source: trunk/source/materials/src/G4IonisParamElm.cc

Last change on this file was 1337, checked in by garnier, 14 years ago

tag geant4.9.4 beta 1 + modifs locales

File size: 9.0 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// $Id: G4IonisParamElm.cc,v 1.16 2010/04/29 11:11:56 vnivanch Exp $
28// GEANT4 tag $Name: geant4-09-04-beta-01 $
29//
30//
31//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
32//
33// 09-07-98, data moved from G4Element. M.Maire
34// 22-11-00, tabulation of ionisation potential from
35//           the ICRU Report N#37. V.Ivanchenko
36// 08-03-01, correct handling of fShellCorrectionVector. M.Maire
37// 17-10-02, Fix excitation energy interpolation. V.Ivanchenko
38// 06-09-04, Update calculated values after any change of ionisation
39//           potential change. V.Ivanchenko
40// 29-04-10, Using G4Pow and mean ionisation energy from NIST V.Ivanchenko
41//
42//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
43
44#include "G4IonisParamElm.hh"
45#include "G4NistManager.hh"
46#include "G4Pow.hh"
47
48//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
49
50G4IonisParamElm::G4IonisParamElm(G4double AtomNumber)
51{
52  G4int Z = G4int(AtomNumber + 0.5);
53  if (Z < 1) {
54    G4Exception (" ERROR! It is not allowed to create an Element with Z < 1" );
55  }
56  G4Pow* g4pow = G4Pow::GetInstance();
57
58  // some basic functions of the atomic number
59  fZ     = Z;
60  fZ3    = g4pow->Z13(Z);
61  fZZ3   = fZ3*g4pow->Z13(Z+1);
62  flogZ3 = g4pow->logZ(Z)/3.;
63   
64  // Parameters for energy loss by ionisation   
65  /*
66  // Mean excitation energy
67  // from "Stopping Powers for Electrons and Positrons"
68  // ICRU Report N#37, 1984  (energy in eV)
69  static double exc[100] = {
70    21.8, 20.9, 13.3, 15.9, 15.2, 13.0, 11.7, 11.9, 11.5, 13.7,
71    13.6, 13.0, 12.8, 12.4, 11.5, 11.3, 10.2, 10.4, 10.0,  9.6,
72    10.3, 10.6, 10.7, 10.7, 10.9, 11.0, 11.0, 11.1, 11.1, 11.0,
73    10.8, 10.4, 10.5, 10.2,  9.8,  9.8,  9.8,  9.6,  9.7,  9.8,
74    10.2, 10.1, 10.0, 10.0, 10.0, 10.2, 10.0,  9.8, 10.0,  9.8,
75     9.5,  9.3,  9.3,  8.9,  8.9,  8.8,  8.8,  8.8,  9.1,  9.1,
76     9.2,  9.3,  9.2,  9.2,  9.4,  9.5,  9.7,  9.7,  9.8,  9.8,
77     9.8,  9.8,  9.8,  9.8,  9.8,  9.8,  9.8, 10.1, 10.0, 10.0,
78    10.0, 10.0,  9.9,  9.9,  9.7,  9.2,  9.5,  9.4,  9.4,  9.4,
79     9.6,  9.7,  9.7,  9.8,  9.8,  9.8,  9.8,  9.9,  9.9,  9.9 };
80
81  fMeanExcitationEnergy = fZ * exc[iz] * eV ;
82  */
83
84  fMeanExcitationEnergy = 
85    G4NistManager::Instance()->GetMeanIonisationEnergy(Z);
86
87  // compute parameters for ion transport
88  // The aproximation from:
89  // J.F.Ziegler, J.P. Biersack, U. Littmark
90  // The Stopping and Range of Ions in Matter,
91  // Vol.1, Pergamon Press, 1985
92  // Fast ions or hadrons
93
94  G4int iz = Z - 1;
95  if(91 < iz) { iz = 91; }
96
97  static G4double vFermi[92] = {
98    1.0309,  0.15976, 0.59782, 1.0781,  1.0486,  1.0,     1.058,   0.93942, 0.74562, 0.3424,
99    0.45259, 0.71074, 0.90519, 0.97411, 0.97184, 0.89852, 0.70827, 0.39816, 0.36552, 0.62712,
100    0.81707, 0.9943,  1.1423,  1.2381,  1.1222,  0.92705, 1.0047,  1.2,     1.0661,  0.97411,
101    0.84912, 0.95,    1.0903,  1.0429,  0.49715, 0.37755, 0.35211, 0.57801, 0.77773, 1.0207,
102    1.029,   1.2542,  1.122,   1.1241,  1.0882,  1.2709,  1.2542,  0.90094, 0.74093, 0.86054,
103    0.93155, 1.0047,  0.55379, 0.43289, 0.32636, 0.5131,  0.695,   0.72591, 0.71202, 0.67413,
104    0.71418, 0.71453, 0.5911,  0.70263, 0.68049, 0.68203, 0.68121, 0.68532, 0.68715, 0.61884,
105    0.71801, 0.83048, 1.1222,  1.2381,  1.045,   1.0733,  1.0953,  1.2381,  1.2879,  0.78654,
106    0.66401, 0.84912, 0.88433, 0.80746, 0.43357, 0.41923, 0.43638, 0.51464, 0.73087, 0.81065,
107    1.9578,  1.0257} ;
108
109  static G4double lFactor[92] = {
110    1.0,  1.0,  1.1,  1.06, 1.01, 1.03, 1.04, 0.99, 0.95, 0.9,
111    0.82, 0.81, 0.83, 0.88, 1.0,  0.95, 0.97, 0.99, 0.98, 0.97,
112    0.98, 0.97, 0.96, 0.93, 0.91, 0.9,  0.88, 0.9,  0.9,  0.9,
113    0.9,  0.85, 0.9,  0.9,  0.91, 0.92, 0.9,  0.9,  0.9,  0.9,
114    0.9,  0.88, 0.9,  0.88, 0.88, 0.9,  0.9,  0.88, 0.9,  0.9,
115    0.9,  0.9,  0.96, 1.2,  0.9,  0.88, 0.88, 0.85, 0.9,  0.9,
116    0.92, 0.95, 0.99, 1.03, 1.05, 1.07, 1.08, 1.1,  1.08, 1.08,
117    1.08, 1.08, 1.09, 1.09, 1.1,  1.11, 1.12, 1.13, 1.14, 1.15,
118    1.17, 1.2,  1.18, 1.17, 1.17, 1.16, 1.16, 1.16, 1.16, 1.16,
119    1.16, 1.16} ;
120
121  fVFermi  = vFermi[iz];
122  fLFactor = lFactor[iz];
123
124  // obsolete parameters for ionisation
125  fTau0 = 0.1*fZ3*MeV/proton_mass_c2;
126  fTaul = 2.*MeV/proton_mass_c2;
127
128  // compute the Bethe-Bloch formula for energy = fTaul*particle mass
129  G4double rate = fMeanExcitationEnergy/electron_mass_c2 ;
130  G4double w = fTaul*(fTaul+2.) ;
131  fBetheBlochLow = (fTaul+1.)*(fTaul+1.)*std::log(2.*w/rate)/w - 1. ;
132  fBetheBlochLow = 2.*fZ*twopi_mc2_rcl2*fBetheBlochLow ; 
133 
134  fClow = std::sqrt(fTaul)*fBetheBlochLow;
135  fAlow = 6.458040 * fClow/fTau0;
136  G4double Taum = 0.035*fZ3*MeV/proton_mass_c2;
137  fBlow =-3.229020*fClow/(fTau0*std::sqrt(Taum));
138
139  // Shell correction parameterization
140  fShellCorrectionVector = new G4double[3]; //[3]
141  rate = 0.001*fMeanExcitationEnergy/eV;
142  G4double rate2 = rate*rate;
143    /*   
144    fShellCorrectionVector[0] = ( 1.10289e5 + 5.14781e8*rate)*rate2 ;
145    fShellCorrectionVector[1] = ( 7.93805e3 - 2.22565e7*rate)*rate2 ;
146    fShellCorrectionVector[2] = (-9.92256e1 + 2.10823e5*rate)*rate2 ;
147    */
148  fShellCorrectionVector[0] = ( 0.422377   + 3.858019*rate)*rate2 ;
149  fShellCorrectionVector[1] = ( 0.0304043  - 0.1667989*rate)*rate2 ;
150  fShellCorrectionVector[2] = (-0.00038106 + 0.00157955*rate)*rate2 ;
151}
152
153//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
154
155// Fake default constructor - sets only member data and allocates memory
156//                            for usage restricted to object persistency
157
158G4IonisParamElm::G4IonisParamElm(__void__&)
159  : fShellCorrectionVector(0)
160{}
161
162//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
163
164G4IonisParamElm::~G4IonisParamElm()
165{
166  if (fShellCorrectionVector) { delete [] fShellCorrectionVector; }
167}
168
169//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
170
171G4IonisParamElm::G4IonisParamElm(G4IonisParamElm& right)
172{
173  *this = right;
174}
175
176//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
177
178const G4IonisParamElm& G4IonisParamElm::operator=(const G4IonisParamElm& right)
179{
180  if (this != &right)
181    {
182      fZ                     = right.fZ;
183      fZ3                    = right.fZ3;
184      fZZ3                   = right.fZZ3;
185      flogZ3                 = right.flogZ3;
186      fTau0                  = right.fTau0;
187      fTaul                  = right.fTaul;
188      fBetheBlochLow         = right.fBetheBlochLow;
189      fAlow                  = right.fAlow;
190      fBlow                  = right.fBlow;
191      fClow                  = right.fClow;
192      fMeanExcitationEnergy  = right.fMeanExcitationEnergy;
193      if (fShellCorrectionVector) { delete [] fShellCorrectionVector; } 
194      fShellCorrectionVector = new G4double[3];           
195      fShellCorrectionVector[0] = right.fShellCorrectionVector[0];
196      fShellCorrectionVector[1] = right.fShellCorrectionVector[1];
197      fShellCorrectionVector[2] = right.fShellCorrectionVector[2];     
198     } 
199  return *this;
200}
201
202//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
203
204G4int G4IonisParamElm::operator==(const G4IonisParamElm& right) const
205{
206  return (this == (G4IonisParamElm *) &right);
207}
208
209//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
210
211G4int G4IonisParamElm::operator!=(const G4IonisParamElm& right) const
212{
213  return (this != (G4IonisParamElm *) &right);
214}
215
216//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
Note: See TracBrowser for help on using the repository browser.