Geant4  10.01.p02
G4VEmProcess.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: G4VEmProcess.hh 85011 2014-10-23 09:41:59Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class header file
31 //
32 //
33 // File name: G4VEmProcess
34 //
35 // Author: Vladimir Ivanchenko
36 //
37 // Creation date: 01.10.2003
38 //
39 // Modifications:
40 // 30-06-04 make destructor virtual (V.Ivanchenko)
41 // 09-08-04 optimise integral option (V.Ivanchenko)
42 // 11-08-04 add protected methods to access cuts (V.Ivanchenko)
43 // 09-09-04 Bug fix for the integral mode with 2 peaks (V.Ivanchneko)
44 // 16-09-04 Add flag for LambdaTable and method RecalculateLambda (VI)
45 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivantchenko)
46 // 08-04-05 Major optimisation of internal interfaces (V.Ivantchenko)
47 // 18-04-05 Use G4ParticleChangeForGamma (V.Ivantchenko)
48 // 09-05-05 Fix problem in logic when path boundary between materials (VI)
49 // 11-01-06 add A to parameters of ComputeCrossSectionPerAtom (VI)
50 // 01-02-06 put default value A=0. to keep compatibility with v5.2 (mma)
51 // 13-05-06 Add method to access model by index (V.Ivanchenko)
52 // 12-09-06 add SetModel() (mma)
53 // 25-09-07 More accurate handling zero xsect in
54 // PostStepGetPhysicalInteractionLength (V.Ivanchenko)
55 // 27-10-07 Virtual functions moved to source (V.Ivanchenko)
56 // 15-07-08 Reorder class members for further multi-thread development (VI)
57 // 17-02-10 Added pointer currentParticle (VI)
58 //
59 // Class Description:
60 //
61 // It is the unified Discrete process
62 
63 // -------------------------------------------------------------------
64 //
65 
66 #ifndef G4VEmProcess_h
67 #define G4VEmProcess_h 1
68 
69 #include <CLHEP/Units/SystemOfUnits.h>
70 
71 #include "G4VDiscreteProcess.hh"
72 #include "globals.hh"
73 #include "G4Material.hh"
74 #include "G4MaterialCutsCouple.hh"
75 #include "G4Track.hh"
76 #include "G4EmModelManager.hh"
77 #include "G4UnitsTable.hh"
78 #include "G4ParticleDefinition.hh"
80 #include "G4EmParameters.hh"
81 
82 class G4Step;
83 class G4VEmModel;
84 class G4DataVector;
85 class G4VParticleChange;
86 class G4PhysicsTable;
87 class G4PhysicsVector;
88 class G4EmBiasingManager;
89 class G4LossTableManager;
90 
91 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
92 
94 {
95 public:
96 
97  G4VEmProcess(const G4String& name,
99 
100  virtual ~G4VEmProcess();
101 
102  //------------------------------------------------------------------------
103  // Virtual methods to be implemented in concrete processes
104  //------------------------------------------------------------------------
105 
106  virtual G4bool IsApplicable(const G4ParticleDefinition& p) = 0;
107 
108  virtual void PrintInfo() = 0;
109 
110 protected:
111 
112  virtual void InitialiseProcess(const G4ParticleDefinition*) = 0;
113 
114  //------------------------------------------------------------------------
115  // Method with standard implementation; may be overwritten if needed
116  //------------------------------------------------------------------------
117 
119  const G4Material*);
120 
121  //------------------------------------------------------------------------
122  // Implementation of virtual methods common to all Discrete processes
123  //------------------------------------------------------------------------
124 
125 public:
126 
127  // Initialise for build of tables
129 
130  // Build physics table during initialisation
132 
133  // Called before tracking of each new G4Track
134  void StartTracking(G4Track*);
135 
136  // implementation of virtual method, specific for G4VEmProcess
138  const G4Track& track,
139  G4double previousStepSize,
141 
142  // implementation of virtual method, specific for G4VEmProcess
143  G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&);
144 
145  // Store PhysicsTable in a file.
146  // Return false in case of failure at I/O
148  const G4String& directory,
149  G4bool ascii = false);
150 
151  // Retrieve Physics from a file.
152  // (return true if the Physics Table can be build by using file)
153  // (return false if the process has no functionality or in case of failure)
154  // File name should is constructed as processName+particleName and the
155  // should be placed under the directory specifed by the argument.
157  const G4String& directory,
158  G4bool ascii);
159 
160  //------------------------------------------------------------------------
161  // Specific methods for Discrete EM post step simulation
162  //------------------------------------------------------------------------
163 
164  // It returns the cross section per volume for energy/ material
166  const G4MaterialCutsCouple* couple);
167 
168  // It returns the cross section of the process per atom
170  G4double Z, G4double A=0.,
171  G4double cut=0.0);
172 
173  G4double MeanFreePath(const G4Track& track);
174 
175  // It returns cross section per volume
176  inline G4double GetLambda(G4double& kinEnergy,
177  const G4MaterialCutsCouple* couple);
178 
179  //------------------------------------------------------------------------
180  // Specific methods to build and access Physics Tables
181  //------------------------------------------------------------------------
182 
183  // Binning for lambda table
184  void SetLambdaBinning(G4int nbins);
185 
186  // Min kinetic energy for tables
187  void SetMinKinEnergy(G4double e);
188 
189  // Min kinetic energy for high energy table
191 
192  // Max kinetic energy for tables
193  void SetMaxKinEnergy(G4double e);
194 
195  // Cross section table pointers
196  inline G4PhysicsTable* LambdaTable() const;
197  inline G4PhysicsTable* LambdaTablePrim() const;
198 
199  //------------------------------------------------------------------------
200  // Define and access particle type
201  //------------------------------------------------------------------------
202 
203  inline const G4ParticleDefinition* Particle() const;
204  inline const G4ParticleDefinition* SecondaryParticle() const;
205 
206  //------------------------------------------------------------------------
207  // Specific methods to set, access, modify models and basic parameters
208  //------------------------------------------------------------------------
209 
210 protected:
211  // Select model in run time
212  inline G4VEmModel* SelectModel(G4double& kinEnergy, size_t index);
213 
214 public:
215  // Select model by energy and region index
216  inline G4VEmModel* SelectModelForMaterial(G4double kinEnergy,
217  size_t& idxRegion) const;
218 
219  // Add model for region, smaller value of order defines which
220  // model will be selected for a given energy interval
221  void AddEmModel(G4int, G4VEmModel*, const G4Region* region = 0);
222 
223  // return the assigned model
224  G4VEmModel* EmModel(G4int index = 1) const;
225 
226  // Assign a model to a process
227  void SetEmModel(G4VEmModel*, G4int index = 1);
228 
229  // Define new energy range for the model identified by the name
230  void UpdateEmModel(const G4String&, G4double, G4double);
231 
232  // Access to models
233  G4VEmModel* GetModelByIndex(G4int idx = 0, G4bool ver = false) const;
234 
235  // access atom on which interaction happens
236  const G4Element* GetCurrentElement() const;
237 
238  // Biasing parameters
239  void SetCrossSectionBiasingFactor(G4double f, G4bool flag = true);
240  inline G4double CrossSectionBiasingFactor() const;
241 
242  // Activate forced interaction
243  void ActivateForcedInteraction(G4double length = 0.0,
244  const G4String& r = "",
245  G4bool flag = true);
246 
247  void ActivateSecondaryBiasing(const G4String& region, G4double factor,
248  G4double energyLimit);
249 
250  inline void SetIntegral(G4bool val);
251 
252  inline void SetBuildTableFlag(G4bool val);
253 
254  //------------------------------------------------------------------------
255  // Other generic methods
256  //------------------------------------------------------------------------
257 
258 protected:
259 
260  G4double GetMeanFreePath(const G4Track& track,
261  G4double previousStepSize,
262  G4ForceCondition* condition);
263 
265 
266  inline G4int LambdaBinning() const;
267 
268  inline G4double MinKinEnergy() const;
269 
270  inline G4double MaxKinEnergy() const;
271 
272  // Single scattering parameters
273  inline G4double PolarAngleLimit() const;
274 
275  inline G4bool IsIntegral() const;
276 
277  inline G4double RecalculateLambda(G4double kinEnergy,
278  const G4MaterialCutsCouple* couple);
279 
281 
282  inline void SetParticle(const G4ParticleDefinition* p);
283 
284  inline void SetSecondaryParticle(const G4ParticleDefinition* p);
285 
286  inline size_t CurrentMaterialCutsCoupleIndex() const;
287 
288  inline G4double GetGammaEnergyCut();
289 
291 
292  inline void SetStartFromNullFlag(G4bool val);
293 
294  inline void SetSplineFlag(G4bool val);
295 
296  inline const G4Element* GetTargetElement() const;
297 
298  inline const G4Isotope* GetTargetIsotope() const;
299 
300 private:
301 
302  void Clear();
303 
304  void BuildLambdaTable();
305 
307 
308  void FindLambdaMax();
309 
310  void PrintWarning(G4String tit, G4double val);
311 
312  inline void DefineMaterial(const G4MaterialCutsCouple* couple);
313 
314  inline void ComputeIntegralLambda(G4double kinEnergy);
315 
316  inline G4double GetLambdaFromTable(G4double kinEnergy);
317 
318  inline G4double GetLambdaFromTablePrim(G4double kinEnergy);
319 
320  inline G4double GetCurrentLambda(G4double kinEnergy);
321 
322  inline G4double ComputeCurrentLambda(G4double kinEnergy);
323 
324  // copy constructor and hide assignment operator
327 
328  // ======== Parameters of the class fixed at construction =========
329 
338 
340 
341  // ======== Parameters of the class fixed at initialisation =======
342 
343  std::vector<G4VEmModel*> emModels;
345 
346  // tables and vectors
349  std::vector<G4double> theEnergyOfCrossSectionMax;
350  std::vector<G4double> theCrossSectionMax;
351 
352  size_t idxLambda;
354 
355  const std::vector<G4double>* theCuts;
356  const std::vector<G4double>* theCutsGamma;
357  const std::vector<G4double>* theCutsElectron;
358  const std::vector<G4double>* theCutsPositron;
359  const std::vector<G4double>* theDensityFactor;
360  const std::vector<G4int>* theDensityIdx;
361 
363 
369 
378 
379  // ======== Cashed values - may be state dependent ================
380 
381 protected:
382 
384 
385 private:
386 
387  std::vector<G4DynamicParticle*> secParticles;
388 
390 
393 
394  // cash
400 
407 
413 };
414 
415 // ======== Run time inline methods ================
416 
418 {
419  return currentCoupleIndex;
420 }
421 
422 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
423 
425 {
426  return (*theCutsGamma)[currentCoupleIndex];
427 }
428 
429 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
430 
432 {
434 }
435 
436 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
437 
439 {
440  if(couple != currentCouple) {
441  currentCouple = couple;
442  currentMaterial = couple->GetMaterial();
444  currentCoupleIndex = couple->GetIndex();
445  basedCoupleIndex = (*theDensityIdx)[currentCoupleIndex];
446  fFactor = biasFactor*(*theDensityFactor)[currentCoupleIndex];
449  idxLambda = idxLambdaPrim = 0;
450  }
451 }
452 
453 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
454 
455 inline
457 {
458  if(1 < numberOfModels) {
459  currentModel = modelManager->SelectModel(kinEnergy, index);
460  }
462  return currentModel;
463 }
464 
465 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
466 
467 inline
469  size_t& idxRegion) const
470 {
471  return modelManager->SelectModel(kinEnergy, idxRegion);
472 }
473 
474 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
475 
477 {
478  return ((*theLambdaTable)[basedCoupleIndex])->Value(e, idxLambda);
479 }
480 
481 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
482 
484 {
485  return ((*theLambdaTablePrim)[basedCoupleIndex])->Value(e, idxLambdaPrim)/e;
486 }
487 
488 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
489 
491 {
494 }
495 
496 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
497 
499 {
500  G4double x;
501  if(e >= minKinEnergyPrim) { x = GetLambdaFromTablePrim(e); }
502  else if(theLambdaTable) { x = GetLambdaFromTable(e); }
503  else { x = ComputeCurrentLambda(e); }
504  return fFactor*x;
505 }
506 
507 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
508 
509 inline G4double
511  const G4MaterialCutsCouple* couple)
512 {
513  DefineMaterial(couple);
514  SelectModel(kinEnergy, currentCoupleIndex);
515  return GetCurrentLambda(kinEnergy);
516 }
517 
518 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
519 
520 inline G4double
522 {
523  DefineMaterial(couple);
525  return fFactor*ComputeCurrentLambda(e);
526 }
527 
528 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
529 
531 {
533  if (e <= mfpKinEnergy) {
535 
536  } else {
538  if(e1 > mfpKinEnergy) {
540  G4double preStepLambda1 = GetCurrentLambda(e1);
541  if(preStepLambda1 > preStepLambda) {
542  mfpKinEnergy = e1;
543  preStepLambda = preStepLambda1;
544  }
545  } else {
547  }
548  }
549 }
550 
551 // ======== Get/Set inline methods used at initialisation ================
552 
553 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
554 
556 {
557  return nLambdaBins;
558 }
559 
560 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
561 
563 {
564  return minKinEnergy;
565 }
566 
567 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
568 
570 {
571  return maxKinEnergy;
572 }
573 
574 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
575 
577 {
578  return theParameters->MscThetaLimit();
579 }
580 
581 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
582 
584 {
585  return biasFactor;
586 }
587 
588 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
589 
591 {
592  return theLambdaTable;
593 }
594 
595 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
596 
598 {
599  return theLambdaTablePrim;
600 }
601 
602 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
603 
605 {
606  return particle;
607 }
608 
609 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
610 
612 {
613  return secondaryParticle;
614 }
615 
616 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
617 
619 {
620  integral = val;
621 }
622 
623 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
624 
626 {
627  return integral;
628 }
629 
630 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
631 
633 {
634  buildLambdaTable = val;
635 }
636 
637 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
638 
640 {
641  return &fParticleChange;
642 }
643 
644 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
645 
647 {
648  particle = p;
649  currentParticle = p;
650 }
651 
652 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
653 
655 {
656  secondaryParticle = p;
657 }
658 
659 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
660 
662 {
663  startFromNull = val;
664 }
665 
666 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
667 
669 {
670  splineFlag = val;
671  actSpline = true;
672 }
673 
674 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
675 
677 {
679 }
680 
681 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
682 
684 {
686 }
687 
688 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
689 
690 #endif
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="", G4bool flag=true)
G4double condition(const G4ErrorSymMatrix &m)
G4PhysicsTable * LambdaTable() const
const std::vector< G4double > * theCutsGamma
std::vector< G4double > theCrossSectionMax
void DefineMaterial(const G4MaterialCutsCouple *couple)
virtual void PrintInfo()=0
G4double GetLambda(G4double &kinEnergy, const G4MaterialCutsCouple *couple)
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:258
G4PhysicsTable * LambdaTablePrim() const
G4VEmModel * SelectModel(G4double &kinEnergy, size_t index)
G4bool weightFlag
G4PhysicsTable * theLambdaTablePrim
G4double GetLambdaFromTable(G4double kinEnergy)
G4VEmProcess(const G4String &name, G4ProcessType type=fElectromagnetic)
Definition: G4VEmProcess.cc:90
G4VEmModel * currentModel
void SetBuildTableFlag(G4bool val)
G4VEmProcess & operator=(const G4VEmProcess &right)
G4double ComputeCrossSectionPerAtom(G4double kineticEnergy, G4double Z, G4double A=0., G4double cut=0.0)
G4VEmModel * EmModel(G4int index=1) const
G4String name
Definition: TRTMaterials.hh:40
G4double MscThetaLimit() const
G4ParticleChangeForGamma fParticleChange
G4double GetElectronEnergyCut()
G4EmParameters * theParameters
void SetSplineFlag(G4bool val)
const std::vector< G4double > * theCutsPositron
virtual void InitialiseProcess(const G4ParticleDefinition *)=0
size_t basedCoupleIndex
G4double MaxKinEnergy() const
virtual ~G4VEmProcess()
const G4ParticleDefinition * SecondaryParticle() const
const G4MaterialCutsCouple * currentCouple
const G4ParticleDefinition * Particle() const
G4double lambdaFactor
void BuildPhysicsTable(const G4ParticleDefinition &)
G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
const G4Material * baseMaterial
G4double fFactor
G4PhysicsTable * theLambdaTable
void SetStartFromNullFlag(G4bool val)
G4double minKinEnergyPrim
int G4int
Definition: G4Types.hh:78
void UpdateEmModel(const G4String &, G4double, G4double)
G4double GetCurrentLambda(G4double kinEnergy)
G4bool applyCuts
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
G4int LambdaBinning() const
void SetMinKinEnergyPrim(G4double e)
G4double mfpKinEnergy
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
virtual G4double MinPrimaryEnergy(const G4ParticleDefinition *, const G4Material *)
G4double MeanFreePath(const G4Track &track)
const G4ParticleDefinition * currentParticle
void SetLambdaBinning(G4int nbins)
G4double minKinEnergy
void SetEmModel(G4VEmModel *, G4int index=1)
G4bool splineFlag
size_t CurrentMaterialCutsCoupleIndex() const
G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii)
G4double PolarAngleLimit() const
G4bool startFromNull
bool G4bool
Definition: G4Types.hh:79
const std::vector< G4int > * theDensityIdx
G4double GetGammaEnergyCut()
void SetParticle(const G4ParticleDefinition *p)
void PrintInfoProcess(const G4ParticleDefinition &)
std::vector< G4DynamicParticle * > secParticles
G4ParticleChangeForGamma * GetParticleChange()
G4double GetMeanFreePath(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
G4bool actSpline
G4LossTableManager * lManager
const G4Element * GetTargetElement() const
const G4ParticleDefinition * secondaryParticle
Definition: G4Step.hh:76
std::vector< G4double > theEnergyOfCrossSectionMax
G4bool actBinning
void SetIntegral(G4bool val)
const G4ParticleDefinition * particle
void PreparePhysicsTable(const G4ParticleDefinition &)
const std::vector< G4double > * theCuts
static const G4double A[nN]
static const G4double e1
std::vector< G4VEmModel * > emModels
G4double ComputeCurrentLambda(G4double kinEnergy)
G4PhysicsVector * LambdaPhysicsVector(const G4MaterialCutsCouple *)
G4VEmModel * SelectModelForMaterial(G4double kinEnergy, size_t &idxRegion) const
const G4Isotope * GetTargetIsotope() const
const G4Material * currentMaterial
void SetMaxKinEnergy(G4double e)
size_t idxLambdaPrim
void FindLambdaMax()
static const G4double factor
const G4ParticleDefinition * theElectron
G4double CrossSectionBiasingFactor() const
void AddEmModel(G4int, G4VEmModel *, const G4Region *region=0)
G4bool actMaxKinEnergy
size_t idxLambda
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:426
void SetSecondaryParticle(const G4ParticleDefinition *p)
G4EmModelManager * modelManager
G4bool IsIntegral() const
const G4Element * GetCurrentElement() const
G4double preStepLambda
void ComputeIntegralLambda(G4double kinEnergy)
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:233
G4EmBiasingManager * biasManager
G4VEmModel * SelectModel(G4double &energy, size_t &index)
G4double GetLambdaFromTablePrim(G4double kinEnergy)
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
G4int numberOfModels
const std::vector< G4double > * theDensityFactor
G4double biasFactor
size_t currentCoupleIndex
G4bool actMinKinEnergy
void SetMinKinEnergy(G4double e)
double G4double
Definition: G4Types.hh:76
const std::vector< G4double > * theCutsElectron
G4ForceCondition
virtual G4bool IsApplicable(const G4ParticleDefinition &p)=0
void StartTracking(G4Track *)
G4double MinKinEnergy() const
#define DBL_MAX
Definition: templates.hh:83
G4double CrossSectionPerVolume(G4double kineticEnergy, const G4MaterialCutsCouple *couple)
void PrintWarning(G4String tit, G4double val)
G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false)
G4double maxKinEnergy
G4double preStepKinEnergy
G4bool buildLambdaTable
G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
const G4ParticleDefinition * theGamma
const G4ParticleDefinition * thePositron
const G4Material * GetMaterial() const
G4int mainSecondaries
void BuildLambdaTable()
G4ProcessType
const G4Isotope * GetCurrentIsotope() const
Definition: G4VEmModel.hh:456
const G4Element * GetCurrentElement() const
Definition: G4VEmModel.hh:449
G4double RecalculateLambda(G4double kinEnergy, const G4MaterialCutsCouple *couple)