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

Last change on this file since 1315 was 819, checked in by garnier, 16 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.