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