Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/utils/src/G4EnergyLossTables.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/utils/src/G4EnergyLossTables.cc (Version 11.3.0) and /processes/electromagnetic/utils/src/G4EnergyLossTables.cc (Version 11.0.p4)


  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 //                                                 25 //
 26 //                                                 26 //
 27 //                                                 27 //
 28 // -------------------------------------------     28 // -------------------------------------------------------------------
 29 // first version created by P.Urban , 06/04/19     29 // first version created by P.Urban , 06/04/1998
 30 // modifications + "precise" functions added b     30 // modifications + "precise" functions added by L.Urban , 27/05/98
 31 // modifications , TOF functions , 26/10/98, L     31 // modifications , TOF functions , 26/10/98, L.Urban
 32 // cache mechanism in order to gain time, 11/0     32 // cache mechanism in order to gain time, 11/02/99, L.Urban
 33 // bug fixed , 12/04/99 , L.Urban                  33 // bug fixed , 12/04/99 , L.Urban
 34 // 10.11.99: moved from RWT hash dictionary to     34 // 10.11.99: moved from RWT hash dictionary to STL map, G.Barrand, M.Maire
 35 // 27.09.01 L.Urban , bug fixed (negative ener     35 // 27.09.01 L.Urban , bug fixed (negative energy deposit)
 36 // 26.10.01 all static functions moved from .i     36 // 26.10.01 all static functions moved from .icc files (mma)
 37 // 15.01.03 Add interfaces required for "cut p     37 // 15.01.03 Add interfaces required for "cut per region" (V.Ivanchenko)
 38 // 12.03.03 Add warnings to obsolete interface     38 // 12.03.03 Add warnings to obsolete interfaces (V.Ivanchenko)
 39 // 10.04.03 Add call to G4LossTableManager is      39 // 10.04.03 Add call to G4LossTableManager is particle is not registered (V.Ivanchenko)
 40 //                                                 40 //
 41 // -------------------------------------------     41 // -------------------------------------------------------------------
 42                                                    42 
 43 #include "G4EnergyLossTables.hh"                   43 #include "G4EnergyLossTables.hh"
 44 #include "G4SystemOfUnits.hh"                      44 #include "G4SystemOfUnits.hh"
 45 #include "G4MaterialCutsCouple.hh"                 45 #include "G4MaterialCutsCouple.hh"
 46 #include "G4RegionStore.hh"                        46 #include "G4RegionStore.hh"
 47 #include "G4LossTableManager.hh"                   47 #include "G4LossTableManager.hh"
 48                                                    48 
 49                                                    49 
 50 //....oooOO0OOooo........oooOO0OOooo........oo     50 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 51                                                    51 
 52 G4EnergyLossTablesHelper *G4EnergyLossTables:: <<  52 G4EnergyLossTablesHelper *G4EnergyLossTables::t = 0;
 53 G4EnergyLossTablesHelper *G4EnergyLossTables:: <<  53 G4EnergyLossTablesHelper *G4EnergyLossTables::null_loss = 0;
 54 G4ParticleDefinition* G4EnergyLossTables::last <<  54 G4ParticleDefinition* G4EnergyLossTables::lastParticle = 0;
 55 G4double G4EnergyLossTables::QQPositron = 1.0;     55 G4double G4EnergyLossTables::QQPositron = 1.0; // e_squared
 56 G4double G4EnergyLossTables::Chargesquare ;        56 G4double G4EnergyLossTables::Chargesquare ;
 57 G4int    G4EnergyLossTables::oldIndex = -1 ;       57 G4int    G4EnergyLossTables::oldIndex = -1 ;
 58 G4double G4EnergyLossTables::rmin = 0. ;           58 G4double G4EnergyLossTables::rmin = 0. ;
 59 G4double G4EnergyLossTables::rmax = 0. ;           59 G4double G4EnergyLossTables::rmax = 0. ;
 60 G4double G4EnergyLossTables::Thigh = 0. ;          60 G4double G4EnergyLossTables::Thigh = 0. ;
 61 G4int    G4EnergyLossTables::let_counter = 0;      61 G4int    G4EnergyLossTables::let_counter = 0;
 62 G4int    G4EnergyLossTables::let_max_num_warni     62 G4int    G4EnergyLossTables::let_max_num_warnings = 100;
 63 G4bool   G4EnergyLossTables::first_loss = true     63 G4bool   G4EnergyLossTables::first_loss = true;
 64                                                    64 
 65 G4EnergyLossTables::helper_map *G4EnergyLossTa <<  65 G4EnergyLossTables::helper_map *G4EnergyLossTables::dict = 0;
 66                                                    66 
 67 //....oooOO0OOooo........oooOO0OOooo........oo     67 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 68                                                    68 
 69 G4EnergyLossTablesHelper::G4EnergyLossTablesHe     69 G4EnergyLossTablesHelper::G4EnergyLossTablesHelper(
 70   const G4PhysicsTable* aDEDXTable,                70   const G4PhysicsTable* aDEDXTable,
 71   const G4PhysicsTable* aRangeTable,               71   const G4PhysicsTable* aRangeTable,
 72   const G4PhysicsTable* anInverseRangeTable,       72   const G4PhysicsTable* anInverseRangeTable,
 73   const G4PhysicsTable* aLabTimeTable,             73   const G4PhysicsTable* aLabTimeTable,
 74   const G4PhysicsTable* aProperTimeTable,          74   const G4PhysicsTable* aProperTimeTable,
 75   G4double aLowestKineticEnergy,                   75   G4double aLowestKineticEnergy,
 76   G4double aHighestKineticEnergy,                  76   G4double aHighestKineticEnergy,
 77   G4double aMassRatio,                             77   G4double aMassRatio,
 78   G4int aNumberOfBins)                             78   G4int aNumberOfBins)
 79   :                                                79   :
 80   theDEDXTable(aDEDXTable), theRangeTable(aRan     80   theDEDXTable(aDEDXTable), theRangeTable(aRangeTable),
 81   theInverseRangeTable(anInverseRangeTable),       81   theInverseRangeTable(anInverseRangeTable),
 82   theLabTimeTable(aLabTimeTable),                  82   theLabTimeTable(aLabTimeTable),
 83   theProperTimeTable(aProperTimeTable),            83   theProperTimeTable(aProperTimeTable),
 84   theLowestKineticEnergy(aLowestKineticEnergy)     84   theLowestKineticEnergy(aLowestKineticEnergy),
 85   theHighestKineticEnergy(aHighestKineticEnerg     85   theHighestKineticEnergy(aHighestKineticEnergy),
 86   theMassRatio(aMassRatio),                        86   theMassRatio(aMassRatio),
 87   theNumberOfBins(aNumberOfBins)                   87   theNumberOfBins(aNumberOfBins)
 88 { }                                                88 { }
 89                                                    89 
 90 //....oooOO0OOooo........oooOO0OOooo........oo     90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 91                                                    91 
 92 G4EnergyLossTablesHelper::G4EnergyLossTablesHe     92 G4EnergyLossTablesHelper::G4EnergyLossTablesHelper()
 93 {                                                  93 { 
 94   theLowestKineticEnergy = 0.0;                    94   theLowestKineticEnergy = 0.0;
 95   theHighestKineticEnergy= 0.0;                    95   theHighestKineticEnergy= 0.0;
 96   theMassRatio = 0.0;                              96   theMassRatio = 0.0;
 97   theNumberOfBins = 0;                             97   theNumberOfBins = 0;
 98   theDEDXTable = theRangeTable = theInverseRan     98   theDEDXTable = theRangeTable = theInverseRangeTable = theLabTimeTable 
 99     = theProperTimeTable = nullptr;                99     = theProperTimeTable = nullptr;
100 }                                                 100 }
101                                                   101 
102 //....oooOO0OOooo........oooOO0OOooo........oo    102 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
103                                                   103 
104 void G4EnergyLossTables::Register(                104 void G4EnergyLossTables::Register(
105   const G4ParticleDefinition* p,                  105   const G4ParticleDefinition* p,
106   const G4PhysicsTable* tDEDX,                    106   const G4PhysicsTable* tDEDX,
107   const G4PhysicsTable* tRange,                   107   const G4PhysicsTable* tRange,
108   const G4PhysicsTable* tInverseRange,            108   const G4PhysicsTable* tInverseRange,
109   const G4PhysicsTable* tLabTime,                 109   const G4PhysicsTable* tLabTime,
110   const G4PhysicsTable* tProperTime,              110   const G4PhysicsTable* tProperTime,
111   G4double lowestKineticEnergy,                   111   G4double lowestKineticEnergy,
112   G4double highestKineticEnergy,                  112   G4double highestKineticEnergy,
113   G4double massRatio,                             113   G4double massRatio,
114   G4int NumberOfBins)                             114   G4int NumberOfBins)
115 {                                                 115 {
116   if (!dict) dict = new G4EnergyLossTables::he    116   if (!dict) dict = new G4EnergyLossTables::helper_map;
117   if (!null_loss) null_loss = new G4EnergyLoss    117   if (!null_loss) null_loss = new G4EnergyLossTablesHelper;
118   if (!t) t = new G4EnergyLossTablesHelper;       118   if (!t) t = new G4EnergyLossTablesHelper;
119                                                   119 
120   (*dict)[p]= G4EnergyLossTablesHelper(tDEDX,     120   (*dict)[p]= G4EnergyLossTablesHelper(tDEDX, tRange,tInverseRange,
121                     tLabTime,tProperTime,lowes    121                     tLabTime,tProperTime,lowestKineticEnergy,
122         highestKineticEnergy, massRatio,Number    122         highestKineticEnergy, massRatio,NumberOfBins);
123                                                   123 
124   *t = GetTables(p) ;    // important for cach    124   *t = GetTables(p) ;    // important for cache !!!!!
125   lastParticle = (G4ParticleDefinition*) p ;      125   lastParticle = (G4ParticleDefinition*) p ;
126   Chargesquare = (p->GetPDGCharge())*(p->GetPD    126   Chargesquare = (p->GetPDGCharge())*(p->GetPDGCharge())/
127                   QQPositron ;                    127                   QQPositron ;
128   if (first_loss ) {                              128   if (first_loss ) {
129     *null_loss = G4EnergyLossTablesHelper(     << 129     *null_loss = G4EnergyLossTablesHelper(0, 0, 0, 0, 0, 0.0, 0.0, 0.0, 0);
130                  nullptr, nullptr, nullptr, nu << 
131     first_loss = false;                           130     first_loss = false;
132   }                                               131   }
133 }                                                 132 }
134                                                   133 
135 //....oooOO0OOooo........oooOO0OOooo........oo    134 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
136                                                   135 
137 const G4PhysicsTable* G4EnergyLossTables::GetD    136 const G4PhysicsTable* G4EnergyLossTables::GetDEDXTable(
138   const G4ParticleDefinition* p)                  137   const G4ParticleDefinition* p)
139 {                                                 138 {
140   if (!dict) dict = new G4EnergyLossTables::he    139   if (!dict) dict = new G4EnergyLossTables::helper_map;
141   helper_map::iterator it;                        140   helper_map::iterator it;
142   if((it=dict->find(p))==dict->end()) return n << 141   if((it=dict->find(p))==dict->end()) return 0;
143   return (*it).second.theDEDXTable;               142   return (*it).second.theDEDXTable;
144 }                                                 143 }
145                                                   144 
146 //....oooOO0OOooo........oooOO0OOooo........oo    145 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
147                                                   146 
148 const G4PhysicsTable* G4EnergyLossTables::GetR    147 const G4PhysicsTable* G4EnergyLossTables::GetRangeTable(
149   const G4ParticleDefinition* p)                  148   const G4ParticleDefinition* p)
150 {                                                 149 {
151   if (!dict) dict = new G4EnergyLossTables::he    150   if (!dict) dict = new G4EnergyLossTables::helper_map;
152   helper_map::iterator it;                        151   helper_map::iterator it;
153   if((it=dict->find(p))==dict->end()) return n << 152   if((it=dict->find(p))==dict->end()) return 0;
154   return (*it).second.theRangeTable;              153   return (*it).second.theRangeTable;
155 }                                                 154 }
156                                                   155 
157 //....oooOO0OOooo........oooOO0OOooo........oo    156 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
158                                                   157 
159 const G4PhysicsTable* G4EnergyLossTables::GetI    158 const G4PhysicsTable* G4EnergyLossTables::GetInverseRangeTable(
160   const G4ParticleDefinition* p)                  159   const G4ParticleDefinition* p)
161 {                                                 160 {
162   if (!dict) dict = new G4EnergyLossTables::he    161   if (!dict) dict = new G4EnergyLossTables::helper_map;
163   helper_map::iterator it;                        162   helper_map::iterator it;
164   if((it=dict->find(p))==dict->end()) return n << 163   if((it=dict->find(p))==dict->end()) return 0;
165   return (*it).second.theInverseRangeTable;       164   return (*it).second.theInverseRangeTable;
166 }                                                 165 }
167                                                   166 
168 //....oooOO0OOooo........oooOO0OOooo........oo    167 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
169                                                   168 
170 const G4PhysicsTable* G4EnergyLossTables::GetL    169 const G4PhysicsTable* G4EnergyLossTables::GetLabTimeTable(
171   const G4ParticleDefinition* p)                  170   const G4ParticleDefinition* p)
172 {                                                 171 {
173   if (!dict) dict = new G4EnergyLossTables::he    172   if (!dict) dict = new G4EnergyLossTables::helper_map;
174   helper_map::iterator it;                        173   helper_map::iterator it;
175   if((it=dict->find(p))==dict->end()) return n << 174   if((it=dict->find(p))==dict->end()) return 0;
176   return (*it).second.theLabTimeTable;            175   return (*it).second.theLabTimeTable;
177 }                                                 176 }
178                                                   177 
179 //....oooOO0OOooo........oooOO0OOooo........oo    178 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
180                                                   179 
181 const G4PhysicsTable* G4EnergyLossTables::GetP    180 const G4PhysicsTable* G4EnergyLossTables::GetProperTimeTable(
182   const G4ParticleDefinition* p)                  181   const G4ParticleDefinition* p)
183 {                                                 182 {
184   if (!dict) dict = new G4EnergyLossTables::he    183   if (!dict) dict = new G4EnergyLossTables::helper_map;
185   helper_map::iterator it;                        184   helper_map::iterator it;
186   if((it=dict->find(p))==dict->end()) return n << 185   if((it=dict->find(p))==dict->end()) return 0;
187   return (*it).second.theProperTimeTable;         186   return (*it).second.theProperTimeTable;
188 }                                                 187 }
189                                                   188 
190 //....oooOO0OOooo........oooOO0OOooo........oo    189 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
191                                                   190 
192 G4EnergyLossTablesHelper G4EnergyLossTables::G    191 G4EnergyLossTablesHelper G4EnergyLossTables::GetTables(
193   const G4ParticleDefinition* p)                  192   const G4ParticleDefinition* p)
194 {                                                 193 {
195   if (!dict) dict = new G4EnergyLossTables::he    194   if (!dict) dict = new G4EnergyLossTables::helper_map;
196   if (!null_loss) null_loss = new G4EnergyLoss    195   if (!null_loss) null_loss = new G4EnergyLossTablesHelper;
197                                                   196 
198   helper_map::iterator it;                        197   helper_map::iterator it;
199   if ((it=dict->find(p))==dict->end()) {          198   if ((it=dict->find(p))==dict->end()) {
200     return *null_loss;                            199     return *null_loss;
201   }                                               200   }
202   return (*it).second;                            201   return (*it).second;
203 }                                                 202 }
204                                                   203 
205 //....oooOO0OOooo........oooOO0OOooo........oo    204 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
206                                                   205 
207 G4double G4EnergyLossTables::GetDEDX(             206 G4double G4EnergyLossTables::GetDEDX(
208     const G4ParticleDefinition *aParticle,        207     const G4ParticleDefinition *aParticle,
209     G4double KineticEnergy,                       208     G4double KineticEnergy,
210     const G4Material *aMaterial)                  209     const G4Material *aMaterial)
211 {                                                 210 {
212   if (!t) t = new G4EnergyLossTablesHelper;       211   if (!t) t = new G4EnergyLossTablesHelper;
213                                                   212 
214   CPRWarning();                                   213   CPRWarning();
215   if(aParticle != (const G4ParticleDefinition*    214   if(aParticle != (const G4ParticleDefinition*) lastParticle)
216   {                                               215   {
217     *t= GetTables(aParticle);                     216     *t= GetTables(aParticle);
218     lastParticle = (G4ParticleDefinition*) aPa    217     lastParticle = (G4ParticleDefinition*) aParticle ;
219     Chargesquare = (aParticle->GetPDGCharge())    218     Chargesquare = (aParticle->GetPDGCharge())*
220                    (aParticle->GetPDGCharge())    219                    (aParticle->GetPDGCharge())/
221                    QQPositron ;                   220                    QQPositron ;
222     oldIndex = -1 ;                               221     oldIndex = -1 ;
223   }                                               222   }
224   const G4PhysicsTable*  dEdxTable= t->theDEDX    223   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
225   if (!dEdxTable) {                               224   if (!dEdxTable) {
226     ParticleHaveNoLoss(aParticle,"dEdx");         225     ParticleHaveNoLoss(aParticle,"dEdx");
227     return 0.0;                                   226     return 0.0;
228   }                                               227   }
229                                                   228 
230   G4int materialIndex = (G4int)aMaterial->GetI << 229   G4int materialIndex = aMaterial->GetIndex();
231   G4double scaledKineticEnergy = KineticEnergy    230   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
232   G4double dEdx;                                  231   G4double dEdx;
233   G4bool isOut;                                   232   G4bool isOut;
234                                                   233 
235   if (scaledKineticEnergy<t->theLowestKineticE    234   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
236                                                   235 
237      dEdx =(*dEdxTable)(materialIndex)->GetVal    236      dEdx =(*dEdxTable)(materialIndex)->GetValue(
238               t->theLowestKineticEnergy,isOut)    237               t->theLowestKineticEnergy,isOut)
239            *std::sqrt(scaledKineticEnergy/t->t    238            *std::sqrt(scaledKineticEnergy/t->theLowestKineticEnergy);
240                                                   239 
241   } else if (scaledKineticEnergy>t->theHighest    240   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
242                                                   241 
243      dEdx = (*dEdxTable)(materialIndex)->GetVa    242      dEdx = (*dEdxTable)(materialIndex)->GetValue(
244         t->theHighestKineticEnergy,isOut);        243         t->theHighestKineticEnergy,isOut);
245                                                   244 
246   } else {                                        245   } else {
247                                                   246 
248     dEdx = (*dEdxTable)(materialIndex)->GetVal    247     dEdx = (*dEdxTable)(materialIndex)->GetValue(
249          scaledKineticEnergy,isOut);              248          scaledKineticEnergy,isOut);
250                                                   249 
251   }                                               250   }
252                                                   251 
253   return dEdx*Chargesquare;                       252   return dEdx*Chargesquare;
254 }                                                 253 }
255                                                   254 
256 //....oooOO0OOooo........oooOO0OOooo........oo    255 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
257                                                   256 
258 G4double G4EnergyLossTables::GetLabTime(          257 G4double G4EnergyLossTables::GetLabTime(
259     const G4ParticleDefinition *aParticle,        258     const G4ParticleDefinition *aParticle,
260     G4double KineticEnergy,                       259     G4double KineticEnergy,
261     const G4Material *aMaterial)                  260     const G4Material *aMaterial)
262 {                                                 261 {
263   if (!t) t = new G4EnergyLossTablesHelper;       262   if (!t) t = new G4EnergyLossTablesHelper;
264                                                   263 
265   CPRWarning();                                   264   CPRWarning();
266   if(aParticle != (const G4ParticleDefinition*    265   if(aParticle != (const G4ParticleDefinition*) lastParticle)
267   {                                               266   {
268     *t= GetTables(aParticle);                     267     *t= GetTables(aParticle);
269     lastParticle = (G4ParticleDefinition*) aPa    268     lastParticle = (G4ParticleDefinition*) aParticle ;
270     oldIndex = -1 ;                               269     oldIndex = -1 ;
271   }                                               270   }
272   const G4PhysicsTable* labtimeTable= t->theLa    271   const G4PhysicsTable* labtimeTable= t->theLabTimeTable;
273   if (!labtimeTable) {                            272   if (!labtimeTable) {
274     ParticleHaveNoLoss(aParticle,"LabTime");      273     ParticleHaveNoLoss(aParticle,"LabTime");
275     return 0.0;                                   274     return 0.0;
276   }                                               275   }
277                                                   276 
278   const G4double parlowen=0.4 , ppar=0.5-parlo    277   const G4double parlowen=0.4 , ppar=0.5-parlowen ;
279   G4int materialIndex = (G4int)aMaterial->GetI << 278   G4int materialIndex = aMaterial->GetIndex();
280   G4double scaledKineticEnergy = KineticEnergy    279   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
281   G4double time;                                  280   G4double time;
282   G4bool isOut;                                   281   G4bool isOut;
283                                                   282 
284   if (scaledKineticEnergy<t->theLowestKineticE    283   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
285                                                   284 
286      time = std::exp(ppar*std::log(scaledKinet    285      time = std::exp(ppar*std::log(scaledKineticEnergy/t->theLowestKineticEnergy))*
287             (*labtimeTable)(materialIndex)->Ge    286             (*labtimeTable)(materialIndex)->GetValue(
288               t->theLowestKineticEnergy,isOut)    287               t->theLowestKineticEnergy,isOut);
289                                                   288 
290                                                   289 
291   } else if (scaledKineticEnergy>t->theHighest    290   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
292                                                   291 
293      time = (*labtimeTable)(materialIndex)->Ge    292      time = (*labtimeTable)(materialIndex)->GetValue(
294               t->theHighestKineticEnergy,isOut    293               t->theHighestKineticEnergy,isOut);
295                                                   294 
296   } else {                                        295   } else {
297                                                   296 
298     time = (*labtimeTable)(materialIndex)->Get    297     time = (*labtimeTable)(materialIndex)->GetValue(
299                scaledKineticEnergy,isOut);        298                scaledKineticEnergy,isOut);
300                                                   299 
301   }                                               300   }
302                                                   301 
303   return time/t->theMassRatio ;                   302   return time/t->theMassRatio ;
304 }                                                 303 }
305                                                   304 
306 //....oooOO0OOooo........oooOO0OOooo........oo    305 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
307                                                   306 
308 G4double G4EnergyLossTables::GetDeltaLabTime(     307 G4double G4EnergyLossTables::GetDeltaLabTime(
309     const G4ParticleDefinition *aParticle,        308     const G4ParticleDefinition *aParticle,
310     G4double KineticEnergyStart,                  309     G4double KineticEnergyStart,
311     G4double KineticEnergyEnd,                    310     G4double KineticEnergyEnd,
312     const G4Material *aMaterial)                  311     const G4Material *aMaterial)
313 {                                                 312 {
314   if (!t) t = new G4EnergyLossTablesHelper;       313   if (!t) t = new G4EnergyLossTablesHelper;
315                                                   314 
316   CPRWarning();                                   315   CPRWarning();
317   if(aParticle != (const G4ParticleDefinition*    316   if(aParticle != (const G4ParticleDefinition*) lastParticle)
318   {                                               317   {
319     *t= GetTables(aParticle);                     318     *t= GetTables(aParticle);
320     lastParticle = (G4ParticleDefinition*) aPa    319     lastParticle = (G4ParticleDefinition*) aParticle ;
321     oldIndex = -1 ;                               320     oldIndex = -1 ;
322   }                                               321   }
323   const G4PhysicsTable* labtimeTable= t->theLa    322   const G4PhysicsTable* labtimeTable= t->theLabTimeTable;
324   if (!labtimeTable) {                            323   if (!labtimeTable) {
325     ParticleHaveNoLoss(aParticle,"LabTime");      324     ParticleHaveNoLoss(aParticle,"LabTime");
326     return 0.0;                                   325     return 0.0;
327   }                                               326   }
328                                                   327 
329   const G4double parlowen=0.4 , ppar=0.5-parlo    328   const G4double parlowen=0.4 , ppar=0.5-parlowen ;
330   const G4double dToverT = 0.05 , facT = 1. -d    329   const G4double dToverT = 0.05 , facT = 1. -dToverT ;
331   G4double timestart,timeend,deltatime,dTT;       330   G4double timestart,timeend,deltatime,dTT;
332   G4bool isOut;                                   331   G4bool isOut;
333                                                   332 
334   G4int materialIndex = (G4int)aMaterial->GetI << 333   G4int materialIndex = aMaterial->GetIndex();
335   G4double scaledKineticEnergy = KineticEnergy    334   G4double scaledKineticEnergy = KineticEnergyStart*t->theMassRatio;
336                                                   335 
337   if (scaledKineticEnergy<t->theLowestKineticE    336   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
338                                                   337 
339      timestart = std::exp(ppar*std::log(scaled    338      timestart = std::exp(ppar*std::log(scaledKineticEnergy/t->theLowestKineticEnergy))*
340                 (*labtimeTable)(materialIndex)    339                 (*labtimeTable)(materialIndex)->GetValue(
341                 t->theLowestKineticEnergy,isOu    340                 t->theLowestKineticEnergy,isOut);
342                                                   341 
343                                                   342 
344   } else if (scaledKineticEnergy>t->theHighest    343   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
345                                                   344 
346      timestart = (*labtimeTable)(materialIndex    345      timestart = (*labtimeTable)(materialIndex)->GetValue(
347                 t->theHighestKineticEnergy,isO    346                 t->theHighestKineticEnergy,isOut);
348                                                   347 
349   } else {                                        348   } else {
350                                                   349 
351     timestart = (*labtimeTable)(materialIndex)    350     timestart = (*labtimeTable)(materialIndex)->GetValue(
352                 scaledKineticEnergy,isOut);       351                 scaledKineticEnergy,isOut);
353                                                   352 
354   }                                               353   }
355                                                   354 
356   dTT = (KineticEnergyStart - KineticEnergyEnd    355   dTT = (KineticEnergyStart - KineticEnergyEnd)/KineticEnergyStart ;
357                                                   356 
358   if( dTT < dToverT )                             357   if( dTT < dToverT )
359     scaledKineticEnergy = facT*KineticEnergySt    358     scaledKineticEnergy = facT*KineticEnergyStart*t->theMassRatio;
360   else                                            359   else
361     scaledKineticEnergy = KineticEnergyEnd*t->    360     scaledKineticEnergy = KineticEnergyEnd*t->theMassRatio;
362                                                   361 
363   if (scaledKineticEnergy<t->theLowestKineticE    362   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
364                                                   363 
365      timeend = std::exp(ppar*std::log(scaledKi    364      timeend = std::exp(ppar*std::log(scaledKineticEnergy/t->theLowestKineticEnergy))*
366                 (*labtimeTable)(materialIndex)    365                 (*labtimeTable)(materialIndex)->GetValue(
367                 t->theLowestKineticEnergy,isOu    366                 t->theLowestKineticEnergy,isOut);
368                                                   367 
369                                                   368 
370   } else if (scaledKineticEnergy>t->theHighest    369   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
371                                                   370 
372      timeend = (*labtimeTable)(materialIndex)-    371      timeend = (*labtimeTable)(materialIndex)->GetValue(
373                 t->theHighestKineticEnergy,isO    372                 t->theHighestKineticEnergy,isOut);
374                                                   373 
375   } else {                                        374   } else {
376                                                   375 
377     timeend = (*labtimeTable)(materialIndex)->    376     timeend = (*labtimeTable)(materialIndex)->GetValue(
378                 scaledKineticEnergy,isOut);       377                 scaledKineticEnergy,isOut);
379                                                   378 
380   }                                               379   }
381                                                   380 
382   deltatime = timestart - timeend ;               381   deltatime = timestart - timeend ;
383                                                   382 
384   if( dTT < dToverT )                             383   if( dTT < dToverT )
385     deltatime *= dTT/dToverT;                     384     deltatime *= dTT/dToverT;
386                                                   385 
387   return deltatime/t->theMassRatio ;              386   return deltatime/t->theMassRatio ;
388 }                                                 387 }
389                                                   388 
390 //....oooOO0OOooo........oooOO0OOooo........oo    389 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
391                                                   390 
392 G4double G4EnergyLossTables::GetProperTime(       391 G4double G4EnergyLossTables::GetProperTime(
393     const G4ParticleDefinition *aParticle,        392     const G4ParticleDefinition *aParticle,
394     G4double KineticEnergy,                       393     G4double KineticEnergy,
395     const G4Material *aMaterial)                  394     const G4Material *aMaterial)
396 {                                                 395 {
397   if (!t) t = new G4EnergyLossTablesHelper;       396   if (!t) t = new G4EnergyLossTablesHelper;
398                                                   397 
399   CPRWarning();                                   398   CPRWarning();
400   if(aParticle != (const G4ParticleDefinition*    399   if(aParticle != (const G4ParticleDefinition*) lastParticle)
401   {                                               400   {
402     *t= GetTables(aParticle);                     401     *t= GetTables(aParticle);
403     lastParticle = (G4ParticleDefinition*) aPa    402     lastParticle = (G4ParticleDefinition*) aParticle ;
404     oldIndex = -1 ;                               403     oldIndex = -1 ;
405   }                                               404   }
406   const G4PhysicsTable* propertimeTable= t->th    405   const G4PhysicsTable* propertimeTable= t->theProperTimeTable;
407   if (!propertimeTable) {                         406   if (!propertimeTable) {
408     ParticleHaveNoLoss(aParticle,"ProperTime")    407     ParticleHaveNoLoss(aParticle,"ProperTime");
409     return 0.0;                                   408     return 0.0;
410   }                                               409   }
411                                                   410 
412   const G4double parlowen=0.4 , ppar=0.5-parlo    411   const G4double parlowen=0.4 , ppar=0.5-parlowen ;
413   G4int materialIndex = (G4int)aMaterial->GetI << 412   G4int materialIndex = aMaterial->GetIndex();
414   G4double scaledKineticEnergy = KineticEnergy    413   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
415   G4double time;                                  414   G4double time;
416   G4bool isOut;                                   415   G4bool isOut;
417                                                   416 
418   if (scaledKineticEnergy<t->theLowestKineticE    417   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
419                                                   418 
420      time = std::exp(ppar*std::log(scaledKinet    419      time = std::exp(ppar*std::log(scaledKineticEnergy/t->theLowestKineticEnergy))*
421             (*propertimeTable)(materialIndex)-    420             (*propertimeTable)(materialIndex)->GetValue(
422               t->theLowestKineticEnergy,isOut)    421               t->theLowestKineticEnergy,isOut);
423                                                   422 
424                                                   423 
425   } else if (scaledKineticEnergy>t->theHighest    424   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
426                                                   425 
427      time = (*propertimeTable)(materialIndex)-    426      time = (*propertimeTable)(materialIndex)->GetValue(
428               t->theHighestKineticEnergy,isOut    427               t->theHighestKineticEnergy,isOut);
429                                                   428 
430   } else {                                        429   } else {
431                                                   430 
432     time = (*propertimeTable)(materialIndex)->    431     time = (*propertimeTable)(materialIndex)->GetValue(
433                scaledKineticEnergy,isOut);        432                scaledKineticEnergy,isOut);
434                                                   433 
435   }                                               434   }
436                                                   435 
437   return time/t->theMassRatio ;                   436   return time/t->theMassRatio ;
438 }                                                 437 }
439                                                   438 
440 //....oooOO0OOooo........oooOO0OOooo........oo    439 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
441                                                   440 
442 G4double G4EnergyLossTables::GetDeltaProperTim    441 G4double G4EnergyLossTables::GetDeltaProperTime(
443     const G4ParticleDefinition *aParticle,        442     const G4ParticleDefinition *aParticle,
444     G4double KineticEnergyStart,                  443     G4double KineticEnergyStart,
445     G4double KineticEnergyEnd,                    444     G4double KineticEnergyEnd,
446     const G4Material *aMaterial)                  445     const G4Material *aMaterial)
447 {                                                 446 {
448   if (!t) t = new G4EnergyLossTablesHelper;       447   if (!t) t = new G4EnergyLossTablesHelper;
449                                                   448 
450   CPRWarning();                                   449   CPRWarning();
451   if(aParticle != (const G4ParticleDefinition*    450   if(aParticle != (const G4ParticleDefinition*) lastParticle)
452   {                                               451   {
453     *t= GetTables(aParticle);                     452     *t= GetTables(aParticle);
454     lastParticle = (G4ParticleDefinition*) aPa    453     lastParticle = (G4ParticleDefinition*) aParticle ;
455     oldIndex = -1 ;                               454     oldIndex = -1 ;
456   }                                               455   }
457   const G4PhysicsTable* propertimeTable= t->th    456   const G4PhysicsTable* propertimeTable= t->theProperTimeTable;
458   if (!propertimeTable) {                         457   if (!propertimeTable) {
459     ParticleHaveNoLoss(aParticle,"ProperTime")    458     ParticleHaveNoLoss(aParticle,"ProperTime");
460     return 0.0;                                   459     return 0.0;
461   }                                               460   }
462                                                   461 
463   const G4double parlowen=0.4 , ppar=0.5-parlo    462   const G4double parlowen=0.4 , ppar=0.5-parlowen ;
464   const G4double dToverT = 0.05 , facT = 1. -d    463   const G4double dToverT = 0.05 , facT = 1. -dToverT ;
465   G4double timestart,timeend,deltatime,dTT;       464   G4double timestart,timeend,deltatime,dTT;
466   G4bool isOut;                                   465   G4bool isOut;
467                                                   466 
468   G4int materialIndex = (G4int)aMaterial->GetI << 467   G4int materialIndex = aMaterial->GetIndex();
469   G4double scaledKineticEnergy = KineticEnergy    468   G4double scaledKineticEnergy = KineticEnergyStart*t->theMassRatio;
470                                                   469 
471   if (scaledKineticEnergy<t->theLowestKineticE    470   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
472                                                   471 
473      timestart = std::exp(ppar*std::log(scaled    472      timestart = std::exp(ppar*std::log(scaledKineticEnergy/t->theLowestKineticEnergy))*
474                 (*propertimeTable)(materialInd    473                 (*propertimeTable)(materialIndex)->GetValue(
475                 t->theLowestKineticEnergy,isOu    474                 t->theLowestKineticEnergy,isOut);
476                                                   475 
477                                                   476 
478   } else if (scaledKineticEnergy>t->theHighest    477   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
479                                                   478 
480      timestart = (*propertimeTable)(materialIn    479      timestart = (*propertimeTable)(materialIndex)->GetValue(
481                 t->theHighestKineticEnergy,isO    480                 t->theHighestKineticEnergy,isOut);
482                                                   481 
483   } else {                                        482   } else {
484                                                   483 
485     timestart = (*propertimeTable)(materialInd    484     timestart = (*propertimeTable)(materialIndex)->GetValue(
486                 scaledKineticEnergy,isOut);       485                 scaledKineticEnergy,isOut);
487                                                   486 
488   }                                               487   }
489                                                   488 
490   dTT = (KineticEnergyStart - KineticEnergyEnd    489   dTT = (KineticEnergyStart - KineticEnergyEnd)/KineticEnergyStart ;
491                                                   490 
492   if( dTT < dToverT )                             491   if( dTT < dToverT )
493     scaledKineticEnergy = facT*KineticEnergySt    492     scaledKineticEnergy = facT*KineticEnergyStart*t->theMassRatio;
494   else                                            493   else
495     scaledKineticEnergy = KineticEnergyEnd*t->    494     scaledKineticEnergy = KineticEnergyEnd*t->theMassRatio;
496                                                   495 
497   if (scaledKineticEnergy<t->theLowestKineticE    496   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
498                                                   497 
499      timeend = std::exp(ppar*std::log(scaledKi    498      timeend = std::exp(ppar*std::log(scaledKineticEnergy/t->theLowestKineticEnergy))*
500                 (*propertimeTable)(materialInd    499                 (*propertimeTable)(materialIndex)->GetValue(
501                 t->theLowestKineticEnergy,isOu    500                 t->theLowestKineticEnergy,isOut);
502                                                   501 
503                                                   502 
504   } else if (scaledKineticEnergy>t->theHighest    503   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
505                                                   504 
506      timeend = (*propertimeTable)(materialInde    505      timeend = (*propertimeTable)(materialIndex)->GetValue(
507                 t->theHighestKineticEnergy,isO    506                 t->theHighestKineticEnergy,isOut);
508                                                   507 
509   } else {                                        508   } else {
510                                                   509 
511     timeend = (*propertimeTable)(materialIndex    510     timeend = (*propertimeTable)(materialIndex)->GetValue(
512                 scaledKineticEnergy,isOut);       511                 scaledKineticEnergy,isOut);
513                                                   512 
514   }                                               513   }
515                                                   514 
516   deltatime = timestart - timeend ;               515   deltatime = timestart - timeend ;
517                                                   516 
518   if( dTT < dToverT )                             517   if( dTT < dToverT )
519     deltatime *= dTT/dToverT ;                    518     deltatime *= dTT/dToverT ;
520                                                   519 
521   return deltatime/t->theMassRatio ;              520   return deltatime/t->theMassRatio ;
522 }                                                 521 }
523                                                   522 
524 //....oooOO0OOooo........oooOO0OOooo........oo    523 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
525                                                   524 
526 G4double G4EnergyLossTables::GetRange(            525 G4double G4EnergyLossTables::GetRange(
527     const G4ParticleDefinition *aParticle,        526     const G4ParticleDefinition *aParticle,
528     G4double KineticEnergy,                       527     G4double KineticEnergy,
529     const G4Material *aMaterial)                  528     const G4Material *aMaterial)
530 {                                                 529 {
531   if (!t) t = new G4EnergyLossTablesHelper;       530   if (!t) t = new G4EnergyLossTablesHelper;
532                                                   531 
533   CPRWarning();                                   532   CPRWarning();
534   if(aParticle != (const G4ParticleDefinition*    533   if(aParticle != (const G4ParticleDefinition*) lastParticle)
535   {                                               534   {
536     *t= GetTables(aParticle);                     535     *t= GetTables(aParticle);
537     lastParticle = (G4ParticleDefinition*) aPa    536     lastParticle = (G4ParticleDefinition*) aParticle ;
538     Chargesquare = (aParticle->GetPDGCharge())    537     Chargesquare = (aParticle->GetPDGCharge())*
539                    (aParticle->GetPDGCharge())    538                    (aParticle->GetPDGCharge())/
540                     QQPositron ;                  539                     QQPositron ;
541     oldIndex = -1 ;                               540     oldIndex = -1 ;
542   }                                               541   }
543   const G4PhysicsTable* rangeTable= t->theRang    542   const G4PhysicsTable* rangeTable= t->theRangeTable;
544   const G4PhysicsTable*  dEdxTable= t->theDEDX    543   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
545   if (!rangeTable) {                              544   if (!rangeTable) {
546     ParticleHaveNoLoss(aParticle,"Range");        545     ParticleHaveNoLoss(aParticle,"Range");
547     return 0.0;                                   546     return 0.0;
548   }                                               547   }
549                                                   548 
550   G4int materialIndex = (G4int)aMaterial->GetI << 549   G4int materialIndex = aMaterial->GetIndex();
551   G4double scaledKineticEnergy = KineticEnergy    550   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
552   G4double Range;                                 551   G4double Range;
553   G4bool isOut;                                   552   G4bool isOut;
554                                                   553 
555   if (scaledKineticEnergy<t->theLowestKineticE    554   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
556                                                   555 
557     Range = std::sqrt(scaledKineticEnergy/t->t    556     Range = std::sqrt(scaledKineticEnergy/t->theLowestKineticEnergy)*
558             (*rangeTable)(materialIndex)->GetV    557             (*rangeTable)(materialIndex)->GetValue(
559               t->theLowestKineticEnergy,isOut)    558               t->theLowestKineticEnergy,isOut);
560                                                   559 
561   } else if (scaledKineticEnergy>t->theHighest    560   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
562                                                   561 
563     Range = (*rangeTable)(materialIndex)->GetV    562     Range = (*rangeTable)(materialIndex)->GetValue(
564         t->theHighestKineticEnergy,isOut)+        563         t->theHighestKineticEnergy,isOut)+
565             (scaledKineticEnergy-t->theHighest    564             (scaledKineticEnergy-t->theHighestKineticEnergy)/
566             (*dEdxTable)(materialIndex)->GetVa    565             (*dEdxTable)(materialIndex)->GetValue(
567               t->theHighestKineticEnergy,isOut    566               t->theHighestKineticEnergy,isOut);
568                                                   567 
569   } else {                                        568   } else {
570                                                   569 
571     Range = (*rangeTable)(materialIndex)->GetV    570     Range = (*rangeTable)(materialIndex)->GetValue(
572          scaledKineticEnergy,isOut);              571          scaledKineticEnergy,isOut);
573                                                   572 
574   }                                               573   }
575                                                   574 
576   return Range/(Chargesquare*t->theMassRatio);    575   return Range/(Chargesquare*t->theMassRatio);
577 }                                                 576 }
578                                                   577 
579 //....oooOO0OOooo........oooOO0OOooo........oo    578 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
580                                                   579 
581 G4double G4EnergyLossTables::GetPreciseEnergyF    580 G4double G4EnergyLossTables::GetPreciseEnergyFromRange(
582                                      const G4P    581                                      const G4ParticleDefinition *aParticle,
583                                            G4d    582                                            G4double range,
584                                      const G4M    583                                      const G4Material *aMaterial)
585 // it returns the value of the kinetic energy     584 // it returns the value of the kinetic energy for a given range
586 {                                                 585 {
587   if (!t) t = new G4EnergyLossTablesHelper;       586   if (!t) t = new G4EnergyLossTablesHelper;
588                                                   587 
589   CPRWarning();                                   588   CPRWarning();
590   if( aParticle != (const G4ParticleDefinition    589   if( aParticle != (const G4ParticleDefinition*) lastParticle)
591   {                                               590   {
592     *t= GetTables(aParticle);                     591     *t= GetTables(aParticle);
593     lastParticle = (G4ParticleDefinition*) aPa    592     lastParticle = (G4ParticleDefinition*) aParticle;
594     Chargesquare = (aParticle->GetPDGCharge())    593     Chargesquare = (aParticle->GetPDGCharge())*
595                    (aParticle->GetPDGCharge())    594                    (aParticle->GetPDGCharge())/
596                     QQPositron ;                  595                     QQPositron ;
597     oldIndex = -1 ;                               596     oldIndex = -1 ;
598   }                                               597   }
599   const G4PhysicsTable*  dEdxTable= t->theDEDX    598   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
600   const G4PhysicsTable*  inverseRangeTable= t-    599   const G4PhysicsTable*  inverseRangeTable= t->theInverseRangeTable;
601   if (!inverseRangeTable) {                       600   if (!inverseRangeTable) {
602     ParticleHaveNoLoss(aParticle,"InverseRange    601     ParticleHaveNoLoss(aParticle,"InverseRange");
603     return 0.0;                                   602     return 0.0;
604   }                                               603   }
605                                                   604 
606   G4double scaledrange,scaledKineticEnergy ;      605   G4double scaledrange,scaledKineticEnergy ;
607   G4bool isOut ;                                  606   G4bool isOut ;
608                                                   607 
609   G4int materialIndex = (G4int)aMaterial->GetI << 608   G4int materialIndex = aMaterial->GetIndex() ;
610                                                   609 
611   if(materialIndex != oldIndex)                   610   if(materialIndex != oldIndex)
612   {                                               611   {
613     oldIndex = materialIndex ;                    612     oldIndex = materialIndex ;
614     rmin = (*inverseRangeTable)(materialIndex)    613     rmin = (*inverseRangeTable)(materialIndex)->
615                               GetLowEdgeEnergy    614                               GetLowEdgeEnergy(0) ;
616     rmax = (*inverseRangeTable)(materialIndex)    615     rmax = (*inverseRangeTable)(materialIndex)->
617                    GetLowEdgeEnergy(t->theNumb    616                    GetLowEdgeEnergy(t->theNumberOfBins-2) ;
618     Thigh = (*inverseRangeTable)(materialIndex    617     Thigh = (*inverseRangeTable)(materialIndex)->
619                               GetValue(rmax,is    618                               GetValue(rmax,isOut) ;
620   }                                               619   }
621                                                   620 
622   scaledrange = range*Chargesquare*t->theMassR    621   scaledrange = range*Chargesquare*t->theMassRatio ;
623                                                   622 
624   if(scaledrange < rmin)                          623   if(scaledrange < rmin)
625   {                                               624   {
626     scaledKineticEnergy = t->theLowestKineticE    625     scaledKineticEnergy = t->theLowestKineticEnergy*
627                    scaledrange*scaledrange/(rm    626                    scaledrange*scaledrange/(rmin*rmin) ;
628   }                                               627   }
629   else                                            628   else
630   {                                               629   {
631     if(scaledrange < rmax)                        630     if(scaledrange < rmax)
632     {                                             631     {
633       scaledKineticEnergy = (*inverseRangeTabl    632       scaledKineticEnergy = (*inverseRangeTable)(materialIndex)->
634                               GetValue( scaled    633                               GetValue( scaledrange,isOut) ;
635     }                                             634     }
636     else                                          635     else
637     {                                             636     {
638       scaledKineticEnergy = Thigh +               637       scaledKineticEnergy = Thigh +
639                       (scaledrange-rmax)*         638                       (scaledrange-rmax)*
640                       (*dEdxTable)(materialInd    639                       (*dEdxTable)(materialIndex)->
641                                  GetValue(Thig    640                                  GetValue(Thigh,isOut) ;
642     }                                             641     }
643   }                                               642   }
644                                                   643 
645   return scaledKineticEnergy/t->theMassRatio ;    644   return scaledKineticEnergy/t->theMassRatio ;
646 }                                                 645 }
647                                                   646 
648 //....oooOO0OOooo........oooOO0OOooo........oo    647 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
649                                                   648 
650  G4double G4EnergyLossTables::GetPreciseDEDX(     649  G4double G4EnergyLossTables::GetPreciseDEDX(
651     const G4ParticleDefinition *aParticle,        650     const G4ParticleDefinition *aParticle,
652     G4double KineticEnergy,                       651     G4double KineticEnergy,
653     const G4Material *aMaterial)                  652     const G4Material *aMaterial)
654 {                                                 653 {
655   if (!t) t = new G4EnergyLossTablesHelper;       654   if (!t) t = new G4EnergyLossTablesHelper;
656                                                   655 
657   CPRWarning();                                   656   CPRWarning();
658   if( aParticle != (const G4ParticleDefinition    657   if( aParticle != (const G4ParticleDefinition*) lastParticle)
659   {                                               658   {
660     *t= GetTables(aParticle);                     659     *t= GetTables(aParticle);
661     lastParticle = (G4ParticleDefinition*) aPa    660     lastParticle = (G4ParticleDefinition*) aParticle;
662     Chargesquare = (aParticle->GetPDGCharge())    661     Chargesquare = (aParticle->GetPDGCharge())*
663                    (aParticle->GetPDGCharge())    662                    (aParticle->GetPDGCharge())/
664                     QQPositron ;                  663                     QQPositron ;
665     oldIndex = -1 ;                               664     oldIndex = -1 ;
666   }                                               665   }
667   const G4PhysicsTable*  dEdxTable= t->theDEDX    666   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
668   if (!dEdxTable) {                               667   if (!dEdxTable) {
669     ParticleHaveNoLoss(aParticle,"dEdx");         668     ParticleHaveNoLoss(aParticle,"dEdx");
670     return 0.0;                                   669     return 0.0;
671   }                                               670   }
672                                                   671 
673   G4int materialIndex = (G4int)aMaterial->GetI << 672   G4int materialIndex = aMaterial->GetIndex();
674   G4double scaledKineticEnergy = KineticEnergy    673   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
675   G4double dEdx;                                  674   G4double dEdx;
676   G4bool isOut;                                   675   G4bool isOut;
677                                                   676 
678   if (scaledKineticEnergy<t->theLowestKineticE    677   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
679                                                   678 
680      dEdx = std::sqrt(scaledKineticEnergy/t->t    679      dEdx = std::sqrt(scaledKineticEnergy/t->theLowestKineticEnergy)
681             *(*dEdxTable)(materialIndex)->GetV    680             *(*dEdxTable)(materialIndex)->GetValue(
682               t->theLowestKineticEnergy,isOut)    681               t->theLowestKineticEnergy,isOut);
683                                                   682 
684   } else if (scaledKineticEnergy>t->theHighest    683   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
685                                                   684 
686      dEdx = (*dEdxTable)(materialIndex)->GetVa    685      dEdx = (*dEdxTable)(materialIndex)->GetValue(
687         t->theHighestKineticEnergy,isOut);        686         t->theHighestKineticEnergy,isOut);
688                                                   687 
689   } else {                                        688   } else {
690                                                   689 
691       dEdx = (*dEdxTable)(materialIndex)->GetV    690       dEdx = (*dEdxTable)(materialIndex)->GetValue(
692                           scaledKineticEnergy,    691                           scaledKineticEnergy,isOut) ;
693                                                   692 
694   }                                               693   }
695                                                   694 
696   return dEdx*Chargesquare;                       695   return dEdx*Chargesquare;
697 }                                                 696 }
698                                                   697 
699 //....oooOO0OOooo........oooOO0OOooo........oo    698 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
700                                                   699 
701  G4double G4EnergyLossTables::GetPreciseRangeF    700  G4double G4EnergyLossTables::GetPreciseRangeFromEnergy(
702     const G4ParticleDefinition *aParticle,        701     const G4ParticleDefinition *aParticle,
703     G4double KineticEnergy,                       702     G4double KineticEnergy,
704     const G4Material *aMaterial)                  703     const G4Material *aMaterial)
705 {                                                 704 {
706   if (!t) t = new G4EnergyLossTablesHelper;       705   if (!t) t = new G4EnergyLossTablesHelper;
707                                                   706 
708   CPRWarning();                                   707   CPRWarning();
709   if( aParticle != (const G4ParticleDefinition    708   if( aParticle != (const G4ParticleDefinition*) lastParticle)
710   {                                               709   {
711     *t= GetTables(aParticle);                     710     *t= GetTables(aParticle);
712     lastParticle = (G4ParticleDefinition*) aPa    711     lastParticle = (G4ParticleDefinition*) aParticle;
713     Chargesquare = (aParticle->GetPDGCharge())    712     Chargesquare = (aParticle->GetPDGCharge())*
714                    (aParticle->GetPDGCharge())    713                    (aParticle->GetPDGCharge())/
715                     QQPositron ;                  714                     QQPositron ;
716     oldIndex = -1 ;                               715     oldIndex = -1 ;
717   }                                               716   }
718   const G4PhysicsTable* rangeTable= t->theRang    717   const G4PhysicsTable* rangeTable= t->theRangeTable;
719   const G4PhysicsTable*  dEdxTable= t->theDEDX    718   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
720   if (!rangeTable) {                              719   if (!rangeTable) {
721     ParticleHaveNoLoss(aParticle,"Range");        720     ParticleHaveNoLoss(aParticle,"Range");
722     return 0.0;                                   721     return 0.0;
723   }                                               722   }
724   G4int materialIndex = (G4int)aMaterial->GetI << 723   G4int materialIndex = aMaterial->GetIndex();
725                                                   724 
726   G4double Thighr = t->theHighestKineticEnergy    725   G4double Thighr = t->theHighestKineticEnergy*t->theLowestKineticEnergy/
727                    (*rangeTable)(materialIndex    726                    (*rangeTable)(materialIndex)->
728                    GetLowEdgeEnergy(1) ;          727                    GetLowEdgeEnergy(1) ;
729                                                   728 
730   G4double scaledKineticEnergy = KineticEnergy    729   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
731   G4double Range;                                 730   G4double Range;
732   G4bool isOut;                                   731   G4bool isOut;
733                                                   732 
734   if (scaledKineticEnergy<t->theLowestKineticE    733   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
735                                                   734 
736     Range = std::sqrt(scaledKineticEnergy/t->t    735     Range = std::sqrt(scaledKineticEnergy/t->theLowestKineticEnergy)*
737             (*rangeTable)(materialIndex)->GetV    736             (*rangeTable)(materialIndex)->GetValue(
738               t->theLowestKineticEnergy,isOut)    737               t->theLowestKineticEnergy,isOut);
739                                                   738 
740   } else if (scaledKineticEnergy>Thighr) {        739   } else if (scaledKineticEnergy>Thighr) {
741                                                   740 
742     Range = (*rangeTable)(materialIndex)->GetV    741     Range = (*rangeTable)(materialIndex)->GetValue(
743         Thighr,isOut)+                            742         Thighr,isOut)+
744             (scaledKineticEnergy-Thighr)/         743             (scaledKineticEnergy-Thighr)/
745             (*dEdxTable)(materialIndex)->GetVa    744             (*dEdxTable)(materialIndex)->GetValue(
746               Thighr,isOut);                      745               Thighr,isOut);
747                                                   746 
748   } else {                                        747   } else {
749                                                   748 
750      Range = (*rangeTable)(materialIndex)->Get    749      Range = (*rangeTable)(materialIndex)->GetValue(
751                        scaledKineticEnergy,isO    750                        scaledKineticEnergy,isOut) ;
752                                                   751 
753   }                                               752   }
754                                                   753 
755   return Range/(Chargesquare*t->theMassRatio);    754   return Range/(Chargesquare*t->theMassRatio);
756 }                                                 755 }
757                                                   756 
758 //....oooOO0OOooo........oooOO0OOooo........oo    757 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
759                                                   758 
760 G4double G4EnergyLossTables::GetDEDX(             759 G4double G4EnergyLossTables::GetDEDX(
761     const G4ParticleDefinition *aParticle,        760     const G4ParticleDefinition *aParticle,
762     G4double KineticEnergy,                       761     G4double KineticEnergy,
763     const G4MaterialCutsCouple *couple,           762     const G4MaterialCutsCouple *couple,
764     G4bool check)                                 763     G4bool check)
765 {                                                 764 {
766   if (!t) t = new G4EnergyLossTablesHelper;       765   if (!t) t = new G4EnergyLossTablesHelper;
767                                                   766 
768   if(aParticle != (const G4ParticleDefinition*    767   if(aParticle != (const G4ParticleDefinition*) lastParticle)
769   {                                               768   {
770     *t= GetTables(aParticle);                     769     *t= GetTables(aParticle);
771     lastParticle = (G4ParticleDefinition*) aPa    770     lastParticle = (G4ParticleDefinition*) aParticle ;
772     Chargesquare = (aParticle->GetPDGCharge())    771     Chargesquare = (aParticle->GetPDGCharge())*
773                    (aParticle->GetPDGCharge())    772                    (aParticle->GetPDGCharge())/
774                    QQPositron ;                   773                    QQPositron ;
775     oldIndex = -1 ;                               774     oldIndex = -1 ;
776   }                                               775   }
777   const G4PhysicsTable*  dEdxTable= t->theDEDX    776   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
778                                                   777   
779   if (!dEdxTable ) {                              778   if (!dEdxTable ) {
780     if (check) return G4LossTableManager::Inst    779     if (check) return G4LossTableManager::Instance()->GetDEDX(aParticle,KineticEnergy,couple);
781     else       ParticleHaveNoLoss(aParticle, "    780     else       ParticleHaveNoLoss(aParticle, "dEdx");
782     return 0.0;                                   781     return 0.0;
783   }                                               782   }
784                                                   783 
785   G4int materialIndex = couple->GetIndex();       784   G4int materialIndex = couple->GetIndex();
786   G4double scaledKineticEnergy = KineticEnergy    785   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
787   G4double dEdx;                                  786   G4double dEdx;
788   G4bool isOut;                                   787   G4bool isOut;
789                                                   788 
790   if (scaledKineticEnergy<t->theLowestKineticE    789   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
791                                                   790 
792      dEdx =(*dEdxTable)(materialIndex)->GetVal    791      dEdx =(*dEdxTable)(materialIndex)->GetValue(
793               t->theLowestKineticEnergy,isOut)    792               t->theLowestKineticEnergy,isOut)
794            *std::sqrt(scaledKineticEnergy/t->t    793            *std::sqrt(scaledKineticEnergy/t->theLowestKineticEnergy);
795                                                   794 
796   } else if (scaledKineticEnergy>t->theHighest    795   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
797                                                   796 
798      dEdx = (*dEdxTable)(materialIndex)->GetVa    797      dEdx = (*dEdxTable)(materialIndex)->GetValue(
799         t->theHighestKineticEnergy,isOut);        798         t->theHighestKineticEnergy,isOut);
800                                                   799 
801   } else {                                        800   } else {
802                                                   801 
803     dEdx = (*dEdxTable)(materialIndex)->GetVal    802     dEdx = (*dEdxTable)(materialIndex)->GetValue(
804          scaledKineticEnergy,isOut);              803          scaledKineticEnergy,isOut);
805                                                   804 
806   }                                               805   }
807                                                   806 
808   return dEdx*Chargesquare;                       807   return dEdx*Chargesquare;
809 }                                                 808 }
810                                                   809 
811 //....oooOO0OOooo........oooOO0OOooo........oo    810 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
812                                                   811 
813 G4double G4EnergyLossTables::GetRange(            812 G4double G4EnergyLossTables::GetRange(
814     const G4ParticleDefinition *aParticle,        813     const G4ParticleDefinition *aParticle,
815     G4double KineticEnergy,                       814     G4double KineticEnergy,
816     const G4MaterialCutsCouple *couple,           815     const G4MaterialCutsCouple *couple,
817     G4bool check)                                 816     G4bool check)
818 {                                                 817 {
819   if (!t) t = new G4EnergyLossTablesHelper;       818   if (!t) t = new G4EnergyLossTablesHelper;
820                                                   819 
821   if(aParticle != (const G4ParticleDefinition*    820   if(aParticle != (const G4ParticleDefinition*) lastParticle)
822   {                                               821   {
823     *t= GetTables(aParticle);                     822     *t= GetTables(aParticle);
824     lastParticle = (G4ParticleDefinition*) aPa    823     lastParticle = (G4ParticleDefinition*) aParticle ;
825     Chargesquare = (aParticle->GetPDGCharge())    824     Chargesquare = (aParticle->GetPDGCharge())*
826                    (aParticle->GetPDGCharge())    825                    (aParticle->GetPDGCharge())/
827                     QQPositron ;                  826                     QQPositron ;
828     oldIndex = -1 ;                               827     oldIndex = -1 ;
829   }                                               828   }
830   const G4PhysicsTable* rangeTable= t->theRang    829   const G4PhysicsTable* rangeTable= t->theRangeTable;
831   const G4PhysicsTable*  dEdxTable= t->theDEDX    830   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
832   if (!rangeTable) {                              831   if (!rangeTable) {
833     if(check) return G4LossTableManager::Insta    832     if(check) return G4LossTableManager::Instance()->GetRange(aParticle,KineticEnergy,couple);
834     else      return DBL_MAX;                     833     else      return DBL_MAX;      
835       //ParticleHaveNoLoss(aParticle,"Range");    834       //ParticleHaveNoLoss(aParticle,"Range");
836   }                                               835   }
837                                                   836 
838   G4int materialIndex = couple->GetIndex();       837   G4int materialIndex = couple->GetIndex();
839   G4double scaledKineticEnergy = KineticEnergy    838   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
840   G4double Range;                                 839   G4double Range;
841   G4bool isOut;                                   840   G4bool isOut;
842                                                   841 
843   if (scaledKineticEnergy<t->theLowestKineticE    842   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
844                                                   843 
845     Range = std::sqrt(scaledKineticEnergy/t->t    844     Range = std::sqrt(scaledKineticEnergy/t->theLowestKineticEnergy)*
846             (*rangeTable)(materialIndex)->GetV    845             (*rangeTable)(materialIndex)->GetValue(
847               t->theLowestKineticEnergy,isOut)    846               t->theLowestKineticEnergy,isOut);
848                                                   847 
849   } else if (scaledKineticEnergy>t->theHighest    848   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
850                                                   849 
851     Range = (*rangeTable)(materialIndex)->GetV    850     Range = (*rangeTable)(materialIndex)->GetValue(
852         t->theHighestKineticEnergy,isOut)+        851         t->theHighestKineticEnergy,isOut)+
853             (scaledKineticEnergy-t->theHighest    852             (scaledKineticEnergy-t->theHighestKineticEnergy)/
854             (*dEdxTable)(materialIndex)->GetVa    853             (*dEdxTable)(materialIndex)->GetValue(
855               t->theHighestKineticEnergy,isOut    854               t->theHighestKineticEnergy,isOut);
856                                                   855 
857   } else {                                        856   } else {
858                                                   857 
859     Range = (*rangeTable)(materialIndex)->GetV    858     Range = (*rangeTable)(materialIndex)->GetValue(
860          scaledKineticEnergy,isOut);              859          scaledKineticEnergy,isOut);
861                                                   860 
862   }                                               861   }
863                                                   862 
864   return Range/(Chargesquare*t->theMassRatio);    863   return Range/(Chargesquare*t->theMassRatio);
865 }                                                 864 }
866                                                   865 
867 //....oooOO0OOooo........oooOO0OOooo........oo    866 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
868                                                   867 
869 G4double G4EnergyLossTables::GetPreciseEnergyF    868 G4double G4EnergyLossTables::GetPreciseEnergyFromRange(
870                                      const G4P    869                                      const G4ParticleDefinition *aParticle,
871                                            G4d    870                                            G4double range,
872                                      const G4M    871                                      const G4MaterialCutsCouple *couple,
873                    G4bool check)                  872                    G4bool check)
874 // it returns the value of the kinetic energy     873 // it returns the value of the kinetic energy for a given range
875 {                                                 874 {
876   if (!t) t = new G4EnergyLossTablesHelper;       875   if (!t) t = new G4EnergyLossTablesHelper;
877                                                   876 
878   if( aParticle != (const G4ParticleDefinition    877   if( aParticle != (const G4ParticleDefinition*) lastParticle)
879   {                                               878   {
880     *t= GetTables(aParticle);                     879     *t= GetTables(aParticle);
881     lastParticle = (G4ParticleDefinition*) aPa    880     lastParticle = (G4ParticleDefinition*) aParticle;
882     Chargesquare = (aParticle->GetPDGCharge())    881     Chargesquare = (aParticle->GetPDGCharge())*
883                    (aParticle->GetPDGCharge())    882                    (aParticle->GetPDGCharge())/
884                     QQPositron ;                  883                     QQPositron ;
885     oldIndex = -1 ;                               884     oldIndex = -1 ;
886   }                                               885   }
887   const G4PhysicsTable*  dEdxTable= t->theDEDX    886   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
888   const G4PhysicsTable*  inverseRangeTable= t-    887   const G4PhysicsTable*  inverseRangeTable= t->theInverseRangeTable;
889                                                   888   
890   if (!inverseRangeTable) {                       889   if (!inverseRangeTable) {
891     if(check) return G4LossTableManager::Insta    890     if(check) return G4LossTableManager::Instance()->GetEnergy(aParticle,range,couple);
892     else      return DBL_MAX;                     891     else      return DBL_MAX;      
893     //    else      ParticleHaveNoLoss(aPartic    892     //    else      ParticleHaveNoLoss(aParticle,"InverseRange");
894   }                                               893   }
895                                                   894 
896   G4double scaledrange,scaledKineticEnergy ;      895   G4double scaledrange,scaledKineticEnergy ;
897   G4bool isOut ;                                  896   G4bool isOut ;
898                                                   897 
899   G4int materialIndex = couple->GetIndex() ;      898   G4int materialIndex = couple->GetIndex() ;
900                                                   899 
901   if(materialIndex != oldIndex)                   900   if(materialIndex != oldIndex)
902   {                                               901   {
903     oldIndex = materialIndex ;                    902     oldIndex = materialIndex ;
904     rmin = (*inverseRangeTable)(materialIndex)    903     rmin = (*inverseRangeTable)(materialIndex)->
905                               GetLowEdgeEnergy    904                               GetLowEdgeEnergy(0) ;
906     rmax = (*inverseRangeTable)(materialIndex)    905     rmax = (*inverseRangeTable)(materialIndex)->
907                    GetLowEdgeEnergy(t->theNumb    906                    GetLowEdgeEnergy(t->theNumberOfBins-2) ;
908     Thigh = (*inverseRangeTable)(materialIndex    907     Thigh = (*inverseRangeTable)(materialIndex)->
909                               GetValue(rmax,is    908                               GetValue(rmax,isOut) ;
910   }                                               909   }
911                                                   910 
912   scaledrange = range*Chargesquare*t->theMassR    911   scaledrange = range*Chargesquare*t->theMassRatio ;
913                                                   912 
914   if(scaledrange < rmin)                          913   if(scaledrange < rmin)
915   {                                               914   {
916     scaledKineticEnergy = t->theLowestKineticE    915     scaledKineticEnergy = t->theLowestKineticEnergy*
917                    scaledrange*scaledrange/(rm    916                    scaledrange*scaledrange/(rmin*rmin) ;
918   }                                               917   }
919   else                                            918   else
920   {                                               919   {
921     if(scaledrange < rmax)                        920     if(scaledrange < rmax)
922     {                                             921     {
923       scaledKineticEnergy = (*inverseRangeTabl    922       scaledKineticEnergy = (*inverseRangeTable)(materialIndex)->
924                               GetValue( scaled    923                               GetValue( scaledrange,isOut) ;
925     }                                             924     }
926     else                                          925     else
927     {                                             926     {
928       scaledKineticEnergy = Thigh +               927       scaledKineticEnergy = Thigh +
929                       (scaledrange-rmax)*         928                       (scaledrange-rmax)*
930                       (*dEdxTable)(materialInd    929                       (*dEdxTable)(materialIndex)->
931                                  GetValue(Thig    930                                  GetValue(Thigh,isOut) ;
932     }                                             931     }
933   }                                               932   }
934                                                   933 
935   return scaledKineticEnergy/t->theMassRatio ;    934   return scaledKineticEnergy/t->theMassRatio ;
936 }                                                 935 }
937                                                   936 
938 //....oooOO0OOooo........oooOO0OOooo........oo    937 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
939                                                   938 
940 G4double G4EnergyLossTables::GetPreciseDEDX(      939 G4double G4EnergyLossTables::GetPreciseDEDX(
941     const G4ParticleDefinition *aParticle,        940     const G4ParticleDefinition *aParticle,
942     G4double KineticEnergy,                       941     G4double KineticEnergy,
943     const G4MaterialCutsCouple *couple)           942     const G4MaterialCutsCouple *couple)
944 {                                                 943 {
945   if (!t) t = new G4EnergyLossTablesHelper;       944   if (!t) t = new G4EnergyLossTablesHelper;
946                                                   945 
947   if( aParticle != (const G4ParticleDefinition    946   if( aParticle != (const G4ParticleDefinition*) lastParticle)
948   {                                               947   {
949     *t= GetTables(aParticle);                     948     *t= GetTables(aParticle);
950     lastParticle = (G4ParticleDefinition*) aPa    949     lastParticle = (G4ParticleDefinition*) aParticle;
951     Chargesquare = (aParticle->GetPDGCharge())    950     Chargesquare = (aParticle->GetPDGCharge())*
952                    (aParticle->GetPDGCharge())    951                    (aParticle->GetPDGCharge())/
953                     QQPositron ;                  952                     QQPositron ;
954     oldIndex = -1 ;                               953     oldIndex = -1 ;
955   }                                               954   }
956   const G4PhysicsTable*  dEdxTable= t->theDEDX    955   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
957   if ( !dEdxTable )                               956   if ( !dEdxTable )
958     return G4LossTableManager::Instance()->Get    957     return G4LossTableManager::Instance()->GetDEDX(aParticle,KineticEnergy,couple);
959                                                   958 
960   G4int materialIndex = couple->GetIndex();       959   G4int materialIndex = couple->GetIndex();
961   G4double scaledKineticEnergy = KineticEnergy    960   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
962   G4double dEdx;                                  961   G4double dEdx;
963   G4bool isOut;                                   962   G4bool isOut;
964                                                   963 
965   if (scaledKineticEnergy<t->theLowestKineticE    964   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
966                                                   965 
967      dEdx = std::sqrt(scaledKineticEnergy/t->t    966      dEdx = std::sqrt(scaledKineticEnergy/t->theLowestKineticEnergy)
968             *(*dEdxTable)(materialIndex)->GetV    967             *(*dEdxTable)(materialIndex)->GetValue(
969               t->theLowestKineticEnergy,isOut)    968               t->theLowestKineticEnergy,isOut);
970                                                   969 
971   } else if (scaledKineticEnergy>t->theHighest    970   } else if (scaledKineticEnergy>t->theHighestKineticEnergy) {
972                                                   971 
973      dEdx = (*dEdxTable)(materialIndex)->GetVa    972      dEdx = (*dEdxTable)(materialIndex)->GetValue(
974         t->theHighestKineticEnergy,isOut);        973         t->theHighestKineticEnergy,isOut);
975                                                   974 
976   } else {                                        975   } else {
977                                                   976 
978       dEdx = (*dEdxTable)(materialIndex)->GetV    977       dEdx = (*dEdxTable)(materialIndex)->GetValue(
979                           scaledKineticEnergy,    978                           scaledKineticEnergy,isOut) ;
980                                                   979 
981   }                                               980   }
982                                                   981 
983   return dEdx*Chargesquare;                       982   return dEdx*Chargesquare;
984 }                                                 983 }
985                                                   984 
986 //....oooOO0OOooo........oooOO0OOooo........oo    985 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
987                                                   986 
988 G4double G4EnergyLossTables::GetPreciseRangeFr    987 G4double G4EnergyLossTables::GetPreciseRangeFromEnergy(
989     const G4ParticleDefinition *aParticle,        988     const G4ParticleDefinition *aParticle,
990     G4double KineticEnergy,                       989     G4double KineticEnergy,
991     const G4MaterialCutsCouple *couple)           990     const G4MaterialCutsCouple *couple)
992 {                                                 991 {
993   if (!t) t = new G4EnergyLossTablesHelper;       992   if (!t) t = new G4EnergyLossTablesHelper;
994                                                   993 
995   if( aParticle != (const G4ParticleDefinition    994   if( aParticle != (const G4ParticleDefinition*) lastParticle)
996   {                                               995   {
997     *t= GetTables(aParticle);                     996     *t= GetTables(aParticle);
998     lastParticle = (G4ParticleDefinition*) aPa    997     lastParticle = (G4ParticleDefinition*) aParticle;
999     Chargesquare = (aParticle->GetPDGCharge())    998     Chargesquare = (aParticle->GetPDGCharge())*
1000                    (aParticle->GetPDGCharge()    999                    (aParticle->GetPDGCharge())/
1001                     QQPositron ;                 1000                     QQPositron ;
1002     oldIndex = -1 ;                              1001     oldIndex = -1 ;
1003   }                                              1002   }
1004   const G4PhysicsTable* rangeTable= t->theRan    1003   const G4PhysicsTable* rangeTable= t->theRangeTable;
1005   const G4PhysicsTable*  dEdxTable= t->theDED    1004   const G4PhysicsTable*  dEdxTable= t->theDEDXTable;
1006   if ( !dEdxTable || !rangeTable)                1005   if ( !dEdxTable || !rangeTable)
1007     return G4LossTableManager::Instance()->Ge    1006     return G4LossTableManager::Instance()->GetDEDX(aParticle,KineticEnergy,couple);
1008                                                  1007 
1009   G4int materialIndex = couple->GetIndex();      1008   G4int materialIndex = couple->GetIndex();
1010                                                  1009 
1011   G4double Thighr = t->theHighestKineticEnerg    1010   G4double Thighr = t->theHighestKineticEnergy*t->theLowestKineticEnergy/
1012                    (*rangeTable)(materialInde    1011                    (*rangeTable)(materialIndex)->
1013                    GetLowEdgeEnergy(1) ;         1012                    GetLowEdgeEnergy(1) ;
1014                                                  1013 
1015   G4double scaledKineticEnergy = KineticEnerg    1014   G4double scaledKineticEnergy = KineticEnergy*t->theMassRatio;
1016   G4double Range;                                1015   G4double Range;
1017   G4bool isOut;                                  1016   G4bool isOut;
1018                                                  1017 
1019   if (scaledKineticEnergy<t->theLowestKinetic    1018   if (scaledKineticEnergy<t->theLowestKineticEnergy) {
1020                                                  1019 
1021     Range = std::sqrt(scaledKineticEnergy/t->    1020     Range = std::sqrt(scaledKineticEnergy/t->theLowestKineticEnergy)*
1022             (*rangeTable)(materialIndex)->Get    1021             (*rangeTable)(materialIndex)->GetValue(
1023               t->theLowestKineticEnergy,isOut    1022               t->theLowestKineticEnergy,isOut);
1024                                                  1023 
1025   } else if (scaledKineticEnergy>Thighr) {       1024   } else if (scaledKineticEnergy>Thighr) {
1026                                                  1025 
1027     Range = (*rangeTable)(materialIndex)->Get    1026     Range = (*rangeTable)(materialIndex)->GetValue(
1028         Thighr,isOut)+                           1027         Thighr,isOut)+
1029             (scaledKineticEnergy-Thighr)/        1028             (scaledKineticEnergy-Thighr)/
1030             (*dEdxTable)(materialIndex)->GetV    1029             (*dEdxTable)(materialIndex)->GetValue(
1031               Thighr,isOut);                     1030               Thighr,isOut);
1032                                                  1031 
1033   } else {                                       1032   } else {
1034                                                  1033 
1035      Range = (*rangeTable)(materialIndex)->Ge    1034      Range = (*rangeTable)(materialIndex)->GetValue(
1036                        scaledKineticEnergy,is    1035                        scaledKineticEnergy,isOut) ;
1037                                                  1036 
1038   }                                              1037   }
1039                                                  1038 
1040   return Range/(Chargesquare*t->theMassRatio)    1039   return Range/(Chargesquare*t->theMassRatio);
1041 }                                                1040 }
1042                                                  1041 
1043 //....oooOO0OOooo........oooOO0OOooo........o    1042 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1044                                                  1043 
1045 void G4EnergyLossTables::CPRWarning()            1044 void G4EnergyLossTables::CPRWarning()
1046 {                                                1045 {
1047   if (let_counter <  let_max_num_warnings) {     1046   if (let_counter <  let_max_num_warnings) {
1048                                                  1047 
1049     G4cout << G4endl;                            1048     G4cout << G4endl;
1050     G4cout << "##### G4EnergyLossTable WARNIN    1049     G4cout << "##### G4EnergyLossTable WARNING: The obsolete interface is used!" << G4endl;
1051     G4cout << "##### RESULTS ARE NOT GARANTEE    1050     G4cout << "##### RESULTS ARE NOT GARANTEED!" << G4endl;
1052     G4cout << "##### Please, substitute G4Mat    1051     G4cout << "##### Please, substitute G4Material by G4MaterialCutsCouple" << G4endl;
1053     G4cout << "##### Obsolete interface will     1052     G4cout << "##### Obsolete interface will be removed soon" << G4endl;
1054     G4cout << G4endl;                            1053     G4cout << G4endl;
1055     let_counter++;                               1054     let_counter++;
1056                                                  1055 
1057   } else if (let_counter == let_max_num_warni    1056   } else if (let_counter == let_max_num_warnings) {
1058                                                  1057 
1059     G4cout << "##### G4EnergyLossTable WARNIN    1058     G4cout << "##### G4EnergyLossTable WARNING closed" << G4endl;
1060     let_counter++;                               1059     let_counter++;
1061   }                                              1060   }
1062 }                                                1061 }
1063                                                  1062 
1064 //....oooOO0OOooo........oooOO0OOooo........o    1063 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1065                                                  1064 
1066 void                                             1065 void 
1067 G4EnergyLossTables::ParticleHaveNoLoss(const     1066 G4EnergyLossTables::ParticleHaveNoLoss(const G4ParticleDefinition*, 
1068                const G4String& /*q*/)            1067                const G4String& /*q*/)
1069 {                                                1068 {
1070 }                                                1069 }
1071                                                  1070 
1072 //....oooOO0OOooo........oooOO0OOooo........o    1071 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1073                                                  1072