source: trunk/source/processes/electromagnetic/lowenergy/src/G4hICRU49He.cc@ 968

Last change on this file since 968 was 819, checked in by garnier, 17 years ago

import all except CVS

File size: 13.3 KB
Line 
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: G4hICRU49He
33//
34// Author: V.Ivanchenko (Vladimir.Ivanchenko@cern.ch)
35//
36// Creation date: 20 July 2000
37//
38// Modifications:
39// 20/07/2000 V.Ivanchenko First implementation
40// 18/09/2000 V.Ivanchenko clean up - all variable are the same as in ICRU
41// 03/10/2000 V.Ivanchenko clean up accoding to CodeWizard
42// 10/05/2001 V.Ivanchenko Clean up againist Linux compilation with -Wall
43// 26/08/2004 V.Ivanchenko Fix a problem of effective charge
44//
45// Class Description:
46//
47// Electronic stopping power parametrised according to
48// ICRU Report N49, 1993. J.F. Ziegler model for He ion.
49//
50// Class Description: End
51//
52// -------------------------------------------------------------------
53//
54//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
55
56#include "G4hICRU49He.hh"
57#include "G4UnitsTable.hh"
58#include "globals.hh"
59#include "G4Material.hh"
60
61//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
62
63G4hICRU49He::G4hICRU49He():G4VhElectronicStoppingPower(),
64 rateMass(4.0026/1.007276),
65 iMolecula(0)
66{;}
67
68//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
69
70G4hICRU49He::~G4hICRU49He()
71{;}
72
73//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
74
75G4bool G4hICRU49He::HasMaterial(const G4Material* material)
76{
77 G4String chFormula = material->GetChemicalFormula() ;
78 G4String myFormula = G4String(" ");
79
80 if (myFormula == chFormula ) {
81 if(1 == (material->GetNumberOfElements())) {return true;}
82 return false ;
83 }
84
85 // ICRU Report N49, 1993. Power's model for He.
86 const size_t numberOfMolecula = 30 ;
87 static G4String name[numberOfMolecula] = {
88 "H_2", "Be-Solid", "C-Solid", "Graphite", "N_2",
89 "O_2", "Al-Solid", "Si-Solid", "Ar-Solid", "Cu-Solid",
90 "Ge", "W-Solid", "Au-Solid", "Pb-Solid", "C_2H_2",
91 "CO_2", "Cellulose-Nitrat", "C_2H_4", "LiF",
92 "CH_4", "Nylon", "Polycarbonate", "(CH_2)_N-Polyetilene", "PMMA",
93 "(C_8H_8)_N", "SiO_2", "CsI", "H_2O", "H_2O-Gas"} ;
94
95 // Special treatment for water in gas state
96
97 myFormula = G4String("H_2O") ;
98 const G4State theState = material->GetState() ;
99 if( theState == kStateGas && myFormula == chFormula) {
100 chFormula = G4String("H_2O-Gas");
101 }
102
103 // Search for the material in the table
104 for (size_t i=0; i<numberOfMolecula; i++) {
105 if (chFormula == name[i]) {
106 SetMoleculaNumber(i) ;
107 return true ;
108 }
109 }
110 return false ;
111}
112
113//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
114
115G4double G4hICRU49He::StoppingPower(const G4Material* material,
116 G4double kineticEnergy)
117{
118 G4double ionloss = 0.0 ;
119
120 // pure material (normally not the case for this function)
121 if(1 == (material->GetNumberOfElements())) {
122 G4double z = material->GetZ() ;
123 ionloss = ElectronicStoppingPower( z, kineticEnergy ) ;
124
125 // The data and the fit from:
126 // ICRU Report N49, 1993. Power's model for He.
127 } else if ( iMolecula < 30 ) {
128
129 // Reduced kinetic energy
130 // in internal units of parametrisation formula (MeV)
131 G4double T = kineticEnergy*rateMass/MeV ;
132
133 static G4double c[30][7] = {
134 {8.0080, 3.6287, 23.0700, 14.9900, 0.8507, 0.60, 2.0
135 },{ 13.3100, 3.7432, 39.4130, 12.1990, 1.0950, 0.38, 1.4
136 },{ 22.7240, 3.6040, 47.1810, 17.5490, 0.9040, 0.40, 1.4
137 },{ 24.4040, 2.4032, 48.9440, 27.9730, 1.2933, 0.40, 1.6
138 },{ 58.4719, 1.5115, 77.6421, 102.490, 1.5811, 0.50, 2.0
139 },{ 60.5408, 1.6297, 91.7601, 94.1260, 1.3662, 0.50, 2.0
140 },{ 48.4480, 6.4323, 59.2890, 18.3810, 0.4937, 0.48, 1.6
141 },{ 59.0346, 5.1305, 47.0866, 30.0857, 0.3500, 0.60, 2.0
142 },{ 71.8691, 2.8250, 51.1658, 57.1235, 0.4477, 0.60, 2.0
143 },{ 78.3520, 4.0961, 136.731, 28.4470, 1.0621, 0.52, 1.2
144 },{ 120.553, 1.5374, 49.8740, 82.2980, 0.8733, 0.45, 1.6
145 },{ 249.896, 0.6996, -37.274, 248.592, 1.1052, 0.50, 1.5
146 },{ 246.698, 0.6219, -58.391, 292.921, 0.8186, 0.56, 1.8
147 },{ 248.563, 0.6235, -36.8968, 306.960, 1.3214, 0.50, 2.0
148 },{ 25.5860, 1.7125, 154.723, 118.620, 2.2580, 0.50, 2.0
149 },{ 138.294, 25.6413, 231.873, 17.3780, 0.3218, 0.58, 1.3
150 },{ 83.2091, 1.1294, 135.7457, 190.865, 2.3461, 0.50, 2.0
151 },{ 263.542, 1.4754, 1541.446, 781.898, 1.9209, 0.40, 2.0
152 },{ 59.5545, 1.5354, 132.1523, 153.3537, 2.0262, 0.50, 2.0
153 },{ 31.7380, 19.820, 125.2100, 6.8910, 0.7242, 0.50, 1.1
154 },{ 31.7549, 1.5682, 97.4777, 106.0774, 2.3204, 0.50, 2.0
155 },{ 230.465, 4.8967, 1845.320, 358.641, 1.0774, 0.46, 1.2
156 },{ 423.444, 5.3761, 1189.114, 319.030, 0.7652, 0.48, 1.5
157 },{ 86.3410, 3.3322, 91.0433, 73.1091, 0.4650, 0.50, 2.0
158 },{ 146.105, 9.4344, 515.1500, 82.8860, 0.6239, 0.55, 1.5
159 },{ 238.050, 5.6901, 372.3575, 146.1835, 0.3992, 0.50, 2.0
160 },{ 124.2338, 2.6730, 133.8175, 99.4109, 0.7776, 0.50, 2.0
161 },{ 221.723, 1.5415, 87.7315, 192.5266, 1.0742, 0.50, 2.0
162 },{ 26.7537, 1.3717, 90.8007, 77.1587, 2.3264, 0.50, 2.0
163 },{ 37.6121, 1.8052, 73.0250, 66.2070, 1.4038, 0.50, 2.0} };
164
165 G4double a1,a2 ;
166
167 // Free electron gas model
168 if ( T < 0.001 ) {
169 G4double T0 = 0.001 ;
170 a1 = 1.0 - std::exp(-c[iMolecula][1]*std::pow(T0,-2.0+c[iMolecula][5])) ;
171 a2 = (c[iMolecula][0]*std::log(T0)/T0 + c[iMolecula][2]/T0) *
172 std::exp(-c[iMolecula][4]*std::pow(T0,-c[iMolecula][6])) +
173 c[iMolecula][3]/(T0*T0) ;
174
175 ionloss *= std::sqrt(T/T0) ;
176
177 // Main parametrisation
178 } else {
179 a1 = 1.0 - std::exp(-c[iMolecula][1]*std::pow(T,-2.0+c[iMolecula][5])) ;
180 a2 = (c[iMolecula][0]*std::log(T)/T + c[iMolecula][2]/T) *
181 std::exp(-c[iMolecula][4]*std::pow(T,-c[iMolecula][6])) +
182 c[iMolecula][3]/(T*T) ;
183 }
184
185 // He effective charge
186 G4double z = (material->GetTotNbOfElectPerVolume()) /
187 (material->GetTotNbOfAtomsPerVolume()) ;
188
189 ionloss = a1*a2 / HeEffChargeSquare(z, kineticEnergy*rateMass) ;
190
191 if ( ionloss < 0.0) ionloss = 0.0 ;
192 }
193
194 return ionloss ;
195}
196
197//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
198
199G4double G4hICRU49He::ElectronicStoppingPower(G4double z,
200 G4double kineticEnergy) const
201{
202 G4double ionloss ;
203 G4int i = G4int(z)-1 ; // index of atom
204 if(i < 0) i = 0 ;
205 if(i > 91) i = 91 ;
206
207 // The data and the fit from:
208 // ICRU Report 49, 1993. Ziegler's type of parametrisations
209 // Reduced kinetic energy
210
211 // He energy in internal units of parametrisation formula (MeV)
212 G4double T = kineticEnergy*rateMass/MeV ;
213
214 static G4double a[92][5] = {
215 {0.35485, 0.6456, 6.01525, 20.8933, 4.3515
216 },{ 0.58, 0.59, 6.3, 130.0, 44.07
217 },{ 1.42, 0.49, 12.25, 32.0, 9.161
218 },{ 2.1895, 0.47183,7.2362, 134.30, 197.96
219 },{ 3.691, 0.4128, 18.48, 50.72, 9.0
220 },{ 3.83523, 0.42993,12.6125, 227.41, 188.97
221 },{ 1.9259, 0.5550, 27.15125, 26.0665, 6.2768
222 },{ 2.81015, 0.4759, 50.0253, 10.556, 1.0382
223 },{ 1.533, 0.531, 40.44, 18.41, 2.718
224 },{ 2.303, 0.4861, 37.01, 37.96, 5.092
225 },{ 9.894, 0.3081, 23.65, 0.384, 92.93
226 },{ 4.3, 0.47, 34.3, 3.3, 12.74
227 },{ 2.5, 0.625, 45.7, 0.1, 4.359
228 },{ 2.1, 0.65, 49.34, 1.788, 4.133
229 },{ 1.729, 0.6562, 53.41, 2.405, 3.845
230 },{ 1.402, 0.6791, 58.98, 3.528, 3.211
231 },{ 1.117, 0.7044, 69.69, 3.705, 2.156
232 },{ 2.291, 0.6284, 73.88, 4.478, 2.066
233 },{ 8.554, 0.3817, 83.61, 11.84, 1.875
234 },{ 6.297, 0.4622, 65.39, 10.14, 5.036
235 },{ 5.307, 0.4918, 61.74, 12.4, 6.665
236 },{ 4.71, 0.5087, 65.28, 8.806, 5.948
237 },{ 6.151, 0.4524, 83.0, 18.31, 2.71
238 },{ 6.57, 0.4322, 84.76, 15.53, 2.779
239 },{ 5.738, 0.4492, 84.6, 14.18, 3.101
240 },{ 5.013, 0.4707, 85.8, 16.55, 3.211
241 },{ 4.32, 0.4947, 76.14, 10.85, 5.441
242 },{ 4.652, 0.4571, 80.73, 22.0, 4.952
243 },{ 3.114, 0.5236, 76.67, 7.62, 6.385
244 },{ 3.114, 0.5236, 76.67, 7.62, 7.502
245 },{ 3.114, 0.5236, 76.67, 7.62, 8.514
246 },{ 5.746, 0.4662, 79.24, 1.185, 7.993
247 },{ 2.792, 0.6346, 106.1, 0.2986, 2.331
248 },{ 4.667, 0.5095, 124.3, 2.102, 1.667
249 },{ 2.44, 0.6346, 105.0, 0.83, 2.851
250 },{ 1.413, 0.7377, 147.9, 1.466, 1.016
251 },{ 11.72, 0.3826, 102.8, 9.231, 4.371
252 },{ 7.126, 0.4804, 119.3, 5.784, 2.454
253 },{ 11.61, 0.3955, 146.7, 7.031, 1.423
254 },{ 10.99, 0.41, 163.9, 7.1, 1.052
255 },{ 9.241, 0.4275, 163.1, 7.954, 1.102
256 },{ 9.276, 0.418, 157.1, 8.038, 1.29
257 },{ 3.999, 0.6152, 97.6, 1.297, 5.792
258 },{ 4.306, 0.5658, 97.99, 5.514, 5.754
259 },{ 3.615, 0.6197, 86.26, 0.333, 8.689
260 },{ 5.8, 0.49, 147.2, 6.903, 1.289
261 },{ 5.6, 0.49, 130.0, 10.0, 2.844
262 },{ 3.55, 0.6068, 124.7, 1.112, 3.119
263 },{ 3.6, 0.62, 105.8, 0.1692, 6.026
264 },{ 5.4, 0.53, 103.1, 3.931, 7.767
265 },{ 3.97, 0.6459, 131.8, 0.2233, 2.723
266 },{ 3.65, 0.64, 126.8, 0.6834, 3.411
267 },{ 3.118, 0.6519, 164.9, 1.208, 1.51
268 },{ 3.949, 0.6209, 200.5, 1.878, 0.9126
269 },{ 14.4, 0.3923, 152.5, 8.354, 2.597
270 },{ 10.99, 0.4599, 138.4, 4.811, 3.726
271 },{ 16.6, 0.3773, 224.1, 6.28, 0.9121
272 },{ 10.54, 0.4533, 159.3, 4.832, 2.529
273 },{ 10.33, 0.4502, 162.0, 5.132, 2.444
274 },{ 10.15, 0.4471, 165.6, 5.378, 2.328
275 },{ 9.976, 0.4439, 168.0, 5.721, 2.258
276 },{ 9.804, 0.4408, 176.2, 5.675, 1.997
277 },{ 14.22, 0.363, 228.4, 7.024, 1.016
278 },{ 9.952, 0.4318, 233.5, 5.065, 0.9244
279 },{ 9.272, 0.4345, 210.0, 4.911, 1.258
280 },{ 10.13, 0.4146, 225.7, 5.525, 1.055
281 },{ 8.949, 0.4304, 213.3, 5.071, 1.221
282 },{ 11.94, 0.3783, 247.2, 6.655, 0.849
283 },{ 8.472, 0.4405, 195.5, 4.051, 1.604
284 },{ 8.301, 0.4399, 203.7, 3.667, 1.459
285 },{ 6.567, 0.4858, 193.0, 2.65, 1.66
286 },{ 5.951, 0.5016, 196.1, 2.662, 1.589
287 },{ 7.495, 0.4523, 251.4, 3.433, 0.8619
288 },{ 6.335, 0.4825, 255.1, 2.834, 0.8228
289 },{ 4.314, 0.5558, 214.8, 2.354, 1.263
290 },{ 4.02, 0.5681, 219.9, 2.402, 1.191
291 },{ 3.836, 0.5765, 210.2, 2.742, 1.305
292 },{ 4.68, 0.5247, 244.7, 2.749, 0.8962
293 },{ 3.223, 0.5883, 232.7, 2.954, 1.05
294 },{ 2.892, 0.6204, 208.6, 2.415, 1.416
295 },{ 4.728, 0.5522, 217.0, 3.091, 1.386
296 },{ 6.18, 0.52, 170.0, 4.0, 3.224
297 },{ 9.0, 0.47, 198.0, 3.8, 2.032
298 },{ 2.324, 0.6997, 216.0, 1.599, 1.399
299 },{ 1.961, 0.7286, 223.0, 1.621, 1.296
300 },{ 1.75, 0.7427, 350.1, 0.9789, 0.5507
301 },{ 10.31, 0.4613, 261.2, 4.738, 0.9899
302 },{ 7.962, 0.519, 235.7, 4.347, 1.313
303 },{ 6.227, 0.5645, 231.9, 3.961, 1.379
304 },{ 5.246, 0.5947, 228.6, 4.027, 1.432
305 },{ 5.408, 0.5811, 235.7, 3.961, 1.358
306 },{ 5.218, 0.5828, 245.0, 3.838, 1.25}
307 };
308
309 // Free electron gas model
310 if ( T < 0.001 ) {
311 G4double slow = a[i][0] ;
312 G4double shigh = std::log( 1.0 + a[i][3]*1000.0 + a[i][4]*0.001 )
313 * a[i][2]*1000.0 ;
314 ionloss = slow*shigh / (slow + shigh) ;
315 ionloss *= std::sqrt(T*1000.0) ;
316
317 // Main parametrisation
318 } else {
319 G4double slow = a[i][0] * std::pow((T*1000.0), a[i][1]) ;
320 G4double shigh = std::log( 1.0 + a[i][3]/T + a[i][4]*T ) * a[i][2]/T ;
321 ionloss = slow*shigh / (slow + shigh) ;
322
323 }
324 if ( ionloss < 0.0) ionloss = 0.0 ;
325
326 // He effective charge
327 ionloss /= HeEffChargeSquare(z, kineticEnergy*rateMass) ;
328
329 return ionloss;
330}
331
332
333
334
335
336
337
338
339
340
341
342
Note: See TracBrowser for help on using the repository browser.