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