Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4VEmProcess.hh
이 파일의 문서화 페이지로 가기
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 109178 2018-04-03 07:13:58Z 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 
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) override = 0;
107 
108  // obsolete
109  virtual void PrintInfo() {};
110 
111  virtual void ProcessDescription(std::ostream& outFile) const override;
112 
113 protected:
114 
115  virtual void StreamProcessInfo(std::ostream&, G4String) const {};
116 
117  virtual void InitialiseProcess(const G4ParticleDefinition*) = 0;
118 
119  //------------------------------------------------------------------------
120  // Method with standard implementation; may be overwritten if needed
121  //------------------------------------------------------------------------
122 
124  const G4Material*);
125 
126  //------------------------------------------------------------------------
127  // Implementation of virtual methods common to all Discrete processes
128  //------------------------------------------------------------------------
129 
130 public:
131 
132  // Initialise for build of tables
133  virtual void PreparePhysicsTable(const G4ParticleDefinition&) override;
134 
135  // Build physics table during initialisation
136  virtual void BuildPhysicsTable(const G4ParticleDefinition&) override;
137 
138  // Called before tracking of each new G4Track
139  virtual void StartTracking(G4Track*) override;
140 
141  // implementation of virtual method, specific for G4VEmProcess
143  const G4Track& track,
144  G4double previousStepSize,
145  G4ForceCondition* condition) override;
146 
147  // implementation of virtual method, specific for G4VEmProcess
148  virtual G4VParticleChange* PostStepDoIt(const G4Track&,
149  const G4Step&) override;
150 
151  // Store PhysicsTable in a file.
152  // Return false in case of failure at I/O
154  const G4String& directory,
155  G4bool ascii = false) override;
156 
157  // Retrieve Physics from a file.
158  // (return true if the Physics Table can be build by using file)
159  // (return false if the process has no functionality or in case of failure)
160  // File name should is constructed as processName+particleName and the
161  // should be placed under the directory specifed by the argument.
163  const G4String& directory,
164  G4bool ascii) override;
165 
166  //------------------------------------------------------------------------
167  // Specific methods for Discrete EM post step simulation
168  //------------------------------------------------------------------------
169 
170  // It returns the cross section per volume for energy/ material
172  const G4MaterialCutsCouple* couple);
173 
174  // It returns the cross section of the process per atom
176  G4double Z, G4double A=0.,
177  G4double cut=0.0);
178 
180 
181  // It returns cross section per volume
182  inline G4double GetLambda(G4double& kinEnergy,
183  const G4MaterialCutsCouple* couple);
184 
185  //------------------------------------------------------------------------
186  // Specific methods to build and access Physics Tables
187  //------------------------------------------------------------------------
188 
189  // Binning for lambda table
190  void SetLambdaBinning(G4int nbins);
191 
192  // Min kinetic energy for tables
193  void SetMinKinEnergy(G4double e);
194 
195  // Min kinetic energy for high energy table
197 
198  // Max kinetic energy for tables
199  void SetMaxKinEnergy(G4double e);
200 
201  // Cross section table pointers
202  inline G4PhysicsTable* LambdaTable() const;
203  inline G4PhysicsTable* LambdaTablePrim() const;
204 
205  //------------------------------------------------------------------------
206  // Define and access particle type
207  //------------------------------------------------------------------------
208 
209  inline const G4ParticleDefinition* Particle() const;
210  inline const G4ParticleDefinition* SecondaryParticle() const;
211 
212  //------------------------------------------------------------------------
213  // Specific methods to set, access, modify models and basic parameters
214  //------------------------------------------------------------------------
215 
216 protected:
217  // Select model in run time
218  inline G4VEmModel* SelectModel(G4double& kinEnergy, size_t index);
219 
220 public:
221  // Select model by energy and region index
222  inline G4VEmModel* SelectModelForMaterial(G4double kinEnergy,
223  size_t& idxRegion) const;
224 
225  // Add model for region, smaller value of order defines which
226  // model will be selected for a given energy interval
227  void AddEmModel(G4int, G4VEmModel*, const G4Region* region = nullptr);
228 
229  // Assign a model to a process local list, to enable the list in run time
230  // the derived process should execute AddEmModel(..) for all such models
231  void SetEmModel(G4VEmModel*, G4int index = 0);
232 
233  // return a model from the local list
234  G4VEmModel* EmModel(size_t index = 0) const;
235 
236  // Define new energy range for the model identified by the name
237  void UpdateEmModel(const G4String&, G4double, G4double);
238 
239  // Access to models
240  G4int GetNumberOfModels() const;
241  G4int GetNumberOfRegionModels(size_t couple_index) const;
242  G4VEmModel* GetRegionModel(G4int idx, size_t couple_index) const;
243  G4VEmModel* GetModelByIndex(G4int idx = 0, G4bool ver = false) const;
244 
245  // Access to active model
246  inline const G4VEmModel* GetCurrentModel() const;
247 
248  // Access to the current G4Element
249  const G4Element* GetCurrentElement() const;
250 
251  // Biasing parameters
252  void SetCrossSectionBiasingFactor(G4double f, G4bool flag = true);
253  inline G4double CrossSectionBiasingFactor() const;
254 
255  // Activate forced interaction
256  void ActivateForcedInteraction(G4double length = 0.0,
257  const G4String& r = "",
258  G4bool flag = true);
259 
260  void ActivateSecondaryBiasing(const G4String& region, G4double factor,
261  G4double energyLimit);
262 
263  inline void SetIntegral(G4bool val);
264 
265  inline void SetBuildTableFlag(G4bool val);
266 
267  //------------------------------------------------------------------------
268  // Other generic methods
269  //------------------------------------------------------------------------
270 
271 protected:
272 
273  virtual G4double GetMeanFreePath(const G4Track& track,
274  G4double previousStepSize,
275  G4ForceCondition* condition) override;
276 
278 
279  inline G4int LambdaBinning() const;
280 
281  inline G4double MinKinEnergy() const;
282 
283  inline G4double MaxKinEnergy() const;
284 
285  // Single scattering parameters
286  inline G4double PolarAngleLimit() const;
287 
288  inline G4bool IsIntegral() const;
289 
290  inline G4double RecalculateLambda(G4double kinEnergy,
291  const G4MaterialCutsCouple* couple);
292 
294 
295  inline void SetParticle(const G4ParticleDefinition* p);
296 
297  inline void SetSecondaryParticle(const G4ParticleDefinition* p);
298 
299  inline size_t CurrentMaterialCutsCoupleIndex() const;
300 
301  inline const G4MaterialCutsCouple* MaterialCutsCouple() const;
302 
303  inline G4bool ApplyCuts() const;
304 
305  inline G4double GetGammaEnergyCut();
306 
308 
309  inline void SetStartFromNullFlag(G4bool val);
310 
311  inline void SetSplineFlag(G4bool val);
312 
313  inline const G4Element* GetTargetElement() const;
314 
315  inline const G4Isotope* GetTargetIsotope() const;
316 
317 private:
318 
319  void Clear();
320 
321  void BuildLambdaTable();
322 
323  void StreamInfo(std::ostream& outFile, const G4ParticleDefinition&,
324  G4String endOfLine=G4String("\n")) const;
325 
326  void FindLambdaMax();
327 
328  void PrintWarning(G4String tit, G4double val);
329 
330  inline void DefineMaterial(const G4MaterialCutsCouple* couple);
331 
332  inline void ComputeIntegralLambda(G4double kinEnergy);
333 
334  inline G4double GetLambdaFromTable(G4double kinEnergy);
335 
336  inline G4double GetLambdaFromTablePrim(G4double kinEnergy);
337 
338  inline G4double GetCurrentLambda(G4double kinEnergy);
339 
340  inline G4double ComputeCurrentLambda(G4double kinEnergy);
341 
342  // hide copy constructor and assignment operator
343  G4VEmProcess(G4VEmProcess &) = delete;
344  G4VEmProcess & operator=(const G4VEmProcess &right) = delete;
345 
346  // ======== Parameters of the class fixed at construction =========
347 
353 
355 
356  // ======== Parameters of the class fixed at initialisation =======
357 
358  std::vector<G4VEmModel*> emModels;
360 
361  // tables and vectors
364  std::vector<G4double> theEnergyOfCrossSectionMax;
365  std::vector<G4double> theCrossSectionMax;
366 
367  size_t idxLambda;
369 
370  const std::vector<G4double>* theCuts;
371  const std::vector<G4double>* theCutsGamma;
372  const std::vector<G4double>* theCutsElectron;
373  const std::vector<G4double>* theCutsPositron;
374  const std::vector<G4double>* theDensityFactor;
375  const std::vector<G4int>* theDensityIdx;
376 
378 
385 
395 
396  // ======== Cashed values - may be state dependent ================
397 
398 protected:
399 
404  std::vector<G4DynamicParticle*> secParticles;
406 
412 
413 private:
414 
416 
419 
420  // cache
425 
432 };
433 
434 // ======== Run time inline methods ================
435 
437 {
438  return applyCuts;
439 }
440 
441 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
442 
444 {
445  return currentCoupleIndex;
446 }
447 
448 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
449 
451 {
452  return currentCouple;
453 }
454 
455 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
456 
458 {
459  return (*theCutsGamma)[currentCoupleIndex];
460 }
461 
462 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
463 
465 {
467 }
468 
469 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
470 
472 {
473  if(couple != currentCouple) {
474  currentCouple = couple;
475  currentMaterial = couple->GetMaterial();
477  currentCoupleIndex = couple->GetIndex();
478  basedCoupleIndex = (*theDensityIdx)[currentCoupleIndex];
479  fFactor = biasFactor*(*theDensityFactor)[currentCoupleIndex];
482  idxLambda = idxLambdaPrim = 0;
483  }
484 }
485 
486 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
487 
488 inline
490 {
491  if(1 < numberOfModels) {
492  currentModel = modelManager->SelectModel(kinEnergy, index);
493  }
495  return currentModel;
496 }
497 
498 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
499 
500 inline
502  size_t& idxRegion) const
503 {
504  return modelManager->SelectModel(kinEnergy, idxRegion);
505 }
506 
507 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
508 
510 {
511  return ((*theLambdaTable)[basedCoupleIndex])->Value(e, idxLambda);
512 }
513 
514 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
515 
517 {
518  return ((*theLambdaTablePrim)[basedCoupleIndex])->Value(e, idxLambdaPrim)/e;
519 }
520 
521 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
522 
524 {
527 }
528 
529 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
530 
532 {
533  G4double x;
534  if(e >= minKinEnergyPrim) { x = GetLambdaFromTablePrim(e); }
535  else if(theLambdaTable) { x = GetLambdaFromTable(e); }
536  else { x = ComputeCurrentLambda(e); }
537  return fFactor*x;
538 }
539 
540 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
541 
542 inline G4double
544  const G4MaterialCutsCouple* couple)
545 {
546  DefineMaterial(couple);
547  SelectModel(kinEnergy, currentCoupleIndex);
548  return GetCurrentLambda(kinEnergy);
549 }
550 
551 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
552 
553 inline G4double
555 {
556  DefineMaterial(couple);
558  return fFactor*ComputeCurrentLambda(e);
559 }
560 
561 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
562 
564 {
566  if (e <= mfpKinEnergy) {
568 
569  } else {
570  G4double e1 = e*lambdaFactor;
571  if(e1 > mfpKinEnergy) {
573  G4double preStepLambda1 = GetCurrentLambda(e1);
574  if(preStepLambda1 > preStepLambda) {
575  mfpKinEnergy = e1;
576  preStepLambda = preStepLambda1;
577  }
578  } else {
580  }
581  }
582 }
583 
584 // ======== Get/Set inline methods used at initialisation ================
585 
586 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
587 
589 {
590  return nLambdaBins;
591 }
592 
593 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
594 
596 {
597  return minKinEnergy;
598 }
599 
600 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
601 
603 {
604  return maxKinEnergy;
605 }
606 
607 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
608 
610 {
611  return theParameters->MscThetaLimit();
612 }
613 
614 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
615 
617 {
618  return biasFactor;
619 }
620 
621 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
622 
624 {
625  return theLambdaTable;
626 }
627 
628 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
629 
631 {
632  return theLambdaTablePrim;
633 }
634 
635 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
636 
638 {
639  return particle;
640 }
641 
642 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
643 
645 {
646  return secondaryParticle;
647 }
648 
649 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
650 
652 {
653  integral = val;
654 }
655 
656 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
657 
659 {
660  return integral;
661 }
662 
663 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
664 
666 {
667  buildLambdaTable = val;
668 }
669 
670 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
671 
673 {
674  return &fParticleChange;
675 }
676 
677 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
678 
680 {
681  particle = p;
682  currentParticle = p;
683 }
684 
685 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
686 
688 {
690 }
691 
692 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
693 
695 {
696  startFromNull = val;
697 }
698 
699 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
700 
702 {
703  splineFlag = val;
704  actSpline = true;
705 }
706 
707 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
708 
710 {
712 }
713 
714 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
715 
717 {
719 }
720 
721 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
722 
724 {
725  return currentModel;
726 }
727 
728 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
729 
730 #endif
virtual void PreparePhysicsTable(const G4ParticleDefinition &) override
Float_t x
Definition: compare.C:6
G4ParticleChangeForGamma * GetParticleChange()
G4double fFactor
void BuildLambdaTable()
const XML_Char * name
Definition: expat.h:151
G4double CrossSectionBiasingFactor() const
void SetMaxKinEnergy(G4double e)
G4double preStepKinEnergy
G4VEmModel * GetRegionModel(G4int idx, size_t couple_index) const
size_t basedCoupleIndex
G4double mfpKinEnergy
void StreamInfo(std::ostream &outFile, const G4ParticleDefinition &, G4String endOfLine=G4String("\n")) const
void SetMinKinEnergy(G4double e)
void ComputeIntegralLambda(G4double kinEnergy)
const G4Element * GetCurrentElement() const
G4bool applyCuts
std::vector< G4double > theCrossSectionMax
G4PhysicsTable * theLambdaTablePrim
const char * p
Definition: xmltok.h:285
void SetBuildTableFlag(G4bool val)
G4VEmProcess(const G4String &name, G4ProcessType type=fElectromagnetic)
Definition: G4VEmProcess.cc:91
G4int mainSecondaries
G4VEmModel * EmModel(size_t index=0) const
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="", G4bool flag=true)
const G4MaterialCutsCouple * MaterialCutsCouple() const
virtual void InitialiseProcess(const G4ParticleDefinition *)=0
virtual void PrintInfo()
void DefineMaterial(const G4MaterialCutsCouple *couple)
const G4ParticleDefinition * SecondaryParticle() const
G4bool splineFlag
G4double condition(const G4ErrorSymMatrix &m)
const G4ParticleDefinition * thePositron
virtual G4double MinPrimaryEnergy(const G4ParticleDefinition *, const G4Material *)
G4bool actMaxKinEnergy
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:455
const G4Material * baseMaterial
void PrintWarning(G4String tit, G4double val)
G4double RecalculateLambda(G4double kinEnergy, const G4MaterialCutsCouple *couple)
virtual G4bool IsApplicable(const G4ParticleDefinition &p) override=0
const G4ParticleDefinition * currentParticle
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
void SetParticle(const G4ParticleDefinition *p)
G4bool actBinning
G4VEmModel * SelectModel(G4double &kinEnergy, size_t index)
const std::vector< G4int > * theDensityIdx
const std::vector< G4double > * theCutsPositron
size_t CurrentMaterialCutsCoupleIndex() const
Float_t Z
std::vector< G4VEmModel * > emModels
G4PhysicsTable * LambdaTablePrim() const
G4double minKinEnergy
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
G4int GetNumberOfModels() const
virtual void ProcessDescription(std::ostream &outFile) const override
G4bool weightFlag
void SetEmModel(G4VEmModel *, G4int index=0)
G4double GetElectronEnergyCut()
const G4Isotope * GetTargetIsotope() const
G4ProcessType
G4PhysicsTable * LambdaTable() const
G4ParticleChangeForGamma fParticleChange
G4double MinKinEnergy() const
const G4ParticleDefinition * theElectron
G4double biasFactor
G4VEmModel * currentModel
G4double GetLambdaFromTable(G4double kinEnergy)
G4int LambdaBinning() const
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
const G4Isotope * GetCurrentIsotope() const
Definition: G4VEmModel.hh:484
G4VEmModel * SelectModel(G4double &energy, size_t &index)
const G4ParticleDefinition * secondaryParticle
G4double ComputeCrossSectionPerAtom(G4double kineticEnergy, G4double Z, G4double A=0., G4double cut=0.0)
double A(double temperature)
G4double GetLambda(G4double &kinEnergy, const G4MaterialCutsCouple *couple)
std::vector< G4DynamicParticle * > secParticles
void SetSplineFlag(G4bool val)
const G4VEmModel * GetCurrentModel() const
const G4ParticleDefinition * Particle() const
G4double maxKinEnergy
G4EmBiasingManager * biasManager
void SetSecondaryParticle(const G4ParticleDefinition *p)
virtual void StartTracking(G4Track *) override
G4EmParameters * theParameters
Definition: G4Step.hh:76
const std::vector< G4double > * theCutsGamma
G4VEmProcess & operator=(const G4VEmProcess &right)=delete
const G4Element * GetTargetElement() const
virtual G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii) override
void SetLambdaBinning(G4int nbins)
std::vector< G4double > theEnergyOfCrossSectionMax
void SetMinKinEnergyPrim(G4double e)
G4bool actSpline
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
G4double MaxKinEnergy() const
G4VEmModel * SelectModelForMaterial(G4double kinEnergy, size_t &idxRegion) const
G4PhysicsTable * theLambdaTable
virtual G4double GetMeanFreePath(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
virtual G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false) override
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:252
int G4int
Definition: G4Types.hh:78
G4double MscThetaLimit() const
G4double PolarAngleLimit() const
G4double minKinEnergyPrim
const std::vector< G4double > * theCutsElectron
G4double GetCurrentLambda(G4double kinEnergy)
const std::vector< G4double > * theDensityFactor
G4PhysicsVector * LambdaPhysicsVector(const G4MaterialCutsCouple *)
G4ForceCondition
G4bool ApplyCuts() const
size_t idxLambda
G4int GetNumberOfRegionModels(size_t couple_index) const
const G4ParticleDefinition * particle
virtual ~G4VEmProcess()
void UpdateEmModel(const G4String &, G4double, G4double)
const std::vector< G4double > * theCuts
void AddEmModel(G4int, G4VEmModel *, const G4Region *region=nullptr)
const G4Element * GetCurrentElement() const
Definition: G4VEmModel.hh:477
virtual void BuildPhysicsTable(const G4ParticleDefinition &) override
const G4Material * GetMaterial() const
G4int numberOfModels
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &) override
G4bool actMinKinEnergy
G4EmModelManager * modelManager
void SetStartFromNullFlag(G4bool val)
G4LossTableManager * lManager
G4double lambdaFactor
G4bool IsIntegral() const
void SetIntegral(G4bool val)
void FindLambdaMax()
size_t currentCoupleIndex
G4bool buildLambdaTable
G4double GetGammaEnergyCut()
G4double ComputeCurrentLambda(G4double kinEnergy)
#define DBL_MAX
Definition: templates.hh:83
G4bool startFromNull
size_t idxLambdaPrim
const G4MaterialCutsCouple * currentCouple
G4double GetLambdaFromTablePrim(G4double kinEnergy)
const G4ParticleDefinition * theGamma
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:234
const G4Material * currentMaterial
G4double CrossSectionPerVolume(G4double kineticEnergy, const G4MaterialCutsCouple *couple)
G4double preStepLambda
G4double MeanFreePath(const G4Track &track)
virtual void StreamProcessInfo(std::ostream &, G4String) const
G4double massRatio