Geant4  10.02.p01
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 95476 2016-02-12 09:39:50Z 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 "G4NistManager.hh"
69 #include "G4VEmProcess.hh"
70 #include "G4VEnergyLossProcess.hh"
71 #include "G4VMultipleScattering.hh"
72 #include "G4Material.hh"
73 #include "G4MaterialCutsCouple.hh"
74 #include "G4ParticleDefinition.hh"
75 #include "G4ParticleTable.hh"
76 #include "G4IonTable.hh"
77 #include "G4PhysicsTable.hh"
78 #include "G4ProductionCutsTable.hh"
79 #include "G4ProcessManager.hh"
80 #include "G4ionEffectiveCharge.hh"
81 #include "G4RegionStore.hh"
82 #include "G4Element.hh"
83 #include "G4EmCorrections.hh"
84 #include "G4GenericIon.hh"
85 #include "G4ProcessVector.hh"
86 #include "G4Gamma.hh"
87 #include "G4Electron.hh"
88 #include "G4Positron.hh"
89 
90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
91 
93 {
97  nLocalMaterials = 0;
98  verbose = 0;
100  currentCouple = nullptr;
101  currentMaterial = cutMaterial = nullptr;
102  currentParticle = nullptr;
103  lambdaParticle = nullptr;
104  baseParticle = nullptr;
105  currentLambda = nullptr;
106  currentModel = nullptr;
107  currentProcess = nullptr;
108  loweModel = nullptr;
109  chargeSquare = 1.0;
110  massRatio = 1.0;
111  mass = 0.0;
112  currentCut = 0.0;
113  cutenergy[0] = cutenergy[1] = cutenergy[2] = DBL_MAX;
116  currentName = "";
117  lambdaName = "";
121  isIon = false;
122  isApplicable = false;
123 }
124 
125 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
126 
128 {
129  delete ionEffCharge;
130  for (G4int i=0; i<nLocalMaterials; ++i) {
131  delete localCouples[i];
132  }
133 }
134 
135 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
136 
138  const G4ParticleDefinition* p,
139  const G4Material* mat,
140  const G4Region* region)
141 {
142  G4double res = 0.0;
143  const G4MaterialCutsCouple* couple = FindCouple(mat, region);
144  if(couple && UpdateParticle(p, kinEnergy) ) {
145  res = manager->GetDEDX(p, kinEnergy, couple);
146 
147  if(isIon) {
148  if(FindEmModel(p, currentProcessName, kinEnergy)) {
149  G4double length = CLHEP::nm;
150  G4double eloss = res*length;
151  //G4cout << "### GetDEDX: E= " << kinEnergy << " dedx0= " << res
152  // << " de= " << eloss << G4endl;;
153  G4double niel = 0.0;
154  dynParticle.SetKineticEnergy(kinEnergy);
155  currentModel->GetChargeSquareRatio(p, mat, kinEnergy);
156  currentModel->CorrectionsAlongStep(couple,&dynParticle,eloss,niel,length);
157  res = eloss/length;
158  //G4cout << " de1= " << eloss << " res1= " << res
159  // << " " << p->GetParticleName() <<G4endl;;
160  }
161  }
162 
163  if(verbose>0) {
164  G4cout << "G4EmCalculator::GetDEDX: E(MeV)= " << kinEnergy/MeV
165  << " DEDX(MeV/mm)= " << res*mm/MeV
166  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
167  << " " << p->GetParticleName()
168  << " in " << mat->GetName()
169  << " isIon= " << isIon
170  << G4endl;
171  }
172  }
173  return res;
174 }
175 
176 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
177 
179  const G4ParticleDefinition* p,
180  const G4Material* mat,
181  const G4Region* region)
182 {
183  G4double res = 0.0;
184  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
185  if(couple && UpdateParticle(p, kinEnergy)) {
186  res = manager->GetRangeFromRestricteDEDX(p, kinEnergy, couple);
187  if(verbose>1) {
188  G4cout << " G4EmCalculator::GetRangeFromRestrictedDEDX: E(MeV)= "
189  << kinEnergy/MeV
190  << " range(mm)= " << res/mm
191  << " " << p->GetParticleName()
192  << " in " << mat->GetName()
193  << G4endl;
194  }
195  }
196  return res;
197 }
198 
199 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
200 
202  const G4ParticleDefinition* p,
203  const G4Material* mat,
204  const G4Region* region)
205 {
206  G4double res = 0.0;
207  if(!manager->BuildCSDARange()) {
209  ed << "G4EmCalculator::GetCSDARange: CSDA table is not built; "
210  << " use UI command: /process/eLoss/CSDARange true";
211  G4Exception("G4EmCalculator::GetCSDARange", "em0077",
212  JustWarning, ed);
213  return res;
214  }
215 
216  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
217  if(couple && UpdateParticle(p, kinEnergy)) {
218  res = manager->GetCSDARange(p, kinEnergy, couple);
219  if(verbose>1) {
220  G4cout << " G4EmCalculator::GetCSDARange: E(MeV)= " << kinEnergy/MeV
221  << " range(mm)= " << res/mm
222  << " " << p->GetParticleName()
223  << " in " << mat->GetName()
224  << G4endl;
225  }
226  }
227  return res;
228 }
229 
230 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
231 
233  const G4ParticleDefinition* p,
234  const G4Material* mat,
235  const G4Region* region)
236 {
237  G4double res = 0.0;
238  if(manager->BuildCSDARange()) {
239  res = GetCSDARange(kinEnergy, p, mat, region);
240  } else {
241  res = GetRangeFromRestricteDEDX(kinEnergy, p, mat, region);
242  }
243  return res;
244 }
245 
246 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
247 
249  const G4ParticleDefinition* p,
250  const G4Material* mat,
251  const G4Region* region)
252 {
253  G4double res = 0.0;
254  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
255  if(couple && UpdateParticle(p, 1.0*GeV)) {
256  res = manager->GetEnergy(p, range, couple);
257  if(verbose>0) {
258  G4cout << "G4EmCalculator::GetKinEnergy: Range(mm)= " << range/mm
259  << " KinE(MeV)= " << res/MeV
260  << " " << p->GetParticleName()
261  << " in " << mat->GetName()
262  << G4endl;
263  }
264  }
265  return res;
266 }
267 
268 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
269 
271  const G4ParticleDefinition* p,
272  const G4String& processName,
273  const G4Material* mat,
274  const G4Region* region)
275 {
276  G4double res = 0.0;
277  const G4MaterialCutsCouple* couple = FindCouple(mat,region);
278 
279  if(couple && UpdateParticle(p, kinEnergy)) {
280  G4int idx = couple->GetIndex();
281  FindLambdaTable(p, processName, kinEnergy);
282 
283  if(currentLambda) {
284  G4double e = kinEnergy*massRatio;
285  res = (((*currentLambda)[idx])->Value(e))*chargeSquare;
286  } else {
287  res = ComputeCrossSectionPerVolume(kinEnergy, p, processName, mat,
288  kinEnergy);
289  }
290  if(verbose>0) {
291  G4cout << "G4EmCalculator::GetXSPerVolume: E(MeV)= " << kinEnergy/MeV
292  << " cross(cm-1)= " << res*cm
293  << " " << p->GetParticleName()
294  << " in " << mat->GetName();
295  if(verbose>1)
296  G4cout << " idx= " << idx << " Escaled((MeV)= "
297  << kinEnergy*massRatio
298  << " q2= " << chargeSquare;
299  G4cout << G4endl;
300  }
301  }
302  return res;
303 }
304 
305 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
306 
308  const G4String& particle,
309  G4int Z,
311  G4double kinEnergy)
312 {
313  G4double res = 0.0;
314  const G4ParticleDefinition* p = FindParticle(particle);
316  if(p && ad) {
317  res = ad->GetShellIonisationCrossSectionPerAtom(p, Z, shell, kinEnergy);
318  }
319  return res;
320 }
321 
322 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
323 
325  const G4ParticleDefinition* p,
326  const G4String& processName,
327  const G4Material* mat,
328  const G4Region* region)
329 {
330  G4double res = DBL_MAX;
331  G4double x = GetCrossSectionPerVolume(kinEnergy,p, processName, mat,region);
332  if(x > 0.0) { res = 1.0/x; }
333  if(verbose>1) {
334  G4cout << "G4EmCalculator::GetMeanFreePath: E(MeV)= " << kinEnergy/MeV
335  << " MFP(mm)= " << res/mm
336  << " " << p->GetParticleName()
337  << " in " << mat->GetName()
338  << G4endl;
339  }
340  return res;
341 }
342 
343 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
344 
346 {
348  G4cout << "##### DEDX Table for " << p->GetParticleName() << G4endl;
349  if(elp) G4cout << *(elp->DEDXTable()) << G4endl;
350 }
351 
352 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
353 
355 {
357  G4cout << "##### Range Table for " << p->GetParticleName() << G4endl;
358  if(elp) G4cout << *(elp->RangeTableForLoss()) << G4endl;
359 }
360 
361 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
362 
364 {
366  G4cout << "### G4EmCalculator: Inverse Range Table for "
367  << p->GetParticleName() << G4endl;
368  if(elp) G4cout << *(elp->InverseRangeTable()) << G4endl;
369 }
370 
371 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
372 
374  const G4ParticleDefinition* p,
375  const G4String& processName,
376  const G4Material* mat,
377  G4double cut)
378 {
379  SetupMaterial(mat);
380  G4double res = 0.0;
381  if(verbose > 1) {
382  G4cout << "### G4EmCalculator::ComputeDEDX: " << p->GetParticleName()
383  << " in " << currentMaterialName
384  << " e(MeV)= " << kinEnergy/MeV << " cut(MeV)= " << cut/MeV
385  << G4endl;
386  }
387  if(UpdateParticle(p, kinEnergy)) {
388  if(FindEmModel(p, processName, kinEnergy)) {
389 
390  // Special case of ICRU'73 model
391  if(currentModel->GetName() == "ParamICRU73") {
392  res = currentModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cut);
393  if(verbose > 1) {
394  G4cout << " ICRU73 ion E(MeV)= " << kinEnergy << " ";
395  G4cout << currentModel->GetName() << ": DEDX(MeV/mm)= " << res*mm/MeV
396  << " DEDX(MeV*cm^2/g)= "
397  << res*gram/(MeV*cm2*mat->GetDensity())
398  << G4endl;
399  }
400  } else {
401 
402  G4double escaled = kinEnergy*massRatio;
403  if(baseParticle) {
405  mat, baseParticle, escaled, cut) * chargeSquare;
406  if(verbose > 1) {
408  << " Escaled(MeV)= " << escaled;
409  }
410  } else {
411  res = currentModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cut);
412  if(verbose > 1) { G4cout << " no basePart E(MeV)= " << kinEnergy << " "; }
413  }
414  if(verbose > 1) {
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 
421  // emulate smoothing procedure
423  // G4cout << "massRatio= " << massRatio << " eth= " << eth << G4endl;
424  if(loweModel) {
425  G4double res0 = 0.0;
426  G4double res1 = 0.0;
427  if(baseParticle) {
428  res1 = currentModel->ComputeDEDXPerVolume(mat, baseParticle, eth, cut)
429  * chargeSquare;
430  res0 = loweModel->ComputeDEDXPerVolume(mat, baseParticle, eth, cut)
431  * chargeSquare;
432  } else {
433  res1 = currentModel->ComputeDEDXPerVolume(mat, p, eth, cut);
434  res0 = loweModel->ComputeDEDXPerVolume(mat, p, eth, cut);
435  }
436  if(verbose > 1) {
437  G4cout << "At boundary energy(MeV)= " << eth/MeV
438  << " DEDX(MeV/mm)= " << res1*mm/MeV
439  << G4endl;
440  }
441 
442  //G4cout << "eth= " << eth << " escaled= " << escaled
443  // << " res0= " << res0 << " res1= "
444  // << res1 << " q2= " << chargeSquare << G4endl;
445 
446  if(res1 > 0.0 && escaled > 0.0) {
447  res *= (1.0 + (res0/res1 - 1.0)*eth/escaled);
448  }
449  }
450 
451  // low energy correction for ions
452  if(isIon) {
453  G4double length = CLHEP::nm;
454  const G4Region* r = 0;
455  const G4MaterialCutsCouple* couple = FindCouple(mat, r);
456  G4double eloss = res*length;
457  G4double niel = 0.0;
458  dynParticle.SetKineticEnergy(kinEnergy);
459  currentModel->GetChargeSquareRatio(p, mat, kinEnergy);
460  currentModel->CorrectionsAlongStep(couple,&dynParticle,eloss,niel,length);
461  res = eloss/length;
462 
463  if(verbose > 1) {
464  G4cout << "After Corrections: DEDX(MeV/mm)= " << res*mm/MeV
465  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
466  << G4endl;
467  }
468  }
469  }
470  }
471  if(verbose > 0) {
472  G4cout << "Sum: E(MeV)= " << kinEnergy/MeV
473  << " DEDX(MeV/mm)= " << res*mm/MeV
474  << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
475  << " cut(MeV)= " << cut/MeV
476  << " " << p->GetParticleName()
477  << " in " << currentMaterialName
478  << " Zi^2= " << chargeSquare
479  << " isIon=" << isIon
480  << G4endl;
481  }
482  }
483  return res;
484 }
485 
486 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
487 
489  const G4ParticleDefinition* part,
490  const G4Material* mat,
491  G4double cut)
492 {
493  SetupMaterial(mat);
494  G4double dedx = 0.0;
495  if(UpdateParticle(part, kinEnergy)) {
496 
498  const std::vector<G4VEnergyLossProcess*> vel =
499  lManager->GetEnergyLossProcessVector();
500  G4int n = vel.size();
501 
502  //G4cout << "ComputeElectronicDEDX for " << part->GetParticleName()
503  // << " n= " << n << G4endl;
504 
505  for(G4int i=0; i<n; ++i) {
506  if(vel[i]) {
507  G4VProcess* p = reinterpret_cast<G4VProcess*>(vel[i]);
508  if(ActiveForParticle(part, p)) {
509  //G4cout << "idx= " << i << " " << (vel[i])->GetProcessName()
510  // << " " << (vel[i])->Particle()->GetParticleName() << G4endl;
511  dedx += ComputeDEDX(kinEnergy,part,(vel[i])->GetProcessName(),mat,cut);
512  }
513  }
514  }
515  }
516  return dedx;
517 }
518 
519 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
520 
522  const G4ParticleDefinition* part,
523  const G4Material* mat,
524  G4double rangecut)
525 {
526  SetupMaterial(mat);
527  G4double dedx = 0.0;
528  if(UpdateParticle(part, kinEnergy)) {
529 
531  const std::vector<G4VEnergyLossProcess*> vel =
532  lManager->GetEnergyLossProcessVector();
533  G4int n = vel.size();
534 
535  if(mat != cutMaterial) {
536  cutMaterial = mat;
540  }
541 
542  //G4cout << "ComputeElectronicDEDX for " << part->GetParticleName()
543  // << " n= " << n << G4endl;
544 
545  for(G4int i=0; i<n; ++i) {
546  if(vel[i]) {
547  G4VProcess* p = reinterpret_cast<G4VProcess*>(vel[i]);
548  if(ActiveForParticle(part, p)) {
549  //G4cout << "idx= " << i << " " << (vel[i])->GetProcessName()
550  // << " " << (vel[i])->Particle()->GetParticleName() << G4endl;
551  const G4ParticleDefinition* sec = (vel[i])->SecondaryParticle();
552  G4int idx = 0;
553  if(sec == G4Electron::Electron()) { idx = 1; }
554  else if(sec == G4Positron::Positron()) { idx = 2; }
555 
556  dedx += ComputeDEDX(kinEnergy,part,(vel[i])->GetProcessName(),
557  mat,cutenergy[idx]);
558  }
559  }
560  }
561  }
562  return dedx;
563 }
564 
565 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
566 
568  const G4ParticleDefinition* part,
569  const G4Material* mat,
570  G4double cut)
571 {
572  G4double dedx = ComputeElectronicDEDX(kinEnergy,part,mat,cut);
573  if(mass > 700.*MeV) { dedx += ComputeNuclearDEDX(kinEnergy,part,mat); }
574  return dedx;
575 }
576 
577 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
578 
580  const G4ParticleDefinition* p,
581  const G4Material* mat)
582 {
583 
584  G4double res = corr->NuclearDEDX(p, mat, kinEnergy, false);
585 
586  if(verbose > 1) {
587  G4cout << p->GetParticleName() << " E(MeV)= " << kinEnergy/MeV
588  << " NuclearDEDX(MeV/mm)= " << res*mm/MeV
589  << " NuclearDEDX(MeV*cm^2/g)= "
590  << res*gram/(MeV*cm2*mat->GetDensity())
591  << G4endl;
592  }
593  return res;
594 }
595 
596 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
597 
599  G4double kinEnergy,
600  const G4ParticleDefinition* p,
601  const G4String& processName,
602  const G4Material* mat,
603  G4double cut)
604 {
605  SetupMaterial(mat);
606  G4double res = 0.0;
607  if(UpdateParticle(p, kinEnergy)) {
608  if(FindEmModel(p, processName, kinEnergy)) {
609  G4double e = kinEnergy;
610  if(baseParticle) {
611  e *= kinEnergy*massRatio;
613  mat, baseParticle, e, cut, e) * chargeSquare;
614  } else {
615  res = currentModel->CrossSectionPerVolume(mat, p, e, cut, e);
616  }
617  if(verbose>0) {
618  G4cout << "G4EmCalculator::ComputeXSPerVolume: E(MeV)= " << kinEnergy/MeV
619  << " cross(cm-1)= " << res*cm
620  << " cut(keV)= " << cut/keV
621  << " " << p->GetParticleName()
622  << " in " << mat->GetName()
623  << G4endl;
624  }
625  }
626  }
627  return res;
628 }
629 
630 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
631 
633  G4double kinEnergy,
634  const G4ParticleDefinition* p,
635  const G4String& processName,
636  G4double Z, G4double A,
637  G4double cut)
638 {
639  G4double res = 0.0;
640  if(UpdateParticle(p, kinEnergy)) {
642  if(FindEmModel(p, processName, kinEnergy)) {
643  G4double e = kinEnergy;
644  if(baseParticle) {
645  e *= kinEnergy*massRatio;
648  baseParticle, e, Z, A, cut) * chargeSquare;
649  } else {
651  res = currentModel->ComputeCrossSectionPerAtom(p, e, Z, A, cut);
652  }
653  if(verbose>0) {
654  G4cout << "E(MeV)= " << kinEnergy/MeV
655  << " cross(barn)= " << res/barn
656  << " " << p->GetParticleName()
657  << " Z= " << Z << " A= " << A/(g/mole) << " g/mole"
658  << G4endl;
659  }
660  }
661  }
662  return res;
663 }
664 
665 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
666 
668  const G4ParticleDefinition* p,
669  const G4String& processName,
670  G4int Z, G4int shellIdx,
671  G4double cut)
672 {
673  G4double res = 0.0;
674  if(UpdateParticle(p, kinEnergy)) {
675  CheckMaterial(Z);
676  if(FindEmModel(p, processName, kinEnergy)) {
677  G4double e = kinEnergy;
678  if(baseParticle) {
679  e *= kinEnergy*massRatio;
682  e, cut) * chargeSquare;
683  } else {
685  res = currentModel->ComputeCrossSectionPerAtom(p, Z, shellIdx, e, cut);
686  }
687  if(verbose>0) {
688  G4cout << "E(MeV)= " << kinEnergy/MeV
689  << " cross(barn)= " << res/barn
690  << " " << p->GetParticleName()
691  << " Z= " << Z << " shellIdx= " << shellIdx << G4endl;
692  }
693  }
694  }
695  return res;
696 }
697 
698 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
699 
700 G4double
702  const G4Material* mat)
703 {
704  G4double res = 0.0;
705  const G4ParticleDefinition* gamma = G4Gamma::Gamma();
706  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "conv", mat, 0.0);
707  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "compt", mat, 0.0);
708  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "phot", mat, 0.0);
709  res += ComputeCrossSectionPerVolume(kinEnergy, gamma, "Rayl", mat, 0.0);
710  if(res > 0.0) { res = 1.0/res; }
711  return res;
712 }
713 
714 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
715 
717  const G4String& particle,
718  G4int Z,
720  G4double kinEnergy,
721  const G4Material* mat)
722 {
723  G4double res = 0.0;
724  const G4ParticleDefinition* p = FindParticle(particle);
726  if(p && ad) {
727  res = ad->ComputeShellIonisationCrossSectionPerAtom(p, Z, shell,
728  kinEnergy, mat);
729  }
730  return res;
731 }
732 
733 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
734 
736  const G4ParticleDefinition* p,
737  const G4String& processName,
738  const G4Material* mat,
739  G4double cut)
740 {
741  G4double mfp = DBL_MAX;
742  G4double x = ComputeCrossSectionPerVolume(kinEnergy, p, processName, mat, cut);
743  if(x > 0.0) { mfp = 1.0/x; }
744  if(verbose>1) {
745  G4cout << "E(MeV)= " << kinEnergy/MeV
746  << " MFP(mm)= " << mfp/mm
747  << " " << p->GetParticleName()
748  << " in " << mat->GetName()
749  << G4endl;
750  }
751  return mfp;
752 }
753 
754 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
755 
757  G4double range,
758  const G4ParticleDefinition* part,
759  const G4Material* mat)
760 {
762  ConvertRangeToEnergy(part, mat, range);
763 }
764 
765 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
766 
768  G4double kinEnergy)
769 {
770  if(p != currentParticle) {
771 
772  // new particle
773  currentParticle = p;
774  dynParticle.SetDefinition(const_cast<G4ParticleDefinition*>(p));
775  dynParticle.SetKineticEnergy(kinEnergy);
776  baseParticle = 0;
778  massRatio = 1.0;
779  mass = p->GetPDGMass();
780  chargeSquare = 1.0;
782  currentProcessName = "";
783  isIon = false;
784 
785  // ionisation process exist
786  if(currentProcess) {
789 
790  // base particle is used
791  if(baseParticle) {
794  chargeSquare = q*q;
795  }
796 
797  if(p->GetParticleType() == "nucleus"
798  && currentParticleName != "deuteron"
799  && currentParticleName != "triton"
800  && currentParticleName != "alpha+"
801  && currentParticleName != "helium"
802  && currentParticleName != "hydrogen"
803  ) {
804  isIon = true;
807  if(verbose>1) {
808  G4cout << "\n G4EmCalculator::UpdateParticle: isIon 1 "
809  << p->GetParticleName()
810  << " in " << currentMaterial->GetName()
811  << " e= " << kinEnergy << G4endl;
812  }
813  }
814  }
815  }
816 
817  // Effective charge for ions
818  if(isIon) {
819  chargeSquare =
822  if(currentProcess) {
824  if(verbose>1) {
825  G4cout <<"\n NewIon: massR= "<< massRatio << " q2= "
826  << chargeSquare << " " << currentProcess << G4endl;
827  }
828  }
829  }
830  return true;
831 }
832 
833 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
834 
836 {
837  const G4ParticleDefinition* p = 0;
838  if(name != currentParticleName) {
840  if(!p) {
841  G4cout << "### WARNING: G4EmCalculator::FindParticle fails to find "
842  << name << G4endl;
843  }
844  } else {
845  p = currentParticle;
846  }
847  return p;
848 }
849 
850 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
851 
853 {
854  const G4ParticleDefinition* p = ionTable->GetIon(Z,A,0);
855  return p;
856 }
857 
858 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
859 
861 {
862  if(name != currentMaterialName) {
864  if(!currentMaterial) {
865  G4cout << "### WARNING: G4EmCalculator::FindMaterial fails to find "
866  << name << G4endl;
867  }
868  }
869  return currentMaterial;
870 }
871 
872 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
873 
875 {
876  const G4Region* r = 0;
877  if(reg != "" && reg != "world") {
879  } else {
880  r = G4RegionStore::GetInstance()->GetRegion("DefaultRegionForTheWorld");
881  }
882  return r;
883 }
884 
885 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
886 
888  const G4Material* material,
889  const G4Region* region)
890 {
891  const G4MaterialCutsCouple* couple = 0;
892  SetupMaterial(material);
893  if(currentMaterial) {
894  // Access to materials
895  const G4ProductionCutsTable* theCoupleTable=
897  const G4Region* r = region;
898  if(r) {
899  couple = theCoupleTable->GetMaterialCutsCouple(material,
900  r->GetProductionCuts());
901  } else {
903  size_t nr = store->size();
904  if(0 < nr) {
905  for(size_t i=0; i<nr; ++i) {
906  couple = theCoupleTable->GetMaterialCutsCouple(
907  material, ((*store)[i])->GetProductionCuts());
908  if(couple) { break; }
909  }
910  }
911  }
912  }
913  if(!couple) {
915  ed << "G4EmCalculator::FindCouple: fail for material <"
916  << currentMaterialName << ">";
917  if(region) { ed << " and region " << region->GetName(); }
918  G4Exception("G4EmCalculator::FindCouple", "em0078",
919  FatalException, ed);
920  }
921  return couple;
922 }
923 
924 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
925 
927 {
928  SetupMaterial(material);
929  if(!currentMaterial) { return false; }
930  for (G4int i=0; i<nLocalMaterials; ++i) {
931  if(material == localMaterials[i] && cut == localCuts[i]) {
934  currentCut = cut;
935  return true;
936  }
937  }
938  const G4MaterialCutsCouple* cc = new G4MaterialCutsCouple(material);
939  localMaterials.push_back(material);
940  localCouples.push_back(cc);
941  localCuts.push_back(cut);
942  nLocalMaterials++;
943  currentCouple = cc;
945  currentCut = cut;
946  return true;
947 }
948 
949 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
950 
952  const G4String& processName,
953  G4double kinEnergy)
954 {
955  // Search for the process
956  if (!currentLambda || p != lambdaParticle || processName != lambdaName) {
957  lambdaName = processName;
958  currentLambda = 0;
959  lambdaParticle = p;
960 
961  const G4ParticleDefinition* part = p;
962  if(isIon) { part = theGenericIon; }
963 
964  // Search for energy loss process
965  currentName = processName;
966  currentModel = 0;
967  loweModel = 0;
968 
969  G4VEnergyLossProcess* elproc = FindEnLossProcess(part, processName);
970  if(elproc) {
971  currentLambda = elproc->LambdaTable();
972  if(currentLambda) {
973  isApplicable = true;
974  if(verbose>1) {
975  G4cout << "G4VEnergyLossProcess is found out: " << currentName
976  << G4endl;
977  }
978  }
979  return;
980  }
981 
982  // Search for discrete process
983  G4VEmProcess* proc = FindDiscreteProcess(part, processName);
984  if(proc) {
985  currentLambda = proc->LambdaTable();
986  if(currentLambda) {
987  isApplicable = true;
988  if(verbose>1) {
989  G4cout << "G4VEmProcess is found out: " << currentName << G4endl;
990  }
991  }
992  return;
993  }
994 
995  // Search for msc process
996  G4VMultipleScattering* msc = FindMscProcess(part, processName);
997  if(msc) {
998  currentModel = msc->SelectModel(kinEnergy,0);
999  if(currentModel) {
1001  if(currentLambda) {
1002  isApplicable = true;
1003  if(verbose>1) {
1004  G4cout << "G4VMultipleScattering is found out: " << currentName
1005  << G4endl;
1006  }
1007  }
1008  }
1009  }
1010  }
1011 }
1012 
1013 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1014 
1016  const G4String& processName,
1017  G4double kinEnergy)
1018 {
1019  isApplicable = false;
1020  if(!p || !currentMaterial) {
1021  G4cout << "G4EmCalculator::FindEmModel WARNING: no particle"
1022  << " or materail defined; particle: " << p << G4endl;
1023  return isApplicable;
1024  }
1025  G4String partname = p->GetParticleName();
1026  const G4ParticleDefinition* part = p;
1027  G4double scaledEnergy = kinEnergy*massRatio;
1028  if(isIon) { part = theGenericIon; }
1029 
1030  if(verbose > 1) {
1031  G4cout << "## G4EmCalculator::FindEmModel for " << partname
1032  << " (type= " << p->GetParticleType()
1033  << ") and " << processName << " at E(MeV)= " << scaledEnergy
1034  << G4endl;
1035  if(p != part) { G4cout << " GenericIon is the base particle" << G4endl; }
1036  }
1037 
1038  // Search for energy loss process
1039  currentName = processName;
1040  currentModel = 0;
1041  loweModel = 0;
1042  size_t idx = 0;
1043 
1044  G4VEnergyLossProcess* elproc = FindEnLossProcess(part, processName);
1045  if(elproc) {
1046  currentModel = elproc->SelectModelForMaterial(scaledEnergy, idx);
1048  currentModel->SetupForMaterial(part, currentMaterial, scaledEnergy);
1050  if(eth > 0.0) {
1051  loweModel = elproc->SelectModelForMaterial(eth - CLHEP::eV, idx);
1052  if(loweModel == currentModel) { loweModel = nullptr; }
1053  else {
1056  }
1057  }
1058  }
1059 
1060  // Search for discrete process
1061  if(!currentModel) {
1062  G4VEmProcess* proc = FindDiscreteProcess(part, processName);
1063  if(proc) {
1064  currentModel = proc->SelectModelForMaterial(kinEnergy, idx);
1066  currentModel->SetupForMaterial(part, currentMaterial, kinEnergy);
1068  if(eth > 0.0) {
1069  loweModel = proc->SelectModelForMaterial(eth - CLHEP::eV, idx);
1070  if(loweModel == currentModel) { loweModel = 0; }
1071  else {
1074  }
1075  }
1076  }
1077  }
1078 
1079  // Search for msc process
1080  if(!currentModel) {
1081  G4VMultipleScattering* proc = FindMscProcess(part, processName);
1082  if(proc) {
1083  currentModel = proc->SelectModel(kinEnergy, idx);
1084  loweModel = 0;
1085  }
1086  }
1087  if(currentModel) {
1088  if(loweModel == currentModel) { loweModel = 0; }
1089  isApplicable = true;
1091  if(loweModel) {
1093  }
1094  if(verbose > 1) {
1095  G4cout << " Model <" << currentModel->GetName()
1096  << "> Emin(MeV)= " << currentModel->LowEnergyLimit()/MeV
1097  << " for " << part->GetParticleName();
1098  if(elproc) {
1099  G4cout << " and " << elproc->GetProcessName() << " " << elproc
1100  << G4endl;
1101  }
1102  if(loweModel) {
1103  G4cout << " LowEnergy model <" << loweModel->GetName() << ">";
1104  }
1105  G4cout << G4endl;
1106  }
1107  }
1108  return isApplicable;
1109 }
1110 
1111 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1112 
1114  const G4ParticleDefinition* p)
1115 {
1116  G4VEnergyLossProcess* elp = 0;
1117  G4String partname = p->GetParticleName();
1118  const G4ParticleDefinition* part = p;
1119 
1120  if(p->GetParticleType() == "nucleus"
1121  && currentParticleName != "deuteron"
1122  && currentParticleName != "triton"
1123  && currentParticleName != "He3"
1124  && currentParticleName != "alpha"
1125  && currentParticleName != "alpha+"
1126  && currentParticleName != "helium"
1127  && currentParticleName != "hydrogen"
1128  ) { part = theGenericIon; }
1129 
1130  elp = manager->GetEnergyLossProcess(part);
1131  /*
1132  G4cout << "\n G4EmCalculator::FindEnergyLossProcess: for " << p->GetParticleName()
1133  << " found " << elp->GetProcessName() << " of "
1134  << elp->Particle()->GetParticleName() << " " << elp << G4endl;
1135  */
1136  return elp;
1137 }
1138 
1139 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1140 
1143  const G4String& processName)
1144 {
1145  G4VEnergyLossProcess* proc = 0;
1146  const std::vector<G4VEnergyLossProcess*> v =
1148  G4int n = v.size();
1149  for(G4int i=0; i<n; ++i) {
1150  if((v[i])->GetProcessName() == processName) {
1151  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1152  if(ActiveForParticle(part, p)) {
1153  proc = v[i];
1154  break;
1155  }
1156  }
1157  }
1158  return proc;
1159 }
1160 
1161 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1162 
1163 G4VEmProcess*
1165  const G4String& processName)
1166 {
1167  G4VEmProcess* proc = 0;
1168  const std::vector<G4VEmProcess*> v =
1170  G4int n = v.size();
1171  for(G4int i=0; i<n; ++i) {
1172  if((v[i])->GetProcessName() == processName) {
1173  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1174  if(ActiveForParticle(part, p)) {
1175  proc = v[i];
1176  break;
1177  }
1178  }
1179  }
1180  return proc;
1181 }
1182 
1183 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1184 
1187  const G4String& processName)
1188 {
1189  G4VMultipleScattering* proc = 0;
1190  const std::vector<G4VMultipleScattering*> v =
1192  G4int n = v.size();
1193  for(G4int i=0; i<n; ++i) {
1194  if((v[i])->GetProcessName() == processName) {
1195  G4VProcess* p = reinterpret_cast<G4VProcess*>(v[i]);
1196  if(ActiveForParticle(part, p)) {
1197  proc = v[i];
1198  break;
1199  }
1200  }
1201  }
1202  return proc;
1203 }
1204 
1205 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1206 
1208  const G4String& processName)
1209 {
1210  G4VProcess* proc = 0;
1211  const G4ProcessManager* procman = part->GetProcessManager();
1212  G4ProcessVector* pv = procman->GetProcessList();
1213  G4int nproc = pv->size();
1214  for(G4int i=0; i<nproc; ++i) {
1215  if(processName == (*pv)[i]->GetProcessName()) {
1216  proc = (*pv)[i];
1217  break;
1218  }
1219  }
1220  return proc;
1221 }
1222 
1223 
1224 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1225 
1227  G4VProcess* proc)
1228 {
1229  G4ProcessManager* pm = part->GetProcessManager();
1230  G4ProcessVector* pv = pm->GetProcessList();
1231  G4int n = pv->size();
1232  G4bool res = false;
1233  for(G4int i=0; i<n; ++i) {
1234  if((*pv)[i] == proc) {
1235  if(pm->GetProcessActivation(i)) { res = true; }
1236  break;
1237  }
1238  }
1239  return res;
1240 }
1241 
1242 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1243 
1245 {
1246  if(mat) {
1247  currentMaterial = mat;
1248  currentMaterialName = mat->GetName();
1249  } else {
1250  currentMaterial = 0;
1251  currentMaterialName = "";
1252  }
1253 }
1254 
1255 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1256 
1258 {
1260 }
1261 
1262 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1263 
1265 {
1266  G4bool isFound = false;
1267  if(currentMaterial) {
1269  for(size_t i=0; i<nn; ++i) {
1270  if(Z == G4lrint(currentMaterial->GetElement(i)->GetZ())) {
1271  isFound = true;
1272  break;
1273  }
1274  }
1275  }
1276  if(!isFound) {
1278  }
1279 }
1280 
1281 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1282 
1284 {
1285  verbose = verb;
1286 }
1287 
1288 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1289 
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
const G4ParticleDefinition * currentParticle
G4PhysicsTable * RangeTableForLoss() const
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4String name
Definition: TRTMaterials.hh:40
G4bool BuildCSDARange() const
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)
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
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]
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
G4NistManager * nist