| 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 //
>> 27 // $Id: G4IonisParamElm.cc 67044 2013-01-30 08:50:06Z gcosmo $
>> 28 //
>> 29 //
>> 30 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 31 //
26 // 09-07-98, data moved from G4Element. M.Mair 32 // 09-07-98, data moved from G4Element. M.Maire
27 // 22-11-00, tabulation of ionisation potentia << 33 // 22-11-00, tabulation of ionisation potential from
28 // the ICRU Report N#37. V.Ivanchenk 34 // the ICRU Report N#37. V.Ivanchenko
29 // 08-03-01, correct handling of fShellCorrect 35 // 08-03-01, correct handling of fShellCorrectionVector. M.Maire
30 // 17-10-02, Fix excitation energy interpolati 36 // 17-10-02, Fix excitation energy interpolation. V.Ivanchenko
31 // 06-09-04, Update calculated values after an << 37 // 06-09-04, Update calculated values after any change of ionisation
32 // potential change. V.Ivanchenko 38 // potential change. V.Ivanchenko
33 // 29-04-10, Using G4Pow and mean ionisation e 39 // 29-04-10, Using G4Pow and mean ionisation energy from NIST V.Ivanchenko
34 // 27.10.11: new scheme for G4Exception (mma) 40 // 27.10.11: new scheme for G4Exception (mma)
>> 41 //
>> 42 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
35 43
36 #include "G4IonisParamElm.hh" 44 #include "G4IonisParamElm.hh"
37 <<
38 #include "G4NistManager.hh" 45 #include "G4NistManager.hh"
39 #include "G4PhysicalConstants.hh" <<
40 #include "G4Pow.hh" 46 #include "G4Pow.hh"
>> 47 #include "G4PhysicalConstants.hh"
41 #include "G4SystemOfUnits.hh" 48 #include "G4SystemOfUnits.hh"
42 49
43 //....oooOO0OOooo........oooOO0OOooo........oo 50 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
44 51
45 G4IonisParamElm::G4IonisParamElm(G4double Atom 52 G4IonisParamElm::G4IonisParamElm(G4double AtomNumber)
46 { 53 {
47 G4int Z = G4lrint(AtomNumber); << 54 G4int Z = G4int(AtomNumber + 0.5);
48 if (Z < 1) { 55 if (Z < 1) {
49 G4Exception("G4IonisParamElm::G4IonisParam << 56 G4Exception("G4IonisParamElm::G4IonisParamElm()", "mat501", FatalException,
50 "It is not allowed to create an Element << 57 "It is not allowed to create an Element with Z<1");
51 } 58 }
52 G4Pow* g4pow = G4Pow::GetInstance(); 59 G4Pow* g4pow = G4Pow::GetInstance();
53 60
54 // some basic functions of the atomic number 61 // some basic functions of the atomic number
55 fZ = Z; << 62 fZ = Z;
56 fZ3 = g4pow->Z13(Z); << 63 fZ3 = g4pow->Z13(Z);
57 fZZ3 = fZ3 * g4pow->Z13(Z + 1); << 64 fZZ3 = fZ3*g4pow->Z13(Z+1);
58 flogZ3 = g4pow->logZ(Z) / 3.; << 65 flogZ3 = g4pow->logZ(Z)/3.;
59 << 66
60 fMeanExcitationEnergy = G4NistManager::Insta << 67 fMeanExcitationEnergy =
>> 68 G4NistManager::Instance()->GetMeanIonisationEnergy(Z);
61 69
62 // compute parameters for ion transport 70 // compute parameters for ion transport
63 // The aproximation from: 71 // The aproximation from:
64 // J.F.Ziegler, J.P. Biersack, U. Littmark 72 // J.F.Ziegler, J.P. Biersack, U. Littmark
65 // The Stopping and Range of Ions in Matter, 73 // The Stopping and Range of Ions in Matter,
66 // Vol.1, Pergamon Press, 1985 74 // Vol.1, Pergamon Press, 1985
67 // Fast ions or hadrons 75 // Fast ions or hadrons
68 76
69 G4int iz = Z - 1; 77 G4int iz = Z - 1;
70 if (91 < iz) { << 78 if(91 < iz) { iz = 91; }
71 iz = 91; <<
72 } <<
73 79
74 // clang-format off << 80 static G4double vFermi[92] = {
75 static const G4double vFermi[92] = { <<
76 1.0309, 0.15976, 0.59782, 1.0781, 1.0486 81 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 82 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 83 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 84 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 85 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 86 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 87 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, 88 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 89 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} ; 90 1.9578, 1.0257} ;
86 91
87 static const G4double lFactor[92] = { << 92 static G4double lFactor[92] = {
88 1.0, 1.0, 1.1, 1.06, 1.01, 1.03, 1.04, 93 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, 94 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, 95 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, 96 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, 97 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, 98 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, 99 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, 100 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, 101 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} ; 102 1.16, 1.16} ;
98 // clang-format on <<
99 103
100 fVFermi = vFermi[iz]; << 104 fVFermi = vFermi[iz];
101 fLFactor = lFactor[iz]; 105 fLFactor = lFactor[iz];
102 106
103 // obsolete parameters for ionisation 107 // obsolete parameters for ionisation
104 fTau0 = 0.1 * fZ3 * MeV / proton_mass_c2; << 108 fTau0 = 0.1*fZ3*MeV/proton_mass_c2;
105 fTaul = 2. * MeV / proton_mass_c2; << 109 fTaul = 2.*MeV/proton_mass_c2;
106 110
107 // compute the Bethe-Bloch formula for energ 111 // compute the Bethe-Bloch formula for energy = fTaul*particle mass
108 G4double rate = fMeanExcitationEnergy / elec << 112 G4double rate = fMeanExcitationEnergy/electron_mass_c2 ;
109 G4double w = fTaul * (fTaul + 2.); << 113 G4double w = fTaul*(fTaul+2.) ;
110 fBetheBlochLow = (fTaul + 1.) * (fTaul + 1.) << 114 fBetheBlochLow = (fTaul+1.)*(fTaul+1.)*std::log(2.*w/rate)/w - 1. ;
111 fBetheBlochLow = 2. * fZ * twopi_mc2_rcl2 * << 115 fBetheBlochLow = 2.*fZ*twopi_mc2_rcl2*fBetheBlochLow ;
112 << 116
113 fClow = std::sqrt(fTaul) * fBetheBlochLow; << 117 fClow = std::sqrt(fTaul)*fBetheBlochLow;
114 fAlow = 6.458040 * fClow / fTau0; << 118 fAlow = 6.458040 * fClow/fTau0;
115 G4double Taum = 0.035 * fZ3 * MeV / proton_m << 119 G4double Taum = 0.035*fZ3*MeV/proton_mass_c2;
116 fBlow = -3.229020 * fClow / (fTau0 * std::sq << 120 fBlow =-3.229020*fClow/(fTau0*std::sqrt(Taum));
117 121
118 // Shell correction parameterization 122 // Shell correction parameterization
119 fShellCorrectionVector = new G4double[3]; / << 123 fShellCorrectionVector = new G4double[3]; //[3]
120 rate = 0.001 * fMeanExcitationEnergy / eV; << 124 rate = 0.001*fMeanExcitationEnergy/eV;
121 G4double rate2 = rate * rate; << 125 G4double rate2 = rate*rate;
122 fShellCorrectionVector[0] = (0.422377 + 3.85 << 126 /*
123 fShellCorrectionVector[1] = (0.0304043 - 0.1 << 127 fShellCorrectionVector[0] = ( 1.10289e5 + 5.14781e8*rate)*rate2 ;
124 fShellCorrectionVector[2] = (-0.00038106 + 0 << 128 fShellCorrectionVector[1] = ( 7.93805e3 - 2.22565e7*rate)*rate2 ;
>> 129 fShellCorrectionVector[2] = (-9.92256e1 + 2.10823e5*rate)*rate2 ;
>> 130 */
>> 131 fShellCorrectionVector[0] = ( 0.422377 + 3.858019*rate)*rate2 ;
>> 132 fShellCorrectionVector[1] = ( 0.0304043 - 0.1667989*rate)*rate2 ;
>> 133 fShellCorrectionVector[2] = (-0.00038106 + 0.00157955*rate)*rate2 ;
>> 134 }
>> 135
>> 136 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 137
>> 138 // Fake default constructor - sets only member data and allocates memory
>> 139 // for usage restricted to object persistency
>> 140
>> 141 G4IonisParamElm::G4IonisParamElm(__void__&)
>> 142 : fShellCorrectionVector(0)
>> 143 {
>> 144 fZ=fZ3=fZZ3=flogZ3=fTau0=fTaul=fBetheBlochLow=fAlow=fBlow=fClow
>> 145 =fMeanExcitationEnergy=fVFermi=fLFactor=0.0;
>> 146 }
>> 147
>> 148 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 149
>> 150 G4IonisParamElm::~G4IonisParamElm()
>> 151 {
>> 152 if (fShellCorrectionVector) { delete [] fShellCorrectionVector; }
125 } 153 }
126 154
127 //....oooOO0OOooo........oooOO0OOooo........oo 155 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
128 156
129 G4IonisParamElm::~G4IonisParamElm() { delete[] << 157 G4IonisParamElm::G4IonisParamElm(G4IonisParamElm& right)
>> 158 {
>> 159 fShellCorrectionVector = 0;
>> 160 *this = right;
>> 161 }
>> 162
>> 163 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 164
>> 165 G4IonisParamElm& G4IonisParamElm::operator=(const G4IonisParamElm& right)
>> 166 {
>> 167 if (this != &right)
>> 168 {
>> 169 fZ = right.fZ;
>> 170 fZ3 = right.fZ3;
>> 171 fZZ3 = right.fZZ3;
>> 172 flogZ3 = right.flogZ3;
>> 173 fTau0 = right.fTau0;
>> 174 fTaul = right.fTaul;
>> 175 fBetheBlochLow = right.fBetheBlochLow;
>> 176 fAlow = right.fAlow;
>> 177 fBlow = right.fBlow;
>> 178 fClow = right.fClow;
>> 179 fMeanExcitationEnergy = right.fMeanExcitationEnergy;
>> 180 fVFermi = right.fVFermi;
>> 181 fLFactor = right.fLFactor;
>> 182 if (fShellCorrectionVector) { delete [] fShellCorrectionVector; }
>> 183 fShellCorrectionVector = new G4double[3];
>> 184 fShellCorrectionVector[0] = right.fShellCorrectionVector[0];
>> 185 fShellCorrectionVector[1] = right.fShellCorrectionVector[1];
>> 186 fShellCorrectionVector[2] = right.fShellCorrectionVector[2];
>> 187 }
>> 188 return *this;
>> 189 }
>> 190
>> 191 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 192
>> 193 G4int G4IonisParamElm::operator==(const G4IonisParamElm& right) const
>> 194 {
>> 195 return (this == (G4IonisParamElm *) &right);
>> 196 }
>> 197
>> 198 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
>> 199
>> 200 G4int G4IonisParamElm::operator!=(const G4IonisParamElm& right) const
>> 201 {
>> 202 return (this != (G4IonisParamElm *) &right);
>> 203 }
>> 204
>> 205 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
130 206