Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/standard/src/G4ionIonisation.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /processes/electromagnetic/standard/src/G4ionIonisation.cc (Version 11.3.0) and /processes/electromagnetic/standard/src/G4ionIonisation.cc (Version 6.2.p2)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                    <<   3 // * DISCLAIMER                                                       *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th <<   5 // * The following disclaimer summarizes all the specific disclaimers *
  6 // * the Geant4 Collaboration.  It is provided <<   6 // * of contributors to this software. The specific disclaimers,which *
  7 // * conditions of the Geant4 Software License <<   7 // * govern, are listed with their locations in:                      *
  8 // * LICENSE and available at  http://cern.ch/ <<   8 // *   http://cern.ch/geant4/license                                  *
  9 // * include a list of copyright holders.      << 
 10 // *                                                9 // *                                                                  *
 11 // * Neither the authors of this software syst     10 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     11 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     12 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     13 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file  <<  14 // * use.                                                             *
 16 // * for the full disclaimer and the limitatio << 
 17 // *                                               15 // *                                                                  *
 18 // * This  code  implementation is the result  <<  16 // * This  code  implementation is the  intellectual property  of the *
 19 // * technical work of the GEANT4 collaboratio <<  17 // * GEANT4 collaboration.                                            *
 20 // * By using,  copying,  modifying or  distri <<  18 // * By copying,  distributing  or modifying the Program (or any work *
 21 // * any work based  on the software)  you  ag <<  19 // * based  on  the Program)  you indicate  your  acceptance of  this *
 22 // * use  in  resulting  scientific  publicati <<  20 // * statement, and all its terms.                                    *
 23 // * acceptance of all terms of the Geant4 Sof << 
 24 // *******************************************     21 // ********************************************************************
 25 //                                                 22 //
                                                   >>  23 // $Id: G4ionIonisation.cc,v 1.23 2004/05/27 17:22:56 vnivanch Exp $
                                                   >>  24 // GEANT4 tag $Name: geant4-06-02 $
 26 //                                                 25 //
 27 // -------------------------------------------     26 // -------------------------------------------------------------------
 28 //                                                 27 //
 29 // GEANT4 Class file                               28 // GEANT4 Class file
 30 //                                                 29 //
 31 //                                                 30 //
 32 // File name:     G4ionIonisation                  31 // File name:     G4ionIonisation
 33 //                                                 32 //
 34 // Author:        Vladimir Ivanchenko              33 // Author:        Vladimir Ivanchenko
 35 //                                                 34 //
 36 // Creation date: 07.05.2002                       35 // Creation date: 07.05.2002
 37 //                                                 36 //
 38 // Modifications:                                  37 // Modifications:
 39 //                                                 38 //
 40 // 23-12-02 Change interface in order to move      39 // 23-12-02 Change interface in order to move to cut per region (V.Ivanchenko)
 41 // 26-12-02 Secondary production moved to deri     40 // 26-12-02 Secondary production moved to derived classes (V.Ivanchenko)
 42 // 13-02-03 SubCutoff regime is assigned to a      41 // 13-02-03 SubCutoff regime is assigned to a region (V.Ivanchenko)
 43 // 18-04-03 Use IonFluctuations (V.Ivanchenko)     42 // 18-04-03 Use IonFluctuations (V.Ivanchenko)
 44 // 03-08-03 Add effective charge (V.Ivanchenko     43 // 03-08-03 Add effective charge (V.Ivanchenko)
 45 // 12-11-03 G4EnergyLossSTD -> G4EnergyLossPro     44 // 12-11-03 G4EnergyLossSTD -> G4EnergyLossProcess (V.Ivanchenko)
 46 // 27-05-04 Set integral to be a default regim <<  45 // 27-05-04 Set integral to be a default regime (V.Ivanchenko) 
 47 // 08-11-04 Migration to new interface of Stor << 
 48 // 08-04-05 Major optimisation of internal int << 
 49 // 10-01-06 SetStepLimits -> SetStepFunction ( << 
 50 // 10-05-06 Add a possibility to download user << 
 51 // 13-05-06 Add data for light ion stopping in << 
 52 // 14-01-07 use SetEmModel() and SetFluctModel << 
 53 // 16-05-07 Add data for light ion stopping on << 
 54 // 07-11-07 Fill non-ionizing energy loss (V.I << 
 55 // 12-09-08 Removed InitialiseMassCharge and C << 
 56 //                                                 46 //
 57 //                                                 47 //
 58 // -------------------------------------------     48 // -------------------------------------------------------------------
 59 //                                                 49 //
 60 //....oooOO0OOooo........oooOO0OOooo........oo     50 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 61 //....oooOO0OOooo........oooOO0OOooo........oo     51 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 62                                                    52 
 63 #include "G4ionIonisation.hh"                      53 #include "G4ionIonisation.hh"
 64 #include "G4PhysicalConstants.hh"              << 
 65 #include "G4SystemOfUnits.hh"                  << 
 66 #include "G4Electron.hh"                           54 #include "G4Electron.hh"
 67 #include "G4GenericIon.hh"                     <<  55 #include "G4Proton.hh"
                                                   >>  56 #include "G4AntiProton.hh"
 68 #include "G4BraggModel.hh"                         57 #include "G4BraggModel.hh"
 69 #include "G4BraggIonModel.hh"                  << 
 70 #include "G4BetheBlochModel.hh"                    58 #include "G4BetheBlochModel.hh"
 71 #include "G4LossTableManager.hh"               <<  59 #include "G4IonFluctuations.hh"
 72 #include "G4EmParameters.hh"                   <<  60 #include "G4UnitsTable.hh"
 73 #include "G4EmStandUtil.hh"                    << 
 74                                                    61 
 75 //....oooOO0OOooo........oooOO0OOooo........oo     62 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 76                                                    63 
 77 G4ionIonisation::G4ionIonisation(const G4Strin     64 G4ionIonisation::G4ionIonisation(const G4String& name)
 78   : G4VEnergyLossProcess(name)                 <<  65   : G4VEnergyLossProcess(name),
                                                   >>  66     theParticle(0),
                                                   >>  67     subCutoff(false)
 79 {                                                  68 {
 80   SetLinearLossLimit(0.02);                    <<  69   InitialiseProcess();
 81   SetProcessSubType(fIonisation);              << 
 82   SetSecondaryParticle(G4Electron::Electron()) << 
 83   eth = 2*CLHEP::MeV;                          << 
 84 }                                                  70 }
 85                                                    71 
 86 //....oooOO0OOooo........oooOO0OOooo........oo     72 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 87                                                    73 
 88 G4bool G4ionIonisation::IsApplicable(const G4P <<  74 G4ionIonisation::~G4ionIonisation()
                                                   >>  75 {}
                                                   >>  76 
                                                   >>  77 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
                                                   >>  78 
                                                   >>  79 void G4ionIonisation::InitialiseProcess()
 89 {                                                  80 {
 90   return true;                                 <<  81   SetSecondaryParticle(G4Electron::Electron());
                                                   >>  82 
                                                   >>  83   SetDEDXBinning(120);
                                                   >>  84   SetLambdaBinning(120);
                                                   >>  85   SetMinKinEnergy(0.1*keV);
                                                   >>  86   SetMaxKinEnergy(100.0*TeV);
                                                   >>  87 
                                                   >>  88   flucModel = new G4IonFluctuations();
                                                   >>  89 
                                                   >>  90   G4VEmModel* em = new G4BraggModel();
                                                   >>  91   em->SetLowEnergyLimit(0.1*keV);
                                                   >>  92   em->SetHighEnergyLimit(2.0*MeV);
                                                   >>  93   AddEmModel(1, em, flucModel);
                                                   >>  94   G4VEmModel* em1 = new G4BetheBlochModel();
                                                   >>  95   em1->SetLowEnergyLimit(2.0*MeV);
                                                   >>  96   em1->SetHighEnergyLimit(100.0*TeV);
                                                   >>  97   AddEmModel(2, em1, flucModel);
                                                   >>  98 
                                                   >>  99   chargeLowLimit = 0.1;
                                                   >> 100   energyLowLimit = 25.*MeV;
                                                   >> 101   SetLinearLossLimit(0.15);
                                                   >> 102   SetStepLimits(0.1, 0.1*mm);
 91 }                                                 103 }
 92                                                   104 
 93 //....oooOO0OOooo........oooOO0OOooo........oo    105 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 94                                                   106 
 95 G4double G4ionIonisation::MinPrimaryEnergy(con << 107 const G4ParticleDefinition* G4ionIonisation::DefineBaseParticle(
 96              const G4Material*,                << 108                       const G4ParticleDefinition* p)
 97              G4double cut)                     << 
 98 {                                                 109 {
 99   return p->GetPDGMass()*(std::sqrt(1. + 0.5*c << 110   if(p) theParticle = p;
                                                   >> 111   theBaseParticle = G4Proton::Proton();
                                                   >> 112   return theBaseParticle;
100 }                                                 113 }
101                                                   114 
102 //....oooOO0OOooo........oooOO0OOooo........oo    115 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
103                                                   116 
104 void G4ionIonisation::InitialiseEnergyLossProc << 117 void G4ionIonisation::PrintInfoDefinition()
105           const G4ParticleDefinition* part,    << 
106           const G4ParticleDefinition* bpart)   << 
107 {                                                 118 {
108   const G4ParticleDefinition* ion = G4GenericI << 119   G4VEnergyLossProcess::PrintInfoDefinition();
109                                                   120 
110   if(!isInitialised) {                         << 121   G4cout << "      Scaling relation is used to proton dE/dx and range"
111     theParticle = part;                        << 122          << G4endl
112                                                << 123          << "      Bether-Bloch model for Escaled > 2 MeV, "
113     // define base particle                    << 124          << "parametrisation of Bragg peak below."
114     const G4ParticleDefinition* theBaseParticl << 125          << G4endl;
115     const G4int pdg = part->GetPDGEncoding();  << 126 }
116                                                   127 
117     if(part == bpart) {                        << 128 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
118       theBaseParticle = nullptr;               << 
119     } else if(nullptr != bpart) {              << 
120       theBaseParticle = bpart;                 << 
121     } else if(part == ion || pdg == 1000020040 << 
122       theBaseParticle = nullptr;               << 
123     } else {                                   << 
124       theBaseParticle = ion;                   << 
125     }                                          << 
126     SetBaseParticle(theBaseParticle);          << 
127                                                   129 
128     // model limit defined for protons         << 130 void G4ionIonisation::SetSubCutoff(G4bool val)
129     eth = 2*CLHEP::MeV*part->GetPDGMass()/CLHE << 131 {
                                                   >> 132   subCutoff = val;
                                                   >> 133 }
130                                                   134 
131     G4EmParameters* param = G4EmParameters::In << 135 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
132     G4double emin = param->MinKinEnergy();     << 
133     G4double emax = param->MaxKinEnergy();     << 
134                                                << 
135     // define model of energy loss fluctuation << 
136     if (nullptr == FluctModel()) {             << 
137       SetFluctModel(G4EmStandUtil::ModelOfFluc << 
138     }                                          << 
139                                                   136 
140     if (nullptr == EmModel(0)) {               << 137 G4double G4ionIonisation::EffectiveCharge(const G4ParticleDefinition* p,
141       if (pdg == 1000020040) {                 << 138                                           const G4Material* material,
142   SetEmModel(new G4BraggIonModel());           << 139                               G4double kineticEnergy)
143       } else {                                 << 140 {
144         SetEmModel(new G4BraggModel());        << 141   G4double mass   = p->GetPDGMass();
                                                   >> 142   G4double charge = p->GetPDGCharge();
                                                   >> 143   G4double Zi     = charge/eplus;
                                                   >> 144 
                                                   >> 145   chargeCorrection = 1.0;
                                                   >> 146 
                                                   >> 147   // The aproximation of ion effective charge from:
                                                   >> 148   // J.F.Ziegler, J.P. Biersack, U. Littmark
                                                   >> 149   // The Stopping and Range of Ions in Matter,
                                                   >> 150   // Vol.1, Pergamon Press, 1985
                                                   >> 151   // Fast ions or hadrons
                                                   >> 152   G4double reducedEnergy = kineticEnergy * proton_mass_c2/mass ;
                                                   >> 153   if( reducedEnergy > energyLowLimit || Zi < 1.5 ) return charge ;
                                                   >> 154 
                                                   >> 155   static G4double vFermi[92] = {
                                                   >> 156     1.0309,  0.15976, 0.59782, 1.0781,  1.0486,  1.0,     1.058,   0.93942, 0.74562, 0.3424,
                                                   >> 157     0.45259, 0.71074, 0.90519, 0.97411, 0.97184, 0.89852, 0.70827, 0.39816, 0.36552, 0.62712,
                                                   >> 158     0.81707, 0.9943,  1.1423,  1.2381,  1.1222,  0.92705, 1.0047,  1.2,     1.0661,  0.97411,
                                                   >> 159     0.84912, 0.95,    1.0903,  1.0429,  0.49715, 0.37755, 0.35211, 0.57801, 0.77773, 1.0207,
                                                   >> 160     1.029,   1.2542,  1.122,   1.1241,  1.0882,  1.2709,  1.2542,  0.90094, 0.74093, 0.86054,
                                                   >> 161     0.93155, 1.0047,  0.55379, 0.43289, 0.32636, 0.5131,  0.695,   0.72591, 0.71202, 0.67413,
                                                   >> 162     0.71418, 0.71453, 0.5911,  0.70263, 0.68049, 0.68203, 0.68121, 0.68532, 0.68715, 0.61884,
                                                   >> 163     0.71801, 0.83048, 1.1222,  1.2381,  1.045,   1.0733,  1.0953,  1.2381,  1.2879,  0.78654,
                                                   >> 164     0.66401, 0.84912, 0.88433, 0.80746, 0.43357, 0.41923, 0.43638, 0.51464, 0.73087, 0.81065,
                                                   >> 165     1.9578,  1.0257} ;
                                                   >> 166 
                                                   >> 167   static G4double lFactor[92] = {
                                                   >> 168     1.0,  1.0,  1.1,  1.06, 1.01, 1.03, 1.04, 0.99, 0.95, 0.9,
                                                   >> 169     0.82, 0.81, 0.83, 0.88, 1.0,  0.95, 0.97, 0.99, 0.98, 0.97,
                                                   >> 170     0.98, 0.97, 0.96, 0.93, 0.91, 0.9,  0.88, 0.9,  0.9,  0.9,
                                                   >> 171     0.9,  0.85, 0.9,  0.9,  0.91, 0.92, 0.9,  0.9,  0.9,  0.9,
                                                   >> 172     0.9,  0.88, 0.9,  0.88, 0.88, 0.9,  0.9,  0.88, 0.9,  0.9,
                                                   >> 173     0.9,  0.9,  0.96, 1.2,  0.9,  0.88, 0.88, 0.85, 0.9,  0.9,
                                                   >> 174     0.92, 0.95, 0.99, 1.03, 1.05, 1.07, 1.08, 1.1,  1.08, 1.08,
                                                   >> 175     1.08, 1.08, 1.09, 1.09, 1.1,  1.11, 1.12, 1.13, 1.14, 1.15,
                                                   >> 176     1.17, 1.2,  1.18, 1.17, 1.17, 1.16, 1.16, 1.16, 1.16, 1.16,
                                                   >> 177     1.16, 1.16} ;
                                                   >> 178 
                                                   >> 179   static G4double c[6] = {0.2865,  0.1266, -0.001429,
                                                   >> 180                           0.02402,-0.01135, 0.001475} ;
                                                   >> 181 
                                                   >> 182   // get elements in the actual material,
                                                   >> 183   const G4ElementVector* theElementVector = material->GetElementVector() ;
                                                   >> 184   const G4double* theAtomicNumDensityVector =
                                                   >> 185                          material->GetAtomicNumDensityVector() ;
                                                   >> 186   const G4int NumberOfElements = material->GetNumberOfElements() ;
                                                   >> 187 
                                                   >> 188   //  loop for the elements in the material
                                                   >> 189   //  to find out average values Z, vF, lF
                                                   >> 190   G4double z = 0.0, vF = 0.0, lF = 0.0, norm = 0.0 ;
                                                   >> 191 
                                                   >> 192   if( 1 == NumberOfElements ) {
                                                   >> 193     z = material->GetZ() ;
                                                   >> 194     G4int iz = G4int(z) - 1 ;
                                                   >> 195     if(iz < 0) iz = 0 ;
                                                   >> 196     else if(iz > 91) iz = 91 ;
                                                   >> 197     vF   = vFermi[iz] ;
                                                   >> 198     lF   = lFactor[iz] ;
                                                   >> 199 
                                                   >> 200   } else {
                                                   >> 201     for (G4int iel=0; iel<NumberOfElements; iel++)
                                                   >> 202       {
                                                   >> 203         const G4Element* element = (*theElementVector)[iel] ;
                                                   >> 204         G4double z2 = element->GetZ() ;
                                                   >> 205         const G4double weight = theAtomicNumDensityVector[iel] ;
                                                   >> 206         norm += weight ;
                                                   >> 207         z    += z2 * weight ;
                                                   >> 208         G4int iz = G4int(z2) - 1 ;
                                                   >> 209         if(iz < 0) iz = 0 ;
                                                   >> 210         else if(iz > 91) iz =91 ;
                                                   >> 211         vF   += vFermi[iz] * weight ;
                                                   >> 212         lF   += lFactor[iz] * weight ;
145       }                                           213       }
                                                   >> 214     z  /= norm ;
                                                   >> 215     vF /= norm ;
                                                   >> 216     lF /= norm ;
                                                   >> 217   }
                                                   >> 218 
                                                   >> 219   G4double q;
                                                   >> 220   // Helium ion case
                                                   >> 221   if( Zi < 2.5 ) {
                                                   >> 222 
                                                   >> 223     // Normalise to He4 mass
                                                   >> 224     G4double e = log(std::max(1.0, kineticEnergy / (keV*4.0026) ) );
                                                   >> 225     G4double x = c[0] ;
                                                   >> 226     G4double y = 1.0 ;
                                                   >> 227     for (G4int i=1; i<6; i++) {
                                                   >> 228       y *= e ;
                                                   >> 229       x += y * c[i] ;
146     }                                             230     }
147     // to compute ranges correctly we have to  << 231     q = 7.6 -  e ;
148     // model even if activation limit is high  << 232     q = 1.0 + ( 0.007 + 0.00005 * z ) * exp( -q*q ) * sqrt(1.0 - exp(-x)) ;
149     EmModel(0)->SetLowEnergyLimit(emin);       << 233     if( q < chargeLowLimit ) q = chargeLowLimit ;
150                                                << 234 
151     // high energy limit may be eth or DBL_MAX << 235     // Heavy ion case
152     G4double emax1 = (EmModel(0)->HighEnergyLi << 236   } else {
153     EmModel(0)->SetHighEnergyLimit(emax1);     << 237     // v1 is ion velocity in vF unit
154     AddEmModel(1, EmModel(0), FluctModel());   << 238     G4double v1 = sqrt( reducedEnergy / (25.0 * keV) )/ vF ;
155                                                << 239     G4double y ;
156     // second model is used if the first does  << 240     G4double z13 = pow(Zi, 0.3333) ;
157     if(emax1 < emax) {                         << 241 
158       if (nullptr == EmModel(1)) { SetEmModel( << 242     // Faster than Fermi velocity
159       EmModel(1)->SetLowEnergyLimit(emax1);    << 243     if ( v1 > 1.0 ) {
160                                                << 244       y = vF * v1 * ( 1.0 + 0.2 / (v1*v1) ) / (z13*z13) ;
161       // for extremely heavy particles upper l << 245 
162       // should be increased                   << 246       // Slower than Fermi velocity
163       emax = std::max(emax, eth*10);           << 247     } else {
164       EmModel(1)->SetHighEnergyLimit(emax);    << 248       y = 0.6923 * vF * (1.0 + 2.0*v1*v1/3.0 + v1*v1*v1*v1/15.0) / (z13*z13) ;
165       AddEmModel(2, EmModel(1), FluctModel()); << 
166     }                                             249     }
167     isInitialised = true;                      << 250 
                                                   >> 251     G4double y3 = pow(y, 0.3) ;
                                                   >> 252     //    G4cout << "y= " << y << " y3= " << y3 << " v1= " << v1 << " vF= " << vF << G4endl; 
                                                   >> 253     q = 1.0 - exp( 0.803*y3 - 1.3167*y3*y3 - 0.38157*y - 0.008983*y*y ) ;
                                                   >> 254 
                                                   >> 255     if( q < chargeLowLimit ) q = chargeLowLimit ;
                                                   >> 256 
                                                   >> 257     G4double s = 7.6 -  log(std::max(1.0, reducedEnergy/keV)) ;
                                                   >> 258     s = 1.0 + ( 0.18 + 0.0015 * z ) * exp( -s*s )/ (Zi*Zi) ;
                                                   >> 259 
                                                   >> 260     // Screen length according to
                                                   >> 261     // J.F.Ziegler and J.M.Manoyan, The stopping of ions in compaunds,
                                                   >> 262     // Nucl. Inst. & Meth. in Phys. Res. B35 (1988) 215-228.
                                                   >> 263 
                                                   >> 264     G4double lambda = 10.0 * vF * pow(1.0-q, 0.6667) / (z13 * (6.0 + q)) ;
                                                   >> 265     chargeCorrection = s * (1.0 + 0.5*(1.0/q - 1.0)*log(1.0 + lambda*lambda)/(vF*vF) );
168   }                                               266   }
                                                   >> 267   //  G4cout << "G4ionIonisation: charge= " << charge << " q= " << q 
                                                   >> 268   //       << " chargeCor= " << chargeCorrection << G4endl;
                                                   >> 269   return charge*q;
169 }                                                 270 }
170                                                   271 
171 //....oooOO0OOooo........oooOO0OOooo........oo    272 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
172                                                   273 
173 void G4ionIonisation::ProcessDescription(std:: << 
174 {                                              << 
175   out << "  Ion ionisation";                   << 
176   G4VEnergyLossProcess::ProcessDescription(out << 
177 }                                              << 
178                                                   274 
179 //....oooOO0OOooo........oooOO0OOooo........oo << 
180                                                   275