| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 << 25 //
26 // 09-07-98, data moved from G4Element. M.Mair << 26 //
27 // 22-11-00, tabulation of ionisation potentia << 27 // $Id: G4IonisParamElm.cc,v 1.15 2008/06/03 14:30:10 vnivanch Exp $
28 // the ICRU Report N#37. V.Ivanchenk << 28 // GEANT4 tag $Name: geant4-09-03-patch-02 $
29 // 08-03-01, correct handling of fShellCorrect << 29 //
30 // 17-10-02, Fix excitation energy interpolati << 30 //
31 // 06-09-04, Update calculated values after an << 31 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 32 //
>> 33 // 06-09-04, Update calculated values after any change of ionisation
32 // potential change. V.Ivanchenko 34 // potential change. V.Ivanchenko
33 // 29-04-10, Using G4Pow and mean ionisation e << 35 // 17-10-02, Fix excitation energy interpolation. V.Ivanchenko
34 // 27.10.11: new scheme for G4Exception (mma) << 36 // 08-03-01, correct handling of fShellCorrectionVector. M.Maire
>> 37 // 22-11-00, tabulation of ionisation potential from
>> 38 // the ICRU Report N#37. V.Ivanchenko
>> 39 // 09-07-98, data moved from G4Element. M.Maire
>> 40 //
>> 41 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
35 42
36 #include "G4IonisParamElm.hh" 43 #include "G4IonisParamElm.hh"
37 44
38 #include "G4NistManager.hh" <<
39 #include "G4PhysicalConstants.hh" <<
40 #include "G4Pow.hh" <<
41 #include "G4SystemOfUnits.hh" <<
42 <<
43 //....oooOO0OOooo........oooOO0OOooo........oo 45 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
44 46
45 G4IonisParamElm::G4IonisParamElm(G4double Atom << 47 G4IonisParamElm::G4IonisParamElm(G4double Z)
46 { 48 {
47 G4int Z = G4lrint(AtomNumber); << 49 if (Z < 1.)
48 if (Z < 1) { << 50 G4Exception (" ERROR! It is not allowed to create an Element with Z < 1" );
49 G4Exception("G4IonisParamElm::G4IonisParam << 51
50 "It is not allowed to create an Element <<
51 } <<
52 G4Pow* g4pow = G4Pow::GetInstance(); <<
53 <<
54 // some basic functions of the atomic number 52 // some basic functions of the atomic number
55 fZ = Z; 53 fZ = Z;
56 fZ3 = g4pow->Z13(Z); << 54 fZ3 = std::pow(fZ, 1./3.);
57 fZZ3 = fZ3 * g4pow->Z13(Z + 1); << 55 fZZ3 = std::pow(fZ*(fZ+1.), 1./3.);
58 flogZ3 = g4pow->logZ(Z) / 3.; << 56 flogZ3 = std::log(fZ)/3.;
59 << 57
60 fMeanExcitationEnergy = G4NistManager::Insta << 58 // Parameters for energy loss by ionisation
>> 59
>> 60 // Mean excitation energy
>> 61 // from "Stopping Powers for Electrons and Positrons"
>> 62 // ICRU Report N#37, 1984 (energy in eV)
>> 63 static double exc[100] = {
>> 64 21.8, 20.9, 13.3, 15.9, 15.2, 13.0, 11.7, 11.9, 11.5, 13.7,
>> 65 13.6, 13.0, 12.8, 12.4, 11.5, 11.3, 10.2, 10.4, 10.0, 9.6,
>> 66 10.3, 10.6, 10.7, 10.7, 10.9, 11.0, 11.0, 11.1, 11.1, 11.0,
>> 67 10.8, 10.4, 10.5, 10.2, 9.8, 9.8, 9.8, 9.6, 9.7, 9.8,
>> 68 10.2, 10.1, 10.0, 10.0, 10.0, 10.2, 10.0, 9.8, 10.0, 9.8,
>> 69 9.5, 9.3, 9.3, 8.9, 8.9, 8.8, 8.8, 8.8, 9.1, 9.1,
>> 70 9.2, 9.3, 9.2, 9.2, 9.4, 9.5, 9.7, 9.7, 9.8, 9.8,
>> 71 9.8, 9.8, 9.8, 9.8, 9.8, 9.8, 9.8, 10.1, 10.0, 10.0,
>> 72 10.0, 10.0, 9.9, 9.9, 9.7, 9.2, 9.5, 9.4, 9.4, 9.4,
>> 73 9.6, 9.7, 9.7, 9.8, 9.8, 9.8, 9.8, 9.9, 9.9, 9.9 };
>> 74
>> 75 G4int iz = (G4int)Z - 1 ;
>> 76 if(0 > iz) iz = 0;
>> 77 else if(99 < iz) iz = 99 ;
>> 78 fMeanExcitationEnergy = fZ * exc[iz] * eV ;
61 79
62 // compute parameters for ion transport 80 // compute parameters for ion transport
63 // The aproximation from: 81 // The aproximation from:
64 // J.F.Ziegler, J.P. Biersack, U. Littmark 82 // J.F.Ziegler, J.P. Biersack, U. Littmark
65 // The Stopping and Range of Ions in Matter, 83 // The Stopping and Range of Ions in Matter,
66 // Vol.1, Pergamon Press, 1985 84 // Vol.1, Pergamon Press, 1985
67 // Fast ions or hadrons 85 // Fast ions or hadrons
68 86
69 G4int iz = Z - 1; << 87 if(91 < iz) iz = 91;
70 if (91 < iz) { <<
71 iz = 91; <<
72 } <<
73 88
74 // clang-format off << 89 static G4double vFermi[92] = {
75 static const G4double vFermi[92] = { <<
76 1.0309, 0.15976, 0.59782, 1.0781, 1.0486 90 1.0309, 0.15976, 0.59782, 1.0781, 1.0486, 1.0, 1.058, 0.93942, 0.74562, 0.3424,
77 0.45259, 0.71074, 0.90519, 0.97411, 0.9718 91 0.45259, 0.71074, 0.90519, 0.97411, 0.97184, 0.89852, 0.70827, 0.39816, 0.36552, 0.62712,
78 0.81707, 0.9943, 1.1423, 1.2381, 1.1222 92 0.81707, 0.9943, 1.1423, 1.2381, 1.1222, 0.92705, 1.0047, 1.2, 1.0661, 0.97411,
79 0.84912, 0.95, 1.0903, 1.0429, 0.4971 93 0.84912, 0.95, 1.0903, 1.0429, 0.49715, 0.37755, 0.35211, 0.57801, 0.77773, 1.0207,
80 1.029, 1.2542, 1.122, 1.1241, 1.0882 94 1.029, 1.2542, 1.122, 1.1241, 1.0882, 1.2709, 1.2542, 0.90094, 0.74093, 0.86054,
81 0.93155, 1.0047, 0.55379, 0.43289, 0.3263 95 0.93155, 1.0047, 0.55379, 0.43289, 0.32636, 0.5131, 0.695, 0.72591, 0.71202, 0.67413,
82 0.71418, 0.71453, 0.5911, 0.70263, 0.6804 96 0.71418, 0.71453, 0.5911, 0.70263, 0.68049, 0.68203, 0.68121, 0.68532, 0.68715, 0.61884,
83 0.71801, 0.83048, 1.1222, 1.2381, 1.045, 97 0.71801, 0.83048, 1.1222, 1.2381, 1.045, 1.0733, 1.0953, 1.2381, 1.2879, 0.78654,
84 0.66401, 0.84912, 0.88433, 0.80746, 0.4335 98 0.66401, 0.84912, 0.88433, 0.80746, 0.43357, 0.41923, 0.43638, 0.51464, 0.73087, 0.81065,
85 1.9578, 1.0257} ; 99 1.9578, 1.0257} ;
86 100
87 static const G4double lFactor[92] = { << 101 static G4double lFactor[92] = {
88 1.0, 1.0, 1.1, 1.06, 1.01, 1.03, 1.04, 102 1.0, 1.0, 1.1, 1.06, 1.01, 1.03, 1.04, 0.99, 0.95, 0.9,
89 0.82, 0.81, 0.83, 0.88, 1.0, 0.95, 0.97, 103 0.82, 0.81, 0.83, 0.88, 1.0, 0.95, 0.97, 0.99, 0.98, 0.97,
90 0.98, 0.97, 0.96, 0.93, 0.91, 0.9, 0.88, 104 0.98, 0.97, 0.96, 0.93, 0.91, 0.9, 0.88, 0.9, 0.9, 0.9,
91 0.9, 0.85, 0.9, 0.9, 0.91, 0.92, 0.9, 105 0.9, 0.85, 0.9, 0.9, 0.91, 0.92, 0.9, 0.9, 0.9, 0.9,
92 0.9, 0.88, 0.9, 0.88, 0.88, 0.9, 0.9, 106 0.9, 0.88, 0.9, 0.88, 0.88, 0.9, 0.9, 0.88, 0.9, 0.9,
93 0.9, 0.9, 0.96, 1.2, 0.9, 0.88, 0.88, 107 0.9, 0.9, 0.96, 1.2, 0.9, 0.88, 0.88, 0.85, 0.9, 0.9,
94 0.92, 0.95, 0.99, 1.03, 1.05, 1.07, 1.08, 108 0.92, 0.95, 0.99, 1.03, 1.05, 1.07, 1.08, 1.1, 1.08, 1.08,
95 1.08, 1.08, 1.09, 1.09, 1.1, 1.11, 1.12, 109 1.08, 1.08, 1.09, 1.09, 1.1, 1.11, 1.12, 1.13, 1.14, 1.15,
96 1.17, 1.2, 1.18, 1.17, 1.17, 1.16, 1.16, 110 1.17, 1.2, 1.18, 1.17, 1.17, 1.16, 1.16, 1.16, 1.16, 1.16,
97 1.16, 1.16} ; 111 1.16, 1.16} ;
98 // clang-format on <<
99 112
100 fVFermi = vFermi[iz]; << 113 fVFermi = vFermi[iz];
101 fLFactor = lFactor[iz]; 114 fLFactor = lFactor[iz];
102 115
103 // obsolete parameters for ionisation 116 // obsolete parameters for ionisation
104 fTau0 = 0.1 * fZ3 * MeV / proton_mass_c2; << 117 fTau0 = 0.1*fZ3*MeV/proton_mass_c2;
105 fTaul = 2. * MeV / proton_mass_c2; << 118 fTaul = 2.*MeV/proton_mass_c2;
106 119
107 // compute the Bethe-Bloch formula for energ 120 // compute the Bethe-Bloch formula for energy = fTaul*particle mass
108 G4double rate = fMeanExcitationEnergy / elec << 121 G4double rate = fMeanExcitationEnergy/electron_mass_c2 ;
109 G4double w = fTaul * (fTaul + 2.); << 122 G4double w = fTaul*(fTaul+2.) ;
110 fBetheBlochLow = (fTaul + 1.) * (fTaul + 1.) << 123 fBetheBlochLow = (fTaul+1.)*(fTaul+1.)*std::log(2.*w/rate)/w - 1. ;
111 fBetheBlochLow = 2. * fZ * twopi_mc2_rcl2 * << 124 fBetheBlochLow = 2.*fZ*twopi_mc2_rcl2*fBetheBlochLow ;
112 << 125
113 fClow = std::sqrt(fTaul) * fBetheBlochLow; << 126 fClow = std::sqrt(fTaul)*fBetheBlochLow;
114 fAlow = 6.458040 * fClow / fTau0; << 127 fAlow = 6.458040 * fClow/fTau0;
115 G4double Taum = 0.035 * fZ3 * MeV / proton_m << 128 G4double Taum = 0.035*fZ3*MeV/proton_mass_c2;
116 fBlow = -3.229020 * fClow / (fTau0 * std::sq << 129 fBlow =-3.229020*fClow/(fTau0*std::sqrt(Taum));
117 130
118 // Shell correction parameterization 131 // Shell correction parameterization
119 fShellCorrectionVector = new G4double[3]; / << 132 fShellCorrectionVector = new G4double[3]; //[3]
120 rate = 0.001 * fMeanExcitationEnergy / eV; << 133 rate = 0.001*fMeanExcitationEnergy/eV;
121 G4double rate2 = rate * rate; << 134 G4double rate2 = rate*rate;
122 fShellCorrectionVector[0] = (0.422377 + 3.85 << 135 /*
123 fShellCorrectionVector[1] = (0.0304043 - 0.1 << 136 fShellCorrectionVector[0] = ( 1.10289e5 + 5.14781e8*rate)*rate2 ;
124 fShellCorrectionVector[2] = (-0.00038106 + 0 << 137 fShellCorrectionVector[1] = ( 7.93805e3 - 2.22565e7*rate)*rate2 ;
>> 138 fShellCorrectionVector[2] = (-9.92256e1 + 2.10823e5*rate)*rate2 ;
>> 139 */
>> 140 fShellCorrectionVector[0] = ( 0.422377 + 3.858019*rate)*rate2 ;
>> 141 fShellCorrectionVector[1] = ( 0.0304043 - 0.1667989*rate)*rate2 ;
>> 142 fShellCorrectionVector[2] = (-0.00038106 + 0.00157955*rate)*rate2 ;
>> 143 }
>> 144
>> 145 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 146
>> 147 // Fake default constructor - sets only member data and allocates memory
>> 148 // for usage restricted to object persistency
>> 149
>> 150 G4IonisParamElm::G4IonisParamElm(__void__&)
>> 151 : fShellCorrectionVector(0)
>> 152 {
>> 153 }
>> 154
>> 155 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 156
>> 157 G4IonisParamElm::~G4IonisParamElm()
>> 158 {
>> 159 if (fShellCorrectionVector) delete [] fShellCorrectionVector;
>> 160 }
>> 161
>> 162 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 163
>> 164 G4IonisParamElm::G4IonisParamElm(G4IonisParamElm& right)
>> 165 {
>> 166 *this = right;
>> 167 }
>> 168
>> 169 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 170
>> 171 const G4IonisParamElm& G4IonisParamElm::operator=(const G4IonisParamElm& right)
>> 172 {
>> 173 if (this != &right)
>> 174 {
>> 175 fZ = right.fZ;
>> 176 fZ3 = right.fZ3;
>> 177 fZZ3 = right.fZZ3;
>> 178 flogZ3 = right.flogZ3;
>> 179 fTau0 = right.fTau0;
>> 180 fTaul = right.fTaul;
>> 181 fBetheBlochLow = right.fBetheBlochLow;
>> 182 fAlow = right.fAlow;
>> 183 fBlow = right.fBlow;
>> 184 fClow = right.fClow;
>> 185 fMeanExcitationEnergy = right.fMeanExcitationEnergy;
>> 186 if (fShellCorrectionVector) delete [] fShellCorrectionVector;
>> 187 fShellCorrectionVector = new G4double[3];
>> 188 fShellCorrectionVector[0] = right.fShellCorrectionVector[0];
>> 189 fShellCorrectionVector[1] = right.fShellCorrectionVector[1];
>> 190 fShellCorrectionVector[2] = right.fShellCorrectionVector[2];
>> 191 }
>> 192 return *this;
125 } 193 }
126 194
127 //....oooOO0OOooo........oooOO0OOooo........oo 195 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
128 196
129 G4IonisParamElm::~G4IonisParamElm() { delete[] << 197 G4int G4IonisParamElm::operator==(const G4IonisParamElm& right) const
>> 198 {
>> 199 return (this == (G4IonisParamElm *) &right);
>> 200 }
>> 201
>> 202 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 203
>> 204 G4int G4IonisParamElm::operator!=(const G4IonisParamElm& right) const
>> 205 {
>> 206 return (this != (G4IonisParamElm *) &right);
>> 207 }
>> 208
>> 209 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
130 210