ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ionEffectiveCharge.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4ionEffectiveCharge.cc
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 file
30 //
31 //
32 // File name: G4ionEffectiveCharge
33 //
34 // Author: Vladimir Ivanchenko
35 //
36 // Creation date: 07.05.2002
37 //
38 // Modifications:
39 // 12.09.2004 Set low energy limit to 1 keV (V.Ivanchenko)
40 // 25.01.2005 Add protection - min Charge 0.1 eplus (V.Ivanchenko)
41 // 28.04.2006 Set upper energy limit to 50 MeV (V.Ivanchenko)
42 // 23.05.2006 Set upper energy limit to Z*10 MeV (V.Ivanchenko)
43 // 15.08.2006 Add protection for not defined material (V.Ivanchenko)
44 // 27-09-2007 Use Fermi energy from material, optimazed formulas (V.Ivanchenko)
45 //
46 
47 // -------------------------------------------------------------------
48 //
49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
50 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
51 
52 #include "G4ionEffectiveCharge.hh"
53 #include "G4PhysicalConstants.hh"
54 #include "G4SystemOfUnits.hh"
55 #include "G4UnitsTable.hh"
56 #include "G4Material.hh"
57 #include "G4NistManager.hh"
58 #include "G4Log.hh"
59 #include "G4Exp.hh"
60 #include "G4Pow.hh"
61 
62 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
63 
64 G4ionEffectiveCharge::G4ionEffectiveCharge()
65 {
66  chargeCorrection = 1.0;
67  energyHighLimit = 20.0*MeV;
68  energyLowLimit = 1.0*keV;
69  energyBohr = 25.*keV;
70  massFactor = amu_c2/(proton_mass_c2*keV);
71  minCharge = 1.0;
72  lastPart = 0;
73  lastMat = 0;
74  lastKinEnergy = 0.0;
75  effCharge = eplus;
76  inveplus = 1.0/CLHEP::eplus;
77  g4calc = G4Pow::GetInstance();
78 }
79 
80 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
81 
82 G4ionEffectiveCharge::~G4ionEffectiveCharge()
83 {}
84 
85 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
86 
87 G4double G4ionEffectiveCharge::EffectiveCharge(const G4ParticleDefinition* p,
88  const G4Material* material,
89  G4double kineticEnergy)
90 {
91  if(p == lastPart && material == lastMat && kineticEnergy == lastKinEnergy)
92  return effCharge;
93 
94  lastPart = p;
95  lastMat = material;
96  lastKinEnergy = kineticEnergy;
97 
98  G4double mass = p->GetPDGMass();
99  G4double charge = p->GetPDGCharge();
100  G4double Zi = charge*inveplus;
101 
102  chargeCorrection = 1.0;
103  effCharge = charge;
104 
105  // The aproximation of ion effective charge from:
106  // J.F.Ziegler, J.P. Biersack, U. Littmark
107  // The Stopping and Range of Ions in Matter,
108  // Vol.1, Pergamon Press, 1985
109  // Fast ions or hadrons
110  G4double reducedEnergy = kineticEnergy * proton_mass_c2/mass ;
111 
112  //G4cout << "e= " << reducedEnergy << " Zi= " << Zi << " "
113  //<< material->GetName() << G4endl;
114 
115  if(Zi < 1.5 || !material || reducedEnergy > Zi*energyHighLimit ) {
116  return charge;
117  }
118  G4double z = material->GetIonisation()->GetZeffective();
119  reducedEnergy = std::max(reducedEnergy,energyLowLimit);
120 
121  // Helium ion case
122  if( Zi < 2.5 ) {
123 
124  static const G4double c[6] =
125  {0.2865,0.1266,-0.001429,0.02402,-0.01135,0.001475};
126 
127  G4double Q = std::max(0.0,G4Log(reducedEnergy*massFactor));
128  G4double x = c[0];
129  G4double y = 1.0;
130  for (G4int i=1; i<6; ++i) {
131  y *= Q;
132  x += y * c[i] ;
133  }
134  G4double ex;
135  if(x < 0.2) { ex = x * (1 - 0.5*x); }
136  else { ex = 1. - G4Exp(-x); }
137 
138  G4double tq = 7.6 - Q;
139  G4double tq2= tq*tq;
140  G4double tt = ( 0.007 + 0.00005 * z );
141  if(tq2 < 0.2) { tt *= (1.0 - tq2 + 0.5*tq2*tq2); }
142  else { tt *= G4Exp(-tq2); }
143 
144  effCharge = charge*(1.0 + tt) * std::sqrt(ex);
145 
146  // Heavy ion case
147  } else {
148 
149  G4double y;
150  G4double zi13 = g4calc->A13(Zi);
151  G4double zi23 = zi13*zi13;
152 
153  // v1 is ion velocity in vF unit
154  G4double eF = material->GetIonisation()->GetFermiEnergy();
155  G4double v1sq = reducedEnergy/eF;
156  G4double vFsq = eF/energyBohr;
157  G4double vF = std::sqrt(eF/energyBohr);
158 
159  // Faster than Fermi velocity
160  if ( v1sq > 1.0 ) {
161  y = vF * std::sqrt(v1sq) * ( 1.0 + 0.2/v1sq ) / zi23 ;
162 
163  // Slower than Fermi velocity
164  } else {
165  y = 0.692308 * vF * (1.0 + 0.666666*v1sq + v1sq*v1sq/15.0) / zi23 ;
166  }
167 
168  G4double q;
169  G4double y3 = std::pow(y, 0.3) ;
170  // G4cout<<"y= "<<y<<" y3= "<<y3<<" v1= "<<v1<<" vF= "<<vF<<G4endl;
171  q = 1.0 - G4Exp( 0.803*y3 - 1.3167*y3*y3 - 0.38157*y - 0.008983*y*y);
172  q = std::max(q, minCharge/Zi);
173 
174  effCharge = q*charge;
175 
176  G4double tq = 7.6 - G4Log(reducedEnergy/keV);
177  G4double tq2= tq*tq;
178  G4double sq = 1.0 + ( 0.18 + 0.0015 * z )*G4Exp(-tq2)/ (Zi*Zi);
179 
180  // G4cout << "sq= " << sq << G4endl;
181 
182  // Screen length according to
183  // J.F.Ziegler and J.M.Manoyan, The stopping of ions in compaunds,
184  // Nucl. Inst. & Meth. in Phys. Res. B35 (1988) 215-228.
185 
186  G4double lambda = 10.0 * vF *g4calc->A23(1.0 - q)/ (zi13 * (6.0 + q));
187 
188  G4double lambda2 = lambda*lambda;
189 
190  G4double xx = (0.5/q - 0.5)*G4Log(1.0 + lambda2)/vFsq;
191 
192  chargeCorrection = sq * (1.0 + xx);
193 
194  }
195  // G4cout << "G4ionEffectiveCharge: charge= " << charge << " q= " << q
196  // << " chargeCor= " << chargeCorrection
197  // << " e(MeV)= " << kineticEnergy/MeV << G4endl;
198  return effCharge;
199 }
200 
201 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
202 
203