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 10.3.p3)


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