Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4EmCalculator.cc
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4EmCalculator.cc 108306 2018-02-02 13:10:06Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4EmCalculator
34 //
35 // Author: Vladimir Ivanchenko
36 //
37 // Creation date: 28.06.2004
38 //
39 // Modifications:
40 // 12.09.2004 Add verbosity (V.Ivanchenko)
41 // 17.11.2004 Change signature of methods, add new methods (V.Ivanchenko)
42 // 08.04.2005 Major optimisation of internal interfaces (V.Ivantchenko)
43 // 08.05.2005 Use updated interfaces (V.Ivantchenko)
44 // 23.10.2005 Fix computations for ions (V.Ivantchenko)
45 // 11.01.2006 Add GetCSDARange (V.Ivantchenko)
46 // 26.01.2006 Rename GetRange -> GetRangeFromRestricteDEDX (V.Ivanchenko)
47 // 14.03.2006 correction in GetCrossSectionPerVolume (mma)
48 // suppress GetCrossSectionPerAtom
49 // elm->GetA() in ComputeCrossSectionPerAtom
50 // 22.03.2006 Add ComputeElectronicDEDX and ComputeTotalDEDX (V.Ivanchenko)
51 // 13.05.2006 Add Corrections for ion stopping (V.Ivanchenko)
52 // 29.09.2006 Uncomment computation of smoothing factor (V.Ivanchenko)
53 // 27.10.2006 Change test energy to access lowEnergy model from
54 // 10 keV to 1 keV (V. Ivanchenko)
55 // 15.03.2007 Add ComputeEnergyCutFromRangeCut methods (V.Ivanchenko)
56 // 21.04.2008 Updated computations for ions (V.Ivanchenko)
57 // 02.02.2018 Fix the MCS (i.e. transport mean free path) case in
58 // GetCrossSectionPerVolume (M. Novak)
59 //
60 // Class Description:
61 //
62 // -------------------------------------------------------------------
63 //
64 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
65 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
66 
67 #include "G4EmCalculator.hh"
68 #include "G4SystemOfUnits.hh"
69 #include "G4LossTableManager.hh"
70 #include "G4EmParameters.hh"
71 #include "G4NistManager.hh"
72 #include "G4VEmProcess.hh"
73 #include "G4VEnergyLossProcess.hh"
74 #include "G4VMultipleScattering.hh"
75 #include "G4Material.hh"
76 #include "G4MaterialCutsCouple.hh"
77 #include "G4ParticleDefinition.hh"
78 #include "G4ParticleTable.hh"
79 #include "G4IonTable.hh"
80 #include "G4PhysicsTable.hh"
81 #include "G4ProductionCutsTable.hh"
82 #include "G4ProcessManager.hh"
83 #include "G4ionEffectiveCharge.hh"
84 #include "G4RegionStore.hh"
85 #include "G4Element.hh"
86 #include "G4EmCorrections.hh"
87 #include "G4GenericIon.hh"
88 #include "G4ProcessVector.hh"
89 #include "G4Gamma.hh"
90 #include "G4Electron.hh"
91 #include "G4Positron.hh"
92 
93 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
94 
96 {
101  nLocalMaterials = 0;
102  verbose = 0;
103  currentCoupleIndex = 0;
104  currentCouple = nullptr;
105  currentMaterial = cutMaterial = nullptr;
106  currentParticle = nullptr;
107  lambdaParticle = nullptr;
108  baseParticle = nullptr;
109  currentLambda = nullptr;
110  currentModel = nullptr;
111  currentProcess = nullptr;
112  curProcess = nullptr;
113  loweModel = nullptr;
114  chargeSquare = 1.0;
115  massRatio = 1.0;
116  mass = 0.0;
117  currentCut = 0.0;
118  cutenergy[0] = cutenergy[1] = cutenergy[2] = DBL_MAX;
121  currentName = "";
122  lambdaName = "";
126  isIon = false;
127  isApplicable = false;
128 }
129 
130 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
131 
133 {
134  delete ionEffCharge;
135  for (G4int i=0; i<nLocalMaterials; ++i) {
136  delete localCouples[i];
137  }
138 }
139 
140 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
141 
143  const G4ParticleDefinition* p,
144  const G4Material* mat,
145  const G4Region* region)
146 {
147  G4double res = 0.0;
148  const G4MaterialCutsCouple* couple = FindCouple(mat, region);
149  if(couple && UpdateParticle(p, kinEnergy) ) {
150  res = manager->GetDEDX(p, kinEnergy, couple);
151 
152  if(isIon) {
153  if(FindEmModel(p, currentProcessName, kinEnergy)) {
154  G4double length = CLHEP::nm;
155  G4double eloss = res*length;
156  //G4cout << "### GetDEDX: E= " << kinEnergy << " dedx0= " << res
157  // << " de= " << eloss << G4endl;;
158  G4double niel = 0.0;
159  dynParticle.SetKineticEnergy(kinEnergy);
160  currentModel->GetChargeSquareRatio(p, mat, kinEnergy);
161  currentModel->CorrectionsAlongStep(couple,&dynParticle,eloss,niel,length);
162  res = eloss/length;
163  //G4cout << " de1= " << eloss << " res1= " << res
164  // << " " << p->GetParticleName() <<G4endl;;
165  }
166  }
167 
168  if(verbose>0) {
169  G4cout << "G4EmCalculator::GetDEDX: E(MeV)= " << kinEnergy/MeV
170  << " DEDX(MeV/mm)= " << res*mm/MeV
171  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
172  << " " << p->GetParticleName()
173  << " in " << mat->GetName()
174  << " isIon= " << isIon
175  << G4endl;
176  }
177  }
178  return res;
179 }
180 
181 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
182 
184  const G4ParticleDefinition* p,
185  const G4Material* mat,
186  const G4Region* region)
187 {
188  G4double res = 0.0;
189  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
190  if(couple && UpdateParticle(p, kinEnergy)) {
191  res = manager->GetRangeFromRestricteDEDX(p, kinEnergy, couple);
192  if(verbose>1) {
193  G4cout << " G4EmCalculator::GetRangeFromRestrictedDEDX: E(MeV)= "
194  << kinEnergy/MeV
195  << " range(mm)= " << res/mm
196  << " " << p->GetParticleName()
197  << " in " << mat->GetName()
198  << G4endl;
199  }
200  }
201  return res;
202 }
203 
204 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
205 
207  const G4ParticleDefinition* p,
208  const G4Material* mat,
209  const G4Region* region)
210 {
211  G4double res = 0.0;
212  if(!theParameters->BuildCSDARange()) {
214  ed << "G4EmCalculator::GetCSDARange: CSDA table is not built; "
215  << " use UI command: /process/eLoss/CSDARange true";
216  G4Exception("G4EmCalculator::GetCSDARange", "em0077",
217  JustWarning, ed);
218  return res;
219  }
220 
221  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
222  if(couple && UpdateParticle(p, kinEnergy)) {
223  res = manager->GetCSDARange(p, kinEnergy, couple);
224  if(verbose>1) {
225  G4cout << " G4EmCalculator::GetCSDARange: E(MeV)= " << kinEnergy/MeV
226  << " range(mm)= " << res/mm
227  << " " << p->GetParticleName()
228  << " in " << mat->GetName()
229  << G4endl;
230  }
231  }
232  return res;
233 }
234 
235 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
236 
238  const G4ParticleDefinition* p,
239  const G4Material* mat,
240  const G4Region* region)
241 {
242  G4double res = 0.0;
244  res = GetCSDARange(kinEnergy, p, mat, region);
245  } else {
246  res = GetRangeFromRestricteDEDX(kinEnergy, p, mat, region);
247  }
248  return res;
249 }
250 
251 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
252 
254  const G4ParticleDefinition* p,
255  const G4Material* mat,
256  const G4Region* region)
257 {
258  G4double res = 0.0;
259  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
260  if(couple && UpdateParticle(p, 1.0*GeV)) {
261  res = manager->GetEnergy(p, range, couple);
262  if(verbose>0) {
263  G4cout << "G4EmCalculator::GetKinEnergy: Range(mm)= " << range/mm
264  << " KinE(MeV)= " << res/MeV
265  << " " << p->GetParticleName()
266  << " in " << mat->GetName()
267  << G4endl;
268  }
269  }
270  return res;
271 }
272 
273 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
274 
276  const G4ParticleDefinition* p,
277  const G4String& processName,
278  const G4Material* mat,
279  const G4Region* region)
280 {
281  G4double res = 0.0;
282  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
283 
284  if(couple && UpdateParticle(p, kinEnergy)) {
285  if(FindEmModel(p, processName, kinEnergy)) {
286  G4int idx = couple->GetIndex();
287  G4int procType = -1;
288  FindLambdaTable(p, processName, kinEnergy, procType);
289 
290  G4VEmProcess* emproc = FindDiscreteProcess(p, processName);
291  if(emproc) {
292  res = emproc->CrossSectionPerVolume(kinEnergy, couple);
293  } else if(currentLambda) {
294  // special tables are built for Msc models (procType is set in FindLambdaTable
295  if(procType==2) {
296  G4VMscModel* mscM = static_cast<G4VMscModel*>(currentModel);
297  mscM->SetCurrentCouple(couple);
298  G4double tr1Mfp = mscM->GetTransportMeanFreePath(p, kinEnergy);
299  if (tr1Mfp<DBL_MAX) {
300  res = 1./tr1Mfp;
301  }
302  } else {
303  G4double e = kinEnergy*massRatio;
304  res = (((*currentLambda)[idx])->Value(e))*chargeSquare;
305  }
306  } else {
307  res = ComputeCrossSectionPerVolume(kinEnergy, p, processName, mat, kinEnergy);
308  }
309  if(verbose>0) {
310  G4cout << "G4EmCalculator::GetXSPerVolume: E(MeV)= " << kinEnergy/MeV
311  << " cross(cm-1)= " << res*cm
312  << " " << p->GetParticleName()
313  << " in " << mat->GetName();
314  if(verbose>1)
315  G4cout << " idx= " << idx << " Escaled((MeV)= "
316  << kinEnergy*massRatio
317  << " q2= " << chargeSquare;
318  G4cout << G4endl;
319  }
320  }
321  }
322  return res;
323 }
324 
325 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
326 
328  const G4String& particle,
329  G4int Z,
331  G4double kinEnergy)
332 {
333  G4double res = 0.0;
334  const G4ParticleDefinition* p = FindParticle(particle);
336  if(p && ad) {
337  res = ad->GetShellIonisationCrossSectionPerAtom(p, Z, shell, kinEnergy);
338  }
339  return res;
340 }
341 
342 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
343 
345  const G4ParticleDefinition* p,
346  const G4String& processName,
347  const G4Material* mat,
348  const G4Region* region)
349 {
350  G4double res = DBL_MAX;
351  G4double x = GetCrossSectionPerVolume(kinEnergy,p, processName, mat,region);
352  if(x > 0.0) { res = 1.0/x; }
353  if(verbose>1) {
354  G4cout << "G4EmCalculator::GetMeanFreePath: E(MeV)= " << kinEnergy/MeV
355  << " MFP(mm)= " << res/mm
356  << " " << p->GetParticleName()
357  << " in " << mat->GetName()
358  << G4endl;
359  }
360  return res;
361 }
362 
363 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
364 
366 {
368  G4cout << "##### DEDX Table for " << p->GetParticleName() << G4endl;
369  if(elp) G4cout << *(elp->DEDXTable()) << G4endl;
370 }
371 
372 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
373 
375 {
377  G4cout << "##### Range Table for " << p->GetParticleName() << G4endl;
378  if(elp) G4cout << *(elp->RangeTableForLoss()) << G4endl;
379 }
380 
381 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
382 
384 {
386  G4cout << "### G4EmCalculator: Inverse Range Table for "
387  << p->GetParticleName() << G4endl;
388  if(elp) G4cout << *(elp->InverseRangeTable()) << G4endl;
389 }
390 
391 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
392 
394  const G4ParticleDefinition* p,
395  const G4String& processName,
396  const G4Material* mat,
397  G4double cut)
398 {
399  SetupMaterial(mat);
400  G4double res = 0.0;
401  if(verbose > 1) {
402  G4cout << "### G4EmCalculator::ComputeDEDX: " << p->GetParticleName()
403  << " in " << currentMaterialName
404  << " e(MeV)= " << kinEnergy/MeV << " cut(MeV)= " << cut/MeV
405  << G4endl;
406  }
407  if(UpdateParticle(p, kinEnergy)) {
408  if(FindEmModel(p, processName, kinEnergy)) {
409 
410  // Special case of ICRU'73 model
411  if(currentModel->GetName() == "ParamICRU73") {
412  res = currentModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cut);
413  if(verbose > 1) {
414  G4cout << " ICRU73 ion E(MeV)= " << kinEnergy << " ";
415  G4cout << currentModel->GetName() << ": DEDX(MeV/mm)= " << res*mm/MeV
416  << " DEDX(MeV*cm^2/g)= "
417  << res*gram/(MeV*cm2*mat->GetDensity())
418  << G4endl;
419  }
420  } else {
421 
422  G4double escaled = kinEnergy*massRatio;
423  if(baseParticle) {
425  mat, baseParticle, escaled, cut) * chargeSquare;
426  if(verbose > 1) {
428  << " Escaled(MeV)= " << escaled;
429  }
430  } else {
431  res = currentModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cut);
432  if(verbose > 1) { G4cout << " no basePart E(MeV)= " << kinEnergy << " "; }
433  }
434  if(verbose > 1) {
435  G4cout << currentModel->GetName() << ": DEDX(MeV/mm)= " << res*mm/MeV
436  << " DEDX(MeV*cm^2/g)= "
437  << res*gram/(MeV*cm2*mat->GetDensity())
438  << G4endl;
439  }
440 
441  // emulate smoothing procedure
443  // G4cout << "massRatio= " << massRatio << " eth= " << eth << G4endl;
444  if(loweModel) {
445  G4double res0 = 0.0;
446  G4double res1 = 0.0;
447  if(baseParticle) {
448  res1 = currentModel->ComputeDEDXPerVolume(mat, baseParticle, eth, cut)
449  * chargeSquare;
450  res0 = loweModel->ComputeDEDXPerVolume(mat, baseParticle, eth, cut)
451  * chargeSquare;
452  } else {
453  res1 = currentModel->ComputeDEDXPerVolume(mat, p, eth, cut);
454  res0 = loweModel->ComputeDEDXPerVolume(mat, p, eth, cut);
455  }
456  if(verbose > 1) {
457  G4cout << "At boundary energy(MeV)= " << eth/MeV
458  << " DEDX(MeV/mm)= " << res1*mm/MeV
459  << G4endl;
460  }
461 
462  //G4cout << "eth= " << eth << " escaled= " << escaled
463  // << " res0= " << res0 << " res1= "
464  // << res1 << " q2= " << chargeSquare << G4endl;
465 
466  if(res1 > 0.0 && escaled > 0.0) {
467  res *= (1.0 + (res0/res1 - 1.0)*eth/escaled);
468  }
469  }
470 
471  // low energy correction for ions
472  if(isIon) {
473  G4double length = CLHEP::nm;
474  const G4Region* r = 0;
475  const G4MaterialCutsCouple* couple = FindCouple(mat, r);
476  G4double eloss = res*length;
477  G4double niel = 0.0;
478  dynParticle.SetKineticEnergy(kinEnergy);
479  currentModel->GetChargeSquareRatio(p, mat, kinEnergy);
480  currentModel->CorrectionsAlongStep(couple,&dynParticle,eloss,niel,length);
481  res = eloss/length;
482 
483  if(verbose > 1) {
484  G4cout << "After Corrections: DEDX(MeV/mm)= " << res*mm/MeV
485  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
486  << G4endl;
487  }
488  }
489  }
490  }
491  if(verbose > 0) {
492  G4cout << "Sum: E(MeV)= " << kinEnergy/MeV
493  << " DEDX(MeV/mm)= " << res*mm/MeV
494  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
495  << " cut(MeV)= " << cut/MeV
496  << " " << p->GetParticleName()
497  << " in " << currentMaterialName
498  << " Zi^2= " << chargeSquare
499  << " isIon=" << isIon
500  << G4endl;
501  }
502  }
503  return res;
504 }
505 
506 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
507 
509  const G4ParticleDefinition* part,
510  const G4Material* mat,
511  G4double cut)
512 {
513  SetupMaterial(mat);
514  G4double dedx = 0.0;
515  if(UpdateParticle(part, kinEnergy)) {
516 
518  const std::vector<G4VEnergyLossProcess*> vel =
519  lManager->GetEnergyLossProcessVector();
520  G4int n = vel.size();
521 
522  //G4cout << "ComputeElectronicDEDX for " << part->GetParticleName()
523  // << " n= " << n << G4endl;
524 
525  for(G4int i=0; i<n; ++i) {
526  if(vel[i]) {
527  G4VProcess* p = reinterpret_cast<G4VProcess*>(vel[i]);
528  if(ActiveForParticle(part, p)) {
529  //G4cout << "idx= " << i << " " << (vel[i])->GetProcessName()
530  // << " " << (vel[i])->Particle()->GetParticleName() << G4endl;
531  dedx += ComputeDEDX(kinEnergy,part,(vel[i])->GetProcessName(),mat,cut);
532  }
533  }
534  }
535  }
536  return dedx;
537 }
538 
539 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
540 
542  const G4ParticleDefinition* part,
543  const G4Material* mat,
544  G4double rangecut)
545 {
546  SetupMaterial(mat);
547  G4double dedx = 0.0;
548  if(UpdateParticle(part, kinEnergy)) {
549 
551  const std::vector<G4VEnergyLossProcess*> vel =
552  lManager->GetEnergyLossProcessVector();
553  G4int n = vel.size();
554 
555  if(mat != cutMaterial) {
556  cutMaterial = mat;
560  }
561 
562  //G4cout << "ComputeElectronicDEDX for " << part->GetParticleName()
563  // << " n= " << n << G4endl;
564 
565  for(G4int i=0; i<n; ++i) {
566  if(vel[i]) {
567  G4VProcess* p = reinterpret_cast<G4VProcess*>(vel[i]);
568  if(ActiveForParticle(part, p)) {
569  //G4cout << "idx= " << i << " " << (vel[i])->GetProcessName()
570  // << " " << (vel[i])->Particle()->GetParticleName() << G4endl;
571  const G4ParticleDefinition* sec = (vel[i])->SecondaryParticle();
572  G4int idx = 0;
573  if(sec == G4Electron::Electron()) { idx = 1; }
574  else if(sec == G4Positron::Positron()) { idx = 2; }
575 
576  dedx += ComputeDEDX(kinEnergy,part,(vel[i])->GetProcessName(),
577  mat,cutenergy[idx]);
578  }
579  }
580  }
581  }
582  return dedx;
583 }
584 
585 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
586 
588  const G4ParticleDefinition* part,
589  const G4Material* mat,
590  G4double cut)
591 {
592  G4double dedx = ComputeElectronicDEDX(kinEnergy,part,mat,cut);
593  if(mass > 700.*MeV) { dedx += ComputeNuclearDEDX(kinEnergy,part,mat); }
594  return dedx;
595 }
596 
597 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
598 
600  const G4ParticleDefinition* p,
601  const G4Material* mat)
602 {
603  G4double res = 0.0;
604  G4VEmProcess* nucst = FindDiscreteProcess(p, "nuclearStopping");
605  if(nucst) {
606  G4VEmModel* mod = nucst->GetModelByIndex();
607  if(mod) {
608  mod->SetFluctuationFlag(false);
609  res = mod->ComputeDEDXPerVolume(mat, p, kinEnergy);
610  }
611  }
612 
613  if(verbose > 1) {
614  G4cout << p->GetParticleName() << " E(MeV)= " << kinEnergy/MeV
615  << " NuclearDEDX(MeV/mm)= " << res*mm/MeV
616  << " NuclearDEDX(MeV*cm^2/g)= "
617  << res*gram/(MeV*cm2*mat->GetDensity())
618  << G4endl;
619  }
620  return res;
621 }
622 
623 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
624 
626  G4double kinEnergy,
627  const G4ParticleDefinition* p,
628  const G4String& processName,
629  const G4Material* mat,
630  G4double cut)
631 {
632  SetupMaterial(mat);
633  G4double res = 0.0;
634  if(UpdateParticle(p, kinEnergy)) {
635  if(FindEmModel(p, processName, kinEnergy)) {
636  G4double e = kinEnergy;
638  if(baseParticle) {
639  e *= kinEnergy*massRatio;
641  mat, baseParticle, e, aCut, e) * chargeSquare;
642  } else {
643  res = currentModel->CrossSectionPerVolume(mat, p, e, aCut, e);
644  }
645  if(verbose>0) {
646  G4cout << "G4EmCalculator::ComputeXSPerVolume: E(MeV)= " << kinEnergy/MeV
647  << " cross(cm-1)= " << res*cm
648  << " cut(keV)= " << aCut/keV
649  << " " << p->GetParticleName()
650  << " in " << mat->GetName()
651  << G4endl;
652  }
653  }
654  }
655  return res;
656 }
657 
658 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
659 
661  G4double kinEnergy,
662  const G4ParticleDefinition* p,
663  const G4String& processName,
664  G4double Z, G4double A,
665  G4double cut)
666 {
667  G4double res = 0.0;
668  if(UpdateParticle(p, kinEnergy)) {
669  G4int iz = G4lrint(Z);
670  CheckMaterial(iz);
671  if(FindEmModel(p, processName, kinEnergy)) {
672  G4double e = kinEnergy;
674  if(baseParticle) {
675  e *= kinEnergy*massRatio;
678  baseParticle, e, Z, A, aCut) * chargeSquare;
679  } else {
681  res = currentModel->ComputeCrossSectionPerAtom(p, e, Z, A, aCut);
682  }
683  if(verbose>0) {
684  G4cout << "E(MeV)= " << kinEnergy/MeV
685  << " cross(barn)= " << res/barn
686  << " " << p->GetParticleName()
687  << " Z= " << Z << " A= " << A/(g/mole) << " g/mole"
688  << " cut(keV)= " << aCut/keV
689  << G4endl;
690  }
691  }
692  }
693  return res;
694 }
695 
696 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
697 
699  const G4ParticleDefinition* p,
700  const G4String& processName,
701  G4int Z, G4int shellIdx,
702  G4double cut)
703 {
704  G4double res = 0.0;
705  if(UpdateParticle(p, kinEnergy)) {
706  CheckMaterial(Z);
707  if(FindEmModel(p, processName, kinEnergy)) {
708  G4double e = kinEnergy;
710  if(baseParticle) {
711  e *= kinEnergy*massRatio;
714  e, aCut) * chargeSquare;
715  } else {
717  res = currentModel->ComputeCrossSectionPerAtom(p, Z, shellIdx, e, aCut);
718  }
719  if(verbose>0) {
720  G4cout << "E(MeV)= " << kinEnergy/MeV
721  << " cross(barn)= " << res/barn
722  << " " << p->GetParticleName()
723  << " Z= " << Z << " shellIdx= " << shellIdx
724  << " cut(keV)= " << aCut/keV
725  << G4endl;
726  }
727  }
728  }
729  return res;
730 }
731 
732 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
733 
734 G4double
736  const G4Material* mat)
737 {
738  G4double res = 0.0;
739  const G4ParticleDefinition* gamma = G4Gamma::Gamma();
740  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "conv", mat, 0.0);
741  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "compt", mat, 0.0);
742  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "phot", mat, 0.0);
743  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "Rayl", mat, 0.0);
744  if(res > 0.0) { res = 1.0/res; }
745  return res;
746 }
747 
748 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
749 
751  const G4String& particle,
752  G4int Z,
754  G4double kinEnergy,
755  const G4Material* mat)
756 {
757  G4double res = 0.0;
758  const G4ParticleDefinition* p = FindParticle(particle);
760  if(p && ad) {
761  res = ad->ComputeShellIonisationCrossSectionPerAtom(p, Z, shell,
762  kinEnergy, mat);
763  }
764  return res;
765 }
766 
767 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
768 
770  const G4ParticleDefinition* p,
771  const G4String& processName,
772  const G4Material* mat,
773  G4double cut)
774 {
775  G4double mfp = DBL_MAX;
776  G4double x = ComputeCrossSectionPerVolume(kinEnergy, p, processName, mat, cut);
777  if(x > 0.0) { mfp = 1.0/x; }
778  if(verbose>1) {
779  G4cout << "E(MeV)= " << kinEnergy/MeV
780  << " MFP(mm)= " << mfp/mm
781  << " " << p->GetParticleName()
782  << " in " << mat->GetName()
783  << G4endl;
784  }
785  return mfp;
786 }
787 
788 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
789 
791  G4double range,
792  const G4ParticleDefinition* part,
793  const G4Material* mat)
794 {
796  ConvertRangeToEnergy(part, mat, range);
797 }
798 
799 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
800 
802  G4double kinEnergy)
803 {
804  if(p != currentParticle) {
805 
806  // new particle
807  currentParticle = p;
808  dynParticle.SetDefinition(const_cast<G4ParticleDefinition*>(p));
809  dynParticle.SetKineticEnergy(kinEnergy);
810  baseParticle = 0;
812  massRatio = 1.0;
813  mass = p->GetPDGMass();
814  chargeSquare = 1.0;
816  currentProcessName = "";
817  isIon = false;
818 
819  // ionisation process exist
820  if(currentProcess) {
823 
824  // base particle is used
825  if(baseParticle) {
828  chargeSquare = q*q;
829  }
830 
831  if(p->GetParticleType() == "nucleus"
832  && currentParticleName != "deuteron"
833  && currentParticleName != "triton"
834  && currentParticleName != "alpha+"
835  && currentParticleName != "helium"
836  && currentParticleName != "hydrogen"
837  ) {
838  isIon = true;
841  if(verbose>1) {
842  G4cout << "\n G4EmCalculator::UpdateParticle: isIon 1 "
843  << p->GetParticleName()
844  << " in " << currentMaterial->GetName()
845  << " e= " << kinEnergy << G4endl;
846  }
847  }
848  }
849  }
850 
851  // Effective charge for ions
852  if(isIon) {
853  chargeSquare =
856  if(currentProcess) {
858  if(verbose>1) {
859  G4cout <<"\n NewIon: massR= "<< massRatio << " q2= "
860  << chargeSquare << " " << currentProcess << G4endl;
861  }
862  }
863  }
864  return true;
865 }
866 
867 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
868 
870 {
871  const G4ParticleDefinition* p = 0;
872  if(name != currentParticleName) {
874  if(!p) {
875  G4cout << "### WARNING: G4EmCalculator::FindParticle fails to find "
876  << name << G4endl;
877  }
878  } else {
879  p = currentParticle;
880  }
881  return p;
882 }
883 
884 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
885 
887 {
888  const G4ParticleDefinition* p = ionTable->GetIon(Z,A,0);
889  return p;
890 }
891 
892 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
893 
895 {
896  if(name != currentMaterialName) {
898  if(!currentMaterial) {
899  G4cout << "### WARNING: G4EmCalculator::FindMaterial fails to find "
900  << name << G4endl;
901  }
902  }
903  return currentMaterial;
904 }
905 
906 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
907 
909 {
910  const G4Region* r = 0;
911  if(reg != "" && reg != "world") {
913  } else {
914  r = G4RegionStore::GetInstance()->GetRegion("DefaultRegionForTheWorld");
915  }
916  return r;
917 }
918 
919 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
920 
922  const G4Material* material,
923  const G4Region* region)
924 {
925  const G4MaterialCutsCouple* couple = nullptr;
926  SetupMaterial(material);
927  if(currentMaterial) {
928  // Access to materials
929  const G4ProductionCutsTable* theCoupleTable=
931  const G4Region* r = region;
932  if(r) {
933  couple = theCoupleTable->GetMaterialCutsCouple(material,
934  r->GetProductionCuts());
935  } else {
937  size_t nr = store->size();
938  if(0 < nr) {
939  for(size_t i=0; i<nr; ++i) {
940  couple = theCoupleTable->GetMaterialCutsCouple(
941  material, ((*store)[i])->GetProductionCuts());
942  if(couple) { break; }
943  }
944  }
945  }
946  }
947  if(!couple) {
949  ed << "G4EmCalculator::FindCouple: fail for material <"
950  << currentMaterialName << ">";
951  if(region) { ed << " and region " << region->GetName(); }
952  G4Exception("G4EmCalculator::FindCouple", "em0078",
953  FatalException, ed);
954  }
955  return couple;
956 }
957 
958 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
959 
961 {
962  SetupMaterial(material);
963  if(!currentMaterial) { return false; }
964  for (G4int i=0; i<nLocalMaterials; ++i) {
965  if(material == localMaterials[i] && cut == localCuts[i]) {
968  currentCut = cut;
969  return true;
970  }
971  }
972  const G4MaterialCutsCouple* cc = new G4MaterialCutsCouple(material);
973  localMaterials.push_back(material);
974  localCouples.push_back(cc);
975  localCuts.push_back(cut);
976  nLocalMaterials++;
977  currentCouple = cc;
979  currentCut = cut;
980  return true;
981 }
982 
983 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
984 
986  const G4String& processName,
987  G4double kinEnergy, G4int& proctype)
988 {
989  // Search for the process
990  if (!currentLambda || p != lambdaParticle || processName != lambdaName) {
991  lambdaName = processName;
992  currentLambda = nullptr;
993  lambdaParticle = p;
994 
995  const G4ParticleDefinition* part = p;
996  if(isIon) { part = theGenericIon; }
997 
998  // Search for energy loss process
999  currentName = processName;
1000  currentModel = nullptr;
1001  loweModel = nullptr;
1002 
1003  G4VEnergyLossProcess* elproc = FindEnLossProcess(part, processName);
1004  if(elproc) {
1005  currentLambda = elproc->LambdaTable();
1006  proctype = 0;
1007  if(currentLambda) {
1008  isApplicable = true;
1009  if(verbose>1) {
1010  G4cout << "G4VEnergyLossProcess is found out: " << currentName
1011  << G4endl;
1012  }
1013  }
1014  curProcess = elproc;
1015  return;
1016  }
1017 
1018  // Search for discrete process
1019  G4VEmProcess* proc = FindDiscreteProcess(part, processName);
1020  if(proc) {
1021  currentLambda = proc->LambdaTable();
1022  proctype = 1;
1023  if(currentLambda) {
1024  isApplicable = true;
1025  if(verbose>1) {
1026  G4cout << "G4VEmProcess is found out: " << currentName << G4endl;
1027  }
1028  }
1029  curProcess = proc;
1030  return;
1031  }
1032 
1033  // Search for msc process
1034  G4VMultipleScattering* msc = FindMscProcess(part, processName);
1035  if(msc) {
1036  currentModel = msc->SelectModel(kinEnergy,0);
1037  proctype = 2;
1038  if(currentModel) {
1040  if(currentLambda) {
1041  isApplicable = true;
1042  if(verbose>1) {
1043  G4cout << "G4VMultipleScattering is found out: " << currentName
1044  << G4endl;
1045  }
1046  }
1047  }
1048  curProcess = msc;
1049  }
1050  }
1051 }
1052 
1053 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1054 
1056  const G4String& processName,
1057  G4double kinEnergy)
1058 {
1059  isApplicable = false;
1060  if(!p || !currentMaterial) {
1061  G4cout << "G4EmCalculator::FindEmModel WARNING: no particle"
1062  << " or materail defined; particle: " << p << G4endl;
1063  return isApplicable;
1064  }
1065  G4String partname = p->GetParticleName();
1066  const G4ParticleDefinition* part = p;
1067  G4double scaledEnergy = kinEnergy*massRatio;
1068  if(isIon) { part = theGenericIon; }
1069 
1070  if(verbose > 1) {
1071  G4cout << "## G4EmCalculator::FindEmModel for " << partname
1072  << " (type= " << p->GetParticleType()
1073  << ") and " << processName << " at E(MeV)= " << scaledEnergy
1074  << G4endl;
1075  if(p != part) { G4cout << " GenericIon is the base particle" << G4endl; }
1076  }
1077 
1078  // Search for energy loss process
1079  currentName = processName;
1080  currentModel = nullptr;
1081  loweModel = nullptr;
1082  size_t idx = 0;
1083 
1084  G4VEnergyLossProcess* elproc = FindEnLossProcess(part, processName);
1085  if(elproc) {
1086  currentModel = elproc->SelectModelForMaterial(scaledEnergy, idx);
1088  currentModel->SetupForMaterial(part, currentMaterial, scaledEnergy);
1090  if(eth > 0.0) {
1091  loweModel = elproc->SelectModelForMaterial(eth - CLHEP::eV, idx);
1092  if(loweModel == currentModel) { loweModel = nullptr; }
1093  else {
1096  }
1097  }
1098  }
1099 
1100  // Search for discrete process
1101  if(!currentModel) {
1102  G4VEmProcess* proc = FindDiscreteProcess(part, processName);
1103  if(proc) {
1104  currentModel = proc->SelectModelForMaterial(kinEnergy, idx);
1106  currentModel->SetupForMaterial(part, currentMaterial, kinEnergy);
1108  if(eth > 0.0) {
1109  loweModel = proc->SelectModelForMaterial(eth - CLHEP::eV, idx);
1110  if(loweModel == currentModel) { loweModel = nullptr; }
1111  else {
1114  }
1115  }
1116  }
1117  }
1118 
1119  // Search for msc process
1120  if(!currentModel) {
1121  G4VMultipleScattering* proc = FindMscProcess(part, processName);
1122  if(proc) {
1123  currentModel = proc->SelectModel(kinEnergy, idx);
1124  loweModel = nullptr;
1125  }
1126  }
1127  if(currentModel) {
1128  if(loweModel == currentModel) { loweModel = nullptr; }
1129  isApplicable = true;
1131  if(loweModel) {
1133  }
1134  if(verbose > 1) {
1135  G4cout << " Model <" << currentModel->GetName()
1136  << "> Emin(MeV)= " << currentModel->LowEnergyLimit()/MeV
1137  << " for " << part->GetParticleName();
1138  if(elproc) {
1139  G4cout << " and " << elproc->GetProcessName() << " " << elproc
1140  << G4endl;
1141  }
1142  if(loweModel) {
1143  G4cout << " LowEnergy model <" << loweModel->GetName() << ">";
1144  }
1145  G4cout << G4endl;
1146  }
1147  }
1148  return isApplicable;
1149 }
1150 
1151 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1152 
1154  const G4ParticleDefinition* p)
1155 {
1156  G4VEnergyLossProcess* elp = nullptr;
1157  G4String partname = p->GetParticleName();
1158  const G4ParticleDefinition* part = p;
1159 
1160  if(p->GetParticleType() == "nucleus"
1161  && currentParticleName != "deuteron"
1162  && currentParticleName != "triton"
1163  && currentParticleName != "He3"
1164  && currentParticleName != "alpha"
1165  && currentParticleName != "alpha+"
1166  && currentParticleName != "helium"
1167  && currentParticleName != "hydrogen"
1168  ) { part = theGenericIon; }
1169 
1170  elp = manager->GetEnergyLossProcess(part);
1171  /*
1172  G4cout << "\n G4EmCalculator::FindEnergyLossProcess: for " << p->GetParticleName()
1173  << " found " << elp->GetProcessName() << " of "
1174  << elp->Particle()->GetParticleName() << " " << elp << G4endl;
1175  */
1176  return elp;
1177 }
1178 
1179 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1180 
1183  const G4String& processName)
1184 {
1185  G4VEnergyLossProcess* proc = 0;
1186  const std::vector<G4VEnergyLossProcess*> v =
1188  G4int n = v.size();
1189  for(G4int i=0; i<n; ++i) {
1190  if((v[i])->GetProcessName() == processName) {
1191  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1192  if(ActiveForParticle(part, p)) {
1193  proc = v[i];
1194  break;
1195  }
1196  }
1197  }
1198  return proc;
1199 }
1200 
1201 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1202 
1203 G4VEmProcess*
1205  const G4String& processName)
1206 {
1207  G4VEmProcess* proc = 0;
1208  const std::vector<G4VEmProcess*> v =
1210  G4int n = v.size();
1211  for(G4int i=0; i<n; ++i) {
1212  if((v[i])->GetProcessName() == processName) {
1213  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1214  if(ActiveForParticle(part, p)) {
1215  proc = v[i];
1216  break;
1217  }
1218  }
1219  }
1220  return proc;
1221 }
1222 
1223 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1224 
1227  const G4String& processName)
1228 {
1229  G4VMultipleScattering* proc = 0;
1230  const std::vector<G4VMultipleScattering*> v =
1232  G4int n = v.size();
1233  for(G4int i=0; i<n; ++i) {
1234  if((v[i])->GetProcessName() == processName) {
1235  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1236  if(ActiveForParticle(part, p)) {
1237  proc = v[i];
1238  break;
1239  }
1240  }
1241  }
1242  return proc;
1243 }
1244 
1245 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1246 
1248  const G4String& processName)
1249 {
1250  G4VProcess* proc = 0;
1251  const G4ProcessManager* procman = part->GetProcessManager();
1252  G4ProcessVector* pv = procman->GetProcessList();
1253  G4int nproc = pv->size();
1254  for(G4int i=0; i<nproc; ++i) {
1255  if(processName == (*pv)[i]->GetProcessName()) {
1256  proc = (*pv)[i];
1257  break;
1258  }
1259  }
1260  return proc;
1261 }
1262 
1263 
1264 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1265 
1267  G4VProcess* proc)
1268 {
1269  G4ProcessManager* pm = part->GetProcessManager();
1270  G4ProcessVector* pv = pm->GetProcessList();
1271  G4int n = pv->size();
1272  G4bool res = false;
1273  for(G4int i=0; i<n; ++i) {
1274  if((*pv)[i] == proc) {
1275  if(pm->GetProcessActivation(i)) { res = true; }
1276  break;
1277  }
1278  }
1279  return res;
1280 }
1281 
1282 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1283 
1285 {
1286  if(mat) {
1287  currentMaterial = mat;
1288  currentMaterialName = mat->GetName();
1289  } else {
1290  currentMaterial = 0;
1291  currentMaterialName = "";
1292  }
1293 }
1294 
1295 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1296 
1298 {
1300 }
1301 
1302 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1303 
1305 {
1306  G4bool isFound = false;
1307  if(currentMaterial) {
1309  for(size_t i=0; i<nn; ++i) {
1310  if(Z == currentMaterial->GetElement(i)->GetZasInt()) {
1311  isFound = true;
1312  break;
1313  }
1314  }
1315  }
1316  if(!isFound) {
1318  }
1319 }
1320 
1321 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1322 
1324 {
1325  verbose = verb;
1326 }
1327 
1328 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1329 
G4double ComputeNuclearDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *)
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
Float_t x
Definition: compare.C:6
G4bool GetProcessActivation(G4VProcess *aProcess) const
G4double GetRangeFromRestricteDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=nullptr)
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:93
void PrintInverseRangeTable(const G4ParticleDefinition *)
const G4PhysicsTable * currentLambda
G4DataVector localCuts
G4double GetMeanFreePath(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, const G4Region *r=nullptr)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
const XML_Char * name
Definition: expat.h:151
virtual void SetupForMaterial(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:447
void FindLambdaTable(const G4ParticleDefinition *, const G4String &processName, G4double kinEnergy, G4int &proctype)
G4PhysicsTable * DEDXTable() const
G4String currentName
const G4Material * currentMaterial
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:45
const G4Region * FindRegion(const G4String &)
void SetupMaterial(const G4Material *)
static G4ParticleTable * GetParticleTable()
G4VAtomDeexcitation * AtomDeexcitation()
static constexpr double MeV
Definition: G4SIunits.hh:214
G4double EffectiveChargeCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
const G4ParticleDefinition * lambdaParticle
void SetKineticEnergy(G4double aEnergy)
static constexpr double keV
Definition: G4SIunits.hh:216
static constexpr double mm
Definition: G4SIunits.hh:115
virtual G4double GetChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:389
const G4String & GetName() const
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
#define G4endl
Definition: G4ios.hh:61
const G4ParticleDefinition * theGenericIon
const char * p
Definition: xmltok.h:285
G4double GetDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=nullptr)
G4bool BuildCSDARange() const
void CheckMaterial(G4int Z)
G4PhysicsTable * RangeTableForLoss() const
G4double ComputeDEDXForCutInRange(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *mat, G4double rangecut=DBL_MAX)
const G4String & GetParticleName() const
const G4String & GetParticleType() const
G4double ComputeElectronicDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *mat, G4double cut=DBL_MAX)
G4IonTable * GetIonTable() const
G4NistManager * nist
G4double GetPDGCharge() const
const G4ParticleDefinition * BaseParticle() const
G4DynamicParticle dynParticle
const G4ParticleDefinition * FindIon(G4int Z, G4int A)
void SetDynamicMassCharge(G4double massratio, G4double charge2ratio)
G4ionEffectiveCharge * ionEffCharge
G4ProcessVector * GetProcessList() const
G4AtomicShellEnumerator
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:455
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:608
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double GetPDGMass() const
G4VProcess * FindProcess(const G4ParticleDefinition *part, const G4String &processName)
const std::vector< G4VEnergyLossProcess * > & GetEnergyLossProcessVector()
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:616
G4double cutenergy[3]
virtual G4double ComputeShellIonisationCrossSectionPerAtom(const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=nullptr)=0
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
const G4String & GetName() const
Definition: G4Material.hh:179
Float_t Z
G4VMultipleScattering * FindMscProcess(const G4ParticleDefinition *, const G4String &processName)
static constexpr double g
Definition: G4SIunits.hh:183
const G4Material * FindMaterial(const G4String &)
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
TString part[npart]
Definition: Style.C:32
G4EmCorrections * EmCorrections()
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:808
G4String currentProcessName
static constexpr double cm2
Definition: G4SIunits.hh:120
const std::vector< G4VEmProcess * > & GetEmProcessVector()
virtual void InitialiseForElement(const G4ParticleDefinition *, G4int Z)
Definition: G4VEmModel.cc:238
G4double GetShellIonisationCrossSectionPerAtom(const G4String &part, G4int Z, G4AtomicShellEnumerator shell, G4double kinEnergy)
G4PhysicsTable * InverseRangeTable() const
const G4ParticleDefinition const G4Material *G4double range
void PrintDEDXTable(const G4ParticleDefinition *)
G4PhysicsTable * LambdaTable() const
G4double GetKinEnergy(G4double range, const G4ParticleDefinition *, const G4Material *, const G4Region *r=nullptr)
G4double ComputeMeanFreePath(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=0.0)
const G4ParticleDefinition * FindParticle(const G4String &)
G4Material * FindOrBuildMaterial(const G4String &name, G4bool isotopes=true, G4bool warning=false)
virtual void CorrectionsAlongStep(const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
Definition: G4VEmModel.cc:406
G4VEmModel * SelectModel(G4double kinEnergy, size_t idx)
G4VEnergyLossProcess * currentProcess
G4VProcess * curProcess
double A(double temperature)
G4double ComputeTotalDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, G4double cut=DBL_MAX)
G4Material * FindOrBuildSimpleMaterial(G4int Z, G4bool warning=false)
const G4String & GetProcessName() const
Definition: G4VProcess.hh:411
static G4Positron * Positron()
Definition: G4Positron.cc:94
G4IonTable * ionTable
static constexpr double eV
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:357
Float_t mat
virtual G4double ComputeCrossSectionPerShell(const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:367
G4VEmModel * currentModel
G4double ComputeGammaAttenuationLength(G4double kinEnergy, const G4Material *)
G4bool UpdateCouple(const G4Material *, G4double cut)
G4double ComputeCrossSectionPerAtom(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, G4double Z, G4double A, G4double cut=0.0)
static G4Electron * Electron()
Definition: G4Electron.cc:94
void PrintRangeTable(const G4ParticleDefinition *)
std::vector< const G4Material * > localMaterials
virtual void InitialiseForMaterial(const G4ParticleDefinition *, const G4Material *)
Definition: G4VEmModel.cc:224
G4ProductionCuts * GetProductionCuts() const
G4int GetZasInt() const
Definition: G4Element.hh:132
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
G4VEmModel * SelectModelForMaterial(G4double kinEnergy, size_t &idxRegion) const
G4VEnergyLossProcess * FindEnLossProcess(const G4ParticleDefinition *, const G4String &processName)
const std::vector< G4VMultipleScattering * > & GetMultipleScatteringVector()
G4double ComputeCrossSectionPerVolume(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=0.0)
static G4ProductionCutsTable * GetProductionCutsTable()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
G4double GetTransportMeanFreePath(const G4ParticleDefinition *part, G4double kinEnergy)
Definition: G4VMscModel.hh:340
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:252
int G4lrint(double ad)
Definition: templates.hh:151
static G4RegionStore * GetInstance()
int G4int
Definition: G4Types.hh:78
static constexpr double nm
Definition: SystemOfUnits.h:92
const G4String & GetName() const
Definition: G4VEmModel.hh:777
G4EmCorrections * corr
G4ProcessManager * GetProcessManager() const
G4double ComputeCrossSectionPerShell(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, G4int Z, G4int shellIdx, G4double cut=0.0)
G4double GetCrossSectionPerVolume(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, const G4Region *r=nullptr)
void SetFluctuationFlag(G4bool val)
Definition: G4VEmModel.hh:686
G4PhysicsTable * LambdaTable() const
static constexpr double barn
Definition: G4SIunits.hh:105
G4double GetCSDARange(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=nullptr)
G4double GetEnergy(const G4ParticleDefinition *aParticle, G4double range, const G4MaterialCutsCouple *couple)
const G4MaterialCutsCouple * currentCouple
G4String currentMaterialName
G4String currentParticleName
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4int size() const
static constexpr double cm
Definition: G4SIunits.hh:119
G4GLOB_DLL std::ostream G4cout
G4LossTableManager * manager
static G4LossTableManager * Instance()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
G4double chargeSquare
Char_t n[5]
G4double GetRangeFromRestricteDEDX(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
G4bool FindEmModel(const G4ParticleDefinition *, const G4String &processName, G4double kinEnergy)
G4VEmProcess * FindDiscreteProcess(const G4ParticleDefinition *, const G4String &processName)
G4VEmModel * SelectModelForMaterial(G4double kinEnergy, size_t &idx) const
G4double GetRange(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=nullptr)
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:504
static const G4double reg
static constexpr double mole
Definition: G4SIunits.hh:286
std::vector< const G4MaterialCutsCouple * > localCouples
const G4MaterialCutsCouple * FindCouple(const G4Material *, const G4Region *r=nullptr)
G4EmParameters * theParameters
void SetVerbose(G4int val)
virtual G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.cc:243
#define DBL_MAX
Definition: templates.hh:83
const G4ParticleDefinition * baseParticle
G4bool UpdateParticle(const G4ParticleDefinition *, G4double kinEnergy)
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:203
G4VEnergyLossProcess * FindEnergyLossProcess(const G4ParticleDefinition *)
const G4ParticleDefinition * currentParticle
virtual G4double GetShellIonisationCrossSectionPerAtom(const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=nullptr)=0
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double LowestElectronEnergy() const
G4double ComputeDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=DBL_MAX)
G4double ComputeShellIonisationCrossSectionPerAtom(const G4String &part, G4int Z, G4AtomicShellEnumerator shell, G4double kinEnergy, const G4Material *mat=nullptr)
const G4Material * cutMaterial
G4double GetCSDARange(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
static G4EmParameters * Instance()
G4double CrossSectionPerVolume(G4double kineticEnergy, const G4MaterialCutsCouple *couple)
size_t GetNumberOfElements() const
Definition: G4Material.hh:187
G4double ComputeEnergyCutFromRangeCut(G4double range, const G4ParticleDefinition *, const G4Material *)
G4bool ActiveForParticle(const G4ParticleDefinition *part, G4VProcess *proc)
G4double GetDEDX(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
static constexpr double gram
Definition: G4SIunits.hh:178
static G4NistManager * Instance()
G4VEmModel * loweModel
G4double GetDensity() const
Definition: G4Material.hh:181