Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4PAIxSection.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 // $Id: G4PAIxSection.hh 96934 2016-05-18 09:10:41Z gcosmo $
28 //
29 //
30 // G4PAIxSection.hh -- header file
31 //
32 // GEANT 4 class header file --- Copyright CERN 1995
33 // CERB Geneva Switzerland
34 //
35 // for information related to this code, please, contact
36 // CERN, CN Division, ASD Group
37 //
38 // Preparation of ionizing collision cross section according to Photo Absorption
39 // Ionization (PAI) model for simulation of ionization energy losses in very thin
40 // absorbers. Author: Vladimir.Grichine@cern.ch
41 //
42 // History:
43 //
44 // 28.10.11, V. Ivanchenko: Migration of exceptions to the new design
45 // 19.10.03, V. Grichine: Integral dEdx was added for G4PAIModel class
46 // 13.05.03, V. Grichine: Numerical instability was fixed in SumOverInterval/Border
47 // functions
48 // 10.02.02, V. Grichine: New functions and arrays/gets for Cerenkov and
49 // plasmon collisions dN/dx
50 // 27.10.99, V. Grichine: Bug fixed in constructors, 3rd constructor and
51 // GetStepEnergyLoss(step) were added, fDelta = 0.005
52 // 30.11.97, V. Grichine: 2nd version
53 // 11.06.97, V. Grichine: 1st version
54 
55 #ifndef G4PAIXSECTION_HH
56 #define G4PAIXSECTION_HH
57 
58 #include "G4ios.hh"
59 #include "globals.hh"
60 #include "Randomize.hh"
61 
62 #include"G4SandiaTable.hh"
63 
65 class G4Sandiatable;
66 
67 
69 {
70 public:
71  // Constructors
72  G4PAIxSection();
74 
75  G4PAIxSection( G4int materialIndex,
76  G4double maxEnergyTransfer );
77 
78  G4PAIxSection( G4int materialIndex, // for proton loss table
79  G4double maxEnergyTransfer,
80  G4double betaGammaSq ,
81  G4double** photoAbsCof, G4int intNumber );
82 
83  G4PAIxSection( G4int materialIndex, // test constructor
84  G4double maxEnergyTransfer,
85  G4double betaGammaSq );
86 
88 
89  void Initialize(const G4Material* material, G4double maxEnergyTransfer,
90  G4double betaGammaSq, G4SandiaTable*);
91 
92  // General control functions
93 
94  void ComputeLowEnergyCof(const G4Material* material);
95  void ComputeLowEnergyCof();
96 
97  void InitPAI();
98 
99  void NormShift( G4double betaGammaSq );
100 
101  void SplainPAI( G4double betaGammaSq );
102 
103  // Physical methods
104 
105  G4double RutherfordIntegral( G4int intervalNumber,
106  G4double limitLow,
107  G4double limitHigh );
108 
109  G4double ImPartDielectricConst( G4int intervalNumber,
110  G4double energy );
111 
112  G4double GetPhotonRange( G4double energy );
114 
116 
117  G4double DifPAIxSection( G4int intervalNumber,
118  G4double betaGammaSq );
119 
120  G4double PAIdNdxCerenkov( G4int intervalNumber,
121  G4double betaGammaSq );
122  G4double PAIdNdxMM( G4int intervalNumber,
123  G4double betaGammaSq );
124 
125  G4double PAIdNdxPlasmon( G4int intervalNumber,
126  G4double betaGammaSq );
127 
128  G4double PAIdNdxResonance( G4int intervalNumber,
129  G4double betaGammaSq );
130 
131 
132  void IntegralPAIxSection();
133  void IntegralCerenkov();
134  void IntegralMM();
135  void IntegralPlasmon();
136  void IntegralResonance();
137 
138  G4double SumOverInterval(G4int intervalNumber);
139  G4double SumOverIntervaldEdx(G4int intervalNumber);
140  G4double SumOverInterCerenkov(G4int intervalNumber);
141  G4double SumOverInterMM(G4int intervalNumber);
142  G4double SumOverInterPlasmon(G4int intervalNumber);
143  G4double SumOverInterResonance(G4int intervalNumber);
144 
145  G4double SumOverBorder( G4int intervalNumber,
146  G4double energy );
147  G4double SumOverBorderdEdx( G4int intervalNumber,
148  G4double energy );
149  G4double SumOverBordCerenkov( G4int intervalNumber,
150  G4double energy );
151  G4double SumOverBordMM( G4int intervalNumber,
152  G4double energy );
153  G4double SumOverBordPlasmon( G4int intervalNumber,
154  G4double energy );
155  G4double SumOverBordResonance( G4int intervalNumber,
156  G4double energy );
157 
163 
170 
171  // Inline access functions
172 
173  G4int GetNumberOfGammas() const { return fNumberOfGammas; }
174 
175  G4int GetSplineSize() const { return fSplineNumber; }
176 
177  G4int GetIntervalNumber() const { return fIntervalNumber; }
178 
179  G4double GetEnergyInterval(G4int i){ return fEnergyInterval[i]; }
180 
181  G4double GetDifPAIxSection(G4int i){ return fDifPAIxSection[i]; }
182  G4double GetPAIdNdxCerenkov(G4int i){ return fdNdxCerenkov[i]; }
183  G4double GetPAIdNdxMM(G4int i){ return fdNdxMM[i]; }
184  G4double GetPAIdNdxPlasmon(G4int i){ return fdNdxPlasmon[i]; }
185  G4double GetPAIdNdxResonance(G4int i){ return fdNdxResonance[i]; }
186 
187  G4double GetMeanEnergyLoss() const {return fIntegralPAIxSection[0]; }
188  G4double GetMeanCerenkovLoss() const {return fIntegralCerenkov[0]; }
189  G4double GetMeanMMLoss() const {return fIntegralMM[0]; }
190  G4double GetMeanPlasmonLoss() const {return fIntegralPlasmon[0]; }
191  G4double GetMeanResonanceLoss() const {return fIntegralResonance[0]; }
192 
193  G4double GetNormalizationCof() const { return fNormalizationCof; }
194 
195  G4double GetLowEnergyCof() const { return fLowEnergyCof; }
196 
197  void SetVerbose(G4int v){fVerbose=v;};
198 
199  inline G4double GetPAItable(G4int i,G4int j) const;
200 
201  inline G4double GetLorentzFactor(G4int i) const;
202 
203  inline G4double GetSplineEnergy(G4int i) const;
204 
205  inline G4double GetIntegralPAIxSection(G4int i) const;
206  inline G4double GetIntegralPAIdEdx(G4int i) const;
207  inline G4double GetIntegralCerenkov(G4int i) const;
208  inline G4double GetIntegralMM(G4int i) const;
209  inline G4double GetIntegralPlasmon(G4int i) const;
210  inline G4double GetIntegralResonance(G4int i) const;
211 
212 private :
213 
214  void CallError(G4int i, const G4String& methodName) const;
215 
216  G4PAIxSection & operator=(const G4PAIxSection &right) = delete;
217  G4PAIxSection(const G4PAIxSection&) = delete;
218 
219  // Local class constants
220 
221  static const G4double fDelta; // energy shift from interval border = 0.001
222  static const G4double fError; // error in lin-log approximation = 0.005
223 
224  static G4int fNumberOfGammas; // = 111;
225  static const G4double fLorentzFactor[112]; // static gamma array
226 
227  static
228  const G4int fRefGammaNumber; //The number of gamma for creation of spline (15)
229 
230  G4int fIntervalNumber ; // The number of energy intervals
231  G4double fNormalizationCof; // Normalization cof for PhotoAbsorptionXsection
232 
233  // G4double fBetaGammaSq; // (beta*gamma)^2
234 
235  G4int fMaterialIndex; // current material index
236  G4double fDensity; // Current density
237  G4double fElectronDensity; // Current electron (number) density
238  G4double fLowEnergyCof; // Correction cof for low energy region
239  G4int fSplineNumber; // Current size of spline
240  G4int fVerbose; // verbose flag
241 
242  // Arrays of Sandia coefficients
243 
244  G4OrderedTable* fMatSandiaMatrix;
245 
246  G4SandiaTable* fSandia;
247 
248  G4DataVector fEnergyInterval;
249  G4DataVector fA1;
250  G4DataVector fA2;
251  G4DataVector fA3;
252  G4DataVector fA4;
253 
254  static
255  const G4int fMaxSplineSize ; // Max size of output splain arrays = 500
256 
257  G4DataVector fSplineEnergy; // energy points of splain
258  G4DataVector fRePartDielectricConst; // Real part of dielectric const
259  G4DataVector fImPartDielectricConst; // Imaginary part of dielectric const
260  G4DataVector fIntegralTerm; // Integral term in PAI cross section
261  G4DataVector fDifPAIxSection; // Differential PAI cross section
262  G4DataVector fdNdxCerenkov; // dNdx of Cerenkov collisions
263  G4DataVector fdNdxPlasmon; // dNdx of Plasmon collisions
264  G4DataVector fdNdxMM; // dNdx of MM-Cerenkov collisions
265  G4DataVector fdNdxResonance; // dNdx of Resonance collisions
266 
267  G4DataVector fIntegralPAIxSection; // Integral PAI cross section ?
268  G4DataVector fIntegralPAIdEdx; // Integral PAI dEdx ?
269  G4DataVector fIntegralCerenkov; // Integral Cerenkov N>omega ?
270  G4DataVector fIntegralPlasmon; // Integral Plasmon N>omega ?
271  G4DataVector fIntegralMM; // Integral MM N>omega ?
272  G4DataVector fIntegralResonance; // Integral resonance N>omega ?
273 
274  G4double fPAItable[500][112]; // Output array
275 
276 };
277 
279 //
280 
282 {
283  return fPAItable[i][j];
284 }
285 
287 {
288  return fLorentzFactor[j];
289 }
290 
292 {
293  if(i < 1 || i > fSplineNumber) { CallError(i, "GetSplineEnergy"); }
294  return fSplineEnergy[i];
295 }
296 
298 {
299  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIxSection"); }
300  return fIntegralPAIxSection[i];
301 }
302 
304 {
305  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIdEdx"); }
306  return fIntegralPAIdEdx[i];
307 }
308 
310 {
311  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralCerenkov"); }
312  return fIntegralCerenkov[i];
313 }
314 
316 {
317  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralMM"); }
318  return fIntegralMM[i];
319 }
320 
322 {
323  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPlasmon"); }
324  return fIntegralPlasmon[i];
325 }
326 
328 {
329  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralResonance"); }
330  return fIntegralResonance[i];
331 }
332 
333 #endif
334 
335 // ----------------- end of G4PAIxSection header file -------------------
void IntegralPAIxSection()
void SplainPAI(G4double betaGammaSq)
G4double GetRutherfordEnergyTransfer()
G4double PAIdNdxResonance(G4int intervalNumber, G4double betaGammaSq)
void IntegralCerenkov()
G4double GetPAIdNdxPlasmon(G4int i)
G4double GetPlasmonEnergyTransfer()
void IntegralPlasmon()
void IntegralResonance()
G4double GetStepResonanceLoss(G4double step)
G4double PAIdNdxCerenkov(G4int intervalNumber, G4double betaGammaSq)
G4double GetDifPAIxSection(G4int i)
G4double SumOverBordPlasmon(G4int intervalNumber, G4double energy)
G4double GetMeanMMLoss() const
G4double SumOverInterMM(G4int intervalNumber)
G4double GetIntegralPlasmon(G4int i) const
void SetVerbose(G4int v)
void NormShift(G4double betaGammaSq)
G4double GetIntegralPAIxSection(G4int i) const
G4double PAIdNdxPlasmon(G4int intervalNumber, G4double betaGammaSq)
G4double SumOverBordResonance(G4int intervalNumber, G4double energy)
G4double SumOverInterPlasmon(G4int intervalNumber)
G4double GetIntegralPAIdEdx(G4int i) const
G4double SumOverInterval(G4int intervalNumber)
G4double RePartDielectricConst(G4double energy)
G4double GetMeanCerenkovLoss() const
void ComputeLowEnergyCof()
int G4int
Definition: G4Types.hh:78
G4double GetNormalizationCof() const
G4double GetLowEnergyCof() const
G4double GetPAIdNdxMM(G4int i)
string material
Definition: eplot.py:19
G4double GetLorentzFactor(G4int i) const
G4double GetResonanceEnergyTransfer()
G4double GetPAIdNdxCerenkov(G4int i)
G4double ImPartDielectricConst(G4int intervalNumber, G4double energy)
G4int GetSplineSize() const
G4double GetIntegralMM(G4int i) const
G4double GetPhotonRange(G4double energy)
G4double SumOverBordMM(G4int intervalNumber, G4double energy)
G4double GetIntegralCerenkov(G4int i) const
G4double GetCerenkovEnergyTransfer()
G4int GetNumberOfGammas() const
G4double SumOverInterResonance(G4int intervalNumber)
tuple v
Definition: test.py:18
G4double GetStepPlasmonLoss(G4double step)
G4double GetPAIdNdxResonance(G4int i)
G4double PAIdNdxMM(G4int intervalNumber, G4double betaGammaSq)
G4double GetStepEnergyLoss(G4double step)
G4int GetIntervalNumber() const
G4double GetMeanResonanceLoss() const
G4double SumOverBordCerenkov(G4int intervalNumber, G4double energy)
G4double GetMMEnergyTransfer()
G4double SumOverBorder(G4int intervalNumber, G4double energy)
G4double energy(const ThreeVector &p, const G4double m)
void Initialize(const G4Material *material, G4double maxEnergyTransfer, G4double betaGammaSq, G4SandiaTable *)
G4double SumOverInterCerenkov(G4int intervalNumber)
G4double GetIntegralResonance(G4int i) const
G4double GetEnergyInterval(G4int i)
G4double GetMeanPlasmonLoss() const
G4double GetPAItable(G4int i, G4int j) const
double G4double
Definition: G4Types.hh:76
G4double GetEnergyTransfer()
G4double GetStepMMLoss(G4double step)
G4double SumOverIntervaldEdx(G4int intervalNumber)
G4double GetStepCerenkovLoss(G4double step)
G4double SumOverBorderdEdx(G4int intervalNumber, G4double energy)
G4double DifPAIxSection(G4int intervalNumber, G4double betaGammaSq)
G4double GetMeanEnergyLoss() const
G4double RutherfordIntegral(G4int intervalNumber, G4double limitLow, G4double limitHigh)
G4double GetSplineEnergy(G4int i) const
G4double GetElectronRange(G4double energy)