Geant4  10.00.p03
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 74376 2013-10-04 08:25:47Z 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 
81 class G4Step;
82 class G4VEmModel;
83 class G4DataVector;
84 class G4VParticleChange;
85 class G4PhysicsTable;
86 class G4PhysicsVector;
87 class G4EmBiasingManager;
88 class G4LossTableManager;
89 
90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
91 
93 {
94 public:
95 
96  G4VEmProcess(const G4String& name,
98 
99  virtual ~G4VEmProcess();
100 
101  //------------------------------------------------------------------------
102  // Virtual methods to be implemented in concrete processes
103  //------------------------------------------------------------------------
104 
105  virtual G4bool IsApplicable(const G4ParticleDefinition& p) = 0;
106 
107  virtual void PrintInfo() = 0;
108 
109 protected:
110 
111  virtual void InitialiseProcess(const G4ParticleDefinition*) = 0;
112 
113  //------------------------------------------------------------------------
114  // Method with standard implementation; may be overwritten if needed
115  //------------------------------------------------------------------------
116 
118  const G4Material*);
119 
120  //------------------------------------------------------------------------
121  // Implementation of virtual methods common to all Discrete processes
122  //------------------------------------------------------------------------
123 
124 public:
125 
126  // Initialise for build of tables
128 
129  // Build physics table during initialisation
131 
132  // Called before tracking of each new G4Track
133  void StartTracking(G4Track*);
134 
135  // implementation of virtual method, specific for G4VEmProcess
137  const G4Track& track,
138  G4double previousStepSize,
140 
141  // implementation of virtual method, specific for G4VEmProcess
142  G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&);
143 
144  // Store PhysicsTable in a file.
145  // Return false in case of failure at I/O
147  const G4String& directory,
148  G4bool ascii = false);
149 
150  // Retrieve Physics from a file.
151  // (return true if the Physics Table can be build by using file)
152  // (return false if the process has no functionality or in case of failure)
153  // File name should is constructed as processName+particleName and the
154  // should be placed under the directory specifed by the argument.
156  const G4String& directory,
157  G4bool ascii);
158 
159  //------------------------------------------------------------------------
160  // Specific methods for Discrete EM post step simulation
161  //------------------------------------------------------------------------
162 
163  // It returns the cross section per volume for energy/ material
165  const G4MaterialCutsCouple* couple);
166 
167  // It returns the cross section of the process per atom
169  G4double Z, G4double A=0.,
170  G4double cut=0.0);
171 
172  G4double MeanFreePath(const G4Track& track);
173 
174  // It returns cross section per volume
175  inline G4double GetLambda(G4double& kinEnergy,
176  const G4MaterialCutsCouple* couple);
177 
178  //------------------------------------------------------------------------
179  // Specific methods to build and access Physics Tables
180  //------------------------------------------------------------------------
181 
182  // Binning for lambda table
183  inline void SetLambdaBinning(G4int nbins);
184  inline G4int LambdaBinning() const;
185 
186  // Min kinetic energy for tables
187  void SetMinKinEnergy(G4double e);
188  inline G4double MinKinEnergy() const;
189 
190  // Max kinetic energy for tables
191  void SetMaxKinEnergy(G4double e);
192  inline G4double MaxKinEnergy() const;
193 
194  // Min kinetic energy for high energy table
195  inline void SetMinKinEnergyPrim(G4double e);
196 
197  // Cross section table pointers
198  inline G4PhysicsTable* LambdaTable() const;
199  inline G4PhysicsTable* LambdaTablePrim() const;
200 
201  //------------------------------------------------------------------------
202  // Define and access particle type
203  //------------------------------------------------------------------------
204 
205  inline const G4ParticleDefinition* Particle() const;
206  inline const G4ParticleDefinition* SecondaryParticle() const;
207 
208  //------------------------------------------------------------------------
209  // Specific methods to set, access, modify models and basic parameters
210  //------------------------------------------------------------------------
211 
212 protected:
213  // Select model in run time
214  inline G4VEmModel* SelectModel(G4double& kinEnergy, size_t index);
215 
216 public:
217  // Select model by energy and region index
218  inline G4VEmModel* SelectModelForMaterial(G4double kinEnergy,
219  size_t& idxRegion) const;
220 
221  // Add model for region, smaller value of order defines which
222  // model will be selected for a given energy interval
223  void AddEmModel(G4int, G4VEmModel*, const G4Region* region = 0);
224 
225  // return the assigned model
226  G4VEmModel* EmModel(G4int index = 1) const;
227 
228  // Assign a model to a process
229  void SetEmModel(G4VEmModel*, G4int index = 1);
230 
231  // Define new energy range for the model identified by the name
232  void UpdateEmModel(const G4String&, G4double, G4double);
233 
234  // Access to models
235  G4VEmModel* GetModelByIndex(G4int idx = 0, G4bool ver = false) const;
236 
237  // access atom on which interaction happens
238  const G4Element* GetCurrentElement() const;
239 
240  // Biasing parameters
241  void SetCrossSectionBiasingFactor(G4double f, G4bool flag = true);
242  inline G4double CrossSectionBiasingFactor() const;
243 
244  // Activate forced interaction
245  void ActivateForcedInteraction(G4double length = 0.0,
246  const G4String& r = "",
247  G4bool flag = true);
248 
249  void ActivateSecondaryBiasing(const G4String& region, G4double factor,
250  G4double energyLimit);
251 
252 
253  // Single scattering parameters
254  inline void SetPolarAngleLimit(G4double a);
255  inline G4double PolarAngleLimit() const;
256 
257  inline void SetLambdaFactor(G4double val);
258 
259  inline void SetIntegral(G4bool val);
260  inline G4bool IsIntegral() const;
261 
262  inline void SetApplyCuts(G4bool val);
263 
264  inline void SetBuildTableFlag(G4bool val);
265 
266  //------------------------------------------------------------------------
267  // Other generic methods
268  //------------------------------------------------------------------------
269 
270 protected:
271 
272  G4double GetMeanFreePath(const G4Track& track,
273  G4double previousStepSize,
274  G4ForceCondition* condition);
275 
277 
278  inline G4double RecalculateLambda(G4double kinEnergy,
279  const G4MaterialCutsCouple* couple);
280 
282 
283  inline void SetParticle(const G4ParticleDefinition* p);
284 
285  inline void SetSecondaryParticle(const G4ParticleDefinition* p);
286 
287  inline size_t CurrentMaterialCutsCoupleIndex() const;
288 
289  inline G4double GetGammaEnergyCut();
290 
292 
293  inline void SetStartFromNullFlag(G4bool val);
294 
295  inline void SetSplineFlag(G4bool val);
296 
297 private:
298 
299  void Clear();
300 
301  void BuildLambdaTable();
302 
304 
305  void FindLambdaMax();
306 
307  inline void DefineMaterial(const G4MaterialCutsCouple* couple);
308 
309  inline void ComputeIntegralLambda(G4double kinEnergy);
310 
311  inline G4double GetLambdaFromTable(G4double kinEnergy);
312 
313  inline G4double GetLambdaFromTablePrim(G4double kinEnergy);
314 
315  inline G4double GetCurrentLambda(G4double kinEnergy);
316 
317  inline G4double ComputeCurrentLambda(G4double kinEnergy);
318 
319  // copy constructor and hide assignment operator
322 
323  // ======== Parameters of the class fixed at construction =========
324 
332 
334 
335  // ======== Parameters of the class fixed at initialisation =======
336 
337  std::vector<G4VEmModel*> emModels;
339 
340  // tables and vectors
343  std::vector<G4double> theEnergyOfCrossSectionMax;
344  std::vector<G4double> theCrossSectionMax;
345 
346  size_t idxLambda;
348 
349  const std::vector<G4double>* theCuts;
350  const std::vector<G4double>* theCutsGamma;
351  const std::vector<G4double>* theCutsElectron;
352  const std::vector<G4double>* theCutsPositron;
353  const std::vector<G4double>* theDensityFactor;
354  const std::vector<G4int>* theDensityIdx;
355 
357 
364 
369 
370  // ======== Cashed values - may be state dependent ================
371 
372 protected:
373 
375 
376 private:
377 
378  std::vector<G4DynamicParticle*> secParticles;
379 
381 
384 
385  // cash
391 
398 
404 };
405 
406 // ======== Run time inline methods ================
407 
409 {
410  return currentCoupleIndex;
411 }
412 
413 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
414 
416 {
417  return (*theCutsGamma)[currentCoupleIndex];
418 }
419 
420 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
421 
423 {
425 }
426 
427 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
428 
430 {
431  if(couple != currentCouple) {
432  currentCouple = couple;
433  currentMaterial = couple->GetMaterial();
435  currentCoupleIndex = couple->GetIndex();
436  basedCoupleIndex = (*theDensityIdx)[currentCoupleIndex];
437  fFactor = biasFactor*(*theDensityFactor)[currentCoupleIndex];
440  idxLambda = idxLambdaPrim = 0;
441  }
442 }
443 
444 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
445 
446 inline
448 {
449  if(1 < numberOfModels) {
450  currentModel = modelManager->SelectModel(kinEnergy, index);
451  }
453  return currentModel;
454 }
455 
456 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
457 
458 inline
460  size_t& idxRegion) const
461 {
462  return modelManager->SelectModel(kinEnergy, idxRegion);
463 }
464 
465 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
466 
468 {
469  return ((*theLambdaTable)[basedCoupleIndex])->Value(e, idxLambda);
470 }
471 
472 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
473 
475 {
476  return ((*theLambdaTablePrim)[basedCoupleIndex])->Value(e, idxLambdaPrim)/e;
477 }
478 
479 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
480 
482 {
485 }
486 
487 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
488 
490 {
491  G4double x;
492  if(e >= minKinEnergyPrim) { x = GetLambdaFromTablePrim(e); }
493  else if(theLambdaTable) { x = GetLambdaFromTable(e); }
494  else { x = ComputeCurrentLambda(e); }
495  return fFactor*x;
496 }
497 
498 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
499 
500 inline G4double
502  const G4MaterialCutsCouple* couple)
503 {
504  DefineMaterial(couple);
505  SelectModel(kinEnergy, currentCoupleIndex);
506  return GetCurrentLambda(kinEnergy);
507 }
508 
509 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
510 
511 inline G4double
513 {
514  DefineMaterial(couple);
516  return fFactor*ComputeCurrentLambda(e);
517 }
518 
519 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
520 
522 {
524  if (e <= mfpKinEnergy) {
526 
527  } else {
529  if(e1 > mfpKinEnergy) {
531  G4double preStepLambda1 = GetCurrentLambda(e1);
532  if(preStepLambda1 > preStepLambda) {
533  mfpKinEnergy = e1;
534  preStepLambda = preStepLambda1;
535  }
536  } else {
538  }
539  }
540 }
541 
542 // ======== Get/Set inline methods used at initialisation ================
543 
545 {
546  nLambdaBins = nbins;
547 }
548 
549 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
550 
552 {
553  return nLambdaBins;
554 }
555 
556 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
557 
559 {
560  return minKinEnergy;
561 }
562 
563 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
564 
566 {
567  return maxKinEnergy;
568 }
569 
570 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
571 
573 {
574  minKinEnergyPrim = e;
575 }
576 
577 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
578 
580 {
581  if(val < 0.0) { polarAngleLimit = 0.0; }
582  else if(val > CLHEP::pi) { polarAngleLimit = CLHEP::pi; }
583  else { polarAngleLimit = val; }
584 }
585 
586 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
587 
589 {
590  return polarAngleLimit;
591 }
592 
593 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
594 
596 {
597  return biasFactor;
598 }
599 
600 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
601 
603 {
604  return theLambdaTable;
605 }
606 
607 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
608 
610 {
611  return theLambdaTablePrim;
612 }
613 
614 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
615 
617 {
618  return particle;
619 }
620 
621 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
622 
624 {
625  return secondaryParticle;
626 }
627 
628 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
629 
631 {
632  if(val > 0.0 && val <= 1.0) { lambdaFactor = val; }
633 }
634 
635 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
636 
638 {
639  if(particle && particle != theGamma) { integral = val; }
640 }
641 
642 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
643 
645 {
646  return integral;
647 }
648 
649 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
650 
652 {
653  applyCuts = val;
654 }
655 
656 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
657 
659 {
660  buildLambdaTable = val;
661 }
662 
663 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
664 
666 {
667  return &fParticleChange;
668 }
669 
670 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
671 
673 {
674  particle = p;
675  currentParticle = p;
676 }
677 
678 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
679 
681 {
682  secondaryParticle = p;
683 }
684 
685 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
686 
688 {
689  startFromNull = val;
690 }
691 
692 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
693 
695 {
696  splineFlag = val;
697 }
698 
699 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
700 
701 #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:245
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
const G4double pi
G4ParticleChangeForGamma fParticleChange
G4double GetElectronEnergyCut()
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 polarAngleLimit
G4double a
Definition: TRTMaterials.hh:39
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 SetApplyCuts(G4bool val)
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)
G4LossTableManager * lManager
const G4ParticleDefinition * secondaryParticle
Definition: G4Step.hh:76
std::vector< G4double > theEnergyOfCrossSectionMax
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 G4Material * currentMaterial
void SetMaxKinEnergy(G4double e)
size_t idxLambdaPrim
void FindLambdaMax()
const G4ParticleDefinition * theElectron
G4double CrossSectionBiasingFactor() const
void AddEmModel(G4int, G4VEmModel *, const G4Region *region=0)
size_t idxLambda
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:419
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:231
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
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)
G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false)
G4double maxKinEnergy
G4double preStepKinEnergy
G4bool buildLambdaTable
void SetPolarAngleLimit(G4double a)
G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
const G4ParticleDefinition * theGamma
const G4ParticleDefinition * thePositron
void SetLambdaFactor(G4double val)
const G4Material * GetMaterial() const
G4int mainSecondaries
void BuildLambdaTable()
G4ProcessType
G4double RecalculateLambda(G4double kinEnergy, const G4MaterialCutsCouple *couple)