Geant4  10.03
G4EmParameters.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: G4EmParameters.hh 66885 2013-01-16 17:37:13Z gunter $
27 //
28 //
29 // -------------------------------------------------------------------
30 //
31 // GEANT4 Class header file
32 //
33 //
34 // File name: G4EmParameters
35 //
36 // Author: Vladimir Ivanchenko for migration to MT
37 //
38 //
39 // Creation date: 17.05.2013
40 //
41 // Modifications:
42 //
43 //
44 // Class Description:
45 //
46 // A utility static class, responsable for keeping parameters
47 // for all EM physics processes and models.
48 //
49 // It is initialized by the master thread but can be updated
50 // at any moment. Parameters may be used in run time or at
51 // initialisation
52 //
53 // -------------------------------------------------------------------
54 //
55 
56 #ifndef G4EmParameters_h
57 #define G4EmParameters_h 1
58 
59 #include "globals.hh"
60 #include "G4ios.hh"
61 #include "G4MscStepLimitType.hh"
63 #include "G4EmSaturation.hh"
64 #include <vector>
65 
68 class G4VEmProcess;
70 class G4StateManager;
71 
73 {
74 public:
75 
76  static G4EmParameters* Instance();
77 
79 
80  void SetDefaults();
81 
82  // printing
83  std::ostream& StreamInfo(std::ostream& os) const;
84  void Dump() const;
85  friend std::ostream& operator<< (std::ostream& os, const G4EmParameters&);
86 
87  // boolean flags
88  void SetLossFluctuations(G4bool val);
89  G4bool LossFluctuation() const;
90 
91  void SetBuildCSDARange(G4bool val);
92  G4bool BuildCSDARange() const;
93 
94  void SetLPM(G4bool val);
95  G4bool LPM() const;
96 
97  void SetSpline(G4bool val);
98  G4bool Spline() const;
99 
100  void SetUseCutAsFinalRange(G4bool val);
101  G4bool UseCutAsFinalRange() const;
102 
103  void SetApplyCuts(G4bool val);
104  G4bool ApplyCuts() const;
105 
106  void SetFluo(G4bool val);
107  G4bool Fluo() const;
108 
109  void SetBeardenFluoDir(G4bool val);
110  G4bool BeardenFluoDir() const;
111 
112  void SetAuger(G4bool val);
113  G4bool Auger() const;
114 
115  void SetAugerCascade(G4bool val);
116  G4bool AugerCascade() const;
117 
118  void SetPixe(G4bool val);
119  G4bool Pixe() const;
120 
123 
124  void SetLateralDisplacement(G4bool val);
125  G4bool LateralDisplacement() const;
126 
129 
132 
135 
136  void SetUseMottCorrection(G4bool val);
137  G4bool UseMottCorrection() const;
138 
139  void SetIntegral(G4bool val);
140  G4bool Integral() const;
141 
142  void SetBirksActive(G4bool val);
143  G4bool BirksActive() const;
144 
147 
148  // double parameters with values
149  void SetMinSubRange(G4double val);
150  G4double MinSubRange() const;
151 
152  void SetMinEnergy(G4double val);
153  G4double MinKinEnergy() const;
154 
155  void SetMaxEnergy(G4double val);
156  G4double MaxKinEnergy() const;
157 
160 
163 
164  void SetLowestMuHadEnergy(G4double val);
165  G4double LowestMuHadEnergy() const;
166 
167  void SetLinearLossLimit(G4double val);
168  G4double LinearLossLimit() const;
169 
170  void SetBremsstrahlungTh(G4double val);
171  G4double BremsstrahlungTh() const;
172 
173  void SetLambdaFactor(G4double val);
174  G4double LambdaFactor() const;
175 
178 
179  void SetMscThetaLimit(G4double val);
180  G4double MscThetaLimit() const;
181 
182  void SetMscRangeFactor(G4double val);
183  G4double MscRangeFactor() const;
184 
187 
188  void SetMscGeomFactor(G4double val);
189  G4double MscGeomFactor() const;
190 
191  void SetMscSkin(G4double val);
192  G4double MscSkin() const;
193 
194  void SetStepFunction(G4double v1, G4double v2);
195 
197 
198  // integer parameters
199  void SetNumberOfBins(G4int val);
200  G4int NumberOfBins() const;
201 
204 
205  void SetVerbose(G4int val);
206  G4int Verbose() const;
207 
208  void SetWorkerVerbose(G4int val);
209  G4int WorkerVerbose() const;
210 
213 
216 
219 
220  // string parameters
221  void SetPIXECrossSectionModel(const G4String&);
223 
226 
227  // parameters per region or per process
228  void AddPAIModel(const G4String& particle,
229  const G4String& region,
230  const G4String& type);
231  const std::vector<G4String>& ParticlesPAI() const;
232  const std::vector<G4String>& RegionsPAI() const;
233  const std::vector<G4String>& TypesPAI() const;
234 
235  void AddMicroElec(const G4String& region);
236  const std::vector<G4String>& RegionsMicroElec() const;
237 
238  void AddDNA(const G4String& region, const G4String& type);
239  const std::vector<G4String>& RegionsDNA() const;
240  const std::vector<G4String>& TypesDNA() const;
241 
242  void AddMsc(const G4String& region, const G4String& type);
243  const std::vector<G4String>& RegionsMsc() const;
244  const std::vector<G4String>& TypesMsc() const;
245 
246  void SetSubCutoff(G4bool val, const G4String& region = "");
247 
248  void SetDeexActiveRegion(const G4String& region, G4bool fdeex,
249  G4bool fauger, G4bool fpixe);
250 
251  void SetProcessBiasingFactor(const G4String& procname,
252  G4double val, G4bool wflag);
253 
254  void ActivateForcedInteraction(const G4String& procname,
255  const G4String& region,
256  G4double length,
257  G4bool wflag);
258 
260  const G4String& region,
261  G4double factor,
262  G4double energyLimit);
263 
264  // initialisation methods
266  G4bool isElectron) const;
267  void DefineRegParamForEM(G4VEmProcess*) const;
269 
270 private:
271 
272  G4EmParameters(G4EmParameters &) = delete;
273  G4EmParameters & operator=(const G4EmParameters &right) = delete;
274 
275  G4EmParameters();
276 
277  void Initialise();
278 
279  G4bool IsLocked() const;
280 
281  G4String CheckRegion(const G4String&) const;
282 
283  void PrintWarning(G4ExceptionDescription& ed) const;
284 
286 
288 
290 
292 
312 
332 
337 
341 
344 
345  std::vector<G4String> m_particlesPAI;
346  std::vector<G4String> m_regnamesPAI;
347  std::vector<G4String> m_typesPAI;
348 
349  std::vector<G4String> m_regnamesME;
350 
351  std::vector<G4String> m_regnamesDNA;
352  std::vector<G4String> m_typesDNA;
353 
354  std::vector<G4String> m_regnamesMsc;
355  std::vector<G4String> m_typesMsc;
356 
357  std::vector<G4String> m_regnamesSubCut;
358  std::vector<G4bool> m_subCuts;
359 
360  std::vector<G4String> m_regnamesDeex;
361  std::vector<G4bool> m_fluo;
362  std::vector<G4bool> m_auger;
363  std::vector<G4bool> m_pixe;
364 
365  std::vector<G4String> m_procBiasedXS;
366  std::vector<G4double> m_factBiasedXS;
367  std::vector<G4bool> m_weightBiasedXS;
368 
369  std::vector<G4String> m_procForced;
370  std::vector<G4String> m_regnamesForced;
371  std::vector<G4double> m_lengthForced;
372  std::vector<G4bool> m_weightForced;
373 
374  std::vector<G4String> m_procBiasedSec;
375  std::vector<G4String> m_regnamesBiasedSec;
376  std::vector<G4double> m_factBiasedSec;
377  std::vector<G4double> m_elimBiasedSec;
378 
379 };
380 
381 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
382 
383 #endif
384 
G4bool UseCutAsFinalRange() const
static G4EmParameters * theInstance
void SetLossFluctuations(G4bool val)
G4int NumberOfBinsPerDecade() const
std::vector< G4String > m_regnamesBiasedSec
void SetApplyCuts(G4bool val)
G4MscStepLimitType mscStepLimit
std::vector< G4String > m_regnamesForced
void SetProcessBiasingFactor(const G4String &procname, G4double val, G4bool wflag)
G4int NumberOfBins() const
void SetVerbose(G4int val)
G4bool Spline() const
G4int WorkerVerbose() const
G4double MaxKinEnergy() const
void SetEmSaturation(G4EmSaturation *)
G4bool isElectron(G4int ityp)
std::vector< G4String > m_typesDNA
void ActivateSecondaryBiasing(const G4String &name, const G4String &region, G4double factor, G4double energyLimit)
void SetDeexcitationIgnoreCut(G4bool val)
void SetUseMottCorrection(G4bool val)
std::vector< G4double > m_lengthForced
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
std::vector< G4String > m_typesPAI
void SetLowestElectronEnergy(G4double val)
void SetLatDisplacementBeyondSafety(G4bool val)
G4bool IsLocked() const
G4MscStepLimitType MscMuHadStepLimitType() const
void SetMscStepLimitType(G4MscStepLimitType val)
G4double linLossLimit
G4bool AugerCascade() const
G4double LowestElectronEnergy() const
G4bool BirksActive() const
void SetBeardenFluoDir(G4bool val)
void SetLinearLossLimit(G4double val)
G4double MscGeomFactor() const
G4double MscMuHadRangeFactor() const
G4double MscThetaLimit() const
std::vector< G4String > m_regnamesME
void SetAuger(G4bool val)
std::vector< G4bool > m_pixe
const std::vector< G4String > & ParticlesPAI() const
G4String nameElectronPIXE
const std::vector< G4String > & RegionsMicroElec() const
G4NuclearFormfactorType nucFormfactor
G4bool Integral() const
G4double maxKinEnergy
G4bool muhadLateralDisplacement
std::ostream & StreamInfo(std::ostream &os) const
void SetNumberOfBins(G4int val)
G4String CheckRegion(const G4String &) const
void SetMinSubRange(G4double val)
G4bool LPM() const
void SetMaxEnergyForCSDARange(G4double val)
const std::vector< G4String > & RegionsPAI() const
const char * name(G4int ptype)
G4bool ApplyCuts() const
void DefineRegParamForLoss(G4VEnergyLossProcess *, G4bool isElectron) const
void SetStepFunctionMuHad(G4double v1, G4double v2)
const std::vector< G4String > & TypesDNA() const
G4bool Fluo() const
void SetPIXEElectronCrossSectionModel(const G4String &)
int G4int
Definition: G4Types.hh:78
void SetDeexActiveRegion(const G4String &region, G4bool fdeex, G4bool fauger, G4bool fpixe)
const std::vector< G4String > & RegionsMsc() const
G4bool useMottCorrection
void SetMaxEnergy(G4double val)
void SetBremsstrahlungTh(G4double val)
std::vector< G4String > m_procForced
G4double MinSubRange() const
G4bool DeexcitationIgnoreCut() const
std::vector< G4String > m_procBiasedXS
G4bool LatDisplacementBeyondSafety() const
std::vector< G4String > m_regnamesMsc
void SetBirksActive(G4bool val)
void AddPAIModel(const G4String &particle, const G4String &region, const G4String &type)
std::vector< G4String > m_regnamesDeex
void SetLateralDisplacement(G4bool val)
G4bool BuildCSDARange() const
G4EmSaturation * GetEmSaturation()
void SetWorkerVerbose(G4int val)
std::vector< G4bool > m_auger
G4bool lateralDisplacement
G4double LinearLossLimit() const
void SetPIXECrossSectionModel(const G4String &)
G4double LambdaFactor() const
std::vector< G4bool > m_subCuts
G4double dRoverRangeMuHad
G4int Verbose() const
G4bool MuHadLateralDisplacement() const
const G4String & PIXECrossSectionModel()
G4double LowestMuHadEnergy() const
G4double lowestMuHadEnergy
bool G4bool
Definition: G4Types.hh:79
void SetMscRangeFactor(G4double val)
friend std::ostream & operator<<(std::ostream &os, const G4EmParameters &)
void SetAugerCascade(G4bool val)
void Dump() const
G4double finalRangeMuHad
void SetNumberOfBinsPerDecade(G4int val)
G4bool LateralDisplacement() const
G4double minSubRange
void SetNuclearFormfactorType(G4NuclearFormfactorType val)
void SetLowestMuHadEnergy(G4double val)
void SetMscGeomFactor(G4double val)
std::vector< G4String > m_regnamesDNA
G4double MinKinEnergy() const
void SetMscMuHadStepLimitType(G4MscStepLimitType val)
std::vector< G4String > m_regnamesPAI
void AddDNA(const G4String &region, const G4String &type)
void SetSubCutoff(G4bool val, const G4String &region="")
std::vector< G4String > m_particlesPAI
G4double lowestElectronEnergy
G4MscStepLimitType mscStepLimitMuHad
G4double lambdaFactor
void PrintWarning(G4ExceptionDescription &ed) const
G4EmSaturation * emSaturation
G4NuclearFormfactorType NuclearFormfactorType() const
void SetBuildCSDARange(G4bool val)
void AddMicroElec(const G4String &region)
G4bool Auger() const
G4bool LossFluctuation() const
const std::vector< G4String > & TypesMsc() const
void SetMscMuHadRangeFactor(G4double val)
void SetPixe(G4bool val)
std::vector< G4bool > m_fluo
void SetSpline(G4bool val)
std::vector< G4String > m_procBiasedSec
void SetMuHadLateralDisplacement(G4bool val)
G4EmParameters & operator=(const G4EmParameters &right)=delete
G4double maxKinEnergyCSDA
void SetMinEnergy(G4double val)
void SetLPM(G4bool val)
G4bool UseAngularGeneratorForIonisation() const
G4NuclearFormfactorType
G4MscStepLimitType
G4double BremsstrahlungTh() const
static G4EmParameters * Instance()
void DefineRegParamForDeex(G4VAtomDeexcitation *) const
std::vector< G4double > m_elimBiasedSec
G4double MscRangeFactor() const
void SetIntegral(G4bool val)
void AddMsc(const G4String &region, const G4String &type)
void ActivateAngularGeneratorForIonisation(G4bool val)
std::vector< G4bool > m_weightBiasedXS
void ActivateForcedInteraction(const G4String &procname, const G4String &region, G4double length, G4bool wflag)
void SetUseCutAsFinalRange(G4bool val)
G4bool latDisplacementBeyondSafety
G4bool UseMottCorrection() const
const std::vector< G4String > & RegionsDNA() const
G4bool Pixe() const
G4double dRoverRange
G4double minKinEnergy
G4MscStepLimitType MscStepLimitType() const
void SetMscThetaLimit(G4double val)
void SetLambdaFactor(G4double val)
void SetFactorForAngleLimit(G4double val)
G4double rangeFactorMuHad
double G4double
Definition: G4Types.hh:76
void DefineRegParamForEM(G4VEmProcess *) const
G4bool useAngGeneratorForIonisation
G4EmParametersMessenger * theMessenger
G4double MaxEnergyForCSDARange() const
std::vector< G4String > m_regnamesSubCut
G4double MscSkin() const
const G4String & PIXEElectronCrossSectionModel()
const std::vector< G4String > & TypesPAI() const
G4StateManager * fStateManager
void SetStepFunction(G4double v1, G4double v2)
std::vector< G4bool > m_weightForced
std::vector< G4double > m_factBiasedXS
void SetFluo(G4bool val)
G4double FactorForAngleLimit() const
G4bool BeardenFluoDir() const
std::vector< G4String > m_typesMsc
void SetMscSkin(G4double val)
G4double rangeFactor
std::vector< G4double > m_factBiasedSec
G4double factorForAngleLimit