| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer << 3 // * DISCLAIMER *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers *
6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which *
7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: *
8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license *
9 // * include a list of copyright holders. <<
10 // * 9 // * *
11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file << 14 // * use. *
16 // * for the full disclaimer and the limitatio <<
17 // * 15 // * *
18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the *
19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. *
20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work *
21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this *
22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. *
23 // * acceptance of all terms of the Geant4 Sof <<
24 // ******************************************* 21 // ********************************************************************
25 // 22 //
26 // G4SPSEneDistribution << 23 ///////////////////////////////////////////////////////////////////////////////
>> 24 //
>> 25 // MODULE: G4SPSEneDistribution.hh
>> 26 //
>> 27 // Version: 1.0
>> 28 // Date: 5/02/04
>> 29 // Author: Fan Lei
>> 30 // Organisation: QinetiQ ltd.
>> 31 // Customer: ESA/ESTEC
>> 32 //
>> 33 ///////////////////////////////////////////////////////////////////////////////
>> 34 //
>> 35 // CHANGE HISTORY
>> 36 // --------------
>> 37 //
>> 38 //
>> 39 // Version 1.0, 05/02/2004, Fan Lei, Created.
>> 40 // Based on the G4GeneralParticleSource class in Geant4 v6.0
>> 41 //
>> 42 ///////////////////////////////////////////////////////////////////////////////
>> 43 //
27 // 44 //
28 // Class Description: 45 // Class Description:
29 // 46 //
30 // To generate the energy of a primary vertex << 47 // To generate the energy of a primary vertex according to the defined distribution
31 // defined distribution. This is a shared clas << 48 //
32 // Only one thread should use the set-methods << 49 ///////////////////////////////////////////////////////////////////////////////
33 // Note that this is exactly what is achieved << 50 //
34 // If you use the set methods to set defaults << 51 // MEMBER FUNCTIONS
35 // care that only one thread is executing them << 52 // ----------------
36 // In addition take care of calling these meth << 53 //
37 // started. Do not use the setters during the << 54 // G4SPSEneDistribution ()
>> 55 // Constructor: Initializes variables
>> 56 //
>> 57 // ~G4SPSEneDistribution ()
>> 58 // Destructor:
>> 59 //
>> 60 // void SetEnergyDisType(G4String)
>> 61 // Allows the user to choose the energy distribution type. The arguments
>> 62 // are Mono (mono-energetic), Lin (linear), Pow (power-law), Exp
>> 63 // (exponential), Gauss (gaussian), Brem (bremsstrahlung), BBody (black-body), Cdg
>> 64 // (cosmic diffuse gamma-ray), User (user-defined), Arb (arbitrary
>> 65 // point-wise), Epn (energy per nucleon).
>> 66 //
>> 67 // void SetEmin(G4double)
>> 68 // Sets the minimum energy.
>> 69 //
>> 70 // void SetEmax(G4double)
>> 71 // Sets the maximum energy.
>> 72 //
>> 73 // void SetMonoEnergy(G4double)
>> 74 // Sets energy for mono-energetic distribution.
>> 75 //
>> 76 // void SetAlpha(G4double)
>> 77 // Sets alpha for a power-law distribution.
>> 78 //
>> 79 // void SetTemp(G4double)
>> 80 // Sets Temperature for a Brem or BBody distributions.
>> 81 //
>> 82 // void SetEzero(G4double)
>> 83 // Sets Ezero for an exponential distribution.
>> 84 //
>> 85 // void SetGradient(G4double)
>> 86 // Sets gradient for a linear distribution.
>> 87 //
>> 88 // void SetInterCept(G4double)
>> 89 // Sets intercept for a linear distribution.
>> 90 //
>> 91 // void UserEnergyHisto(G4ThreeVector)
>> 92 // Allows user to defined a histogram for the energy distribution.
>> 93 //
>> 94 // void ArbEnergyHisto(G4ThreeVector)
>> 95 // Allows the user to define an Arbitrary set of points for the
>> 96 // energy distribution.
>> 97 //
>> 98 // void EpnEnergyHisto(G4ThreeVector)
>> 99 // Allows the user to define an Energy per nucleon histogram.
>> 100 //
>> 101 // void Calculate()
>> 102 // Controls the calculation of Integral PDF for the Cdg and BBody
>> 103 // distributions.
>> 104 //
>> 105 // void InputEnergySpectra(G4bool)
>> 106 // Allows the user to choose between momentum and energy histograms
>> 107 // for user-defined histograms and arbitrary point-wise spectr.
>> 108 // The default is true (energy).
>> 109 //
>> 110 // void InputDifferentialSpectra(G4bool)
>> 111 // Allows the user to choose between integral and differential
>> 112 // distributions when using the arbitrary point-wise option.
>> 113 //
>> 114 // void ArbInterpolate(G4String)
>> 115 // ArbInterpolate allows the user to specify the type of function to
>> 116 // interpolate the Arbitrary points spectrum with.
>> 117 //
>> 118 // void SetBiasRndm (G4SPSRandomGenerator* a)
>> 119 // Sets the biased random number generator
>> 120 //
>> 121 // G4double GenerateOne(G4ParticleDefinition*);
>> 122 // Generate one random energy for the specified particle
>> 123 //
>> 124 // void ReSetHist(G4String);
>> 125 // Re-sets the histogram for user defined distribution
>> 126 //
>> 127 // void SetVerbosity(G4int)
>> 128 // Sets the verbosity level.
>> 129 //
>> 130 ///////////////////////////////////////////////////////////////////////////////
38 131
39 // Author: Fan Lei, QinetiQ ltd. << 132 #ifndef G4SPSEneDistribution_h
40 // Customer: ESA/ESTEC << 133 #define G4SPSEneDistribution_h 1
41 // History: <<
42 // - 05/02/2004, Fan Lei - Created. <<
43 // Based on the G4GeneralParticleSource cl <<
44 // - 26/03/2014, Andrew Green. <<
45 // Modification to use STL vectors instead <<
46 // Also moved to dynamically allocated mem <<
47 // ExpInterpolation() and LogInterpolation <<
48 // - 06/06/2014, Andrea Dotti. <<
49 // For thread safety: this is a shared obj <<
50 // Added mutex to control access to shared <<
51 // in Getters and Setters, mutex is NOT us <<
52 // is assumed that properties are not chan <<
53 // - 24/11/2017, Fan Lei <<
54 // Added cutoff power-law distribution opti <<
55 // to that of the BlackBody one. <<
56 // ------------------------------------------- <<
57 #ifndef G4SPSEneDistribution_hh <<
58 #define G4SPSEneDistribution_hh 1 <<
59 134
60 #include "G4PhysicsFreeVector.hh" << 135 #include "G4PhysicsOrderedFreeVector.hh"
61 #include "G4ParticleMomentum.hh" 136 #include "G4ParticleMomentum.hh"
62 #include "G4ParticleDefinition.hh" 137 #include "G4ParticleDefinition.hh"
63 #include "G4DataInterpolation.hh" 138 #include "G4DataInterpolation.hh"
64 #include "G4Threading.hh" <<
65 #include "G4Cache.hh" <<
66 #include <vector> <<
67 139
>> 140 //
68 #include "G4SPSRandomGenerator.hh" 141 #include "G4SPSRandomGenerator.hh"
69 142
70 class G4SPSEneDistribution 143 class G4SPSEneDistribution
71 { 144 {
72 public: << 145 public:
>> 146 G4SPSEneDistribution ();
>> 147 ~G4SPSEneDistribution ();
>> 148
>> 149 void SetEnergyDisType(G4String);
>> 150 inline G4String GetEnergyDisType() {return EnergyDisType;};
>> 151 void SetEmin(G4double);
>> 152 inline G4double GetEmin() {return Emin;} ;
>> 153 inline G4double GetArbEmin() {return ArbEmin;} ;
>> 154 void SetEmax(G4double);
>> 155 inline G4double GetEmax() {return Emax;} ;
>> 156 inline G4double GetArbEmax() {return ArbEmax;};
>> 157 void SetMonoEnergy(G4double);
>> 158 void SetAlpha(G4double);
>> 159 void SetTemp(G4double);
>> 160 void SetBeamSigmaInE(G4double);
>> 161 void SetEzero(G4double);
>> 162 void SetGradient(G4double);
>> 163 void SetInterCept(G4double);
>> 164 void UserEnergyHisto(G4ThreeVector);
>> 165 void ArbEnergyHisto(G4ThreeVector);
>> 166 void EpnEnergyHisto(G4ThreeVector);
>> 167
>> 168 void InputEnergySpectra(G4bool);
>> 169 void InputDifferentialSpectra(G4bool);
>> 170 void ArbInterpolate(G4String);
>> 171 inline G4String GetIntType() {return IntType;};
>> 172 void Calculate();
>> 173 //
>> 174 void SetBiasRndm(G4SPSRandomGenerator* a) {eneRndm = a; };
>> 175 // method to re-set the histograms
>> 176 void ReSetHist(G4String);
>> 177 // Set the verbosity level.
>> 178 void SetVerbosity(G4int a) {verbosityLevel = a; } ;
>> 179 //x
>> 180 G4double GenerateOne(G4ParticleDefinition*);
>> 181
>> 182 private:
>> 183 void LinearInterpolation();
>> 184 void LogInterpolation();
>> 185 void ExpInterpolation();
>> 186 void SplineInterpolation();
>> 187 void CalculateCdgSpectrum();
>> 188 void CalculateBbodySpectrum();
>> 189
>> 190 // The following methods generate energies according to the spectral
>> 191 // parameters defined above.
>> 192 void GenerateMonoEnergetic();
>> 193 void GenerateLinearEnergies(G4bool);
>> 194 void GeneratePowEnergies(G4bool);
>> 195 void GenerateExpEnergies(G4bool );
>> 196 void GenerateGaussEnergies();
>> 197 void GenerateBremEnergies();
>> 198 void GenerateBbodyEnergies();
>> 199 void GenerateCdgEnergies();
>> 200 void GenUserHistEnergies();
>> 201 void GenEpnHistEnergies();
>> 202 void GenArbPointEnergies();
>> 203 // converts energy per nucleon to energy.
>> 204 void ConvertEPNToEnergy();
>> 205
>> 206 private:
>> 207
>> 208 G4String EnergyDisType; // energy dis type Variable - Mono,Lin,Exp,etc
>> 209 G4double MonoEnergy; //Mono-energteic energy
>> 210 G4double SE ; // Standard deviation for Gaussion distrbution in energy
>> 211 G4double Emin, Emax; // emin and emax
>> 212 G4double alpha, Ezero, Temp; // alpha (pow), E0 (exp) and Temp (bbody,brem)
>> 213 G4double grad, cept; // gradient and intercept for linear spectra
>> 214 G4bool EnergySpec; // true - energy spectra, false - momentum spectra
>> 215 G4bool DiffSpec; // true - differential spec, false integral spec
>> 216 G4bool ApplyRig; // false no rigidity cutoff, true then apply one
>> 217 G4double ERig; // energy of rigidity cutoff
>> 218 G4PhysicsOrderedFreeVector UDefEnergyH; // energy hist data
>> 219 G4PhysicsOrderedFreeVector IPDFEnergyH;
>> 220 G4bool IPDFEnergyExist, IPDFArbExist, Epnflag;
>> 221 G4PhysicsOrderedFreeVector ArbEnergyH; // Arb x,y histogram
>> 222 G4PhysicsOrderedFreeVector IPDFArbEnergyH; // IPDF for Arb
>> 223 G4PhysicsOrderedFreeVector EpnEnergyH;
>> 224 G4double CDGhist[3]; // cumulative histo for cdg
>> 225 G4double BBHist[10001], Bbody_x[10001];
>> 226 G4String IntType; // Interpolation type
>> 227 G4double Arb_grad[1024], Arb_cept[1024]; // grad and cept for 1024 segments
>> 228 G4double Arb_alpha[1024], Arb_Const[1024]; // alpha and constants
>> 229 G4double Arb_ezero[1024]; // ezero
>> 230 G4double ArbEmin, ArbEmax; // Emin and Emax for the whole arb distribution used primarily for debug.
>> 231
>> 232 G4double particle_energy;
>> 233 G4ParticleDefinition* particle_definition;
>> 234
>> 235 G4SPSRandomGenerator* eneRndm;
>> 236
>> 237 // Verbosity
>> 238 G4int verbosityLevel;
>> 239
>> 240 G4PhysicsOrderedFreeVector ZeroPhysVector ; // for re-set only
>> 241
>> 242 G4DataInterpolation *SplineInt; // holds Spline stuff
73 243
74 G4SPSEneDistribution(); <<
75 // Constructor: initializes variables <<
76 ~G4SPSEneDistribution(); <<
77 // Destructor <<
78 <<
79 void SetEnergyDisType(const G4String&); <<
80 // Allows the user to choose the energy <<
81 // The arguments are: Mono (mono-energet <<
82 // Pow (power-law), Exp (exponential), G <<
83 // Brem (bremsstrahlung), BBody (black-b <<
84 // Cdg (cosmic diffuse gamma-ray), User <<
85 // Arb (arbitrary point-wise), Epn (ener <<
86 <<
87 const G4String& GetEnergyDisType(); <<
88 <<
89 void SetEmin(G4double); <<
90 // Sets the minimum energy <<
91 <<
92 G4double GetEmin() const; <<
93 G4double GetArbEmin(); <<
94 <<
95 void SetEmax(G4double); <<
96 // Sets the maximum energy <<
97 <<
98 G4double GetEmax() const; <<
99 G4double GetArbEmax(); <<
100 <<
101 void SetMonoEnergy(G4double); <<
102 // Sets energy for mono-energetic distri <<
103 <<
104 void SetAlpha(G4double); <<
105 // Sets alpha for a power-law distributi <<
106 <<
107 void SetBiasAlpha(G4double); <<
108 <<
109 void SetTemp(G4double); <<
110 // Sets Temperature for a Brem or BBody <<
111 <<
112 void SetBeamSigmaInE(G4double); <<
113 <<
114 void SetEzero(G4double); <<
115 // Sets Ezero for an exponential distrib <<
116 <<
117 void SetGradient(G4double); <<
118 // Sets gradient for a linear distributi <<
119 <<
120 void SetInterCept(G4double); <<
121 // Sets intercept for a linear distribut <<
122 <<
123 void UserEnergyHisto(const G4ThreeVector&) <<
124 // Allows user to defined a histogram fo <<
125 <<
126 void ArbEnergyHisto(const G4ThreeVector&); <<
127 // Allows the user to define an Arbitrar <<
128 // energy distribution <<
129 <<
130 void ArbEnergyHistoFile(const G4String&); <<
131 <<
132 void EpnEnergyHisto(const G4ThreeVector&); <<
133 // Allows the user to define an Energy p <<
134 <<
135 void InputEnergySpectra(G4bool); <<
136 // Allows the user to choose between mom <<
137 // for user-defined histograms and arbit <<
138 // The default is true (energy) <<
139 <<
140 void InputDifferentialSpectra(G4bool); <<
141 // Allows the user to choose between int <<
142 // distributions when using the arbitrar <<
143 <<
144 void ArbInterpolate(const G4String&); <<
145 // Allows the user to specify the type o <<
146 // interpolate the Arbitrary points spec <<
147 <<
148 const G4String& GetIntType(); <<
149 <<
150 void Calculate(); <<
151 // Controls the calculation of Integral <<
152 // distributions <<
153 <<
154 void SetBiasRndm(G4SPSRandomGenerator* a); <<
155 // Sets the biased random number generat <<
156 <<
157 void ReSetHist(const G4String&); <<
158 // Resets the histogram for user defined <<
159 <<
160 void SetVerbosity(G4int a); <<
161 // Sets the verbosity level <<
162 <<
163 G4double GetWeight() const; <<
164 <<
165 G4double GetMonoEnergy(); <<
166 // Mono-energetic energy <<
167 <<
168 G4double GetSE(); <<
169 // Standard deviation for Gaussian distr <<
170 <<
171 G4double Getalpha() const; <<
172 // Alpha (pow) <<
173 <<
174 G4double GetEzero() const; <<
175 // E0 (exp) <<
176 <<
177 G4double GetTemp(); <<
178 // Temp (bbody,brem) <<
179 <<
180 G4double Getgrad() const; <<
181 // Gradient and intercept for linear spe <<
182 <<
183 G4double Getcept() const; <<
184 <<
185 G4PhysicsFreeVector GetUserDefinedEnergyHi <<
186 <<
187 G4PhysicsFreeVector GetArbEnergyHisto(); <<
188 <<
189 G4double GenerateOne(G4ParticleDefinition* <<
190 // Generate one random energy for the s <<
191 <<
192 G4double GetProbability (G4double); <<
193 <<
194 G4double GetArbEneWeight(G4double); <<
195 <<
196 inline void ApplyEnergyWeight(G4bool val) <<
197 inline G4bool IfApplyEnergyWeight() const <<
198 <<
199 private: <<
200 <<
201 void LinearInterpolation(); <<
202 void LogInterpolation(); <<
203 void ExpInterpolation(); <<
204 void SplineInterpolation(); <<
205 void CalculateCdgSpectrum(); <<
206 void CalculateBbodySpectrum(); <<
207 void CalculateCPowSpectrum(); <<
208 <<
209 // The following methods generate energies <<
210 // to the spectral parameters defined abov <<
211 <<
212 void GenerateMonoEnergetic(); <<
213 void GenerateBiasPowEnergies(); <<
214 void GenerateGaussEnergies(); <<
215 void GenerateBremEnergies(); <<
216 void GenerateBbodyEnergies(); <<
217 void GenerateCdgEnergies(); <<
218 void GenUserHistEnergies(); <<
219 void GenEpnHistEnergies(); <<
220 void GenArbPointEnergies(); // NOTE: REQUI <<
221 void GenerateExpEnergies(G4bool); <<
222 void GenerateLinearEnergies(G4bool); <<
223 void GeneratePowEnergies(G4bool); <<
224 void GenerateCPowEnergies(); <<
225 <<
226 void ConvertEPNToEnergy(); <<
227 // Converts energy per nucleon to energy <<
228 <<
229 void BBInitHists(); <<
230 void CPInitHists(); <<
231 <<
232 private: // Non invariant data members beco <<
233 <<
234 G4String EnergyDisType; // energy dis type <<
235 G4double weight; // particle weight //// N <<
236 G4double MonoEnergy; //Mono-energteic ener <<
237 G4double SE; // Standard deviation for Gau <<
238 <<
239 G4double Emin, Emax; // emin and emax //// <<
240 G4double alpha, Ezero;// alpha (pow), E0 ( <<
241 G4double Temp; // Temp (bbody,brem) <<
242 G4double biasalpha; // biased power index <<
243 G4double grad, cept; // gradient and inter <<
244 G4double prob_norm; // normalisation facto <<
245 G4bool Biased = false; // biased to power- <<
246 G4bool EnergySpec = true; // energy spectr <<
247 G4bool DiffSpec = true; // differential sp <<
248 <<
249 G4PhysicsFreeVector UDefEnergyH; // energy <<
250 G4PhysicsFreeVector IPDFEnergyH; <<
251 G4bool IPDFEnergyExist = false, IPDFArbExi <<
252 G4PhysicsFreeVector ArbEnergyH; // Arb x,y <<
253 G4PhysicsFreeVector IPDFArbEnergyH; // IPD <<
254 G4PhysicsFreeVector EpnEnergyH; <<
255 G4double CDGhist[3]; // cumulative histo f <<
256 <<
257 std::vector<G4double>* BBHist = nullptr; <<
258 std::vector<G4double>* Bbody_x = nullptr; <<
259 G4bool BBhistInit = false; <<
260 G4bool BBhistCalcd = false; <<
261 <<
262 // For cutoff power-law <<
263 // <<
264 std::vector<G4double>* CPHist = nullptr; <<
265 std::vector<G4double>* CP_x = nullptr; <<
266 G4bool CPhistInit = false; <<
267 G4bool CPhistCalcd = false; <<
268 <<
269 G4String IntType; // Interpolation type <<
270 G4double* Arb_grad = nullptr; <<
271 G4double* Arb_cept = nullptr; <<
272 G4bool Arb_grad_cept_flag = false; <<
273 G4double* Arb_alpha = nullptr; <<
274 G4double* Arb_Const = nullptr; <<
275 G4bool Arb_alpha_Const_flag = false; <<
276 G4double* Arb_ezero = nullptr; <<
277 G4bool Arb_ezero_flag = false; <<
278 <<
279 G4bool applyEvergyWeight = false; <<
280 <<
281 G4double ArbEmin, ArbEmax; <<
282 // Emin and Emax for the whole arb distr <<
283 <<
284 G4double particle_energy; <<
285 <<
286 G4SPSRandomGenerator* eneRndm = nullptr; <<
287 <<
288 G4int verbosityLevel; <<
289 <<
290 G4PhysicsFreeVector ZeroPhysVector; // for <<
291 <<
292 std::vector<G4DataInterpolation*> SplineIn <<
293 // Holds Spline stuff required for sampl <<
294 G4DataInterpolation* Splinetemp = nullptr; <<
295 // Holds a temp Spline used for calculat <<
296 <<
297 G4Mutex mutex; // protect access to shared <<
298 <<
299 // Thread local data (non-invariant during <<
300 // These are copied from master one at the <<
301 // generation of each event <<
302 // <<
303 struct threadLocal_t <<
304 { <<
305 G4double Emin; <<
306 G4double Emax; <<
307 G4double alpha; <<
308 G4double Ezero; <<
309 G4double grad; <<
310 G4double cept; <<
311 G4ParticleDefinition* particle_definitio <<
312 G4double weight; <<
313 G4double particle_energy; <<
314 }; <<
315 G4Cache<threadLocal_t> threadLocalData; <<
316 }; 244 };
317 245
>> 246
318 #endif 247 #endif
>> 248
>> 249
>> 250
>> 251
319 252