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.0.p2)


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