Geant4  10.01.p03
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 75582 2013-11-04 12:13:01Z 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 
87  // G4PAIxSection(const G4PAIxSection& right);
88 
89  // Destructor
90 
92 
93  void Initialize(const G4Material* material,
94  G4double maxEnergyTransfer, G4double betaGammaSq, G4SandiaTable*);
95 
96  // Operators
97  // G4PAIxSection& operator=(const G4PAIxSection& right);
98  // G4int operator==(const G4PAIxSection& right)const;
99  // G4int operator!=(const G4PAIxSection& right)const;
100 
101  // Methods
102 
103  // General control functions
104 
105  void ComputeLowEnergyCof(const G4Material* material);
106  void ComputeLowEnergyCof();
107 
108  void InitPAI();
109 
110  void NormShift( G4double betaGammaSq );
111 
112  void SplainPAI( G4double betaGammaSq );
113 
114  // Physical methods
115 
116 
117  G4double RutherfordIntegral( G4int intervalNumber,
118  G4double limitLow,
119  G4double limitHigh );
120 
121  G4double ImPartDielectricConst( G4int intervalNumber,
122  G4double energy );
123 
124  G4double GetPhotonRange( G4double energy );
126 
128 
129  G4double DifPAIxSection( G4int intervalNumber,
130  G4double betaGammaSq );
131 
132  G4double PAIdNdxCerenkov( G4int intervalNumber,
133  G4double betaGammaSq );
134  G4double PAIdNdxMM( G4int intervalNumber,
135  G4double betaGammaSq );
136 
137  G4double PAIdNdxPlasmon( G4int intervalNumber,
138  G4double betaGammaSq );
139 
140  G4double PAIdNdxResonance( G4int intervalNumber,
141  G4double betaGammaSq );
142 
143 
144  void IntegralPAIxSection();
145  void IntegralCerenkov();
146  void IntegralMM();
147  void IntegralPlasmon();
148  void IntegralResonance();
149 
150  G4double SumOverInterval(G4int intervalNumber);
151  G4double SumOverIntervaldEdx(G4int intervalNumber);
152  G4double SumOverInterCerenkov(G4int intervalNumber);
153  G4double SumOverInterMM(G4int intervalNumber);
154  G4double SumOverInterPlasmon(G4int intervalNumber);
155  G4double SumOverInterResonance(G4int intervalNumber);
156 
157  G4double SumOverBorder( G4int intervalNumber,
158  G4double energy );
159  G4double SumOverBorderdEdx( G4int intervalNumber,
160  G4double energy );
161  G4double SumOverBordCerenkov( G4int intervalNumber,
162  G4double energy );
163  G4double SumOverBordMM( G4int intervalNumber,
164  G4double energy );
165  G4double SumOverBordPlasmon( G4int intervalNumber,
166  G4double energy );
167  G4double SumOverBordResonance( G4int intervalNumber,
168  G4double energy );
169 
175 
182 
183  // Inline access functions
184 
186 
187  G4int GetSplineSize() const { return fSplineNumber; }
188 
190 
192 
195  G4double GetPAIdNdxMM(G4int i){ return fdNdxMM[i]; }
198 
201  G4double GetMeanMMLoss() const {return fIntegralMM[0]; }
204 
206 
208 
209  void SetVerbose(G4int v){fVerbose=v;};
210 
211  inline G4double GetPAItable(G4int i,G4int j) const;
212 
213  inline G4double GetLorentzFactor(G4int i) const;
214 
215  inline G4double GetSplineEnergy(G4int i) const;
216 
217  inline G4double GetIntegralPAIxSection(G4int i) const;
218  inline G4double GetIntegralPAIdEdx(G4int i) const;
219  inline G4double GetIntegralCerenkov(G4int i) const;
220  inline G4double GetIntegralMM(G4int i) const;
221  inline G4double GetIntegralPlasmon(G4int i) const;
222  inline G4double GetIntegralResonance(G4int i) const;
223 
224 private :
225 
226  void CallError(G4int i, const G4String& methodName) const;
227 
230 
231 // Local class constants
232 
233 static const G4double fDelta; // energy shift from interval border = 0.001
234 static const G4double fError; // error in lin-log approximation = 0.005
235 
236 static G4int fNumberOfGammas; // = 111;
237 static const G4double fLorentzFactor[112]; // static gamma array
238 
239 static
240 const G4int fRefGammaNumber ; // The number of gamma for creation of spline (15)
241 
242 G4int fIntervalNumber ; // The number of energy intervals
243 G4double fNormalizationCof; // Normalization cof for PhotoAbsorptionXsection
244 
245 // G4double fBetaGammaSq; // (beta*gamma)^2
246 
247  G4int fMaterialIndex; // current material index
248  G4double fDensity; // Current density
249  G4double fElectronDensity; // Current electron (number) density
250  G4double fLowEnergyCof; // Correction cof for low energy region
251  G4int fSplineNumber; // Current size of spline
252  G4int fVerbose; // verbose flag
253 
254 // Arrays of Sandia coefficients
255 
257 
259 
265 
266 static
267 const G4int fMaxSplineSize ; // Max size of output splain arrays = 500
268 
269 
270  G4DataVector fSplineEnergy; // energy points of splain
271  G4DataVector fRePartDielectricConst; // Real part of dielectric const
272  G4DataVector fImPartDielectricConst; // Imaginary part of dielectric const
273  G4DataVector fIntegralTerm; // Integral term in PAI cross section
274  G4DataVector fDifPAIxSection; // Differential PAI cross section
275  G4DataVector fdNdxCerenkov; // dNdx of Cerenkov collisions
276  G4DataVector fdNdxPlasmon; // dNdx of Plasmon collisions
277  G4DataVector fdNdxMM; // dNdx of MM-Cerenkov collisions
278  G4DataVector fdNdxResonance; // dNdx of Resonance collisions
279 
280  G4DataVector fIntegralPAIxSection; // Integral PAI cross section ?
281  G4DataVector fIntegralPAIdEdx; // Integral PAI dEdx ?
282  G4DataVector fIntegralCerenkov; // Integral Cerenkov N>omega ?
283  G4DataVector fIntegralPlasmon; // Integral Plasmon N>omega ?
284  G4DataVector fIntegralMM; // Integral MM N>omega ?
285  G4DataVector fIntegralResonance; // Integral resonance N>omega ?
286 
287 
288 
289 G4double fPAItable[500][112]; // Output array
290 
291 };
292 
294 //
295 
297 {
298  return fPAItable[i][j];
299 }
300 
302 {
303  return fLorentzFactor[j];
304 }
305 
307 {
308  if(i < 1 || i > fSplineNumber) { CallError(i, "GetSplineEnergy"); }
309  return fSplineEnergy[i];
310 }
311 
313 {
314  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIxSection"); }
315  return fIntegralPAIxSection[i];
316 }
317 
319 {
320  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIdEdx"); }
321  return fIntegralPAIdEdx[i];
322 }
323 
325 {
326  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralCerenkov"); }
327  return fIntegralCerenkov[i];
328 }
329 
331 {
332  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralMM"); }
333  return fIntegralMM[i];
334 }
335 
337 {
338  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPlasmon"); }
339  return fIntegralPlasmon[i];
340 }
341 
343 {
344  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralResonance"); }
345  return fIntegralResonance[i];
346 }
347 
348 #endif
349 
350 // ----------------- end of G4PAIxSection header file -------------------
G4DataVector fA1
G4DataVector fIntegralTerm
G4DataVector fIntegralResonance
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 fNormalizationCof
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)
static G4int fNumberOfGammas
G4DataVector fIntegralMM
G4DataVector fdNdxMM
G4double SumOverInterPlasmon(G4int intervalNumber)
G4double GetIntegralPAIdEdx(G4int i) const
G4double SumOverInterval(G4int intervalNumber)
G4double RePartDielectricConst(G4double energy)
G4double fDensity
G4double GetMeanCerenkovLoss() const
G4DataVector fdNdxCerenkov
void ComputeLowEnergyCof()
int G4int
Definition: G4Types.hh:78
void CallError(G4int i, const G4String &methodName) const
static const G4int fMaxSplineSize
G4DataVector fIntegralCerenkov
G4double GetNormalizationCof() const
G4double GetLowEnergyCof() const
G4double fPAItable[500][112]
G4double GetPAIdNdxMM(G4int i)
G4double GetLorentzFactor(G4int i) const
G4double GetResonanceEnergyTransfer()
G4DataVector fRePartDielectricConst
G4DataVector fA3
G4DataVector fdNdxPlasmon
G4double GetPAIdNdxCerenkov(G4int i)
G4double ImPartDielectricConst(G4int intervalNumber, G4double energy)
G4int GetSplineSize() const
G4DataVector fDifPAIxSection
G4double GetIntegralMM(G4int i) const
G4double GetPhotonRange(G4double energy)
G4double SumOverBordMM(G4int intervalNumber, G4double energy)
G4double GetIntegralCerenkov(G4int i) const
G4DataVector fIntegralPAIxSection
G4double GetCerenkovEnergyTransfer()
static const G4int fRefGammaNumber
static const G4double fError
G4OrderedTable * fMatSandiaMatrix
G4int GetNumberOfGammas() const
G4double SumOverInterResonance(G4int intervalNumber)
G4DataVector fA4
G4DataVector fdNdxResonance
G4double GetStepPlasmonLoss(G4double step)
G4double GetPAIdNdxResonance(G4int i)
G4double PAIdNdxMM(G4int intervalNumber, G4double betaGammaSq)
G4DataVector fIntegralPlasmon
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)
static const G4double fDelta
void Initialize(const G4Material *material, G4double maxEnergyTransfer, G4double betaGammaSq, G4SandiaTable *)
G4double SumOverInterCerenkov(G4int intervalNumber)
G4DataVector fEnergyInterval
G4double fLowEnergyCof
G4DataVector fA2
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)
G4PAIxSection & operator=(const G4PAIxSection &right)
G4DataVector fIntegralPAIdEdx
G4DataVector fImPartDielectricConst
G4double GetStepCerenkovLoss(G4double step)
G4double SumOverBorderdEdx(G4int intervalNumber, G4double energy)
G4double DifPAIxSection(G4int intervalNumber, G4double betaGammaSq)
G4DataVector fSplineEnergy
G4double fElectronDensity
G4double GetMeanEnergyLoss() const
static const G4double fLorentzFactor[112]
G4double RutherfordIntegral(G4int intervalNumber, G4double limitLow, G4double limitHigh)
G4double GetSplineEnergy(G4int i) const
G4SandiaTable * fSandia
G4double GetElectronRange(G4double energy)