Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/materials/src/G4IonisParamMat.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 /materials/src/G4IonisParamMat.cc (Version 11.3.0) and /materials/src/G4IonisParamMat.cc (Version 10.5)


  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 //
                                                   >>  26 //
                                                   >>  27 // 
                                                   >>  28 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
 25                                                    29 
 26 // 09-07-98, data moved from G4Material, M.Mai     30 // 09-07-98, data moved from G4Material, M.Maire
 27 // 18-07-98, bug corrected in ComputeDensityEf     31 // 18-07-98, bug corrected in ComputeDensityEffect() for gas
 28 // 16-01-01, bug corrected in ComputeDensityEf     32 // 16-01-01, bug corrected in ComputeDensityEffect() E100eV (L.Urban)
 29 // 08-02-01, fShellCorrectionVector correctly      33 // 08-02-01, fShellCorrectionVector correctly handled (mma)
 30 // 28-10-02, add setMeanExcitationEnergy (V.Iv     34 // 28-10-02, add setMeanExcitationEnergy (V.Ivanchenko)
 31 // 06-09-04, factor 2 to shell correction term <<  35 // 06-09-04, factor 2 to shell correction term (V.Ivanchenko) 
 32 // 10-05-05, add a missing coma in FindMeanExc     36 // 10-05-05, add a missing coma in FindMeanExcitationEnergy() - Bug#746 (mma)
 33 // 27-09-07, add computation of parameters for     37 // 27-09-07, add computation of parameters for ions (V.Ivanchenko)
 34 // 04-03-08, remove reference to G4NistManager     38 // 04-03-08, remove reference to G4NistManager. Add fBirks constant (mma)
 35 // 30-10-09, add G4DensityEffectData class and     39 // 30-10-09, add G4DensityEffectData class and density effect computation (VI)
 36                                                    40 
 37 #include "G4IonisParamMat.hh"                  <<  41 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
 38                                                    42 
 39 #include "G4AtomicShells.hh"                   <<  43 #include "G4IonisParamMat.hh"
 40 #include "G4AutoLock.hh"                       << 
 41 #include "G4DensityEffectData.hh"              << 
 42 #include "G4Exp.hh"                            << 
 43 #include "G4Log.hh"                            << 
 44 #include "G4Material.hh"                           44 #include "G4Material.hh"
                                                   >>  45 #include "G4DensityEffectData.hh"
 45 #include "G4NistManager.hh"                        46 #include "G4NistManager.hh"
 46 #include "G4PhysicalConstants.hh"              << 
 47 #include "G4Pow.hh"                                47 #include "G4Pow.hh"
                                                   >>  48 #include "G4PhysicalConstants.hh"
 48 #include "G4SystemOfUnits.hh"                      49 #include "G4SystemOfUnits.hh"
 49                                                    50 
 50 G4DensityEffectData* G4IonisParamMat::fDensity     51 G4DensityEffectData* G4IonisParamMat::fDensityData = nullptr;
 51                                                    52 
 52 namespace                                      <<  53 #ifdef G4MULTITHREADED
 53 {                                              <<  54   G4Mutex G4IonisParamMat::ionisMutex = G4MUTEX_INITIALIZER;
 54   G4Mutex ionisMutex = G4MUTEX_INITIALIZER;    <<  55 #endif
 55 }                                              << 
 56                                                    56 
 57 //....oooOO0OOooo........oooOO0OOooo........oo     57 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
 58                                                    58 
 59 G4IonisParamMat::G4IonisParamMat(const G4Mater <<  59 G4IonisParamMat::G4IonisParamMat(const G4Material* material)
                                                   >>  60   : fMaterial(material)
 60 {                                                  61 {
 61   fBirks = 0.;                                     62   fBirks = 0.;
 62   fMeanEnergyPerIon = 0.0;                         63   fMeanEnergyPerIon = 0.0;
 63   twoln10 = 2. * G4Pow::GetInstance()->logZ(10 <<  64   twoln10 = 2.*G4Pow::GetInstance()->logZ(10);
 64                                                    65 
 65   // minimal set of default parameters for den     66   // minimal set of default parameters for density effect
 66   fCdensity = 0.0;                                 67   fCdensity = 0.0;
 67   fD0density = 0.0;                                68   fD0density = 0.0;
 68   fAdjustmentFactor = 1.0;                         69   fAdjustmentFactor = 1.0;
 69   if (fDensityData == nullptr) {               <<  70   if(fDensityData == nullptr) { fDensityData = new G4DensityEffectData(); }
 70     fDensityData = new G4DensityEffectData();  << 
 71   }                                            << 
 72   fDensityEffectCalc = nullptr;                << 
 73                                                    71 
 74   // compute parameters                            72   // compute parameters
 75   ComputeMeanParameters();                         73   ComputeMeanParameters();
 76   ComputeDensityEffectParameters();            <<  74   ComputeDensityEffect();
 77   ComputeFluctModel();                             75   ComputeFluctModel();
 78   ComputeIonParameters();                          76   ComputeIonParameters();
 79 }                                                  77 }
 80                                                    78 
 81 //....oooOO0OOooo........oooOO0OOooo........oo     79 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
 82                                                    80 
 83 G4IonisParamMat::~G4IonisParamMat()            <<  81 // Fake default constructor - sets only member data and allocates memory
                                                   >>  82 //                            for usage restricted to object persistency
                                                   >>  83 
                                                   >>  84 G4IonisParamMat::G4IonisParamMat(__void__&)
                                                   >>  85   : fMaterial(nullptr), fShellCorrectionVector(nullptr)
 84 {                                                  86 {
 85   delete fDensityEffectCalc;                   <<  87   fMeanExcitationEnergy = 0.0;
 86   delete[] fShellCorrectionVector;             <<  88   fLogMeanExcEnergy = 0.0;
 87   delete fDensityData;                         <<  89   fTaul = 0.0;
 88   fDensityData = nullptr;                      <<  90   fCdensity = 0.0;
 89   fShellCorrectionVector = nullptr;            <<  91   fMdensity = 0.0;
 90   fDensityEffectCalc = nullptr;                <<  92   fAdensity = 0.0;
                                                   >>  93   fX0density = 0.0;
                                                   >>  94   fX1density = 0.0;
                                                   >>  95   fD0density = 0.0;
                                                   >>  96   fPlasmaEnergy = 0.0;
                                                   >>  97   fAdjustmentFactor = 0.0;
                                                   >>  98   fF1fluct = 0.0;          
                                                   >>  99   fF2fluct = 0.0;                       
                                                   >> 100   fEnergy1fluct = 0.0;
                                                   >> 101   fLogEnergy1fluct = 0.0;
                                                   >> 102   fEnergy2fluct = 0.0;
                                                   >> 103   fLogEnergy2fluct = 0.0;
                                                   >> 104   fEnergy0fluct = 0.0;
                                                   >> 105   fRateionexcfluct = 0.0;
                                                   >> 106   fZeff = 0.0;
                                                   >> 107   fFermiEnergy = 0.0;
                                                   >> 108   fLfactor = 0.0;
                                                   >> 109   fInvA23 = 0.0;
                                                   >> 110   fBirks = 0.0;
                                                   >> 111   fMeanEnergyPerIon = 0.0;
                                                   >> 112   twoln10 = 2.*G4Pow::GetInstance()->logZ(10);
                                                   >> 113 
                                                   >> 114   if(fDensityData == nullptr) { fDensityData = new G4DensityEffectData(); }
 91 }                                                 115 }
 92                                                   116 
 93 //....oooOO0OOooo........oooOO0OOooo........oo    117 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
 94                                                   118 
 95 G4double G4IonisParamMat::GetDensityCorrection << 119 G4IonisParamMat::~G4IonisParamMat()
 96 {                                                 120 {
 97   // x = log10(beta*gamma)                     << 121   delete [] fShellCorrectionVector; 
 98   G4double y = 0.0;                            << 122   delete fDensityData;
 99   if (x < fX0density) {                        << 123   fDensityData = nullptr;
100     if (fD0density > 0.0) {                    << 124   fShellCorrectionVector = nullptr;
101       y = fD0density * G4Exp(twoln10 * (x - fX << 
102     }                                          << 
103   }                                            << 
104   else if (x >= fX1density) {                  << 
105     y = twoln10 * x - fCdensity;               << 
106   }                                            << 
107   else {                                       << 
108     y = twoln10 * x - fCdensity + fAdensity *  << 
109   }                                            << 
110   return y;                                    << 
111 }                                                 125 }
112                                                   126 
113 //....oooOO0OOooo........oooOO0OOooo........oo    127 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
114                                                   128 
115 void G4IonisParamMat::ComputeMeanParameters()     129 void G4IonisParamMat::ComputeMeanParameters()
116 {                                                 130 {
117   // compute mean excitation energy and shell     131   // compute mean excitation energy and shell correction vector
118   fTaul = (*(fMaterial->GetElementVector()))[0    132   fTaul = (*(fMaterial->GetElementVector()))[0]->GetIonisation()->GetTaul();
119                                                   133 
120   std::size_t nElements = fMaterial->GetNumber << 134   fMeanExcitationEnergy = 0.;
                                                   >> 135   fLogMeanExcEnergy = 0.;
                                                   >> 136 
                                                   >> 137   size_t nElements = fMaterial->GetNumberOfElements();
121   const G4ElementVector* elmVector = fMaterial    138   const G4ElementVector* elmVector = fMaterial->GetElementVector();
122   const G4double* nAtomsPerVolume = fMaterial-    139   const G4double* nAtomsPerVolume = fMaterial->GetVecNbOfAtomsPerVolume();
123                                                << 140  
124   fMeanExcitationEnergy = FindMeanExcitationEn    141   fMeanExcitationEnergy = FindMeanExcitationEnergy(fMaterial);
125   fLogMeanExcEnergy = 0.;                      << 
126                                                   142 
127   // Chemical formula defines mean excitation     143   // Chemical formula defines mean excitation energy
128   if (fMeanExcitationEnergy > 0.0) {           << 144   if(fMeanExcitationEnergy > 0.0) {
129     fLogMeanExcEnergy = G4Log(fMeanExcitationE    145     fLogMeanExcEnergy = G4Log(fMeanExcitationEnergy);
130                                                   146 
131     // Compute average                         << 147     // Compute average 
132   }                                            << 148   } else {
133   else {                                       << 149     for (size_t i=0; i < nElements; i++) {
134     for (std::size_t i = 0; i < nElements; ++i << 
135       const G4Element* elm = (*elmVector)[i];     150       const G4Element* elm = (*elmVector)[i];
136       fLogMeanExcEnergy +=                     << 151       fLogMeanExcEnergy += nAtomsPerVolume[i]*elm->GetZ()
137         nAtomsPerVolume[i] * elm->GetZ() * G4L << 152   *G4Log(elm->GetIonisation()->GetMeanExcitationEnergy());
138     }                                             153     }
139     fLogMeanExcEnergy /= fMaterial->GetTotNbOf    154     fLogMeanExcEnergy /= fMaterial->GetTotNbOfElectPerVolume();
140     fMeanExcitationEnergy = G4Exp(fLogMeanExcE    155     fMeanExcitationEnergy = G4Exp(fLogMeanExcEnergy);
141   }                                               156   }
142                                                   157 
143   fShellCorrectionVector = new G4double[3];    << 158   fShellCorrectionVector = new G4double[3]; 
144                                                   159 
145   for (G4int j = 0; j <= 2; ++j) {             << 160   for (G4int j=0; j<=2; j++)
                                                   >> 161   {
146     fShellCorrectionVector[j] = 0.;               162     fShellCorrectionVector[j] = 0.;
147                                                   163 
148     for (std::size_t k = 0; k < nElements; ++k << 164     for (size_t k=0; k<nElements; k++) {
149       fShellCorrectionVector[j] +=             << 165       fShellCorrectionVector[j] += nAtomsPerVolume[k]
150         nAtomsPerVolume[k] * (((*elmVector)[k] << 166   *(((*elmVector)[k])->GetIonisation()->GetShellCorrectionVector())[j];
151     }                                             167     }
152     fShellCorrectionVector[j] *= 2.0 / fMateri << 168     fShellCorrectionVector[j] *= 2.0/fMaterial->GetTotNbOfElectPerVolume();
153   }                                            << 169   } 
154 }                                                 170 }
155                                                   171 
156 //....oooOO0OOooo........oooOO0OOooo........oo    172 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
157                                                   173 
158 G4DensityEffectData* G4IonisParamMat::GetDensi << 174 G4DensityEffectData* G4IonisParamMat::GetDensityEffectData()
                                                   >> 175 {
159   return fDensityData;                            176   return fDensityData;
160 }                                                 177 }
161                                                   178 
162 //....oooOO0OOooo........oooOO0OOooo........oo    179 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
163                                                << 180                     
164 void G4IonisParamMat::ComputeDensityEffectPara << 181 void G4IonisParamMat::ComputeDensityEffect()
165 {                                                 182 {
166   G4State State = fMaterial->GetState();          183   G4State State = fMaterial->GetState();
167   G4double density = fMaterial->GetDensity();  << 
168                                                   184 
169   // Check if density effect data exist in the    185   // Check if density effect data exist in the table
170   // R.M. Sternheimer, Atomic Data and Nuclear    186   // R.M. Sternheimer, Atomic Data and Nuclear Data Tables, 30: 261 (1984)
171   // or is assign to one of data set in this t << 
172   G4int idx = fDensityData->GetIndex(fMaterial    187   G4int idx = fDensityData->GetIndex(fMaterial->GetName());
173   auto nelm = (G4int)fMaterial->GetNumberOfEle << 188   G4int nelm= fMaterial->GetNumberOfElements();
174   G4int Z0 = ((*(fMaterial->GetElementVector() << 189   G4int Z0  = ((*(fMaterial->GetElementVector()))[0])->GetZasInt();
175   const G4Material* bmat = fMaterial->GetBaseM << 
176   G4NistManager* nist = G4NistManager::Instanc << 
177                                                << 
178   // arbitrary empirical limits                << 
179   // parameterisation with very different dens << 
180   static const G4double corrmax = 1.;          << 
181   static const G4double massfracmax = 0.9;     << 
182                                                   190 
183   // for simple non-NIST materials             << 191   // for simple non-NIST materials 
184   G4double corr = 0.0;                            192   G4double corr = 0.0;
185                                                << 193   if(idx < 0 && 1 == nelm) {
186   if (idx < 0 && 1 == nelm) {                  << 194     idx = fDensityData->GetElementIndex(Z0, fMaterial->GetState());
187     G4int z = (1 == Z0 && State == kStateLiqui << 
188     idx = fDensityData->GetElementIndex(z);    << 
189                                                   195 
190     // Correction for base material or for non    196     // Correction for base material or for non-nominal density
191     // Except cases of very different density     197     // Except cases of very different density defined in user code
192     if (idx >= 0 && 0 < z) {                   << 198     if(idx >= 0) {
193       G4double dens = nist->GetNominalDensity( << 199       const G4Material* bmat = fMaterial->GetBaseMaterial();
194       if (dens <= 0.0) {                       << 200       if(bmat) {
195         idx = -1;                              << 201   corr = G4Log(bmat->GetDensity()/fMaterial->GetDensity());
196       }                                        << 202       } else {
197       else {                                   << 203   G4double dens = G4NistManager::Instance()->GetNominalDensity(Z0);
198         corr = G4Log(dens / density);          << 204   if(dens <= 0.0) { idx = -1; }
199         if (std::abs(corr) > corrmax) {        << 205   else {
200           idx = -1;                            << 206     corr = G4Log(dens/fMaterial->GetDensity());
201         }                                      << 207   }
202       }                                        << 208       }
203     }                                          << 209       // 1.0 is an arbitrary empirical limit 
204   }                                            << 210       // parameterisation with very different density is not applicable
205   // for base material case                    << 211       if(std::abs(corr) > 1.0) { idx = -1; }
206   if (idx < 0 && nullptr != bmat) {            << 
207     idx = fDensityData->GetIndex(bmat->GetName << 
208     if (idx >= 0) {                            << 
209       corr = G4Log(bmat->GetDensity() / densit << 
210       if (std::abs(corr) > corrmax) {          << 
211         idx = -1;                              << 
212       }                                        << 
213     }                                             212     }
214   }                                               213   }
215                                                   214 
216   // for compound non-NIST materials with one  << 215   //G4cout<<"DensityEffect for "<<fMaterial->GetName()<<"  "<< idx << G4endl; 
217   if (idx < 0 && 1 < nelm) {                   << 216 
218     const G4double tot = fMaterial->GetTotNbOf << 217   if(idx >= 0) {
219     for (G4int i = 0; i < nelm; ++i) {         << 
220       const G4double frac = fMaterial->GetVecN << 
221       if (frac > massfracmax) {                << 
222         Z0 = ((*(fMaterial->GetElementVector() << 
223         idx = fDensityData->GetElementIndex(Z0 << 
224         G4double dens = nist->GetNominalDensit << 
225         if (idx >= 0 && dens > 0.0) {          << 
226           corr = G4Log(dens / density);        << 
227           if (std::abs(corr) > corrmax) {      << 
228             idx = -1;                          << 
229           }                                    << 
230           else {                               << 
231             break;                             << 
232           }                                    << 
233         }                                      << 
234       }                                        << 
235     }                                          << 
236   }                                            << 
237                                                   218 
238   if (idx >= 0) {                              << 
239     // Take parameters for the density effect     219     // Take parameters for the density effect correction from
240     // R.M. Sternheimer et al. Density Effect  << 220     // R.M. Sternheimer et al. Density Effect For The Ionization Loss 
241     // of Charged Particles in Various Substan << 221     // of Charged Particles in Various Substances. 
242     // Atom. Data Nucl. Data Tabl. 30 (1984) 2 << 222     // Atom. Data Nucl. Data Tabl. 30 (1984) 261-271. 
243                                                << 223 
244     fCdensity = fDensityData->GetCdensity(idx) << 224     fCdensity  = fDensityData->GetCdensity(idx); 
245     fMdensity = fDensityData->GetMdensity(idx) << 225     fMdensity  = fDensityData->GetMdensity(idx);
246     fAdensity = fDensityData->GetAdensity(idx) << 226     fAdensity  = fDensityData->GetAdensity(idx);
247     fX0density = fDensityData->GetX0density(id    227     fX0density = fDensityData->GetX0density(idx);
248     fX1density = fDensityData->GetX1density(id    228     fX1density = fDensityData->GetX1density(idx);
249     fD0density = fDensityData->GetDelta0densit    229     fD0density = fDensityData->GetDelta0density(idx);
250     fPlasmaEnergy = fDensityData->GetPlasmaEne    230     fPlasmaEnergy = fDensityData->GetPlasmaEnergy(idx);
251     fAdjustmentFactor = fDensityData->GetAdjus    231     fAdjustmentFactor = fDensityData->GetAdjustmentFactor(idx);
252                                                   232 
253     // parameter C is computed and not taken f    233     // parameter C is computed and not taken from Sternheimer tables
254     // fCdensity = 1. + 2*G4Log(fMeanExcitatio << 234     //fCdensity = 1. + 2*G4Log(fMeanExcitationEnergy/fPlasmaEnergy);
255     // G4cout << "IonisParamMat: " << fMateria << 235     //G4cout << "IonisParamMat: " << fMaterial->GetName() 
256     //     << "  Cst= " << Cdensity << " C= "     236     //     << "  Cst= " << Cdensity << " C= " << fCdensity << G4endl;
257                                                   237 
258     // correction on nominal density              238     // correction on nominal density
259     fCdensity += corr;                         << 239     fCdensity  += corr;
260     fX0density += corr / twoln10;              << 240     fX0density += corr/twoln10;
261     fX1density += corr / twoln10;              << 241     fX1density += corr/twoln10;
262   }                                            << 242 
263   else {                                       << 243   } else {
264     static const G4double Cd2 = 4 * CLHEP::pi  << 244 
265     fPlasmaEnergy = std::sqrt(Cd2 * fMaterial- << 245     static const G4double Cd2 = 4*pi*hbarc_squared*classic_electr_radius;
                                                   >> 246     fPlasmaEnergy = std::sqrt(Cd2*fMaterial->GetTotNbOfElectPerVolume());
266                                                   247 
267     // Compute parameters for the density effe    248     // Compute parameters for the density effect correction in DE/Dx formula.
268     // The parametrization is from R.M. Sternh    249     // The parametrization is from R.M. Sternheimer, Phys. Rev.B,3:3681 (1971)
269     G4int icase;                                  250     G4int icase;
270                                                << 251     
271     fCdensity = 1. + 2 * G4Log(fMeanExcitation << 252     fCdensity = 1. + 2*G4Log(fMeanExcitationEnergy/fPlasmaEnergy);
272     //                                            253     //
273     // condensed materials                        254     // condensed materials
274     //                                         << 255     //  
275     if ((State == kStateSolid) || (State == kS << 256     if ((State == kStateSolid)||(State == kStateLiquid)) {
276       static const G4double E100eV = 100. * CL << 
277       static const G4double ClimiS[] = {3.681, << 
278       static const G4double X0valS[] = {1.0, 1 << 
279       static const G4double X1valS[] = {2.0, 3 << 
280                                                << 
281       if (fMeanExcitationEnergy < E100eV) {    << 
282         icase = 0;                             << 
283       }                                        << 
284       else {                                   << 
285         icase = 1;                             << 
286       }                                        << 
287                                                << 
288       if (fCdensity < ClimiS[icase]) {         << 
289         fX0density = 0.2;                      << 
290       }                                        << 
291       else {                                   << 
292         fX0density = 0.326 * fCdensity - X0val << 
293       }                                        << 
294                                                << 
295       fX1density = X1valS[icase];              << 
296       fMdensity = 3.0;                         << 
297                                                   257 
298       // special: Hydrogen                     << 258       static const G4double E100eV  = 100.*eV; 
                                                   >> 259       static const G4double ClimiS[] = {3.681 , 5.215 };
                                                   >> 260       static const G4double X0valS[] = {1.0   , 1.5   };
                                                   >> 261       static const G4double X1valS[] = {2.0   , 3.0   };
                                                   >> 262                                 
                                                   >> 263       if(fMeanExcitationEnergy < E100eV) { icase = 0; }
                                                   >> 264       else                               { icase = 1; } 
                                                   >> 265 
                                                   >> 266       if(fCdensity < ClimiS[icase])    { fX0density = 0.2; }
                                                   >> 267       else { fX0density = 0.326*fCdensity - X0valS[icase]; }
                                                   >> 268 
                                                   >> 269       fX1density = X1valS[icase]; fMdensity = 3.0;
                                                   >> 270       
                                                   >> 271       //special: Hydrogen
299       if (1 == nelm && 1 == Z0) {                 272       if (1 == nelm && 1 == Z0) {
300         fX0density = 0.425;                    << 273          fX0density = 0.425; fX1density = 2.0; fMdensity = 5.949;
301         fX1density = 2.0;                      << 
302         fMdensity = 5.949;                     << 
303       }                                           274       }
304     }                                          << 275     } else {
305     else {                                     << 
306       //                                          276       //
307       // gases                                    277       // gases
308       //                                          278       //
309       fMdensity = 3.;                             279       fMdensity = 3.;
310       fX1density = 4.0;                           280       fX1density = 4.0;
311                                                << 281       //static const G4double ClimiG[] = {10.,10.5,11.,11.5,12.25,13.804};
312       if (fCdensity <= 10.) {                  << 282       //static const G4double X0valG[] = {1.6,1.7,1.8,1.9,2.0,2.0};
313         fX0density = 1.6;                      << 283       //static const G4double X1valG[] = {4.0,4.0,4.0,4.0,4.0,5.0};
314       }                                        << 284 
315       else if (fCdensity <= 10.5) {            << 285       if(fCdensity < 10.) {
316         fX0density = 1.7;                      << 286   fX0density = 1.6; 
317       }                                        << 287       } else if(fCdensity < 11.5) { 
318       else if (fCdensity <= 11.0) {            << 288   fX0density = 1.6 + 0.2*(fCdensity - 10.); 
319         fX0density = 1.8;                      << 289       } else if(fCdensity < 12.25) { 
320       }                                        << 290   fX0density = 1.9 + (fCdensity - 11.5)/7.5; 
321       else if (fCdensity <= 11.5) {            << 291       } else if(fCdensity < 13.804) { 
322         fX0density = 1.9;                      << 292   fX0density = 2.0; 
323       }                                        << 293   fX1density = 4.0 + (fCdensity - 12.25)/1.554;
324       else if (fCdensity <= 12.25) {           << 294       } else {
325         fX0density = 2.0;                      << 295   fX0density = 0.326*fCdensity-2.5; fX1density = 5.0; 
326       }                                        << 
327       else if (fCdensity <= 13.804) {          << 
328         fX0density = 2.0;                      << 
329         fX1density = 5.0;                      << 
330       }                                           296       }
331       else {                                   << 297       
332         fX0density = 0.326 * fCdensity - 2.5;  << 298       //special: Hydrogen
333         fX1density = 5.0;                      << 
334       }                                        << 
335                                                << 
336       // special: Hydrogen                     << 
337       if (1 == nelm && 1 == Z0) {                 299       if (1 == nelm && 1 == Z0) {
338         fX0density = 1.837;                    << 300          fX0density = 1.837; fX1density = 3.0; fMdensity = 4.754;
339         fX1density = 3.0;                      << 
340         fMdensity = 4.754;                     << 
341       }                                           301       }
342                                                << 302       
343       // special: Helium                       << 303       //special: Helium
344       if (1 == nelm && 2 == Z0) {                 304       if (1 == nelm && 2 == Z0) {
345         fX0density = 2.191;                    << 305          fX0density = 2.191; fX1density = 3.0; fMdensity = 3.297;
346         fX1density = 3.0;                      << 
347         fMdensity = 3.297;                     << 
348       }                                           306       }
349     }                                             307     }
350   }                                               308   }
351                                                   309 
352   // change parameters if the gas is not in ST    310   // change parameters if the gas is not in STP.
353   // For the correction the density(STP) is ne << 311   // For the correction the density(STP) is needed. 
354   // Density(STP) is calculated here :         << 312   // Density(STP) is calculated here : 
355                                                << 313   
356   if (State == kStateGas) {                    << 314     
                                                   >> 315   if (State == kStateGas) { 
                                                   >> 316     G4double Density  = fMaterial->GetDensity();
357     G4double Pressure = fMaterial->GetPressure    317     G4double Pressure = fMaterial->GetPressure();
358     G4double Temp = fMaterial->GetTemperature( << 318     G4double Temp     = fMaterial->GetTemperature();
359                                                << 319       
360     G4double DensitySTP = density * STP_Pressu << 320     G4double DensitySTP = Density*STP_Pressure*Temp/(Pressure*NTP_Temperature);
361                                                << 321 
362     G4double ParCorr = G4Log(density / Density << 322     G4double ParCorr = G4Log(Density/DensitySTP);
363                                                << 323   
364     fCdensity -= ParCorr;                      << 324     fCdensity  -= ParCorr;
365     fX0density -= ParCorr / twoln10;           << 325     fX0density -= ParCorr/twoln10;
366     fX1density -= ParCorr / twoln10;           << 326     fX1density -= ParCorr/twoln10;
367   }                                            << 327   }
368                                                << 328 
369   // fAdensity parameter can be fixed for not  << 329   // fAdensity parameter can be fixed for not conductive materials 
370   if (0.0 == fD0density) {                     << 330   if(0.0 == fD0density) {
371     G4double Xa = fCdensity / twoln10;         << 331     G4double Xa = fCdensity/twoln10;
372     fAdensity = twoln10 * (Xa - fX0density) /  << 332     fAdensity = twoln10*(Xa-fX0density)
373   }                                            << 333       /std::pow((fX1density-fX0density),fMdensity);
                                                   >> 334   }
                                                   >> 335   /*  
                                                   >> 336   G4cout << "G4IonisParamMat: density effect data for <" 
                                                   >> 337          << fMaterial->GetName() 
                                                   >> 338    << "> " << G4endl;
                                                   >> 339   G4cout << "Eplasma(eV)= " << fPlasmaEnergy/eV
                                                   >> 340    << " rho= " << fAdjustmentFactor
                                                   >> 341    << " -C= " << fCdensity 
                                                   >> 342    << " x0= " << fX0density
                                                   >> 343    << " x1= " << fX1density
                                                   >> 344    << " a= " << fAdensity
                                                   >> 345    << " m= " << fMdensity
                                                   >> 346    << G4endl;
                                                   >> 347   */
374 }                                                 348 }
375                                                   349 
376 //....oooOO0OOooo........oooOO0OOooo........oo    350 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
377                                                   351 
378 void G4IonisParamMat::ComputeFluctModel()         352 void G4IonisParamMat::ComputeFluctModel()
379 {                                                 353 {
380   // compute parameters for the energy loss fl    354   // compute parameters for the energy loss fluctuation model
381   // needs an 'effective Z'                    << 355   // needs an 'effective Z' 
382   G4double Zeff = 0.;                             356   G4double Zeff = 0.;
383   for (std::size_t i = 0; i < fMaterial->GetNu << 357   for (size_t i=0;i<fMaterial->GetNumberOfElements();i++) {
384     Zeff += (fMaterial->GetFractionVector())[i << 358      Zeff += (fMaterial->GetFractionVector())[i]
385   }                                            << 359              *((*(fMaterial->GetElementVector()))[i]->GetZ());
386   if (Zeff > 2.1) {                            << 360   }
387     fF2fluct = 2.0 / Zeff;                     << 361   if (Zeff > 2.) { fF2fluct = 2./Zeff; }
388     fF1fluct = 1. - fF2fluct;                  << 362   else           { fF2fluct = 0.; }
389     fEnergy2fluct = 10. * Zeff * Zeff * CLHEP: << 363 
390     fLogEnergy2fluct = G4Log(fEnergy2fluct);   << 364   fF1fluct         = 1. - fF2fluct;
391     fLogEnergy1fluct = (fLogMeanExcEnergy - fF << 365   fEnergy2fluct    = 10.*Zeff*Zeff*eV;
392   } else if (Zeff > 1.1) {                     << 366   fLogEnergy2fluct = G4Log(fEnergy2fluct);
393     fF2fluct = 0.0;                            << 367   fLogEnergy1fluct = (fLogMeanExcEnergy - fF2fluct*fLogEnergy2fluct)
394     fF1fluct = 1.0;                            << 368                      /fF1fluct;
395     fEnergy2fluct = 40. * CLHEP::eV;           << 369   fEnergy1fluct    = G4Exp(fLogEnergy1fluct);
396     fLogEnergy2fluct = G4Log(fEnergy2fluct);   << 370   fEnergy0fluct    = 10.*eV;
397     fLogEnergy1fluct = fLogMeanExcEnergy;      << 
398   } else {                                     << 
399     fF2fluct = 0.0;                            << 
400     fF1fluct = 1.0;                            << 
401     fEnergy2fluct = 10. * CLHEP::eV;           << 
402     fLogEnergy2fluct = G4Log(fEnergy2fluct);   << 
403     fLogEnergy1fluct = fLogMeanExcEnergy;      << 
404   }                                            << 
405   fEnergy1fluct = G4Exp(fLogEnergy1fluct);     << 
406   fEnergy0fluct = 10. * CLHEP::eV;             << 
407   fRateionexcfluct = 0.4;                         371   fRateionexcfluct = 0.4;
408 }                                                 372 }
409                                                   373 
410 //....oooOO0OOooo........oooOO0OOooo........oo    374 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
411                                                   375 
412 void G4IonisParamMat::ComputeIonParameters()      376 void G4IonisParamMat::ComputeIonParameters()
413 {                                                 377 {
414   // get elements in the actual material,         378   // get elements in the actual material,
415   const G4ElementVector* theElementVector = fM << 379   const G4ElementVector* theElementVector = fMaterial->GetElementVector() ;
416   const G4double* theAtomicNumDensityVector =  << 380   const G4double* theAtomicNumDensityVector =
417   const auto NumberOfElements = (G4int)fMateri << 381                          fMaterial->GetAtomicNumDensityVector() ;
                                                   >> 382   const G4int NumberOfElements = fMaterial->GetNumberOfElements() ;
418                                                   383 
419   //  loop for the elements in the material       384   //  loop for the elements in the material
420   //  to find out average values Z, vF, lF        385   //  to find out average values Z, vF, lF
421   G4double z(0.0), vF(0.0), lF(0.0), a23(0.0); << 386   G4double z(0.0), vF(0.0), lF(0.0), norm(0.0), a23(0.0);
422                                                   387 
423   G4Pow* g4pow = G4Pow::GetInstance();            388   G4Pow* g4pow = G4Pow::GetInstance();
424   if (1 == NumberOfElements) {                 << 389   if( 1 == NumberOfElements ) {
425     const G4Element* element = (*theElementVec    390     const G4Element* element = (*theElementVector)[0];
426     z = element->GetZ();                          391     z = element->GetZ();
427     vF = element->GetIonisation()->GetFermiVel << 392     vF= element->GetIonisation()->GetFermiVelocity();
428     lF = element->GetIonisation()->GetLFactor( << 393     lF= element->GetIonisation()->GetLFactor();
429     a23 = 1.0 / g4pow->A23(element->GetN());   << 394     a23 = 1.0/g4pow->A23(element->GetN());
430   }                                            << 395 
431   else {                                       << 396   } else {
432     G4double norm(0.0);                        << 397     for (G4int iel=0; iel<NumberOfElements; iel++)
433     for (G4int iel = 0; iel < NumberOfElements << 398       {
434       const G4Element* element = (*theElementV << 399         const G4Element* element = (*theElementVector)[iel];
435       const G4double weight = theAtomicNumDens << 400         const G4double weight = theAtomicNumDensityVector[iel];
436       norm += weight;                          << 401         norm += weight ;
437       z += element->GetZ() * weight;           << 402         z    += element->GetZ() * weight;
438       vF += element->GetIonisation()->GetFermi << 403         vF   += element->GetIonisation()->GetFermiVelocity() * weight;
439       lF += element->GetIonisation()->GetLFact << 404         lF   += element->GetIonisation()->GetLFactor() * weight;
440       a23 += weight / g4pow->A23(element->GetN << 405   a23  += weight/g4pow->A23(element->GetN());
441     }                                          << 406       }
442     if (norm > 0.0) { norm = 1.0/norm; }       << 407     z  /= norm;
443     z *= norm;                                 << 408     vF /= norm;
444     vF *= norm;                                << 409     lF /= norm;
445     lF *= norm;                                << 410     a23 /= norm;
446     a23 *= norm;                               << 411   }  
447   }                                            << 412   fZeff        = z;
448   fZeff = z;                                   << 413   fLfactor     = lF;
449   fLfactor = lF;                               << 414   fFermiEnergy = 25.*keV*vF*vF;
450   fFermiEnergy = 25. * CLHEP::keV * vF * vF;   << 415   fInvA23      = a23;
451   fInvA23 = a23;                               << 
452 }                                                 416 }
453                                                   417 
454 //....oooOO0OOooo........oooOO0OOooo........oo    418 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
455                                                   419 
456 void G4IonisParamMat::SetMeanExcitationEnergy(    420 void G4IonisParamMat::SetMeanExcitationEnergy(G4double value)
457 {                                                 421 {
458   if (value == fMeanExcitationEnergy || value  << 422   if(value == fMeanExcitationEnergy || value <= 0.0) { return; }
459     return;                                    << 
460   }                                            << 
461   if (G4NistManager::Instance()->GetVerbose()     423   if (G4NistManager::Instance()->GetVerbose() > 1) {
462     G4cout << "G4Material: Mean excitation ene << 424     G4cout << "G4Material: Mean excitation energy is changed for "
463            << " Iold= " << fMeanExcitationEner << 425            << fMaterial->GetName()
                                                   >> 426            << " Iold= " << fMeanExcitationEnergy/eV
                                                   >> 427            << "eV; Inew= " << value/eV << " eV;"
464            << G4endl;                             428            << G4endl;
465   }                                               429   }
466                                                << 430   
467   fMeanExcitationEnergy = value;                  431   fMeanExcitationEnergy = value;
468                                                   432 
469   // add corrections to density effect            433   // add corrections to density effect
470   G4double newlog = G4Log(value);                 434   G4double newlog = G4Log(value);
471   G4double corr = 2 * (newlog - fLogMeanExcEne << 435   G4double corr = 2*(newlog - fLogMeanExcEnergy);
472   fCdensity += corr;                           << 436   fCdensity  += corr;
473   fX0density += corr / twoln10;                << 437   fX0density += corr/twoln10;
474   fX1density += corr / twoln10;                << 438   fX1density += corr/twoln10;
475                                                   439 
476   // recompute parameters of fluctuation model    440   // recompute parameters of fluctuation model
477   fLogMeanExcEnergy = newlog;                     441   fLogMeanExcEnergy = newlog;
478   ComputeFluctModel();                            442   ComputeFluctModel();
479 }                                                 443 }
480                                                   444 
481 //....oooOO0OOooo........oooOO0OOooo........oo    445 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
482                                                   446 
483 void G4IonisParamMat::SetDensityEffectParamete << 447 void G4IonisParamMat::SetDensityEffectParameters(G4double cd, G4double md, 
484   G4double cd, G4double md, G4double ad, G4dou << 448              G4double ad, G4double x0, 
485 {                                              << 449              G4double x1, G4double d0)
486   // no check on consistence of user parameter << 450 {
487   G4AutoLock l(&ionisMutex);                   << 451   // no check on consistence of user parameters  
488   fCdensity = cd;                              << 452 #ifdef G4MULTITHREADED
489   fMdensity = md;                              << 453   G4MUTEXLOCK(&ionisMutex);
490   fAdensity = ad;                              << 454 #endif
                                                   >> 455   fCdensity  = cd; 
                                                   >> 456   fMdensity  = md;
                                                   >> 457   fAdensity  = ad;
491   fX0density = x0;                                458   fX0density = x0;
492   fX1density = x1;                                459   fX1density = x1;
493   fD0density = d0;                                460   fD0density = d0;
494   l.unlock();                                  << 461 #ifdef G4MULTITHREADED
                                                   >> 462   G4MUTEXUNLOCK(&ionisMutex);
                                                   >> 463 #endif
495 }                                                 464 }
496                                                   465 
497 //....oooOO0OOooo........oooOO0OOooo........oo    466 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
498                                                   467 
499 void G4IonisParamMat::SetDensityEffectParamete    468 void G4IonisParamMat::SetDensityEffectParameters(const G4Material* bmat)
500 {                                                 469 {
501   G4AutoLock l(&ionisMutex);                   << 470 #ifdef G4MULTITHREADED
                                                   >> 471   G4MUTEXLOCK(&ionisMutex);
                                                   >> 472 #endif
502   const G4IonisParamMat* ipm = bmat->GetIonisa    473   const G4IonisParamMat* ipm = bmat->GetIonisation();
503   fCdensity = ipm->GetCdensity();              << 474   fCdensity  = ipm->GetCdensity(); 
504   fMdensity = ipm->GetMdensity();              << 475   fMdensity  = ipm->GetMdensity();
505   fAdensity = ipm->GetAdensity();              << 476   fAdensity  = ipm->GetAdensity();
506   fX0density = ipm->GetX0density();               477   fX0density = ipm->GetX0density();
507   fX1density = ipm->GetX1density();               478   fX1density = ipm->GetX1density();
508   fD0density = ipm->GetD0density();               479   fD0density = ipm->GetD0density();
509                                                   480 
510   G4double corr = G4Log(bmat->GetDensity() / f << 481   G4double corr = G4Log(bmat->GetDensity()/fMaterial->GetDensity());
511   fCdensity += corr;                           << 482   fCdensity  += corr;
512   fX0density += corr / twoln10;                << 483   fX0density += corr/twoln10;
513   fX1density += corr / twoln10;                << 484   fX1density += corr/twoln10;
514   l.unlock();                                  << 485 #ifdef G4MULTITHREADED
515 }                                              << 486   G4MUTEXUNLOCK(&ionisMutex);
516                                                << 487 #endif
517 //....oooOO0OOooo........oooOO0OOooo........oo << 
518                                                << 
519 void G4IonisParamMat::ComputeDensityEffectOnFl << 
520 {                                              << 
521   if (val) {                                   << 
522     if (nullptr == fDensityEffectCalc) {       << 
523       G4int n = 0;                             << 
524       for (std::size_t i = 0; i < fMaterial->G << 
525         const G4int Z = fMaterial->GetElement( << 
526         n += G4AtomicShells::GetNumberOfShells << 
527       }                                        << 
528       // The last level is the conduction leve << 
529       // make a dummy conductor level with zer << 
530       fDensityEffectCalc = new G4DensityEffect << 
531     }                                          << 
532   }                                            << 
533   else {                                       << 
534     delete fDensityEffectCalc;                 << 
535     fDensityEffectCalc = nullptr;              << 
536   }                                            << 
537 }                                                 488 }
538                                                   489 
539 //....oooOO0OOooo........oooOO0OOooo........oo    490 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
540                                                   491 
541 G4double G4IonisParamMat::FindMeanExcitationEn    492 G4double G4IonisParamMat::FindMeanExcitationEnergy(const G4Material* mat) const
542 {                                                 493 {
543   G4double res = 0.0;                             494   G4double res = 0.0;
544   // data from density effect data                495   // data from density effect data
545   if (fDensityData != nullptr) {               << 496   if(fDensityData) {
546     G4int idx = fDensityData->GetIndex(mat->Ge    497     G4int idx = fDensityData->GetIndex(mat->GetName());
547     if (idx >= 0) {                            << 498     if(idx >= 0) {
548       res = fDensityData->GetMeanIonisationPot    499       res = fDensityData->GetMeanIonisationPotential(idx);
549     }                                             500     }
550   }                                               501   }
551                                                   502 
552   // The data on mean excitation energy for co    503   // The data on mean excitation energy for compaunds
553   // from "Stopping Powers for Electrons and P    504   // from "Stopping Powers for Electrons and Positrons"
554   // ICRU Report N#37, 1984  (energy in eV)       505   // ICRU Report N#37, 1984  (energy in eV)
555   // this value overwrites Density effect data    506   // this value overwrites Density effect data
556   G4String chFormula = mat->GetChemicalFormula    507   G4String chFormula = mat->GetChemicalFormula();
557   if (! chFormula.empty()) {                   << 508   if(chFormula != "") { 
558     static const size_t numberOfMolecula = 54; << 509 
559     // clang-format off                        << 510     static const size_t numberOfMolecula = 54; 
560     static const G4String name[numberOfMolecul    511     static const G4String name[numberOfMolecula] = {
561       // gas 0 - 12                               512       // gas 0 - 12
562       "NH_3",       "C_4H_10",     "CO_2",        513       "NH_3",       "C_4H_10",     "CO_2",      "C_2H_6",      "C_7H_16-Gas",
563       // "G4_AMMONIA", "G4_BUTANE","G4_CARBON_    514       // "G4_AMMONIA", "G4_BUTANE","G4_CARBON_DIOXIDE","G4_ETHANE", "G4_N-HEPTANE"
564       "C_6H_14-Gas",   "CH_4",     "NO",          515       "C_6H_14-Gas",   "CH_4",     "NO",        "N_2O",        "C_8H_18-Gas",      
565       // "G4_N-HEXANE" , "G4_METHANE", "x", "G    516       // "G4_N-HEXANE" , "G4_METHANE", "x", "G4_NITROUS_OXIDE", "G4_OCTANE"
566       "C_5H_12-Gas",   "C_3H_8",   "H_2O-Gas",    517       "C_5H_12-Gas",   "C_3H_8",   "H_2O-Gas",                        
567       // "G4_N-PENTANE", "G4_PROPANE", "G4_WAT    518       // "G4_N-PENTANE", "G4_PROPANE", "G4_WATER_VAPOR"
568                                                   519 
569       // liquid 13 - 39                           520       // liquid 13 - 39
570       "C_3H_6O",    "C_6H_5NH_2",  "C_6H_6",      521       "C_3H_6O",    "C_6H_5NH_2",  "C_6H_6",    "C_4H_9OH",    "CCl_4", 
571       //"G4_ACETONE","G4_ANILINE","G4_BENZENE" << 522       //"G4_ACETONE","G4_ANILINE","G4_BENZENE","G4_N-BUTYL_ALCOHOL","G4_CARBON_TETRACHLORIDE"  
572       "C_6H_5Cl",   "CHCl_3",      "C_6H_12",     523       "C_6H_5Cl",   "CHCl_3",      "C_6H_12",   "C_6H_4Cl_2",  "C_4Cl_2H_8O",
573       //"G4_CHLOROBENZENE","G4_CHLOROFORM","G4    524       //"G4_CHLOROBENZENE","G4_CHLOROFORM","G4_CYCLOHEXANE","G4_1,2-DICHLOROBENZENE",
574       //"G4_DICHLORODIETHYL_ETHER"                525       //"G4_DICHLORODIETHYL_ETHER"
575       "C_2Cl_2H_4", "(C_2H_5)_2O", "C_2H_5OH",    526       "C_2Cl_2H_4", "(C_2H_5)_2O", "C_2H_5OH",  "C_3H_5(OH)_3","C_7H_16",
576       //"G4_1,2-DICHLOROETHANE","G4_DIETHYL_ET    527       //"G4_1,2-DICHLOROETHANE","G4_DIETHYL_ETHER","G4_ETHYL_ALCOHOL","G4_GLYCEROL","G4_N-HEPTANE"
577       "C_6H_14",    "CH_3OH",      "C_6H_5NO_2    528       "C_6H_14",    "CH_3OH",      "C_6H_5NO_2","C_5H_12",     "C_3H_7OH",
578       //"G4_N-HEXANE","G4_METHANOL","G4_NITROB << 529       //"G4_N-HEXANE","G4_METHANOL","G4_NITROBENZENE","G4_N-PENTANE","G4_N-PROPYL_ALCOHOL",  
579       "C_5H_5N",    "C_8H_8",      "C_2Cl_4",     530       "C_5H_5N",    "C_8H_8",      "C_2Cl_4",   "C_7H_8",      "C_2Cl_3H",
580       //"G4_PYRIDINE","G4_POLYSTYRENE","G4_TET    531       //"G4_PYRIDINE","G4_POLYSTYRENE","G4_TETRACHLOROETHYLENE","G4_TOLUENE","G4_TRICHLOROETHYLENE"
581       "H_2O",       "C_8H_10",                    532       "H_2O",       "C_8H_10",             
582       // "G4_WATER", "G4_XYLENE"                  533       // "G4_WATER", "G4_XYLENE"
583                                                   534 
584       // solid 40 - 53                            535       // solid 40 - 53
585       "C_5H_5N_5",  "C_5H_5N_5O",  "(C_6H_11NO    536       "C_5H_5N_5",  "C_5H_5N_5O",  "(C_6H_11NO)-nylon",  "C_25H_52",
586       // "G4_ADENINE", "G4_GUANINE", "G4_NYLON    537       // "G4_ADENINE", "G4_GUANINE", "G4_NYLON-6-6", "G4_PARAFFIN"
587       "(C_2H_4)-Polyethylene",     "(C_5H_8O_2    538       "(C_2H_4)-Polyethylene",     "(C_5H_8O_2)-Polymethil_Methacrylate",
588       // "G4_ETHYLENE", "G4_PLEXIGLASS"           539       // "G4_ETHYLENE", "G4_PLEXIGLASS"
589       "(C_8H_8)-Polystyrene",      "A-150-tiss    540       "(C_8H_8)-Polystyrene",      "A-150-tissue",       "Al_2O_3",  "CaF_2",
590       // "G4_POLYSTYRENE", "G4_A-150_TISSUE",     541       // "G4_POLYSTYRENE", "G4_A-150_TISSUE", "G4_ALUMINUM_OXIDE", "G4_CALCIUM_FLUORIDE"
591       "LiF",        "Photo_Emulsion",  "(C_2F_    542       "LiF",        "Photo_Emulsion",  "(C_2F_4)-Teflon",  "SiO_2"
592       // "G4_LITHIUM_FLUORIDE", "G4_PHOTO_EMUL    543       // "G4_LITHIUM_FLUORIDE", "G4_PHOTO_EMULSION", "G4_TEFLON", "G4_SILICON_DIOXIDE"
593     } ;                                           544     } ;
594                                                << 545     
595     static const G4double meanExcitation[numbe    546     static const G4double meanExcitation[numberOfMolecula] = {
596                                                   547 
597       53.7,   48.3,  85.0,  45.4,  49.2,          548       53.7,   48.3,  85.0,  45.4,  49.2,
598       49.1,   41.7,  87.8,  84.9,  49.5,          549       49.1,   41.7,  87.8,  84.9,  49.5,
599       48.2,   47.1,  71.6,                        550       48.2,   47.1,  71.6,
600                                                   551 
601       64.2,   66.2,  63.4,  59.9,  166.3,         552       64.2,   66.2,  63.4,  59.9,  166.3,
602       89.1,  156.0,  56.4, 106.5,  103.3,         553       89.1,  156.0,  56.4, 106.5,  103.3, 
603       111.9,   60.0,  62.9,  72.6,   54.4,        554       111.9,   60.0,  62.9,  72.6,   54.4,  
604       54.0,  67.6,   75.8,  53.6,   61.1,         555       54.0,  67.6,   75.8,  53.6,   61.1,  
605       66.2,  64.0,  159.2,  62.5,  148.1,         556       66.2,  64.0,  159.2,  62.5,  148.1,  
606       75.0,  61.8,                                557       75.0,  61.8,
607                                                   558 
608       71.4,  75.0,   63.9,  48.3,   57.4,         559       71.4,  75.0,   63.9,  48.3,   57.4,
609       74.0,  68.7,   65.1, 145.2,  166.,          560       74.0,  68.7,   65.1, 145.2,  166.,
610       94.0, 331.0,   99.1, 139.2                  561       94.0, 331.0,   99.1, 139.2 
611     };                                            562     };
612     // clang-format on                         << 
613                                                   563 
614     for (std::size_t i = 0; i < numberOfMolecu << 564     for(size_t i=0; i<numberOfMolecula; i++) {
615       if (chFormula == name[i]) {              << 565       if(chFormula == name[i]) {
616         res = meanExcitation[i] * CLHEP::eV;   << 566   res = meanExcitation[i]*eV;
617         break;                                 << 567   break;
618       }                                           568       }
619     }                                             569     }
620   }                                               570   }
621   return res;                                     571   return res;
622 }                                                 572 }
                                                   >> 573 
                                                   >> 574 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
                                                   >> 575 
623                                                   576