Geant4  10.02.p02
G4EmCalculator.cc
Go to the documentation of this file.
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 97616 2016-06-06 12:37:57Z 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 //
58 // Class Description:
59 //
60 // -------------------------------------------------------------------
61 //
62 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
63 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
64 
65 #include "G4EmCalculator.hh"
66 #include "G4SystemOfUnits.hh"
67 #include "G4LossTableManager.hh"
68 #include "G4EmParameters.hh"
69 #include "G4NistManager.hh"
70 #include "G4VEmProcess.hh"
71 #include "G4VEnergyLossProcess.hh"
72 #include "G4VMultipleScattering.hh"
73 #include "G4Material.hh"
74 #include "G4MaterialCutsCouple.hh"
75 #include "G4ParticleDefinition.hh"
76 #include "G4ParticleTable.hh"
77 #include "G4IonTable.hh"
78 #include "G4PhysicsTable.hh"
79 #include "G4ProductionCutsTable.hh"
80 #include "G4ProcessManager.hh"
81 #include "G4ionEffectiveCharge.hh"
82 #include "G4RegionStore.hh"
83 #include "G4Element.hh"
84 #include "G4EmCorrections.hh"
85 #include "G4GenericIon.hh"
86 #include "G4ProcessVector.hh"
87 #include "G4Gamma.hh"
88 #include "G4Electron.hh"
89 #include "G4Positron.hh"
90 
91 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
92 
94 {
99  nLocalMaterials = 0;
100  verbose = 0;
101  currentCoupleIndex = 0;
102  currentCouple = nullptr;
103  currentMaterial = cutMaterial = nullptr;
104  currentParticle = nullptr;
105  lambdaParticle = nullptr;
106  baseParticle = nullptr;
107  currentLambda = nullptr;
108  currentModel = nullptr;
109  currentProcess = nullptr;
110  loweModel = nullptr;
111  chargeSquare = 1.0;
112  massRatio = 1.0;
113  mass = 0.0;
114  currentCut = 0.0;
115  cutenergy[0] = cutenergy[1] = cutenergy[2] = DBL_MAX;
118  currentName = "";
119  lambdaName = "";
123  isIon = false;
124  isApplicable = false;
125 }
126 
127 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
128 
130 {
131  delete ionEffCharge;
132  for (G4int i=0; i<nLocalMaterials; ++i) {
133  delete localCouples[i];
134  }
135 }
136 
137 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
138 
140  const G4ParticleDefinition* p,
141  const G4Material* mat,
142  const G4Region* region)
143 {
144  G4double res = 0.0;
145  const G4MaterialCutsCouple* couple = FindCouple(mat, region);
146  if(couple && UpdateParticle(p, kinEnergy) ) {
147  res = manager->GetDEDX(p, kinEnergy, couple);
148 
149  if(isIon) {
150  if(FindEmModel(p, currentProcessName, kinEnergy)) {
151  G4double length = CLHEP::nm;
152  G4double eloss = res*length;
153  //G4cout << "### GetDEDX: E= " << kinEnergy << " dedx0= " << res
154  // << " de= " << eloss << G4endl;;
155  G4double niel = 0.0;
156  dynParticle.SetKineticEnergy(kinEnergy);
157  currentModel->GetChargeSquareRatio(p, mat, kinEnergy);
158  currentModel->CorrectionsAlongStep(couple,&dynParticle,eloss,niel,length);
159  res = eloss/length;
160  //G4cout << " de1= " << eloss << " res1= " << res
161  // << " " << p->GetParticleName() <<G4endl;;
162  }
163  }
164 
165  if(verbose>0) {
166  G4cout << "G4EmCalculator::GetDEDX: E(MeV)= " << kinEnergy/MeV
167  << " DEDX(MeV/mm)= " << res*mm/MeV
168  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
169  << " " << p->GetParticleName()
170  << " in " << mat->GetName()
171  << " isIon= " << isIon
172  << G4endl;
173  }
174  }
175  return res;
176 }
177 
178 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
179 
181  const G4ParticleDefinition* p,
182  const G4Material* mat,
183  const G4Region* region)
184 {
185  G4double res = 0.0;
186  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
187  if(couple && UpdateParticle(p, kinEnergy)) {
188  res = manager->GetRangeFromRestricteDEDX(p, kinEnergy, couple);
189  if(verbose>1) {
190  G4cout << " G4EmCalculator::GetRangeFromRestrictedDEDX: E(MeV)= "
191  << kinEnergy/MeV
192  << " range(mm)= " << res/mm
193  << " " << p->GetParticleName()
194  << " in " << mat->GetName()
195  << G4endl;
196  }
197  }
198  return res;
199 }
200 
201 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
202 
204  const G4ParticleDefinition* p,
205  const G4Material* mat,
206  const G4Region* region)
207 {
208  G4double res = 0.0;
209  if(!theParameters->BuildCSDARange()) {
211  ed << "G4EmCalculator::GetCSDARange: CSDA table is not built; "
212  << " use UI command: /process/eLoss/CSDARange true";
213  G4Exception("G4EmCalculator::GetCSDARange", "em0077",
214  JustWarning, ed);
215  return res;
216  }
217 
218  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
219  if(couple && UpdateParticle(p, kinEnergy)) {
220  res = manager->GetCSDARange(p, kinEnergy, couple);
221  if(verbose>1) {
222  G4cout << " G4EmCalculator::GetCSDARange: E(MeV)= " << kinEnergy/MeV
223  << " range(mm)= " << res/mm
224  << " " << p->GetParticleName()
225  << " in " << mat->GetName()
226  << G4endl;
227  }
228  }
229  return res;
230 }
231 
232 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
233 
235  const G4ParticleDefinition* p,
236  const G4Material* mat,
237  const G4Region* region)
238 {
239  G4double res = 0.0;
241  res = GetCSDARange(kinEnergy, p, mat, region);
242  } else {
243  res = GetRangeFromRestricteDEDX(kinEnergy, p, mat, region);
244  }
245  return res;
246 }
247 
248 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
249 
251  const G4ParticleDefinition* p,
252  const G4Material* mat,
253  const G4Region* region)
254 {
255  G4double res = 0.0;
256  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
257  if(couple && UpdateParticle(p, 1.0*GeV)) {
258  res = manager->GetEnergy(p, range, couple);
259  if(verbose>0) {
260  G4cout << "G4EmCalculator::GetKinEnergy: Range(mm)= " << range/mm
261  << " KinE(MeV)= " << res/MeV
262  << " " << p->GetParticleName()
263  << " in " << mat->GetName()
264  << G4endl;
265  }
266  }
267  return res;
268 }
269 
270 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
271 
273  const G4ParticleDefinition* p,
274  const G4String& processName,
275  const G4Material* mat,
276  const G4Region* region)
277 {
278  G4double res = 0.0;
279  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
280 
281  if(couple && UpdateParticle(p, kinEnergy)) {
282  G4int idx = couple->GetIndex();
283  FindLambdaTable(p, processName, kinEnergy);
284 
285  if(currentLambda) {
286  G4double e = kinEnergy*massRatio;
287  res = (((*currentLambda)[idx])->Value(e))*chargeSquare;
288  } else {
289  res = ComputeCrossSectionPerVolume(kinEnergy, p, processName, mat,
290  kinEnergy);
291  }
292  if(verbose>0) {
293  G4cout << "G4EmCalculator::GetXSPerVolume: E(MeV)= " << kinEnergy/MeV
294  << " cross(cm-1)= " << res*cm
295  << " " << p->GetParticleName()
296  << " in " << mat->GetName();
297  if(verbose>1)
298  G4cout << " idx= " << idx << " Escaled((MeV)= "
299  << kinEnergy*massRatio
300  << " q2= " << chargeSquare;
301  G4cout << G4endl;
302  }
303  }
304  return res;
305 }
306 
307 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
308 
310  const G4String& particle,
311  G4int Z,
313  G4double kinEnergy)
314 {
315  G4double res = 0.0;
316  const G4ParticleDefinition* p = FindParticle(particle);
318  if(p && ad) {
319  res = ad->GetShellIonisationCrossSectionPerAtom(p, Z, shell, kinEnergy);
320  }
321  return res;
322 }
323 
324 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
325 
327  const G4ParticleDefinition* p,
328  const G4String& processName,
329  const G4Material* mat,
330  const G4Region* region)
331 {
332  G4double res = DBL_MAX;
333  G4double x = GetCrossSectionPerVolume(kinEnergy,p, processName, mat,region);
334  if(x > 0.0) { res = 1.0/x; }
335  if(verbose>1) {
336  G4cout << "G4EmCalculator::GetMeanFreePath: E(MeV)= " << kinEnergy/MeV
337  << " MFP(mm)= " << res/mm
338  << " " << p->GetParticleName()
339  << " in " << mat->GetName()
340  << G4endl;
341  }
342  return res;
343 }
344 
345 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
346 
348 {
350  G4cout << "##### DEDX Table for " << p->GetParticleName() << G4endl;
351  if(elp) G4cout << *(elp->DEDXTable()) << G4endl;
352 }
353 
354 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
355 
357 {
359  G4cout << "##### Range Table for " << p->GetParticleName() << G4endl;
360  if(elp) G4cout << *(elp->RangeTableForLoss()) << G4endl;
361 }
362 
363 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
364 
366 {
368  G4cout << "### G4EmCalculator: Inverse Range Table for "
369  << p->GetParticleName() << G4endl;
370  if(elp) G4cout << *(elp->InverseRangeTable()) << G4endl;
371 }
372 
373 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
374 
376  const G4ParticleDefinition* p,
377  const G4String& processName,
378  const G4Material* mat,
379  G4double cut)
380 {
381  SetupMaterial(mat);
382  G4double res = 0.0;
383  if(verbose > 1) {
384  G4cout << "### G4EmCalculator::ComputeDEDX: " << p->GetParticleName()
385  << " in " << currentMaterialName
386  << " e(MeV)= " << kinEnergy/MeV << " cut(MeV)= " << cut/MeV
387  << G4endl;
388  }
389  if(UpdateParticle(p, kinEnergy)) {
390  if(FindEmModel(p, processName, kinEnergy)) {
391 
392  // Special case of ICRU'73 model
393  if(currentModel->GetName() == "ParamICRU73") {
394  res = currentModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cut);
395  if(verbose > 1) {
396  G4cout << " ICRU73 ion E(MeV)= " << kinEnergy << " ";
397  G4cout << currentModel->GetName() << ": DEDX(MeV/mm)= " << res*mm/MeV
398  << " DEDX(MeV*cm^2/g)= "
399  << res*gram/(MeV*cm2*mat->GetDensity())
400  << G4endl;
401  }
402  } else {
403 
404  G4double escaled = kinEnergy*massRatio;
405  if(baseParticle) {
407  mat, baseParticle, escaled, cut) * chargeSquare;
408  if(verbose > 1) {
410  << " Escaled(MeV)= " << escaled;
411  }
412  } else {
413  res = currentModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cut);
414  if(verbose > 1) { G4cout << " no basePart E(MeV)= " << kinEnergy << " "; }
415  }
416  if(verbose > 1) {
417  G4cout << currentModel->GetName() << ": DEDX(MeV/mm)= " << res*mm/MeV
418  << " DEDX(MeV*cm^2/g)= "
419  << res*gram/(MeV*cm2*mat->GetDensity())
420  << G4endl;
421  }
422 
423  // emulate smoothing procedure
425  // G4cout << "massRatio= " << massRatio << " eth= " << eth << G4endl;
426  if(loweModel) {
427  G4double res0 = 0.0;
428  G4double res1 = 0.0;
429  if(baseParticle) {
430  res1 = currentModel->ComputeDEDXPerVolume(mat, baseParticle, eth, cut)
431  * chargeSquare;
432  res0 = loweModel->ComputeDEDXPerVolume(mat, baseParticle, eth, cut)
433  * chargeSquare;
434  } else {
435  res1 = currentModel->ComputeDEDXPerVolume(mat, p, eth, cut);
436  res0 = loweModel->ComputeDEDXPerVolume(mat, p, eth, cut);
437  }
438  if(verbose > 1) {
439  G4cout << "At boundary energy(MeV)= " << eth/MeV
440  << " DEDX(MeV/mm)= " << res1*mm/MeV
441  << G4endl;
442  }
443 
444  //G4cout << "eth= " << eth << " escaled= " << escaled
445  // << " res0= " << res0 << " res1= "
446  // << res1 << " q2= " << chargeSquare << G4endl;
447 
448  if(res1 > 0.0 && escaled > 0.0) {
449  res *= (1.0 + (res0/res1 - 1.0)*eth/escaled);
450  }
451  }
452 
453  // low energy correction for ions
454  if(isIon) {
455  G4double length = CLHEP::nm;
456  const G4Region* r = 0;
457  const G4MaterialCutsCouple* couple = FindCouple(mat, r);
458  G4double eloss = res*length;
459  G4double niel = 0.0;
460  dynParticle.SetKineticEnergy(kinEnergy);
461  currentModel->GetChargeSquareRatio(p, mat, kinEnergy);
462  currentModel->CorrectionsAlongStep(couple,&dynParticle,eloss,niel,length);
463  res = eloss/length;
464 
465  if(verbose > 1) {
466  G4cout << "After Corrections: DEDX(MeV/mm)= " << res*mm/MeV
467  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
468  << G4endl;
469  }
470  }
471  }
472  }
473  if(verbose > 0) {
474  G4cout << "Sum: E(MeV)= " << kinEnergy/MeV
475  << " DEDX(MeV/mm)= " << res*mm/MeV
476  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
477  << " cut(MeV)= " << cut/MeV
478  << " " << p->GetParticleName()
479  << " in " << currentMaterialName
480  << " Zi^2= " << chargeSquare
481  << " isIon=" << isIon
482  << G4endl;
483  }
484  }
485  return res;
486 }
487 
488 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
489 
491  const G4ParticleDefinition* part,
492  const G4Material* mat,
493  G4double cut)
494 {
495  SetupMaterial(mat);
496  G4double dedx = 0.0;
497  if(UpdateParticle(part, kinEnergy)) {
498 
500  const std::vector<G4VEnergyLossProcess*> vel =
501  lManager->GetEnergyLossProcessVector();
502  G4int n = vel.size();
503 
504  //G4cout << "ComputeElectronicDEDX for " << part->GetParticleName()
505  // << " n= " << n << G4endl;
506 
507  for(G4int i=0; i<n; ++i) {
508  if(vel[i]) {
509  G4VProcess* p = reinterpret_cast<G4VProcess*>(vel[i]);
510  if(ActiveForParticle(part, p)) {
511  //G4cout << "idx= " << i << " " << (vel[i])->GetProcessName()
512  // << " " << (vel[i])->Particle()->GetParticleName() << G4endl;
513  dedx += ComputeDEDX(kinEnergy,part,(vel[i])->GetProcessName(),mat,cut);
514  }
515  }
516  }
517  }
518  return dedx;
519 }
520 
521 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
522 
524  const G4ParticleDefinition* part,
525  const G4Material* mat,
526  G4double rangecut)
527 {
528  SetupMaterial(mat);
529  G4double dedx = 0.0;
530  if(UpdateParticle(part, kinEnergy)) {
531 
533  const std::vector<G4VEnergyLossProcess*> vel =
534  lManager->GetEnergyLossProcessVector();
535  G4int n = vel.size();
536 
537  if(mat != cutMaterial) {
538  cutMaterial = mat;
542  }
543 
544  //G4cout << "ComputeElectronicDEDX for " << part->GetParticleName()
545  // << " n= " << n << G4endl;
546 
547  for(G4int i=0; i<n; ++i) {
548  if(vel[i]) {
549  G4VProcess* p = reinterpret_cast<G4VProcess*>(vel[i]);
550  if(ActiveForParticle(part, p)) {
551  //G4cout << "idx= " << i << " " << (vel[i])->GetProcessName()
552  // << " " << (vel[i])->Particle()->GetParticleName() << G4endl;
553  const G4ParticleDefinition* sec = (vel[i])->SecondaryParticle();
554  G4int idx = 0;
555  if(sec == G4Electron::Electron()) { idx = 1; }
556  else if(sec == G4Positron::Positron()) { idx = 2; }
557 
558  dedx += ComputeDEDX(kinEnergy,part,(vel[i])->GetProcessName(),
559  mat,cutenergy[idx]);
560  }
561  }
562  }
563  }
564  return dedx;
565 }
566 
567 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
568 
570  const G4ParticleDefinition* part,
571  const G4Material* mat,
572  G4double cut)
573 {
574  G4double dedx = ComputeElectronicDEDX(kinEnergy,part,mat,cut);
575  if(mass > 700.*MeV) { dedx += ComputeNuclearDEDX(kinEnergy,part,mat); }
576  return dedx;
577 }
578 
579 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
580 
582  const G4ParticleDefinition* p,
583  const G4Material* mat)
584 {
585 
586  G4double res = corr->NuclearDEDX(p, mat, kinEnergy, false);
587 
588  if(verbose > 1) {
589  G4cout << p->GetParticleName() << " E(MeV)= " << kinEnergy/MeV
590  << " NuclearDEDX(MeV/mm)= " << res*mm/MeV
591  << " NuclearDEDX(MeV*cm^2/g)= "
592  << res*gram/(MeV*cm2*mat->GetDensity())
593  << G4endl;
594  }
595  return res;
596 }
597 
598 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
599 
601  G4double kinEnergy,
602  const G4ParticleDefinition* p,
603  const G4String& processName,
604  const G4Material* mat,
605  G4double cut)
606 {
607  SetupMaterial(mat);
608  G4double res = 0.0;
609  if(UpdateParticle(p, kinEnergy)) {
610  if(FindEmModel(p, processName, kinEnergy)) {
611  G4double e = kinEnergy;
613  if(baseParticle) {
614  e *= kinEnergy*massRatio;
616  mat, baseParticle, e, aCut, e) * chargeSquare;
617  } else {
618  res = currentModel->CrossSectionPerVolume(mat, p, e, aCut, e);
619  }
620  if(verbose>0) {
621  G4cout << "G4EmCalculator::ComputeXSPerVolume: E(MeV)= " << kinEnergy/MeV
622  << " cross(cm-1)= " << res*cm
623  << " cut(keV)= " << aCut/keV
624  << " " << p->GetParticleName()
625  << " in " << mat->GetName()
626  << G4endl;
627  }
628  }
629  }
630  return res;
631 }
632 
633 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
634 
636  G4double kinEnergy,
637  const G4ParticleDefinition* p,
638  const G4String& processName,
639  G4double Z, G4double A,
640  G4double cut)
641 {
642  G4double res = 0.0;
643  if(UpdateParticle(p, kinEnergy)) {
645  if(FindEmModel(p, processName, kinEnergy)) {
646  G4double e = kinEnergy;
648  if(baseParticle) {
649  e *= kinEnergy*massRatio;
652  baseParticle, e, Z, A, aCut) * chargeSquare;
653  } else {
655  res = currentModel->ComputeCrossSectionPerAtom(p, e, Z, A, aCut);
656  }
657  if(verbose>0) {
658  G4cout << "E(MeV)= " << kinEnergy/MeV
659  << " cross(barn)= " << res/barn
660  << " " << p->GetParticleName()
661  << " Z= " << Z << " A= " << A/(g/mole) << " g/mole"
662  << " cut(keV)= " << aCut/keV
663  << G4endl;
664  }
665  }
666  }
667  return res;
668 }
669 
670 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
671 
673  const G4ParticleDefinition* p,
674  const G4String& processName,
675  G4int Z, G4int shellIdx,
676  G4double cut)
677 {
678  G4double res = 0.0;
679  if(UpdateParticle(p, kinEnergy)) {
680  CheckMaterial(Z);
681  if(FindEmModel(p, processName, kinEnergy)) {
682  G4double e = kinEnergy;
684  if(baseParticle) {
685  e *= kinEnergy*massRatio;
688  e, aCut) * chargeSquare;
689  } else {
691  res = currentModel->ComputeCrossSectionPerAtom(p, Z, shellIdx, e, aCut);
692  }
693  if(verbose>0) {
694  G4cout << "E(MeV)= " << kinEnergy/MeV
695  << " cross(barn)= " << res/barn
696  << " " << p->GetParticleName()
697  << " Z= " << Z << " shellIdx= " << shellIdx
698  << " cut(keV)= " << aCut/keV
699  << G4endl;
700  }
701  }
702  }
703  return res;
704 }
705 
706 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
707 
708 G4double
710  const G4Material* mat)
711 {
712  G4double res = 0.0;
713  const G4ParticleDefinition* gamma = G4Gamma::Gamma();
714  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "conv", mat, 0.0);
715  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "compt", mat, 0.0);
716  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "phot", mat, 0.0);
717  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "Rayl", mat, 0.0);
718  if(res > 0.0) { res = 1.0/res; }
719  return res;
720 }
721 
722 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
723 
725  const G4String& particle,
726  G4int Z,
728  G4double kinEnergy,
729  const G4Material* mat)
730 {
731  G4double res = 0.0;
732  const G4ParticleDefinition* p = FindParticle(particle);
734  if(p && ad) {
735  res = ad->ComputeShellIonisationCrossSectionPerAtom(p, Z, shell,
736  kinEnergy, mat);
737  }
738  return res;
739 }
740 
741 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
742 
744  const G4ParticleDefinition* p,
745  const G4String& processName,
746  const G4Material* mat,
747  G4double cut)
748 {
749  G4double mfp = DBL_MAX;
750  G4double x = ComputeCrossSectionPerVolume(kinEnergy, p, processName, mat, cut);
751  if(x > 0.0) { mfp = 1.0/x; }
752  if(verbose>1) {
753  G4cout << "E(MeV)= " << kinEnergy/MeV
754  << " MFP(mm)= " << mfp/mm
755  << " " << p->GetParticleName()
756  << " in " << mat->GetName()
757  << G4endl;
758  }
759  return mfp;
760 }
761 
762 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
763 
765  G4double range,
766  const G4ParticleDefinition* part,
767  const G4Material* mat)
768 {
770  ConvertRangeToEnergy(part, mat, range);
771 }
772 
773 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
774 
776  G4double kinEnergy)
777 {
778  if(p != currentParticle) {
779 
780  // new particle
781  currentParticle = p;
782  dynParticle.SetDefinition(const_cast<G4ParticleDefinition*>(p));
783  dynParticle.SetKineticEnergy(kinEnergy);
784  baseParticle = 0;
786  massRatio = 1.0;
787  mass = p->GetPDGMass();
788  chargeSquare = 1.0;
790  currentProcessName = "";
791  isIon = false;
792 
793  // ionisation process exist
794  if(currentProcess) {
797 
798  // base particle is used
799  if(baseParticle) {
802  chargeSquare = q*q;
803  }
804 
805  if(p->GetParticleType() == "nucleus"
806  && currentParticleName != "deuteron"
807  && currentParticleName != "triton"
808  && currentParticleName != "alpha+"
809  && currentParticleName != "helium"
810  && currentParticleName != "hydrogen"
811  ) {
812  isIon = true;
815  if(verbose>1) {
816  G4cout << "\n G4EmCalculator::UpdateParticle: isIon 1 "
817  << p->GetParticleName()
818  << " in " << currentMaterial->GetName()
819  << " e= " << kinEnergy << G4endl;
820  }
821  }
822  }
823  }
824 
825  // Effective charge for ions
826  if(isIon) {
827  chargeSquare =
830  if(currentProcess) {
832  if(verbose>1) {
833  G4cout <<"\n NewIon: massR= "<< massRatio << " q2= "
834  << chargeSquare << " " << currentProcess << G4endl;
835  }
836  }
837  }
838  return true;
839 }
840 
841 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
842 
844 {
845  const G4ParticleDefinition* p = 0;
846  if(name != currentParticleName) {
848  if(!p) {
849  G4cout << "### WARNING: G4EmCalculator::FindParticle fails to find "
850  << name << G4endl;
851  }
852  } else {
853  p = currentParticle;
854  }
855  return p;
856 }
857 
858 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
859 
861 {
862  const G4ParticleDefinition* p = ionTable->GetIon(Z,A,0);
863  return p;
864 }
865 
866 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
867 
869 {
870  if(name != currentMaterialName) {
872  if(!currentMaterial) {
873  G4cout << "### WARNING: G4EmCalculator::FindMaterial fails to find "
874  << name << G4endl;
875  }
876  }
877  return currentMaterial;
878 }
879 
880 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
881 
883 {
884  const G4Region* r = 0;
885  if(reg != "" && reg != "world") {
887  } else {
888  r = G4RegionStore::GetInstance()->GetRegion("DefaultRegionForTheWorld");
889  }
890  return r;
891 }
892 
893 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
894 
896  const G4Material* material,
897  const G4Region* region)
898 {
899  const G4MaterialCutsCouple* couple = 0;
900  SetupMaterial(material);
901  if(currentMaterial) {
902  // Access to materials
903  const G4ProductionCutsTable* theCoupleTable=
905  const G4Region* r = region;
906  if(r) {
907  couple = theCoupleTable->GetMaterialCutsCouple(material,
908  r->GetProductionCuts());
909  } else {
911  size_t nr = store->size();
912  if(0 < nr) {
913  for(size_t i=0; i<nr; ++i) {
914  couple = theCoupleTable->GetMaterialCutsCouple(
915  material, ((*store)[i])->GetProductionCuts());
916  if(couple) { break; }
917  }
918  }
919  }
920  }
921  if(!couple) {
923  ed << "G4EmCalculator::FindCouple: fail for material <"
924  << currentMaterialName << ">";
925  if(region) { ed << " and region " << region->GetName(); }
926  G4Exception("G4EmCalculator::FindCouple", "em0078",
927  FatalException, ed);
928  }
929  return couple;
930 }
931 
932 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
933 
935 {
936  SetupMaterial(material);
937  if(!currentMaterial) { return false; }
938  for (G4int i=0; i<nLocalMaterials; ++i) {
939  if(material == localMaterials[i] && cut == localCuts[i]) {
942  currentCut = cut;
943  return true;
944  }
945  }
946  const G4MaterialCutsCouple* cc = new G4MaterialCutsCouple(material);
947  localMaterials.push_back(material);
948  localCouples.push_back(cc);
949  localCuts.push_back(cut);
950  nLocalMaterials++;
951  currentCouple = cc;
953  currentCut = cut;
954  return true;
955 }
956 
957 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
958 
960  const G4String& processName,
961  G4double kinEnergy)
962 {
963  // Search for the process
964  if (!currentLambda || p != lambdaParticle || processName != lambdaName) {
965  lambdaName = processName;
966  currentLambda = 0;
967  lambdaParticle = p;
968 
969  const G4ParticleDefinition* part = p;
970  if(isIon) { part = theGenericIon; }
971 
972  // Search for energy loss process
973  currentName = processName;
974  currentModel = 0;
975  loweModel = 0;
976 
977  G4VEnergyLossProcess* elproc = FindEnLossProcess(part, processName);
978  if(elproc) {
979  currentLambda = elproc->LambdaTable();
980  if(currentLambda) {
981  isApplicable = true;
982  if(verbose>1) {
983  G4cout << "G4VEnergyLossProcess is found out: " << currentName
984  << G4endl;
985  }
986  }
987  return;
988  }
989 
990  // Search for discrete process
991  G4VEmProcess* proc = FindDiscreteProcess(part, processName);
992  if(proc) {
993  currentLambda = proc->LambdaTable();
994  if(currentLambda) {
995  isApplicable = true;
996  if(verbose>1) {
997  G4cout << "G4VEmProcess is found out: " << currentName << G4endl;
998  }
999  }
1000  return;
1001  }
1002 
1003  // Search for msc process
1004  G4VMultipleScattering* msc = FindMscProcess(part, processName);
1005  if(msc) {
1006  currentModel = msc->SelectModel(kinEnergy,0);
1007  if(currentModel) {
1009  if(currentLambda) {
1010  isApplicable = true;
1011  if(verbose>1) {
1012  G4cout << "G4VMultipleScattering is found out: " << currentName
1013  << G4endl;
1014  }
1015  }
1016  }
1017  }
1018  }
1019 }
1020 
1021 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1022 
1024  const G4String& processName,
1025  G4double kinEnergy)
1026 {
1027  isApplicable = false;
1028  if(!p || !currentMaterial) {
1029  G4cout << "G4EmCalculator::FindEmModel WARNING: no particle"
1030  << " or materail defined; particle: " << p << G4endl;
1031  return isApplicable;
1032  }
1033  G4String partname = p->GetParticleName();
1034  const G4ParticleDefinition* part = p;
1035  G4double scaledEnergy = kinEnergy*massRatio;
1036  if(isIon) { part = theGenericIon; }
1037 
1038  if(verbose > 1) {
1039  G4cout << "## G4EmCalculator::FindEmModel for " << partname
1040  << " (type= " << p->GetParticleType()
1041  << ") and " << processName << " at E(MeV)= " << scaledEnergy
1042  << G4endl;
1043  if(p != part) { G4cout << " GenericIon is the base particle" << G4endl; }
1044  }
1045 
1046  // Search for energy loss process
1047  currentName = processName;
1048  currentModel = 0;
1049  loweModel = 0;
1050  size_t idx = 0;
1051 
1052  G4VEnergyLossProcess* elproc = FindEnLossProcess(part, processName);
1053  if(elproc) {
1054  currentModel = elproc->SelectModelForMaterial(scaledEnergy, idx);
1056  currentModel->SetupForMaterial(part, currentMaterial, scaledEnergy);
1058  if(eth > 0.0) {
1059  loweModel = elproc->SelectModelForMaterial(eth - CLHEP::eV, idx);
1060  if(loweModel == currentModel) { loweModel = nullptr; }
1061  else {
1064  }
1065  }
1066  }
1067 
1068  // Search for discrete process
1069  if(!currentModel) {
1070  G4VEmProcess* proc = FindDiscreteProcess(part, processName);
1071  if(proc) {
1072  currentModel = proc->SelectModelForMaterial(kinEnergy, idx);
1074  currentModel->SetupForMaterial(part, currentMaterial, kinEnergy);
1076  if(eth > 0.0) {
1077  loweModel = proc->SelectModelForMaterial(eth - CLHEP::eV, idx);
1078  if(loweModel == currentModel) { loweModel = 0; }
1079  else {
1082  }
1083  }
1084  }
1085  }
1086 
1087  // Search for msc process
1088  if(!currentModel) {
1089  G4VMultipleScattering* proc = FindMscProcess(part, processName);
1090  if(proc) {
1091  currentModel = proc->SelectModel(kinEnergy, idx);
1092  loweModel = 0;
1093  }
1094  }
1095  if(currentModel) {
1096  if(loweModel == currentModel) { loweModel = 0; }
1097  isApplicable = true;
1099  if(loweModel) {
1101  }
1102  if(verbose > 1) {
1103  G4cout << " Model <" << currentModel->GetName()
1104  << "> Emin(MeV)= " << currentModel->LowEnergyLimit()/MeV
1105  << " for " << part->GetParticleName();
1106  if(elproc) {
1107  G4cout << " and " << elproc->GetProcessName() << " " << elproc
1108  << G4endl;
1109  }
1110  if(loweModel) {
1111  G4cout << " LowEnergy model <" << loweModel->GetName() << ">";
1112  }
1113  G4cout << G4endl;
1114  }
1115  }
1116  return isApplicable;
1117 }
1118 
1119 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1120 
1122  const G4ParticleDefinition* p)
1123 {
1124  G4VEnergyLossProcess* elp = 0;
1125  G4String partname = p->GetParticleName();
1126  const G4ParticleDefinition* part = p;
1127 
1128  if(p->GetParticleType() == "nucleus"
1129  && currentParticleName != "deuteron"
1130  && currentParticleName != "triton"
1131  && currentParticleName != "He3"
1132  && currentParticleName != "alpha"
1133  && currentParticleName != "alpha+"
1134  && currentParticleName != "helium"
1135  && currentParticleName != "hydrogen"
1136  ) { part = theGenericIon; }
1137 
1138  elp = manager->GetEnergyLossProcess(part);
1139  /*
1140  G4cout << "\n G4EmCalculator::FindEnergyLossProcess: for " << p->GetParticleName()
1141  << " found " << elp->GetProcessName() << " of "
1142  << elp->Particle()->GetParticleName() << " " << elp << G4endl;
1143  */
1144  return elp;
1145 }
1146 
1147 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1148 
1151  const G4String& processName)
1152 {
1153  G4VEnergyLossProcess* proc = 0;
1154  const std::vector<G4VEnergyLossProcess*> v =
1156  G4int n = v.size();
1157  for(G4int i=0; i<n; ++i) {
1158  if((v[i])->GetProcessName() == processName) {
1159  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1160  if(ActiveForParticle(part, p)) {
1161  proc = v[i];
1162  break;
1163  }
1164  }
1165  }
1166  return proc;
1167 }
1168 
1169 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1170 
1171 G4VEmProcess*
1173  const G4String& processName)
1174 {
1175  G4VEmProcess* proc = 0;
1176  const std::vector<G4VEmProcess*> v =
1178  G4int n = v.size();
1179  for(G4int i=0; i<n; ++i) {
1180  if((v[i])->GetProcessName() == processName) {
1181  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1182  if(ActiveForParticle(part, p)) {
1183  proc = v[i];
1184  break;
1185  }
1186  }
1187  }
1188  return proc;
1189 }
1190 
1191 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1192 
1195  const G4String& processName)
1196 {
1197  G4VMultipleScattering* proc = 0;
1198  const std::vector<G4VMultipleScattering*> v =
1200  G4int n = v.size();
1201  for(G4int i=0; i<n; ++i) {
1202  if((v[i])->GetProcessName() == processName) {
1203  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1204  if(ActiveForParticle(part, p)) {
1205  proc = v[i];
1206  break;
1207  }
1208  }
1209  }
1210  return proc;
1211 }
1212 
1213 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1214 
1216  const G4String& processName)
1217 {
1218  G4VProcess* proc = 0;
1219  const G4ProcessManager* procman = part->GetProcessManager();
1220  G4ProcessVector* pv = procman->GetProcessList();
1221  G4int nproc = pv->size();
1222  for(G4int i=0; i<nproc; ++i) {
1223  if(processName == (*pv)[i]->GetProcessName()) {
1224  proc = (*pv)[i];
1225  break;
1226  }
1227  }
1228  return proc;
1229 }
1230 
1231 
1232 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1233 
1235  G4VProcess* proc)
1236 {
1237  G4ProcessManager* pm = part->GetProcessManager();
1238  G4ProcessVector* pv = pm->GetProcessList();
1239  G4int n = pv->size();
1240  G4bool res = false;
1241  for(G4int i=0; i<n; ++i) {
1242  if((*pv)[i] == proc) {
1243  if(pm->GetProcessActivation(i)) { res = true; }
1244  break;
1245  }
1246  }
1247  return res;
1248 }
1249 
1250 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1251 
1253 {
1254  if(mat) {
1255  currentMaterial = mat;
1256  currentMaterialName = mat->GetName();
1257  } else {
1258  currentMaterial = 0;
1259  currentMaterialName = "";
1260  }
1261 }
1262 
1263 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1264 
1266 {
1268 }
1269 
1270 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1271 
1273 {
1274  G4bool isFound = false;
1275  if(currentMaterial) {
1277  for(size_t i=0; i<nn; ++i) {
1278  if(Z == G4lrint(currentMaterial->GetElement(i)->GetZ())) {
1279  isFound = true;
1280  break;
1281  }
1282  }
1283  }
1284  if(!isFound) {
1286  }
1287 }
1288 
1289 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1290 
1292 {
1293  verbose = verb;
1294 }
1295 
1296 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1297 
static const double cm
Definition: G4SIunits.hh:118
G4PhysicsTable * LambdaTable() const
static const double gram
Definition: G4SIunits.hh:175
G4VEmModel * loweModel
G4VEmProcess * FindDiscreteProcess(const G4ParticleDefinition *, const G4String &processName)
G4ProductionCuts * GetProductionCuts() const
const G4PhysicsTable * currentLambda
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:641
void FindLambdaTable(const G4ParticleDefinition *, const G4String &processName, G4double kinEnergy)
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:258
G4Material * FindOrBuildMaterial(const G4String &name, G4bool isotopes=true, G4bool warning=false)
const G4ParticleDefinition * theGenericIon
static const double MeV
Definition: G4SIunits.hh:211
G4VEmModel * SelectModel(G4double kinEnergy, size_t idx)
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
virtual void InitialiseForElement(const G4ParticleDefinition *, G4int Z)
Definition: G4VEmModel.cc:244
const G4String & GetName() const
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
static G4LossTableManager * Instance()
static const double cm2
Definition: G4SIunits.hh:119
virtual void CorrectionsAlongStep(const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
Definition: G4VEmModel.cc:362
void SetDynamicMassCharge(G4double massratio, G4double charge2ratio)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4double GetDEDX(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
const G4Material * FindMaterial(const G4String &)
G4Material * FindOrBuildSimpleMaterial(G4int Z, G4bool warning=false)
const G4Material * cutMaterial
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:604
std::vector< const G4MaterialCutsCouple * > localCouples
G4double LowestElectronEnergy() const
const G4ParticleDefinition * currentParticle
G4PhysicsTable * RangeTableForLoss() const
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4String name
Definition: TRTMaterials.hh:40
const G4String & GetName() const
Definition: G4Material.hh:178
void PrintRangeTable(const G4ParticleDefinition *)
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:491
virtual G4double GetShellIonisationCrossSectionPerAtom(const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=0)=0
virtual void SetupForMaterial(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:403
G4double GetZ() const
Definition: G4Element.hh:131
G4bool UpdateCouple(const G4Material *, G4double cut)
G4double GetDensity() const
Definition: G4Material.hh:180
G4double ComputeTotalDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, G4double cut=DBL_MAX)
G4double ComputeElectronicDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *mat, G4double cut=DBL_MAX)
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:826
G4double EffectiveChargeCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4EmCorrections * corr
G4ionEffectiveCharge * ionEffCharge
const G4MaterialCutsCouple * currentCouple
G4double GetCrossSectionPerVolume(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, const G4Region *r=0)
G4double GetMeanFreePath(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, const G4Region *r=0)
void SetupMaterial(const G4Material *)
void PrintInverseRangeTable(const G4ParticleDefinition *)
G4double GetRange(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=0)
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:202
G4ProcessManager * GetProcessManager() const
const std::vector< G4VEmProcess * > & GetEmProcessVector()
int G4int
Definition: G4Types.hh:78
G4DynamicParticle dynParticle
static G4NistManager * Instance()
const G4String & GetParticleName() const
G4double ComputeShellIonisationCrossSectionPerAtom(const G4String &part, G4int Z, G4AtomicShellEnumerator shell, G4double kinEnergy, const G4Material *mat=0)
G4double GetCSDARange(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
G4double GetDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=0)
static const G4double reg
static G4RegionStore * GetInstance()
G4double GetCSDARange(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=0)
G4bool BuildCSDARange() const
G4double cutenergy[3]
G4double GetShellIonisationCrossSectionPerAtom(const G4String &part, G4int Z, G4AtomicShellEnumerator shell, G4double kinEnergy)
G4PhysicsTable * LambdaTable() const
G4VEnergyLossProcess * FindEnergyLossProcess(const G4ParticleDefinition *)
G4IonTable * GetIonTable() const
virtual void InitialiseForMaterial(const G4ParticleDefinition *, const G4Material *)
Definition: G4VEmModel.cc:225
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
G4VEnergyLossProcess * currentProcess
const G4MaterialCutsCouple * FindCouple(const G4Material *, const G4Region *r=0)
G4PhysicsTable * DEDXTable() const
G4double ComputeDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=DBL_MAX)
G4String currentName
G4double NuclearDEDX(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy, G4bool fluct=true)
const G4ParticleDefinition * baseParticle
bool G4bool
Definition: G4Types.hh:79
G4EmCorrections * EmCorrections()
G4bool FindEmModel(const G4ParticleDefinition *, const G4String &processName, G4double kinEnergy)
static const double nm
Definition: G4SIunits.hh:111
G4double ComputeNuclearDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *)
G4double GetRangeFromRestricteDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=0)
G4DataVector localCuts
static const double GeV
Definition: G4SIunits.hh:214
const G4String & GetParticleType() const
const G4ParticleDefinition * BaseParticle() const
const G4int n
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4double GetEnergy(const G4ParticleDefinition *aParticle, G4double range, const G4MaterialCutsCouple *couple)
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void CheckMaterial(G4int Z)
G4bool ActiveForParticle(const G4ParticleDefinition *part, G4VProcess *proc)
const std::vector< G4VEnergyLossProcess * > & GetEnergyLossProcessVector()
G4VMultipleScattering * FindMscProcess(const G4ParticleDefinition *, const G4String &processName)
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:313
void SetKineticEnergy(G4double aEnergy)
virtual G4double ComputeShellIonisationCrossSectionPerAtom(const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=0)=0
G4double GetRangeFromRestricteDEDX(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
std::vector< const G4Material * > localMaterials
G4VEnergyLossProcess * FindEnLossProcess(const G4ParticleDefinition *, const G4String &processName)
G4int size() const
virtual G4double GetChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:345
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
static G4ProductionCutsTable * GetProductionCutsTable()
static const double eV
Definition: G4SIunits.hh:212
static G4Positron * Positron()
Definition: G4Positron.cc:94
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:93
virtual G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.cc:249
G4VEmModel * SelectModelForMaterial(G4double kinEnergy, size_t &idxRegion) const
const G4ParticleDefinition * FindParticle(const G4String &)
G4PhysicsTable * InverseRangeTable() const
G4double GetPDGMass() const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
G4double chargeSquare
static G4ParticleTable * GetParticleTable()
int G4lrint(double ad)
Definition: templates.hh:163
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4double ComputeMeanFreePath(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=0.0)
G4double ComputeCrossSectionPerShell(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, G4int Z, G4int shellIdx, G4double cut=0.0)
G4double GetKinEnergy(G4double range, const G4ParticleDefinition *, const G4Material *, const G4Region *r=0)
static const double g
Definition: G4SIunits.hh:180
const G4Region * FindRegion(const G4String &)
G4IonTable * ionTable
G4String currentMaterialName
const G4ParticleDefinition * FindIon(G4int Z, G4int A)
G4String currentProcessName
G4String currentParticleName
const G4double x[NPOINTSGL]
static G4EmParameters * Instance()
G4double ComputeCrossSectionPerVolume(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=0.0)
virtual G4double ComputeCrossSectionPerShell(const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:323
G4double ComputeGammaAttenuationLength(G4double kinEnergy, const G4Material *)
static const double mole
Definition: G4SIunits.hh:283
G4bool GetProcessActivation(G4VProcess *aProcess) const
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
G4double ComputeCrossSectionPerAtom(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, G4double Z, G4double A, G4double cut=0.0)
const G4String & GetName() const
Definition: G4VEmModel.hh:795
G4VProcess * FindProcess(const G4ParticleDefinition *part, const G4String &processName)
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
void SetVerbose(G4int val)
void PrintDEDXTable(const G4ParticleDefinition *)
G4VAtomDeexcitation * AtomDeexcitation()
static const double barn
Definition: G4SIunits.hh:104
static const double keV
Definition: G4SIunits.hh:213
const G4ParticleDefinition * lambdaParticle
double G4double
Definition: G4Types.hh:76
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double ComputeEnergyCutFromRangeCut(G4double range, const G4ParticleDefinition *, const G4Material *)
const G4Material * currentMaterial
G4VEmModel * currentModel
G4double GetPDGCharge() const
G4bool UpdateParticle(const G4ParticleDefinition *, G4double kinEnergy)
#define DBL_MAX
Definition: templates.hh:83
static const double mm
Definition: G4SIunits.hh:114
G4LossTableManager * manager
G4double ComputeDEDXForCutInRange(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *mat, G4double rangecut=DBL_MAX)
G4AtomicShellEnumerator
const std::vector< G4VMultipleScattering * > & GetMultipleScatteringVector()
G4VEmModel * SelectModelForMaterial(G4double kinEnergy, size_t &idx) const
G4ProcessVector * GetProcessList() const
G4EmParameters * theParameters
G4NistManager * nist