Geant4_10
G4SPSEneDistribution.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 //
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 //
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 //
46 //
47 //
48 // Class Description:
49 //
50 // To generate the energy of a primary vertex according to the defined distribution
51 //
53 //
54 // MEMBER FUNCTIONS
55 // ----------------
56 //
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 //
134 
135 #ifndef G4SPSEneDistribution_h
136 #define G4SPSEneDistribution_h 1
137 
139 #include "G4ParticleMomentum.hh"
140 #include "G4ParticleDefinition.hh"
141 #include "G4DataInterpolation.hh"
142 
143 //
144 #include "G4SPSRandomGenerator.hh"
145 
147 public:
150 
153  return EnergyDisType;
154  }
155  ;
156  void SetEmin(G4double);
157  inline G4double GetEmin() {
158  return Emin;
159  }
160  ;
161  inline G4double GetArbEmin() {
162  return ArbEmin;
163  }
164  ;
165  void SetEmax(G4double);
166  inline G4double GetEmax() {
167  return Emax;
168  }
169  ;
170  inline G4double GetArbEmax() {
171  return ArbEmax;
172  }
173  ;
174  void SetMonoEnergy(G4double);
175  void SetAlpha(G4double);
176  void SetBiasAlpha(G4double);
177  void SetTemp(G4double);
179  void SetEzero(G4double);
180  void SetGradient(G4double);
181  void SetInterCept(G4double);
186 
189  void ArbInterpolate(G4String);
190  inline G4String GetIntType() {
191  return IntType;
192  }
193  ;
194  void Calculate();
195  //
197  eneRndm = a;
198  }
199  ;
200  // method to re-set the histograms
201  void ReSetHist(G4String);
202  // Set the verbosity level.
204  verbosityLevel = a;
205  }
206  ;
207  //x
209  return weight;
210  }
211 
213  return MonoEnergy;
214  }
215  ; //Mono-energteic energy
217  return SE;
218  }
219  ; // Standard deviation for Gaussion distrbution in energy
221  return alpha;
222  }
223  ; // alpha (pow)
225  return Ezero;
226  }
227  ; // E0 (exp)
229  return Temp;
230  }
231  ; // Temp (bbody,brem)
233  return grad;
234  }
235  ; // gradient and intercept for linear spectra
237  return cept;
238  }
239  ;
240 
242  return UDefEnergyH;
243  }
244  ;
246  return ArbEnergyH;
247  }
248  ;
249 
252 
253 
254 private:
255  void LinearInterpolation();
256  void LogInterpolation();
257  void ExpInterpolation();
258  void SplineInterpolation();
259  void CalculateCdgSpectrum();
260  void CalculateBbodySpectrum();
261 
262  // The following methods generate energies according to the spectral
263  // parameters defined above.
264  void GenerateMonoEnergetic();
265  void GenerateLinearEnergies(G4bool);
266  void GeneratePowEnergies(G4bool);
267  void GenerateBiasPowEnergies();
268  void GenerateExpEnergies(G4bool);
269  void GenerateGaussEnergies();
270  void GenerateBremEnergies();
271  void GenerateBbodyEnergies();
272  void GenerateCdgEnergies();
273  void GenUserHistEnergies();
274  void GenEpnHistEnergies();
275  void GenArbPointEnergies();
276  // converts energy per nucleon to energy.
277  void ConvertEPNToEnergy();
278 
279 
280 private:
281 
282  G4String EnergyDisType; // energy dis type Variable - Mono,Lin,Exp,etc
283  G4double weight; // particle weight
284  G4double MonoEnergy; //Mono-energteic energy
285  G4double SE; // Standard deviation for Gaussion distrbution in energy
286  G4double Emin, Emax; // emin and emax
287  G4double alpha, Ezero, Temp; // alpha (pow), E0 (exp) and Temp (bbody,brem)
288  G4double biasalpha; // biased power index
289  G4double grad, cept; // gradient and intercept for linear spectra
290  G4double prob_norm; // normalisation factor use in calculate the probability
291  G4bool Biased; // true - biased to power-law
292  G4bool EnergySpec; // true - energy spectra, false - momentum spectra
293  G4bool DiffSpec; // true - differential spec, false integral spec
294  //G4bool ApplyRig; // false no rigidity cutoff, true then apply one
295  //G4double ERig; // energy of rigidity cutoff
296  G4PhysicsOrderedFreeVector UDefEnergyH; // energy hist data
297  G4PhysicsOrderedFreeVector IPDFEnergyH;
298  G4bool IPDFEnergyExist, IPDFArbExist, Epnflag;
299  G4PhysicsOrderedFreeVector ArbEnergyH; // Arb x,y histogram
300  G4PhysicsOrderedFreeVector IPDFArbEnergyH; // IPDF for Arb
301  G4PhysicsOrderedFreeVector EpnEnergyH;
302  G4double CDGhist[3]; // cumulative histo for cdg
303  G4double BBHist[10001], Bbody_x[10001];
304  G4String IntType; // Interpolation type
305  G4double Arb_grad[1024], Arb_cept[1024]; // grad and cept for 1024 segments
306  G4double Arb_alpha[1024], Arb_Const[1024]; // alpha and constants
307  G4double Arb_ezero[1024]; // ezero
308  G4double ArbEmin, ArbEmax; // Emin and Emax for the whole arb distribution used primarily for debug.
309 
310  G4double particle_energy;
311  G4ParticleDefinition* particle_definition;
312 
313  G4SPSRandomGenerator* eneRndm;
314 
315  // Verbosity
316  G4int verbosityLevel;
317 
318  G4PhysicsOrderedFreeVector ZeroPhysVector; // for re-set only
319 
320  G4DataInterpolation *SplineInt[1024]; // holds Spline stuff required for sampling
321  G4DataInterpolation *Splinetemp; // holds a temp Spline used for calculating area
322 
323 };
324 
325 #endif
326 
tuple a
Definition: test.py:11
void ArbEnergyHisto(G4ThreeVector)
G4PhysicsOrderedFreeVector GetUserDefinedEnergyHisto()
int G4int
Definition: G4Types.hh:78
G4PhysicsOrderedFreeVector GetArbEnergyHisto()
void SetEnergyDisType(G4String)
bool G4bool
Definition: G4Types.hh:79
G4double GetProbability(G4double)
G4double GenerateOne(G4ParticleDefinition *)
void SetBiasRndm(G4SPSRandomGenerator *a)
void EpnEnergyHisto(G4ThreeVector)
void UserEnergyHisto(G4ThreeVector)
void InputDifferentialSpectra(G4bool)
void ArbEnergyHistoFile(G4String)
double G4double
Definition: G4Types.hh:76