Geant4-11
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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// -------------------------------------------------------------------
27//
28// GEANT4 Class header file
29//
30//
31// File name: G4VEmProcess
32//
33// Author: Vladimir Ivanchenko
34//
35// Creation date: 01.10.2003
36//
37// Modifications: Vladimir Ivanchenko
38//
39// Class Description:
40//
41// It is the base class - EM discrete and rest/discrete process
42
43// -------------------------------------------------------------------
44//
45
46#ifndef G4VEmProcess_h
47#define G4VEmProcess_h 1
48
50
51#include "G4VDiscreteProcess.hh"
52#include "globals.hh"
53#include "G4Material.hh"
55#include "G4Track.hh"
56#include "G4UnitsTable.hh"
59#include "G4EmParameters.hh"
60#include "G4EmDataHandler.hh"
61#include "G4EmTableType.hh"
62#include "G4EmModelManager.hh"
64
65class G4Step;
66class G4VEmModel;
67class G4DataVector;
69class G4PhysicsTable;
70class G4PhysicsVector;
73
74//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
75
77{
78public:
79
81
82 ~G4VEmProcess() override;
83
84 //------------------------------------------------------------------------
85 // Virtual methods to be implemented in concrete processes
86 //------------------------------------------------------------------------
87
88 virtual G4bool IsApplicable(const G4ParticleDefinition& p) override = 0;
89
90 void ProcessDescription(std::ostream& outFile) const override;
91
92protected:
93
94 virtual void StreamProcessInfo(std::ostream&) const {};
95
96 virtual void InitialiseProcess(const G4ParticleDefinition*) = 0;
97
98 //------------------------------------------------------------------------
99 // Method with standard implementation; may be overwritten if needed
100 //------------------------------------------------------------------------
101
103 const G4Material*);
104
105 //------------------------------------------------------------------------
106 // Implementation of virtual methods common to all Discrete processes
107 //------------------------------------------------------------------------
108
109public:
110
111 // Initialise for build of tables
112 void PreparePhysicsTable(const G4ParticleDefinition&) override;
113
114 // Build physics table during initialisation
115 void BuildPhysicsTable(const G4ParticleDefinition&) override;
116
117 // Called before tracking of each new G4Track
118 void StartTracking(G4Track*) override;
119
120 // implementation of virtual method, specific for G4VEmProcess
122 const G4Track& track,
123 G4double previousStepSize,
124 G4ForceCondition* condition) override;
125
126 // implementation of virtual method, specific for G4VEmProcess
127 G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&) override;
128
129 // Store PhysicsTable in a file.
130 // Return false in case of failure at I/O
132 const G4String& directory,
133 G4bool ascii = false) override;
134
135 // Retrieve Physics from a file.
136 // (return true if the Physics Table can be build by using file)
137 // (return false if the process has no functionality or in case of failure)
138 // File name should is constructed as processName+particleName and the
139 // should be placed under the directory specified by the argument.
141 const G4String& directory,
142 G4bool ascii) override;
143
144 // allowing check process name
145 virtual G4VEmProcess* GetEmProcess(const G4String& name);
146
147 //------------------------------------------------------------------------
148 // Specific methods for Discrete EM post step simulation
149 //------------------------------------------------------------------------
150
151 // It returns the cross section per volume for energy/ material
153 const G4MaterialCutsCouple* couple,
154 G4double logKinEnergy = DBL_MAX);
155
156 // It returns the cross section of the process per atom
158 G4double Z, G4double A=0.,
159 G4double cut=0.0);
160
161 G4double MeanFreePath(const G4Track& track);
162
163 // Obsolete method to access cross section per volume
164 G4double GetLambda(G4double kinEnergy, const G4MaterialCutsCouple* couple);
165
166 // The main method to access cross section per volume
167 inline G4double GetLambda(G4double kinEnergy,
168 const G4MaterialCutsCouple* couple,
169 G4double logKinEnergy);
170
171 //------------------------------------------------------------------------
172 // Specific methods to build and access Physics Tables
173 //------------------------------------------------------------------------
174
175 // Binning for lambda table
176 void SetLambdaBinning(G4int nbins);
177
178 // Min kinetic energy for tables
180
181 // Min kinetic energy for high energy table
183
184 // Max kinetic energy for tables
186
187 // for cross section with one peak
188 void SetEnergyOfCrossSectionMax(std::vector<G4double>*);
189
190 // Cross section table pointers
191 inline G4PhysicsTable* LambdaTable() const;
192 inline G4PhysicsTable* LambdaTablePrim() const;
193 inline std::vector<G4double>* EnergyOfCrossSectionMax() const;
194 //------------------------------------------------------------------------
195 // Define and access particle type
196 //------------------------------------------------------------------------
197
198 inline const G4ParticleDefinition* Particle() const;
199 inline const G4ParticleDefinition* SecondaryParticle() const;
200
201 //------------------------------------------------------------------------
202 // Specific methods to set, access, modify models and basic parameters
203 //------------------------------------------------------------------------
204
205protected:
206
207 // Select model in run time
208 inline G4VEmModel* SelectModel(G4double kinEnergy, size_t);
209
210public:
211
212 // Select model by energy and couple index
213 inline G4VEmModel* SelectModelForMaterial(G4double kinEnergy,
214 size_t idxCouple) const;
215
216 // Add model for region, smaller value of order defines which
217 // model will be selected for a given energy interval
218 void AddEmModel(G4int, G4VEmModel*, const G4Region* region = nullptr);
219
220 // Assign a model to a process local list, to enable the list in run time
221 // the derived process should execute AddEmModel(..) for all such models
222 void SetEmModel(G4VEmModel*, G4int index = 0);
223
224 inline G4int NumberOfModels() const;
225
226 // return a model from the local list
227 inline G4VEmModel* EmModel(size_t index = 0) const;
228
229 // Access to active model
230 inline const G4VEmModel* GetCurrentModel() const;
231
232 // Access to models
233 G4VEmModel* GetModelByIndex(G4int idx = 0, G4bool ver = false) const;
234
235 // Access to the current G4Element
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 std::vector<G4double>* FindLambdaMax();
251
252 inline void SetEmMasterProcess(const G4VEmProcess*);
253
255
256 inline void SetBuildTableFlag(G4bool val);
257
259
260 inline G4bool UseBaseMaterial() const;
261
262 // hide copy constructor and assignment operator
264 G4VEmProcess & operator=(const G4VEmProcess &right) = delete;
265
266 //------------------------------------------------------------------------
267 // Other generic methods
268 //------------------------------------------------------------------------
269
270protected:
271
272 G4double GetMeanFreePath(const G4Track& track,
273 G4double previousStepSize,
274 G4ForceCondition* condition) override;
275
277
278 inline void DefineMaterial(const G4MaterialCutsCouple* couple);
279
280 inline G4int LambdaBinning() const;
281
282 inline G4double MinKinEnergy() const;
283
284 inline G4double MaxKinEnergy() const;
285
286 // Single scattering parameters
287 inline G4double PolarAngleLimit() const;
288
290
291 inline G4double RecalculateLambda(G4double kinEnergy,
292 const G4MaterialCutsCouple* couple);
293
295
296 inline void SetParticle(const G4ParticleDefinition* p);
297
298 inline void SetSecondaryParticle(const G4ParticleDefinition* p);
299
300 inline size_t CurrentMaterialCutsCoupleIndex() const;
301
302 inline const G4MaterialCutsCouple* MaterialCutsCouple() const;
303
304 inline G4bool ApplyCuts() const;
305
307
309
310 inline void SetStartFromNullFlag(G4bool val);
311
312 inline void SetSplineFlag(G4bool val);
313
314 inline const G4Element* GetTargetElement() const;
315
316 inline const G4Isotope* GetTargetIsotope() const;
317
318 // these two methods assume that vectors are initilized
319 // and idx is within vector length
320 inline G4int DensityIndex(G4int idx) const;
321 inline G4double DensityFactor(G4int idx) const;
322
323private:
324
325 void Clear();
326
327 void BuildLambdaTable();
328
329 void StreamInfo(std::ostream& outFile, const G4ParticleDefinition&,
330 G4bool rst=false) const;
331
332 void PrintWarning(G4String tit, G4double val);
333
334 void ComputeIntegralLambda(G4double kinEnergy, G4double logKinEnergy);
335
336 inline G4double GetLambdaFromTable(G4double kinEnergy);
337
338 inline G4double GetLambdaFromTable(G4double kinEnergy, G4double logKinEnergy);
339
340 inline G4double GetLambdaFromTablePrim(G4double kinEnergy);
341
342 inline G4double GetLambdaFromTablePrim(G4double kinEnergy, G4double logKinEnergy);
343
344 inline G4double GetCurrentLambda(G4double kinEnergy);
345
346 inline G4double GetCurrentLambda(G4double kinEnergy, G4double logKinEnergy);
347
348 inline G4double ComputeCurrentLambda(G4double kinEnergy);
349
350 // ======== pointers =========
351
359 const G4VEmProcess* masterProc = nullptr;
364 const G4Material* baseMaterial = nullptr;
365
366 // ======== tables and vectors ========
369
370 const std::vector<G4double>* theCuts = nullptr;
371 const std::vector<G4double>* theCutsGamma = nullptr;
372 const std::vector<G4double>* theCutsElectron = nullptr;
373 const std::vector<G4double>* theCutsPositron = nullptr;
374
375protected:
376
377 // ======== pointers =========
378
380 const G4Material* currentMaterial = nullptr;
382 std::vector<G4double>* theEnergyOfCrossSectionMax = nullptr;
383
384private:
385
386 const std::vector<G4double>* theDensityFactor = nullptr;
387 const std::vector<G4int>* theDensityIdx = nullptr;
388
389 // ======== parameters =========
400
401protected:
402
407
408private:
409
411
414
415protected:
416
425 size_t coupleIdxLambda = 0;
426 size_t idxLambda = 0;
427
430
431private:
432
440 G4bool isIon = false;
443
444protected:
445
446 // ======== particle change =========
447 std::vector<G4DynamicParticle*> secParticles;
449
450private:
451
452 // ======== local vectors =========
453 std::vector<G4VEmModel*> emModels;
454
455};
456
457// ======== Run time inline methods ================
458
460{
461 return applyCuts;
462}
463
464//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
465
467{
468 return currentCoupleIndex;
469}
470
471//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
472
474{
475 return currentCouple;
476}
477
478//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
479
481{
483}
484
485//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
486
488{
490}
491
492//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
493
495{
496 if(couple != currentCouple) {
497 currentCouple = couple;
502 if(baseMat) {
503 basedCoupleIndex = (*theDensityIdx)[currentCoupleIndex];
504 if(nullptr != currentMaterial->GetBaseMaterial())
506 fFactor *= (*theDensityFactor)[currentCoupleIndex];
507 }
508 }
509}
510
511//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
512
513inline
515{
516 if(1 < numberOfModels) {
518 }
520 return currentModel;
521}
522
523//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
524
525inline
527 size_t idxCouple) const
528{
529 return modelManager->SelectModel(kinEnergy, idxCouple);
530}
531
532//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
533
535{
536 return ((*theLambdaTable)[basedCoupleIndex])->Value(e, idxLambda);
537}
538
539//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
540
542{
543 return ((*theLambdaTable)[basedCoupleIndex])->LogVectorValue(e, loge);
544}
545
546//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
547
549{
550 return ((*theLambdaTablePrim)[basedCoupleIndex])->Value(e, idxLambda)/e;
551}
552
553//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
554
556{
557 return ((*theLambdaTablePrim)[basedCoupleIndex])->LogVectorValue(e, loge)/e;
558}
559
560//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
561
563{
565}
566
567//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
568
570{
573 fLambdaEnergy = e;
575 else if(nullptr != theLambdaTable) { fLambda = GetLambdaFromTable(e); }
576 else { fLambda = ComputeCurrentLambda(e); }
577 fLambda *= fFactor;
578 }
579 return fLambda;
580}
581
582//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
583
585{
588 fLambdaEnergy = e;
589 if(e >= minKinEnergyPrim) { fLambda = GetLambdaFromTablePrim(e, loge); }
590 else if(nullptr != theLambdaTable) { fLambda = GetLambdaFromTable(e, loge); }
591 else { fLambda = ComputeCurrentLambda(e); }
592 fLambda *= fFactor;
593 }
594 return fLambda;
595}
596
597//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
598
599inline void
601{
602 DefineMaterial(couple);
604}
605
606//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
607
608inline G4double
610 G4double logKinEnergy)
611{
612 CurrentSetup(couple, kinEnergy);
613 return GetCurrentLambda(kinEnergy, logKinEnergy);
614}
615
616//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
617
618inline G4double
620{
621 CurrentSetup(couple, e);
623}
624
625// ======== Get/Set inline methods used at initialisation ================
626
627//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
628
630{
631 return nLambdaBins;
632}
633
634//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
635
637{
638 return minKinEnergy;
639}
640
641//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
642
644{
645 return maxKinEnergy;
646}
647
648//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
649
651{
652 return biasFactor;
653}
654
655//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
656
658{
659 return theLambdaTable;
660}
661
662//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
663
665{
666 return theLambdaTablePrim;
667}
668
669//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
670
671inline std::vector<G4double>* G4VEmProcess::EnergyOfCrossSectionMax() const
672{
674}
675
676//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
677
679{
680 return particle;
681}
682
683//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
684
686{
687 return secondaryParticle;
688}
689
690//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
691
693{
694 fXSType = val;
695}
696
697//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
698
700{
701 return fXSType;
702}
703
704//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
705
707{
708 buildLambdaTable = val;
709}
710
711//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
712
714{
715 return &fParticleChange;
716}
717
718//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
719
721{
722 particle = p;
723 currentParticle = p;
724}
725
726//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
727
729{
731}
732
733//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
734
736{
737 startFromNull = val;
738}
739
740//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
741
743{
744 splineFlag = val;
745}
746
747//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
748
750{
752}
753
754//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
755
757{
759}
760
761//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
762
764{
765 return (*theDensityIdx)[idx];
766}
767
768//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
769
771{
772 return (*theDensityFactor)[idx];
773}
774
775//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
776
778{
779 return baseMat;
780}
781
782//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
783
785{
786 return currentModel;
787}
788
789//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
790
792{
793 masterProc = ptr;
794}
795
796//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
797
799{
800 return numberOfModels;
801}
802
803//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
804
805inline G4VEmModel* G4VEmProcess::EmModel(size_t index) const
806{
807 return (index < emModels.size()) ? emModels[index] : nullptr;
808}
809
810//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
811
812#endif
G4CrossSectionType
@ fEmNoIntegral
G4double condition(const G4ErrorSymMatrix &m)
G4ForceCondition
G4ProcessType
@ fElectromagnetic
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
int G4int
Definition: G4Types.hh:85
const G4int Z[17]
const G4double A[17]
G4VEmModel * SelectModel(G4double energy, size_t index)
const G4Material * GetMaterial() const
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:229
Definition: G4Step.hh:62
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:237
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:472
const G4Element * GetCurrentElement() const
Definition: G4VEmModel.hh:505
const G4Isotope * GetCurrentIsotope() const
Definition: G4VEmModel.hh:512
G4CrossSectionType CrossSectionType() const
const std::vector< G4double > * theCutsElectron
void DefineMaterial(const G4MaterialCutsCouple *couple)
G4double RecalculateLambda(G4double kinEnergy, const G4MaterialCutsCouple *couple)
G4double fLambdaEnergy
G4bool ApplyCuts() const
G4double GetGammaEnergyCut()
G4double maxKinEnergy
void ComputeIntegralLambda(G4double kinEnergy, G4double logKinEnergy)
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
G4double MeanFreePath(const G4Track &track)
G4bool applyCuts
G4int mainSecondaries
G4PhysicsVector * LambdaPhysicsVector(const G4MaterialCutsCouple *)
void CurrentSetup(const G4MaterialCutsCouple *, G4double energy)
G4double CrossSectionPerVolume(G4double kineticEnergy, const G4MaterialCutsCouple *couple, G4double logKinEnergy=DBL_MAX)
const G4ParticleDefinition * Particle() const
virtual void StreamProcessInfo(std::ostream &) const
Definition: G4VEmProcess.hh:94
G4double massRatio
G4VEmProcess(const G4String &name, G4ProcessType type=fElectromagnetic)
Definition: G4VEmProcess.cc:79
G4LossTableManager * lManager
G4VEmModel * SelectModel(G4double kinEnergy, size_t)
G4EmBiasingManager * biasManager
const G4ParticleDefinition * secondaryParticle
G4bool weightFlag
void BuildPhysicsTable(const G4ParticleDefinition &) override
G4double preStepLambda
G4VEmModel * EmModel(size_t index=0) const
void SetBuildTableFlag(G4bool val)
G4double preStepLogKinEnergy
virtual G4double MinPrimaryEnergy(const G4ParticleDefinition *, const G4Material *)
const std::vector< G4double > * theCutsGamma
G4double ComputeCrossSectionPerAtom(G4double kineticEnergy, G4double Z, G4double A=0., G4double cut=0.0)
G4VEmModel * currentModel
G4VEmProcess(G4VEmProcess &)=delete
void SetEnergyOfCrossSectionMax(std::vector< G4double > *)
G4bool splineFlag
G4PhysicsTable * LambdaTable() const
std::vector< G4double > * theEnergyOfCrossSectionMax
virtual G4bool IsApplicable(const G4ParticleDefinition &p) override=0
G4bool actMaxKinEnergy
void SetMinKinEnergy(G4double e)
G4int numberOfModels
size_t basedCoupleIndex
G4double GetElectronEnergyCut()
const G4ParticleDefinition * currentParticle
G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &) override
G4int DensityIndex(G4int idx) const
std::vector< G4VEmModel * > emModels
size_t idxLambda
void StartTracking(G4Track *) override
G4double CrossSectionBiasingFactor() const
void AddEmModel(G4int, G4VEmModel *, const G4Region *region=nullptr)
G4bool actMinKinEnergy
void SetEmModel(G4VEmModel *, G4int index=0)
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
G4double mfpKinEnergy
void SetEmMasterProcess(const G4VEmProcess *)
virtual void InitialiseProcess(const G4ParticleDefinition *)=0
void StreamInfo(std::ostream &outFile, const G4ParticleDefinition &, G4bool rst=false) const
size_t coupleIdxLambda
G4double GetLambdaFromTable(G4double kinEnergy)
~G4VEmProcess() override
void SetLambdaBinning(G4int nbins)
const std::vector< G4double > * theDensityFactor
void PrintWarning(G4String tit, G4double val)
G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii) override
void SetSecondaryParticle(const G4ParticleDefinition *p)
const G4VEmModel * GetCurrentModel() const
void SetSplineFlag(G4bool val)
G4VEmModel * SelectModelForMaterial(G4double kinEnergy, size_t idxCouple) const
std::vector< G4double > * EnergyOfCrossSectionMax() const
G4bool UseBaseMaterial() const
G4CrossSectionType fXSType
G4int NumberOfModels() const
void SetCrossSectionType(G4CrossSectionType val)
void ProcessDescription(std::ostream &outFile) const override
G4double GetCurrentLambda(G4double kinEnergy)
G4PhysicsTable * theLambdaTablePrim
G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
G4VEmProcess & operator=(const G4VEmProcess &right)=delete
G4double fFactor
const G4Element * GetTargetElement() const
G4double DensityFactor(G4int idx) const
G4bool actBinning
G4double lambdaFactor
virtual G4VEmProcess * GetEmProcess(const G4String &name)
G4double MaxKinEnergy() const
G4double GetLambdaFromTablePrim(G4double kinEnergy)
G4double logLambdaFactor
G4ParticleChangeForGamma * GetParticleChange()
G4bool startFromNull
const G4Isotope * GetTargetIsotope() const
const G4ParticleDefinition * thePositron
G4bool isTheMaster
G4double MinKinEnergy() const
G4double biasFactor
std::vector< G4DynamicParticle * > secParticles
std::vector< G4double > * FindLambdaMax()
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="", G4bool flag=true)
const G4ParticleDefinition * theElectron
G4EmParameters * theParameters
const G4MaterialCutsCouple * currentCouple
const G4ParticleDefinition * theGamma
G4double minKinEnergyPrim
G4bool buildLambdaTable
void SetStartFromNullFlag(G4bool val)
G4PhysicsTable * LambdaTablePrim() const
G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false) override
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
const G4MaterialCutsCouple * MaterialCutsCouple() const
G4double preStepKinEnergy
G4double PolarAngleLimit() const
const G4VEmProcess * masterProc
const std::vector< G4double > * theCuts
size_t currentCoupleIndex
G4double GetLambda(G4double kinEnergy, const G4MaterialCutsCouple *couple)
G4ParticleChangeForGamma fParticleChange
G4PhysicsTable * theLambdaTable
const std::vector< G4int > * theDensityIdx
const std::vector< G4double > * theCutsPositron
G4double minKinEnergy
G4int LambdaBinning() const
void SetParticle(const G4ParticleDefinition *p)
const G4ParticleDefinition * particle
G4double fLambda
void SetMinKinEnergyPrim(G4double e)
const G4Material * baseMaterial
void PreparePhysicsTable(const G4ParticleDefinition &) override
void SetMaxKinEnergy(G4double e)
G4EmModelManager * modelManager
G4EmDataHandler * theData
const G4Material * currentMaterial
size_t CurrentMaterialCutsCoupleIndex() const
G4double ComputeCurrentLambda(G4double kinEnergy)
const G4Element * GetCurrentElement() const
G4double GetMeanFreePath(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
void BuildLambdaTable()
const G4ParticleDefinition * SecondaryParticle() const
G4double energy(const ThreeVector &p, const G4double m)
const char * name(G4int ptype)
#define LOG_EKIN_MIN
Definition: templates.hh:98
#define DBL_MAX
Definition: templates.hh:62