Geant4  10.02.p01
G4EmCalculator.hh
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.hh 90333 2015-05-26 08:28:09Z gcosmo $
27 //
28 //
29 // -------------------------------------------------------------------
30 //
31 // GEANT4 Class header file
32 //
33 //
34 // File name: G4EmCalculator
35 //
36 // Author: Vladimir Ivanchenko
37 //
38 // Creation date: 27.06.2004
39 //
40 // Modifications:
41 // 17.11.2004 Change signature of methods, add new methods (V.Ivanchenko)
42 // 11.01.2006 Add GetCSDARange (V.Ivanchenko)
43 // 26.01.2006 Rename GetRange -> GetRangeFromRestricteDEDX (V.Ivanchenko)
44 // 22.03.2006 Add ComputeElectronicDEDX and ComputeTotalDEDX (V.Ivanchenko)
45 // 29.09.2006 Add member loweModel (V.Ivanchenko)
46 // 15.03.2007 Add ComputeEnergyCutFromRangeCut methods (V.Ivanchenko)
47 //
48 // Class Description:
49 //
50 // Provide access to dE/dx and cross sections
51 
52 // -------------------------------------------------------------------
53 //
54 
55 #ifndef G4EmCalculator_h
56 #define G4EmCalculator_h 1
57 
58 #include <vector>
59 #include "globals.hh"
60 #include "G4DataVector.hh"
61 #include "G4DynamicParticle.hh"
62 #include "G4VAtomDeexcitation.hh"
63 
64 class G4LossTableManager;
65 class G4NistManager;
66 class G4Material;
69 class G4PhysicsTable;
70 class G4VEmModel;
72 class G4VEmProcess;
74 class G4VProcess;
76 class G4Region;
77 class G4Element;
78 class G4EmCorrections;
79 class G4IonTable;
80 
82 {
83 
84 public:
85 
87 
89 
90  //===========================================================================
91  // Methods to access precalculated dE/dx and cross sections
92  // Materials should exist in the list of the G4MaterialCutsCouple
93  //===========================================================================
94 
95  G4double GetDEDX(G4double kinEnergy, const G4ParticleDefinition*,
96  const G4Material*,
97  const G4Region* r = 0);
98  inline G4double GetDEDX(G4double kinEnergy, const G4String& part,
99  const G4String& mat,
100  const G4String& s = "world");
101 
103  const G4ParticleDefinition*,
104  const G4Material*,
105  const G4Region* r = 0);
106  inline G4double GetRangeFromRestricteDEDX(G4double kinEnergy,
107  const G4String& part,
108  const G4String& mat,
109  const G4String& s = "world");
110 
112  const G4Material*,
113  const G4Region* r = 0);
114  inline G4double GetCSDARange(G4double kinEnergy, const G4String& part,
115  const G4String& mat,
116  const G4String& s = "world");
117 
118  G4double GetRange(G4double kinEnergy, const G4ParticleDefinition*,
119  const G4Material*,
120  const G4Region* r = 0);
121  inline G4double GetRange(G4double kinEnergy, const G4String& part,
122  const G4String& mat,
123  const G4String& s = "world");
124 
126  const G4Material*,
127  const G4Region* r = 0);
128  inline G4double GetKinEnergy(G4double range, const G4String& part,
129  const G4String& mat,
130  const G4String& s = "world");
131 
133  G4double kinEnergy, const G4ParticleDefinition*,
134  const G4String& processName, const G4Material*,
135  const G4Region* r = 0);
137  G4double kinEnergy, const G4String& part, const G4String& proc,
138  const G4String& mat, const G4String& s = "world");
139 
141  const G4String& part, G4int Z,
143  G4double kinEnergy);
144 
146  const G4String& processName, const G4Material*,
147  const G4Region* r = 0);
148  inline G4double GetMeanFreePath(G4double kinEnergy, const G4String& part,
149  const G4String& proc, const G4String& mat,
150  const G4String& s = "world");
151 
153 
155 
157 
158  //===========================================================================
159  // Methods to calculate dE/dx and cross sections "on fly"
160  // Existing tables and G4MaterialCutsCouples are not used
161  //===========================================================================
162 
164  const G4String& processName, const G4Material*,
165  G4double cut = DBL_MAX);
166  inline G4double ComputeDEDX(G4double kinEnergy, const G4String& part,
167  const G4String& proc,
168  const G4String& mat, G4double cut = DBL_MAX);
169 
171  const G4ParticleDefinition*,
172  const G4Material* mat, G4double cut = DBL_MAX);
173  inline G4double ComputeElectronicDEDX(G4double kinEnergy, const G4String& part,
174  const G4String& mat, G4double cut = DBL_MAX);
175 
177  const G4ParticleDefinition*,
178  const G4Material* mat, G4double rangecut = DBL_MAX);
179  inline G4double ComputeDEDXForCutInRange(G4double kinEnergy, const G4String& part,
180  const G4String& mat,
181  G4double rangecut = DBL_MAX);
182 
184  const G4Material*);
185  inline G4double ComputeNuclearDEDX(G4double kinEnergy, const G4String& part,
186  const G4String& mat);
187 
189  const G4Material*, G4double cut = DBL_MAX);
190  inline G4double ComputeTotalDEDX(G4double kinEnergy, const G4String& part,
191  const G4String& mat, G4double cut = DBL_MAX);
192 
194  G4double kinEnergy, const G4ParticleDefinition*,
195  const G4String& processName, const G4Material*,
196  G4double cut = 0.0);
198  G4double kinEnergy, const G4String& part,
199  const G4String& proc,
200  const G4String& mat, G4double cut = 0.0);
201 
203  G4double kinEnergy, const G4ParticleDefinition*,
204  const G4String& processName, G4double Z, G4double A,
205  G4double cut = 0.0);
207  G4double kinEnergy, const G4String& part,
208  const G4String& processName, const G4Element*,
209  G4double cut = 0.0);
210 
212  G4double kinEnergy, const G4ParticleDefinition*,
213  const G4String& processName, G4int Z, G4int shellIdx,
214  G4double cut = 0.0);
216  G4double kinEnergy, const G4String& part,
217  const G4String& processName, const G4Element*,
218  G4int shellIdx,
219  G4double cut = 0.0);
220 
222  const G4Material*);
223 
225  const G4String& part, G4int Z,
227  G4double kinEnergy,
228  const G4Material* mat = 0);
229 
231  G4double kinEnergy, const G4ParticleDefinition*,
232  const G4String& processName, const G4Material*,
233  G4double cut = 0.0);
235  G4double kinEnergy, const G4String&, const G4String&,
236  const G4String& processName, G4double cut = 0.0);
237 
239  G4double range, const G4ParticleDefinition*,
240  const G4Material*);
242  G4double range, const G4String&,
243  const G4String&);
244 
245  //===========================================================================
246  // Methods to access particles, materials, regions, processes
247  //===========================================================================
248 
250 
252 
253  const G4Material* FindMaterial(const G4String&);
254 
255  const G4Region* FindRegion(const G4String&);
256 
258  const G4Region* r = 0);
259 
261  const G4String& processName);
262 
263  void SetupMaterial(const G4Material*);
264 
265  void SetupMaterial(const G4String&);
266 
267  void SetVerbose(G4int val);
268 
269  //===========================================================================
270  // Private methods
271  //===========================================================================
272 
273 private:
274 
276 
277  G4bool UpdateCouple(const G4Material*, G4double cut);
278 
280  const G4String& processName,
281  G4double kinEnergy);
282 
284  const G4String& processName,
285  G4double kinEnergy);
286 
288 
290  const G4String& processName);
291 
293  const G4String& processName);
294 
296  const G4String& processName);
297 
299  G4VProcess* proc);
300 
301  void CheckMaterial(G4int Z);
302 
305 
306  std::vector<const G4Material*> localMaterials;
307  std::vector<const G4MaterialCutsCouple*> localCouples;
308 
315 
317 
318  // cash
330 
334 
344 
348 };
349 
350 //....oooOO0OOooo.......oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
351 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
352 
353 inline
355  const G4String& material, const G4String& reg)
356 {
357  return GetDEDX(kinEnergy,FindParticle(particle),
358  FindMaterial(material),FindRegion(reg));
359 }
360 
361 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
362 
363 inline
365  const G4String& particle,
366  const G4String& material,
367  const G4String& reg)
368 {
369  return GetRangeFromRestricteDEDX(kinEnergy,FindParticle(particle),
370  FindMaterial(material),FindRegion(reg));
371 }
372 
373 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
374 
375 inline
377  const G4String& particle,
378  const G4String& material,
379  const G4String& reg)
380 {
381  return GetCSDARange(kinEnergy,FindParticle(particle),
382  FindMaterial(material),FindRegion(reg));
383 }
384 
385 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
386 
387 inline
389  const G4String& particle,
390  const G4String& material,
391  const G4String& reg)
392 {
393  return GetRange(kinEnergy,FindParticle(particle),
394  FindMaterial(material),FindRegion(reg));
395 }
396 
397 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
398 
399 inline
401  const G4String& material, const G4String& reg)
402 {
403  return GetKinEnergy(range,FindParticle(particle),
404  FindMaterial(material),FindRegion(reg));
405 }
406 
407 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
408 
409 inline
411  const G4String& particle,
412  const G4String& processName,
413  const G4String& material,
414  const G4String& reg)
415 {
416  return GetCrossSectionPerVolume(kinEnergy,FindParticle(particle),processName,
417  FindMaterial(material),FindRegion(reg));
418 }
419 
420 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
421 
422 inline
424  const G4String& particle,
425  const G4String& processName,
426  const G4String& material,
427  const G4String& reg)
428 {
429  return GetMeanFreePath(kinEnergy,FindParticle(particle),processName,
430  FindMaterial(material),FindRegion(reg));
431 }
432 
433 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
434 
435 inline G4double
437  const G4String& mat, G4double cut)
438 {
439  return
440  ComputeElectronicDEDX(kinEnergy,FindParticle(part),FindMaterial(mat),cut);
441 }
442 
443 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
444 
445 inline G4double
447  const G4String& part,
448  const G4String& mat,
449  G4double rangecut)
450 {
451  return ComputeDEDXForCutInRange(kinEnergy,FindParticle(part),
452  FindMaterial(mat), rangecut);
453 }
454 
455 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
456 
457 inline
459  const G4String& part,
460  const G4String& mat,
461  G4double cut)
462 {
463  return ComputeTotalDEDX(kinEnergy,FindParticle(part),FindMaterial(mat),cut);
464 }
465 
466 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
467 
468 inline
470  const G4String& particle,
471  const G4String& processName,
472  const G4String& material,
473  G4double cut)
474 {
475  return ComputeDEDX(kinEnergy,FindParticle(particle),processName,
476  FindMaterial(material),cut);
477 }
478 
479 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
480 
481 inline
483  const G4String& particle,
484  const G4String& material)
485 {
486  return ComputeNuclearDEDX(kinEnergy,FindParticle(particle),
487  FindMaterial(material));
488 }
489 
490 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
491 
492 inline
494  G4double kinEnergy,
495  const G4String& particle,
496  const G4String& processName,
497  const G4String& material,
498  G4double cut)
499 {
500  return ComputeCrossSectionPerVolume(kinEnergy,FindParticle(particle),
501  processName,
502  FindMaterial(material),cut);
503 }
504 
505 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
506 
507 inline
509  const G4String& particle,
510  const G4String& processName,
511  const G4Element* elm,
512  G4double cut)
513 {
514  return ComputeCrossSectionPerAtom(kinEnergy,FindParticle(particle),
515  processName,
516  elm->GetZ(),elm->GetN(),cut);
517 }
518 
519 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
520 
522  G4double kinEnergy, const G4String& part,
523  const G4String& processName, const G4Element* elm,
524  G4int shellIdx, G4double cut)
525 {
526  return ComputeCrossSectionPerShell(kinEnergy, FindParticle(part),
527  processName, G4lrint(elm->GetZ()),
528  shellIdx, cut);
529 }
530 
531 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
532 
533 inline
535  G4double range,
536  const G4String& particle,
537  const G4String& material)
538 {
539  return ComputeEnergyCutFromRangeCut(range,FindParticle(particle),
540  FindMaterial(material));
541 }
542 
543 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
544 
545 inline
547  const G4String& particle,
548  const G4String& processName,
549  const G4String& material,
550  G4double cut)
551 {
552  return ComputeMeanFreePath(kinEnergy,FindParticle(particle),processName,
553  FindMaterial(material),cut);
554 }
555 
556 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
557 
558 #endif
G4VEmModel * loweModel
G4VEmProcess * FindDiscreteProcess(const G4ParticleDefinition *, const G4String &processName)
const G4PhysicsTable * currentLambda
void FindLambdaTable(const G4ParticleDefinition *, const G4String &processName, G4double kinEnergy)
const G4ParticleDefinition * theGenericIon
const G4Material * FindMaterial(const G4String &)
const G4Material * cutMaterial
G4double GetN() const
Definition: G4Element.hh:134
std::vector< const G4MaterialCutsCouple * > localCouples
const G4ParticleDefinition * currentParticle
void PrintRangeTable(const G4ParticleDefinition *)
G4double GetZ() const
Definition: G4Element.hh:131
G4bool UpdateCouple(const G4Material *, G4double cut)
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)
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)
int G4int
Definition: G4Types.hh:78
G4DynamicParticle dynParticle
G4double ComputeShellIonisationCrossSectionPerAtom(const G4String &part, G4int Z, G4AtomicShellEnumerator shell, G4double kinEnergy, const G4Material *mat=0)
G4double GetDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=0)
static const G4double reg
static const double s
Definition: G4SIunits.hh:168
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)
G4VEnergyLossProcess * FindEnergyLossProcess(const G4ParticleDefinition *)
double A(double temperature)
G4VEnergyLossProcess * currentProcess
const G4MaterialCutsCouple * FindCouple(const G4Material *, const G4Region *r=0)
G4double ComputeDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=DBL_MAX)
G4String currentName
const G4ParticleDefinition * baseParticle
bool G4bool
Definition: G4Types.hh:79
G4bool FindEmModel(const G4ParticleDefinition *, const G4String &processName, G4double kinEnergy)
G4double ComputeNuclearDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *)
G4double GetRangeFromRestricteDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=0)
G4DataVector localCuts
G4EmCalculator & operator=(const G4EmCalculator &right)
void CheckMaterial(G4int Z)
G4bool ActiveForParticle(const G4ParticleDefinition *part, G4VProcess *proc)
G4VMultipleScattering * FindMscProcess(const G4ParticleDefinition *, const G4String &processName)
std::vector< const G4Material * > localMaterials
G4VEnergyLossProcess * FindEnLossProcess(const G4ParticleDefinition *, const G4String &processName)
const G4ParticleDefinition * FindParticle(const G4String &)
G4double chargeSquare
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)
const G4Region * FindRegion(const G4String &)
G4IonTable * ionTable
G4String currentMaterialName
const G4ParticleDefinition * FindIon(G4int Z, G4int A)
G4String currentProcessName
G4String currentParticleName
G4double ComputeCrossSectionPerVolume(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=0.0)
G4double ComputeGammaAttenuationLength(G4double kinEnergy, const G4Material *)
G4double ComputeCrossSectionPerAtom(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, G4double Z, G4double A, G4double cut=0.0)
G4VProcess * FindProcess(const G4ParticleDefinition *part, const G4String &processName)
void SetVerbose(G4int val)
void PrintDEDXTable(const G4ParticleDefinition *)
const G4ParticleDefinition * lambdaParticle
double G4double
Definition: G4Types.hh:76
G4double ComputeEnergyCutFromRangeCut(G4double range, const G4ParticleDefinition *, const G4Material *)
const G4Material * currentMaterial
G4VEmModel * currentModel
G4bool UpdateParticle(const G4ParticleDefinition *, G4double kinEnergy)
#define DBL_MAX
Definition: templates.hh:83
G4LossTableManager * manager
G4double ComputeDEDXForCutInRange(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *mat, G4double rangecut=DBL_MAX)
G4AtomicShellEnumerator
G4NistManager * nist