Geant4-11
G4EmCorrections.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// -------------------------------------------------------------------
28//
29// GEANT4 Class header file
30//
31//
32// File name: G4EmCorrections
33//
34// Author: Vladimir Ivanchenko
35//
36// Creation date: 13.01.2005
37//
38// Modifications:
39// 28.04.2006 General cleanup, add finite size corrections (V.Ivanchenko)
40// 13.05.2006 Add corrections for ion stopping (V.Ivanhcenko)
41// 20.05.2008 Removed Finite Size correction (V.Ivanchenko)
42// 12.09.2008 Added inlined interfaces to effective charge (V.Ivanchenko)
43// 19.04.2012 Fix reproducibility problem (A.Ribon)
44//
45// Class Description:
46//
47// This class provides calculation of EM corrections to ionisation
48//
49
50// -------------------------------------------------------------------
51//
52
53#ifndef G4EmCorrections_h
54#define G4EmCorrections_h 1
55
57
58#include "globals.hh"
60#include "G4Material.hh"
62#include "G4Threading.hh"
63
64class G4VEmModel;
65class G4PhysicsVector;
66class G4IonTable;
69class G4Pow;
70
72{
73
74public:
75
76 explicit G4EmCorrections(G4int verb);
77
79
81 const G4Material*,
82 G4double kineticEnergy,
83 G4double cutEnergy);
84
87 G4double kineticEnergy);
88
90 const G4Material*,
91 G4double kineticEnergy);
92
94 const G4Material*,
95 G4double kineticEnergy);
96
98 const G4Material*,
99 G4double kineticEnergy);
100
102 const G4Material*,
103 G4double kineticEnergy);
104
106 const G4Material*,
107 G4double kineticEnergy);
108
110 const G4Material*,
111 G4double kineticEnergy);
112
114 const G4Material*,
115 G4double kineticEnergy);
116
118 const G4Material*,
119 G4double kineticEnergy);
120
122 const G4Material*,
123 G4double kineticEnergy);
124
126 const G4Material*,
127 G4double kineticEnergy);
128
130 const G4Material*,
131 G4double kineticEnergy);
132
134 const G4Material*,
135 G4double kineticEnergy);
136
138 G4PhysicsVector* dVector);
139
140 void InitialiseForNewRun();
141
142 // effective charge correction using stopping power data
144 const G4Material*,
145 G4double kineticEnergy);
146
147 // effective charge of an ion
149 const G4Material*,
150 G4double kineticEnergy);
151
152 inline
154 const G4Material*,
155 G4double kineticEnergy);
156
157 // ionisation models for ions
158 inline void SetIonisationModels(G4VEmModel* mod1 = nullptr,
159 G4VEmModel* mod2 = nullptr);
160
161 inline G4int GetNumberOfStoppingVectors() const;
162
163 inline void SetVerbose(G4int verb);
164
165 // hide assignment operator
166 G4EmCorrections & operator=(const G4EmCorrections &right) = delete;
168
169private:
170
171 void Initialise();
172
174
176 const G4Material*,
177 G4double kineticEnergy);
178
179 G4double KShell(G4double theta, G4double eta);
180
181 G4double LShell(G4double theta, G4double eta);
182
183 G4int Index(G4double x, const G4double* y, G4int n) const;
184
186 G4double y1, G4double y2) const;
187
189 G4double y1, G4double y2, G4double z11, G4double z21,
190 G4double z12, G4double z22) const;
191
194
197 const G4Material* material = nullptr;
198 const G4Material* curMaterial = nullptr;
200 const G4double* atomDensity = nullptr;
201
203
206
222
224
225 size_t ncouples = 0;
227
232
235
236 // Ion stopping data
240
241 std::vector<G4int> Zion;
242 std::vector<G4int> Aion;
243 std::vector<G4String> materialName;
244
245 std::vector<const G4ParticleDefinition*> ionList;
246
247 std::map< G4int, std::vector<G4double> > thcorr;
248
249 std::vector<const G4Material*> currmat;
250 std::vector<const G4Material*> materialList;
251 std::vector<G4PhysicsVector*> stopData;
252
254
255 static const G4double ZD[11];
256 static const G4double UK[20];
257 static const G4double VK[20];
258 static G4double ZK[20];
259 static const G4double Eta[29];
260 static G4double CK[20][29];
261 static G4double CL[26][28];
262 static const G4double UL[26];
263 static G4double VL[26];
264
268
269#ifdef G4MULTITHREADED
270 static G4Mutex theCorrMutex;
271#endif
272};
273
274inline G4int
276{
277 G4int iddd = n-1;
278 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
279 do {--iddd;} while (iddd>0 && x<y[iddd]);
280 return iddd;
281}
282
284 G4double y1, G4double y2) const
285{
286 return y1 + (y2 - y1)*(xv - x1)/(x2 - x1);
287}
288
290 G4double x1, G4double x2,
291 G4double y1, G4double y2,
292 G4double z11, G4double z21,
293 G4double z12, G4double z22) const
294{
295 return ( z11*(x2-xv)*(y2-yv) + z22*(xv-x1)*(yv-y1) +
296 z12*(x2-xv)*(yv-y1) + z21*(xv-x1)*(y2-yv) )
297 / ((x2-x1)*(y2-y1));
298}
299
300inline void
302{
303 if(nullptr != mod1) { ionLEModel = mod1; }
304 if(nullptr != mod2) { ionHEModel = mod2; }
305}
306
308{
309 return nIons;
310}
311
312inline G4double
314 const G4Material* mat,
315 G4double kineticEnergy)
316{
317 return effCharge.EffectiveCharge(p,mat,kineticEnergy);
318}
319
320inline G4double
322 const G4Material* mat,
323 G4double kineticEnergy)
324{
325 return effCharge.EffectiveChargeSquareRatio(p,mat,kineticEnergy);
326}
327
329{
330 verbose = verb;
331}
332
333//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
334
335#endif
std::vector< const G4Element * > G4ElementVector
std::mutex G4Mutex
Definition: G4Threading.hh:81
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
int G4int
Definition: G4Types.hh:85
const G4int Z[17]
const G4double A[17]
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
void SetIonisationModels(G4VEmModel *mod1=nullptr, G4VEmModel *mod2=nullptr)
static G4double VL[26]
G4double LShellCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4VEmModel * ionLEModel
static G4double CK[20][29]
void AddStoppingData(G4int Z, G4int A, const G4String &materialName, G4PhysicsVector *dVector)
std::vector< const G4ParticleDefinition * > ionList
G4VEmModel * ionHEModel
const G4Material * material
static G4double ZK[20]
G4double MottCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
const G4ParticleDefinition * curParticle
G4double KShellCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
const G4Material * curMaterial
G4PhysicsVector * curVector
G4double IonHighOrderCorrections(const G4ParticleDefinition *, const G4MaterialCutsCouple *, G4double kineticEnergy)
void BuildCorrectionVector()
G4double BarkasCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double Value(G4double xv, G4double x1, G4double x2, G4double y1, G4double y2) const
std::vector< G4String > materialName
G4double ComputeIonCorrections(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double IonBarkasCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4IonTable * ionTable
G4double Value2(G4double xv, G4double yv, G4double x1, G4double x2, G4double y1, G4double y2, G4double z11, G4double z21, G4double z12, G4double z22) const
G4double SpinCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
static const G4double ZD[11]
G4double HighOrderCorrections(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy, G4double cutEnergy)
G4EmCorrections & operator=(const G4EmCorrections &right)=delete
G4double KShell(G4double theta, G4double eta)
std::map< G4int, std::vector< G4double > > thcorr
static const G4double Eta[29]
std::vector< G4PhysicsVector * > stopData
static const G4double UL[26]
void SetVerbose(G4int verb)
static const G4double UK[20]
G4int GetNumberOfStoppingVectors() const
G4int Index(G4double x, const G4double *y, G4int n) const
G4double ShellCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
const G4ElementVector * theElementVector
G4double Bethe(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
const G4ParticleDefinition * particle
G4double EffectiveChargeCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
std::vector< G4int > Zion
G4double ShellCorrectionSTD(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
static G4double CL[26][28]
const G4double * atomDensity
static G4PhysicsFreeVector * sBarkasCorr
void SetupKinematics(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double LShell(G4double theta, G4double eta)
static G4PhysicsFreeVector * sThetaK
G4ionEffectiveCharge effCharge
G4double DensityCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
std::vector< const G4Material * > materialList
G4EmCorrections(const G4EmCorrections &)=delete
G4EmCorrections(G4int verb)
static const G4double VK[20]
static G4PhysicsFreeVector * sThetaL
G4double BlochCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
std::vector< G4int > Aion
std::vector< const G4Material * > currmat
Definition: G4Pow.hh:49
G4double EffectiveCharge(const G4ParticleDefinition *p, const G4Material *material, G4double kineticEnergy)
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *p, const G4Material *material, G4double kineticEnergy)