| 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 // 26 //
27 // ------------------------------------------- 27 // -------------------------------------------------------------------
28 // 28 //
29 // GEANT4 Class file 29 // GEANT4 Class file
30 // 30 //
31 // 31 //
32 // File name: G4hParametrisedLossModel 32 // File name: G4hParametrisedLossModel
33 // 33 //
34 // Author: V.Ivanchenko (Vladimir.Ivanc 34 // Author: V.Ivanchenko (Vladimir.Ivanchenko@cern.ch)
35 // << 35 //
36 // Creation date: 20 July 2000 36 // Creation date: 20 July 2000
37 // 37 //
38 // Modifications: << 38 // Modifications:
39 // 20/07/2000 V.Ivanchenko First implementati 39 // 20/07/2000 V.Ivanchenko First implementation
40 // 18/08/2000 V.Ivanchenko TRIM85 model is ad 40 // 18/08/2000 V.Ivanchenko TRIM85 model is added
41 // 03/10/2000 V.Ivanchenko CodeWizard clean u 41 // 03/10/2000 V.Ivanchenko CodeWizard clean up
42 // 10/05/2001 V.Ivanchenko Clean up againist 42 // 10/05/2001 V.Ivanchenko Clean up againist Linux compilation with -Wall
43 // 30/12/2003 V.Ivanchenko SRIM2003 model is 43 // 30/12/2003 V.Ivanchenko SRIM2003 model is added
44 // 07/05/2004 V.Ivanchenko Fix Graphite probl 44 // 07/05/2004 V.Ivanchenko Fix Graphite problem, add QAO model
45 // 45 //
46 // Class Description: 46 // Class Description:
47 // 47 //
48 // Low energy protons/ions electronic stopping 48 // Low energy protons/ions electronic stopping power parametrisation
49 // 49 //
50 // Class Description: End 50 // Class Description: End
51 // 51 //
52 // ------------------------------------------- 52 // -------------------------------------------------------------------
53 // 53 //
54 //....oooOO0OOooo........oooOO0OOooo........oo 54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
55 //....oooOO0OOooo........oooOO0OOooo........oo 55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
56 56
57 #include "G4hParametrisedLossModel.hh" 57 #include "G4hParametrisedLossModel.hh"
58 <<
59 #include "globals.hh" <<
60 #include "G4PhysicalConstants.hh" <<
61 #include "G4SystemOfUnits.hh" <<
62 #include "G4UnitsTable.hh" 58 #include "G4UnitsTable.hh"
>> 59 #include "globals.hh"
>> 60 #include "G4hZiegler1977p.hh"
>> 61 #include "G4hZiegler1977He.hh"
63 #include "G4hZiegler1985p.hh" 62 #include "G4hZiegler1985p.hh"
>> 63 #include "G4hSRIM2000p.hh"
>> 64 //#include "G4hQAOModel.hh"
64 #include "G4hICRU49p.hh" 65 #include "G4hICRU49p.hh"
65 #include "G4hICRU49He.hh" 66 #include "G4hICRU49He.hh"
66 #include "G4DynamicParticle.hh" 67 #include "G4DynamicParticle.hh"
67 #include "G4ParticleDefinition.hh" 68 #include "G4ParticleDefinition.hh"
68 #include "G4ElementVector.hh" 69 #include "G4ElementVector.hh"
69 #include "G4Material.hh" 70 #include "G4Material.hh"
70 #include "G4Exp.hh" <<
71 71
72 //....oooOO0OOooo........oooOO0OOooo........oo 72 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
73 73
74 G4hParametrisedLossModel::G4hParametrisedLossM 74 G4hParametrisedLossModel::G4hParametrisedLossModel(const G4String& name)
75 :G4VLowEnergyModel(name), modelName(name) 75 :G4VLowEnergyModel(name), modelName(name)
76 { 76 {
77 InitializeMe(); 77 InitializeMe();
78 } 78 }
79 79
80 //....oooOO0OOooo........oooOO0OOooo........oo 80 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
81 81
82 void G4hParametrisedLossModel::InitializeMe() 82 void G4hParametrisedLossModel::InitializeMe()
83 { 83 {
84 expStopPower125 = 0.; <<
85 <<
86 theZieglerFactor = eV*cm2*1.0e-15 ; 84 theZieglerFactor = eV*cm2*1.0e-15 ;
87 85
88 // Registration of parametrisation models 86 // Registration of parametrisation models
89 const G4String& blank(" "); << 87 G4String blank = G4String(" ") ;
90 const G4String& ir49p("ICRU_R49p"); << 88 G4String zi77p = G4String("Ziegler1977p") ;
91 const G4String& ir49He("ICRU_R49He"); << 89 G4String zi77He = G4String("Ziegler1977He") ;
92 const G4String& zi85p("Ziegler1985p"); << 90 G4String ir49p = G4String("ICRU_R49p") ;
93 if(zi85p == modelName) { << 91 G4String ir49He = G4String("ICRU_R49He") ;
>> 92 G4String zi85p = G4String("Ziegler1985p") ;
>> 93 G4String zi00p = G4String("SRIM2000p") ;
>> 94 G4String qao = G4String("QAO") ;
>> 95 if(zi77p == modelName) {
>> 96 eStopingPowerTable = new G4hZiegler1977p();
>> 97 highEnergyLimit = 100.0*MeV;
>> 98 lowEnergyLimit = 1.0*keV;
>> 99
>> 100 } else if(zi77He == modelName) {
>> 101 eStopingPowerTable = new G4hZiegler1977He();
>> 102 highEnergyLimit = 10.0*MeV/4.0;
>> 103 lowEnergyLimit = 1.0*keV/4.0;
>> 104
>> 105 } else if(zi85p == modelName) {
94 eStopingPowerTable = new G4hZiegler1985p 106 eStopingPowerTable = new G4hZiegler1985p();
95 highEnergyLimit = 100.0*MeV; 107 highEnergyLimit = 100.0*MeV;
96 lowEnergyLimit = 1.0*keV; 108 lowEnergyLimit = 1.0*keV;
97 109
>> 110 } else if(zi00p == modelName ) {
>> 111 eStopingPowerTable = new G4hSRIM2000p();
>> 112 highEnergyLimit = 100.0*MeV;
>> 113 lowEnergyLimit = 1.0*keV;
>> 114
98 } else if(ir49p == modelName || blank == mod 115 } else if(ir49p == modelName || blank == modelName) {
99 eStopingPowerTable = new G4hICRU49p(); 116 eStopingPowerTable = new G4hICRU49p();
100 highEnergyLimit = 2.0*MeV; 117 highEnergyLimit = 2.0*MeV;
101 lowEnergyLimit = 1.0*keV; 118 lowEnergyLimit = 1.0*keV;
102 119
103 } else if(ir49He == modelName) { 120 } else if(ir49He == modelName) {
104 eStopingPowerTable = new G4hICRU49He(); 121 eStopingPowerTable = new G4hICRU49He();
105 highEnergyLimit = 10.0*MeV/4.0; 122 highEnergyLimit = 10.0*MeV/4.0;
106 lowEnergyLimit = 1.0*keV/4.0; 123 lowEnergyLimit = 1.0*keV/4.0;
107 << 124 /*
>> 125 } else if(qao == modelName) {
>> 126 eStopingPowerTable = new G4hQAOModel();
>> 127 highEnergyLimit = 2.0*MeV;
>> 128 lowEnergyLimit = 5.0*keV;
>> 129 */
108 } else { 130 } else {
109 eStopingPowerTable = new G4hICRU49p(); 131 eStopingPowerTable = new G4hICRU49p();
110 highEnergyLimit = 2.0*MeV; 132 highEnergyLimit = 2.0*MeV;
111 lowEnergyLimit = 1.0*keV; 133 lowEnergyLimit = 1.0*keV;
112 G4cout << "G4hParametrisedLossModel Warn << 134 G4cout << "G4hParametrisedLossModel Warning: <" << modelName
113 << "> is unknown - default <" 135 << "> is unknown - default <"
114 << ir49p << ">" << " is used for 136 << ir49p << ">" << " is used for Electronic Stopping"
115 << G4endl; 137 << G4endl;
116 modelName = ir49p; 138 modelName = ir49p;
117 } 139 }
118 /* 140 /*
119 G4cout << "G4hParametrisedLossModel: the 141 G4cout << "G4hParametrisedLossModel: the model <"
120 << modelName << ">" << " is used 142 << modelName << ">" << " is used for Electronic Stopping"
121 << G4endl; 143 << G4endl;
122 */ 144 */
123 } 145 }
124 146
125 //....oooOO0OOooo........oooOO0OOooo........oo 147 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
126 148
127 G4hParametrisedLossModel::~G4hParametrisedLoss 149 G4hParametrisedLossModel::~G4hParametrisedLossModel()
128 { 150 {
129 delete eStopingPowerTable; 151 delete eStopingPowerTable;
130 } 152 }
131 153
132 //....oooOO0OOooo........oooOO0OOooo........oo 154 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
133 155
134 G4double G4hParametrisedLossModel::TheValue(co 156 G4double G4hParametrisedLossModel::TheValue(const G4DynamicParticle* particle,
135 const G4Material* material) 157 const G4Material* material)
136 { 158 {
137 G4double scaledEnergy = (particle->GetKineti 159 G4double scaledEnergy = (particle->GetKineticEnergy())
138 * proton_mass_c2/(part 160 * proton_mass_c2/(particle->GetMass());
139 G4double factor = theZieglerFactor; 161 G4double factor = theZieglerFactor;
140 if (scaledEnergy < lowEnergyLimit) { 162 if (scaledEnergy < lowEnergyLimit) {
141 if (modelName != "QAO") factor *= std::sqr 163 if (modelName != "QAO") factor *= std::sqrt(scaledEnergy/lowEnergyLimit);
142 scaledEnergy = lowEnergyLimit; 164 scaledEnergy = lowEnergyLimit;
143 } 165 }
144 G4double eloss = StoppingPower(material,scal 166 G4double eloss = StoppingPower(material,scaledEnergy) * factor;
145 167
146 return eloss; 168 return eloss;
147 } 169 }
148 170
149 //....oooOO0OOooo........oooOO0OOooo........oo 171 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
150 172
151 G4double G4hParametrisedLossModel::TheValue(co 173 G4double G4hParametrisedLossModel::TheValue(const G4ParticleDefinition* aParticle,
152 const G4Material* material, 174 const G4Material* material,
153 G4double kineticEnergy) 175 G4double kineticEnergy)
154 { 176 {
155 G4double scaledEnergy = kineticEnergy 177 G4double scaledEnergy = kineticEnergy
156 * proton_mass_c2/(aPar 178 * proton_mass_c2/(aParticle->GetPDGMass());
157 179
158 G4double factor = theZieglerFactor; 180 G4double factor = theZieglerFactor;
159 if (scaledEnergy < lowEnergyLimit) { 181 if (scaledEnergy < lowEnergyLimit) {
160 if (modelName != "QAO") factor *= std::sqr 182 if (modelName != "QAO") factor *= std::sqrt(scaledEnergy/lowEnergyLimit);
161 scaledEnergy = lowEnergyLimit; 183 scaledEnergy = lowEnergyLimit;
162 } 184 }
163 G4double eloss = StoppingPower(material,scal 185 G4double eloss = StoppingPower(material,scaledEnergy) * factor;
164 186
165 return eloss; 187 return eloss;
166 } 188 }
167 189
168 //....oooOO0OOooo........oooOO0OOooo........oo 190 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
169 191
170 G4double G4hParametrisedLossModel::LowEnergyLi 192 G4double G4hParametrisedLossModel::LowEnergyLimit(const G4ParticleDefinition* ,
171 const G4Material*) const 193 const G4Material*) const
172 { 194 {
173 return lowEnergyLimit; 195 return lowEnergyLimit;
174 } 196 }
175 197
176 //....oooOO0OOooo........oooOO0OOooo........oo 198 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
177 199
178 G4double G4hParametrisedLossModel::HighEnergyL 200 G4double G4hParametrisedLossModel::HighEnergyLimit(const G4ParticleDefinition* ,
179 const G4Material*) const 201 const G4Material*) const
180 { 202 {
181 return highEnergyLimit; 203 return highEnergyLimit;
182 } 204 }
183 //....oooOO0OOooo........oooOO0OOooo........oo 205 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
184 206
185 G4double G4hParametrisedLossModel::LowEnergyLi 207 G4double G4hParametrisedLossModel::LowEnergyLimit(const G4ParticleDefinition* ) const
186 { 208 {
187 return lowEnergyLimit; 209 return lowEnergyLimit;
188 } 210 }
189 211
190 //....oooOO0OOooo........oooOO0OOooo........oo 212 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
191 213
192 G4double G4hParametrisedLossModel::HighEnergyL 214 G4double G4hParametrisedLossModel::HighEnergyLimit(const G4ParticleDefinition* ) const
193 { 215 {
194 return highEnergyLimit; 216 return highEnergyLimit;
195 } 217 }
196 218
197 //....oooOO0OOooo........oooOO0OOooo........oo 219 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
198 220
199 G4bool G4hParametrisedLossModel::IsInCharge(co 221 G4bool G4hParametrisedLossModel::IsInCharge(const G4DynamicParticle* ,
200 const G4Material*) const 222 const G4Material*) const
201 { 223 {
202 return true; 224 return true;
203 } 225 }
204 226
205 //....oooOO0OOooo........oooOO0OOooo........oo 227 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
206 228
207 G4bool G4hParametrisedLossModel::IsInCharge(co 229 G4bool G4hParametrisedLossModel::IsInCharge(const G4ParticleDefinition* ,
208 const G4Material*) const 230 const G4Material*) const
209 { 231 {
210 return true; 232 return true;
211 } 233 }
212 234
213 //....oooOO0OOooo........oooOO0OOooo........oo 235 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
214 236
215 G4double G4hParametrisedLossModel::StoppingPow 237 G4double G4hParametrisedLossModel::StoppingPower(const G4Material* material,
216 G4double kineticEnergy) 238 G4double kineticEnergy)
217 { 239 {
218 G4double eloss = 0.0; 240 G4double eloss = 0.0;
219 241
220 const std::size_t numberOfElements = materia << 242 const G4int numberOfElements = material->GetNumberOfElements() ;
221 const G4double* theAtomicNumDensityVector = 243 const G4double* theAtomicNumDensityVector =
222 material->GetAtomicNumDensityVector() ; 244 material->GetAtomicNumDensityVector() ;
223 245
224 246
225 // compound material with parametrisation 247 // compound material with parametrisation
226 if( (eStopingPowerTable->HasMaterial(materia 248 if( (eStopingPowerTable->HasMaterial(material)) ) {
227 249
228 eloss = eStopingPowerTable->StoppingPower( 250 eloss = eStopingPowerTable->StoppingPower(material, kineticEnergy);
229 if ("QAO" != modelName) { 251 if ("QAO" != modelName) {
230 eloss *= material->GetTotNbOfAtomsPerVo 252 eloss *= material->GetTotNbOfAtomsPerVolume();
231 if(1 < numberOfElements) { 253 if(1 < numberOfElements) {
232 G4int nAtoms = 0; 254 G4int nAtoms = 0;
233 255
234 const G4int* theAtomsVector = material 256 const G4int* theAtomsVector = material->GetAtomsVector();
235 for (std::size_t iel=0; iel<numberOfEl << 257 for (G4int iel=0; iel<numberOfElements; iel++) {
236 nAtoms += theAtomsVector[iel]; 258 nAtoms += theAtomsVector[iel];
237 } 259 }
238 eloss /= nAtoms; 260 eloss /= nAtoms;
239 } 261 }
240 } 262 }
241 263
242 // pure material 264 // pure material
243 } else if(1 == numberOfElements) { 265 } else if(1 == numberOfElements) {
244 266
245 G4double z = material->GetZ(); 267 G4double z = material->GetZ();
246 eloss = (eStopingPowerTable->ElectronicSto 268 eloss = (eStopingPowerTable->ElectronicStoppingPower(z, kineticEnergy))
247 * (material->Ge 269 * (material->GetTotNbOfAtomsPerVolume()) ;
248 270
249 // Experimental data exist only for kinetic 271 // Experimental data exist only for kinetic energy 125 keV
250 } else if( MolecIsInZiegler1988(material)) { 272 } else if( MolecIsInZiegler1988(material)) {
251 273
252 // Cycle over elements - calculation based o 274 // Cycle over elements - calculation based on Bragg's rule
253 G4double eloss125 = 0.0 ; 275 G4double eloss125 = 0.0 ;
254 const G4ElementVector* theElementVector = 276 const G4ElementVector* theElementVector =
255 material->GetElemen 277 material->GetElementVector() ;
256 278
257 279
258 // loop for the elements in the material 280 // loop for the elements in the material
259 for (std::size_t i=0; i<numberOfElements; << 281 for (G4int i=0; i<numberOfElements; i++) {
260 const G4Element* element = (*theElementV 282 const G4Element* element = (*theElementVector)[i] ;
261 G4double z = element->GetZ() ; 283 G4double z = element->GetZ() ;
262 eloss +=(eStopingPowerTable->Electron 284 eloss +=(eStopingPowerTable->ElectronicStoppingPower(z,kineticEnergy))
263 * theAtomi 285 * theAtomicNumDensityVector[i] ;
264 eloss125 +=(eStopingPowerTable->Electron 286 eloss125 +=(eStopingPowerTable->ElectronicStoppingPower(z,125.0*keV))
265 * theAtomi 287 * theAtomicNumDensityVector[i] ;
266 } 288 }
267 289
268 // Chemical factor is taken into account 290 // Chemical factor is taken into account
269 if (eloss125 > 0.0) { << 291 eloss *= ChemicalFactor(kineticEnergy, eloss125) ;
270 eloss *= ChemicalFactor(kineticEnergy, e <<
271 } <<
272 292
273 // Brugg's rule calculation 293 // Brugg's rule calculation
274 } else { 294 } else {
275 const G4ElementVector* theElementVector = 295 const G4ElementVector* theElementVector =
276 material->GetElemen 296 material->GetElementVector() ;
277 297
278 // loop for the elements in the material 298 // loop for the elements in the material
279 for (std::size_t i=0; i<numberOfElements; << 299 for (G4int i=0; i<numberOfElements; i++)
280 { 300 {
281 const G4Element* element = (*theElementV 301 const G4Element* element = (*theElementVector)[i] ;
282 G4double z = element->GetZ() ; 302 G4double z = element->GetZ() ;
283 eloss += (eStopingPowerTable->Electron 303 eloss += (eStopingPowerTable->ElectronicStoppingPower(z,kineticEnergy))
284 * theAtomic 304 * theAtomicNumDensityVector[i];
285 } 305 }
286 } 306 }
287 return eloss; 307 return eloss;
288 } 308 }
289 309
290 //....oooOO0OOooo........oooOO0OOooo........oo 310 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
291 311
292 G4bool G4hParametrisedLossModel::MolecIsInZieg 312 G4bool G4hParametrisedLossModel::MolecIsInZiegler1988(
293 const G4Mater 313 const G4Material* material)
294 { 314 {
295 // The list of molecules from 315 // The list of molecules from
296 // J.F.Ziegler and J.M.Manoyan, The stopping 316 // J.F.Ziegler and J.M.Manoyan, The stopping of ions in compaunds,
297 // Nucl. Inst. & Meth. in Phys. Res. B35 (19 317 // Nucl. Inst. & Meth. in Phys. Res. B35 (1988) 215-228.
298 318
299 G4String myFormula = G4String(" ") ; 319 G4String myFormula = G4String(" ") ;
300 const G4String chFormula = material->GetChem 320 const G4String chFormula = material->GetChemicalFormula() ;
301 if (myFormula == chFormula ) return false ; 321 if (myFormula == chFormula ) return false ;
302 << 322
303 // There are no evidence for difference of 323 // There are no evidence for difference of stopping power depended on
304 // phase of the compound except for water. << 324 // phase of the compound except for water. The stopping power of the
305 // water in gas phase can be predicted usin 325 // water in gas phase can be predicted using Bragg's rule.
306 // << 326 //
307 // No chemical factor for water-gas << 327 // No chemical factor for water-gas
308 << 328
309 myFormula = G4String("H_2O") ; 329 myFormula = G4String("H_2O") ;
310 const G4State theState = material->GetState( 330 const G4State theState = material->GetState() ;
311 if( theState == kStateGas && myFormula == ch 331 if( theState == kStateGas && myFormula == chFormula) return false ;
312 << 332
313 const std::size_t numberOfMolecula = 53 ; << 333 const size_t numberOfMolecula = 53 ;
314 334
315 // The coffecient from Table.4 of Ziegler & 335 // The coffecient from Table.4 of Ziegler & Manoyan
316 static const G4double HeEff = 2.8735 ; << 336 const G4double HeEff = 2.8735 ;
317 << 337
318 static const G4String name[numberOfMolecula] << 338 static G4String name[numberOfMolecula] = {
319 "H_2O", "C_2H_4O", "C_3H_6O", "C_2 << 339 "H_2O", "C_2H_4O", "C_3H_6O", "C_2H_2", "C_H_3OH",
320 "C_2H_5OH", "C_3H_7OH", "C_3H_4", "NH 340 "C_2H_5OH", "C_3H_7OH", "C_3H_4", "NH_3", "C_14H_10",
321 "C_6H_6", "C_4H_10", "C_4H_6", "C_ 341 "C_6H_6", "C_4H_10", "C_4H_6", "C_4H_8O", "CCl_4",
322 "CF_4", "C_6H_8", "C_6H_12", "C_ 342 "CF_4", "C_6H_8", "C_6H_12", "C_6H_10O", "C_6H_10",
323 "C_8H_16", "C_5H_10", "C_5H_8", "C_ 343 "C_8H_16", "C_5H_10", "C_5H_8", "C_3H_6-Cyclopropane","C_2H_4F_2",
324 "C_2H_2F_2", "C_4H_8O_2", "C_2H_6", "C_ 344 "C_2H_2F_2", "C_4H_8O_2", "C_2H_6", "C_2F_6", "C_2H_6O",
325 "C_3H_6O", "C_4H_10O", "C_2H_4", "C_ 345 "C_3H_6O", "C_4H_10O", "C_2H_4", "C_2H_4O", "C_2H_4S",
326 "SH_2", "CH_4", "CCLF_3", "CC 346 "SH_2", "CH_4", "CCLF_3", "CCl_2F_2", "CHCl_2F",
327 "(CH_3)_2S", "N_2O", "C_5H_10O", "C_ << 347 "(CH_3)_2S", "N_2O", "C_5H_10O", "C_8H_6", "(CH_2)_N",
328 "(C_3H_6)_N","(C_8H_8)_N", "C_3H_8", "C_ 348 "(C_3H_6)_N","(C_8H_8)_N", "C_3H_8", "C_3H_6-Propylene", "C_3H_6O",
329 "C_3H_6S", "C_4H_4S", "C_7H_8" 349 "C_3H_6S", "C_4H_4S", "C_7H_8"
330 }; << 350 } ;
331 << 351
332 static const G4double expStopping[numberOfMo << 352 static G4double expStopping[numberOfMolecula] = {
333 66.1, 190.4, 258.7, 42.2, 141.5, << 353 66.1, 190.4, 258.7, 42.2, 141.5,
334 210.9, 279.6, 198.8, 31.0, 267.5, 354 210.9, 279.6, 198.8, 31.0, 267.5,
335 122.8, 311.4, 260.3, 328.9, 391.3, 355 122.8, 311.4, 260.3, 328.9, 391.3,
336 206.6, 374.0, 422.0, 432.0, 398.0, 356 206.6, 374.0, 422.0, 432.0, 398.0,
337 554.0, 353.0, 326.0, 74.6, 220.5, 357 554.0, 353.0, 326.0, 74.6, 220.5,
338 197.4, 362.0, 170.0, 330.5, 211.3, 358 197.4, 362.0, 170.0, 330.5, 211.3,
339 262.3, 349.6, 51.3, 187.0, 236.9, 359 262.3, 349.6, 51.3, 187.0, 236.9,
340 121.9, 35.8, 247.0, 292.6, 268.0, 360 121.9, 35.8, 247.0, 292.6, 268.0,
341 262.3, 49.0, 398.9, 444.0, 22.91, 361 262.3, 49.0, 398.9, 444.0, 22.91,
342 68.0, 155.0, 84.0, 74.2, 254.7, 362 68.0, 155.0, 84.0, 74.2, 254.7,
343 306.8, 324.4, 420.0 363 306.8, 324.4, 420.0
344 } ; 364 } ;
345 365
346 static const G4double expCharge[numberOfMole << 366 static G4double expCharge[numberOfMolecula] = {
347 HeEff, HeEff, HeEff, 1.0, HeEff, << 367 HeEff, HeEff, HeEff, 1.0, HeEff,
348 HeEff, HeEff, HeEff, 1.0, 1.0, 368 HeEff, HeEff, HeEff, 1.0, 1.0,
349 1.0, HeEff, HeEff, HeEff, HeEff, 369 1.0, HeEff, HeEff, HeEff, HeEff,
350 HeEff, HeEff, HeEff, HeEff, HeEff, 370 HeEff, HeEff, HeEff, HeEff, HeEff,
351 HeEff, HeEff, HeEff, 1.0, HeEff, 371 HeEff, HeEff, HeEff, 1.0, HeEff,
352 HeEff, HeEff, HeEff, HeEff, HeEff, 372 HeEff, HeEff, HeEff, HeEff, HeEff,
353 HeEff, HeEff, 1.0, HeEff, HeEff, 373 HeEff, HeEff, 1.0, HeEff, HeEff,
354 HeEff, 1.0, HeEff, HeEff, HeEff, 374 HeEff, 1.0, HeEff, HeEff, HeEff,
355 HeEff, 1.0, HeEff, HeEff, 1.0, 375 HeEff, 1.0, HeEff, HeEff, 1.0,
356 1.0, 1.0, 1.0, 1.0, HeEff, 376 1.0, 1.0, 1.0, 1.0, HeEff,
357 HeEff, HeEff, HeEff 377 HeEff, HeEff, HeEff
358 } ; 378 } ;
359 379
360 static const G4double numberOfAtomsPerMolecu << 380 static G4double numberOfAtomsPerMolecula[numberOfMolecula] = {
361 3.0, 7.0, 10.0, 4.0, 6.0, << 381 3.0, 7.0, 10.0, 4.0, 6.0,
362 9.0, 12.0, 7.0, 4.0, 24.0, 382 9.0, 12.0, 7.0, 4.0, 24.0,
363 12.0, 14.0, 10.0, 13.0, 5.0, 383 12.0, 14.0, 10.0, 13.0, 5.0,
364 5.0, 14.0, 18.0, 17.0, 17.0, 384 5.0, 14.0, 18.0, 17.0, 17.0,
365 24.0, 15.0, 13.0, 9.0, 8.0, 385 24.0, 15.0, 13.0, 9.0, 8.0,
366 6.0, 14.0, 8.0, 8.0, 9.0, 386 6.0, 14.0, 8.0, 8.0, 9.0,
367 10.0, 15.0, 6.0, 7.0, 7.0, 387 10.0, 15.0, 6.0, 7.0, 7.0,
368 3.0, 5.0, 5.0, 5.0, 5.0, 388 3.0, 5.0, 5.0, 5.0, 5.0,
369 9.0, 3.0, 16.0, 14.0, 3.0, 389 9.0, 3.0, 16.0, 14.0, 3.0,
370 9.0, 16.0, 11.0, 9.0, 10.0, 390 9.0, 16.0, 11.0, 9.0, 10.0,
371 10.0, 9.0, 15.0 391 10.0, 9.0, 15.0
372 } ; 392 } ;
373 393
374 // Search for the compaund in the table 394 // Search for the compaund in the table
375 for (std::size_t i=0; i<numberOfMolecula; ++ << 395 for (size_t i=0; i<numberOfMolecula; i++)
376 { << 396 {
377 if(chFormula == name[i]) { << 397 if(chFormula == name[i]) {
378 G4double exp125 = expStopping[i] * << 398 G4double exp125 = expStopping[i] *
379 (material->GetTotNbOfAtoms 399 (material->GetTotNbOfAtomsPerVolume()) /
380 (expCharge[i] * numberOfAt 400 (expCharge[i] * numberOfAtomsPerMolecula[i]) ;
381 SetExpStopPower125(exp125) ; 401 SetExpStopPower125(exp125) ;
382 return true ; 402 return true ;
383 } 403 }
384 } 404 }
385 << 405
386 return false ; 406 return false ;
387 } 407 }
388 408
389 //....oooOO0OOooo........oooOO0OOooo........oo 409 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
390 410
391 G4double G4hParametrisedLossModel::ChemicalFac 411 G4double G4hParametrisedLossModel::ChemicalFactor(
392 G4double kineticEn 412 G4double kineticEnergy, G4double eloss125) const
393 { 413 {
394 // Approximation of Chemical Factor accordin 414 // Approximation of Chemical Factor according to
395 // J.F.Ziegler and J.M.Manoyan, The stopping 415 // J.F.Ziegler and J.M.Manoyan, The stopping of ions in compaunds,
396 // Nucl. Inst. & Meth. in Phys. Res. B35 (19 416 // Nucl. Inst. & Meth. in Phys. Res. B35 (1988) 215-228.
397 << 417
398 G4double gamma = 1.0 + kineticEnergy/prot << 418 G4double gamma = 1.0 + kineticEnergy/proton_mass_c2 ;
399 G4double gamma25 = 1.0 + 25.0*keV /proton_m 419 G4double gamma25 = 1.0 + 25.0*keV /proton_mass_c2 ;
400 G4double gamma125 = 1.0 + 125.0*keV/proton_m 420 G4double gamma125 = 1.0 + 125.0*keV/proton_mass_c2 ;
401 G4double beta = std::sqrt(1.0 - 1.0/(gam 421 G4double beta = std::sqrt(1.0 - 1.0/(gamma*gamma)) ;
402 G4double beta25 = std::sqrt(1.0 - 1.0/(gam 422 G4double beta25 = std::sqrt(1.0 - 1.0/(gamma25*gamma25)) ;
403 G4double beta125 = std::sqrt(1.0 - 1.0/(gam 423 G4double beta125 = std::sqrt(1.0 - 1.0/(gamma125*gamma125)) ;
404 << 424
405 G4double factor = 1.0 + (expStopPower125/elo 425 G4double factor = 1.0 + (expStopPower125/eloss125 - 1.0) *
406 (1.0 + G4Exp( 1.48 * ( beta << 426 (1.0 + std::exp( 1.48 * ( beta125/beta25 - 7.0 ) ) ) /
407 (1.0 + G4Exp( 1.48 * ( beta << 427 (1.0 + std::exp( 1.48 * ( beta/beta25 - 7.0 ) ) ) ;
408 << 428
409 return factor ; 429 return factor ;
410 } 430 }
411 431
412 //....oooOO0OOooo........oooOO0OOooo........oo 432 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
413 433