Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/adjoint/src/G4ContinuousGainOfEnergy.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/adjoint/src/G4ContinuousGainOfEnergy.cc (Version 11.3.0) and /processes/electromagnetic/adjoint/src/G4ContinuousGainOfEnergy.cc (Version 9.5)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                         3 // * License and Disclaimer                                           *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th      5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided      6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License      7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/      8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.           9 // * include a list of copyright holders.                             *
 10 // *                                               10 // *                                                                  *
 11 // * Neither the authors of this software syst     11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file      15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitatio     16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                               17 // *                                                                  *
 18 // * This  code  implementation is the result      18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboratio     19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distri     20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  ag     21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publicati     22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Sof     23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // *******************************************     24 // ********************************************************************
 25 //                                                 25 //
                                                   >>  26 // $Id: G4ContinuousGainOfEnergy.cc,v 1.5 2010-11-11 11:51:56 ldesorgh Exp $
                                                   >>  27 // GEANT4 tag $Name: not supported by cvs2svn $
 26 //                                                 28 //
 27                                                << 
 28 #include "G4ContinuousGainOfEnergy.hh"             29 #include "G4ContinuousGainOfEnergy.hh"
 29                                                << 
 30 #include "G4EmCorrections.hh"                  << 
 31 #include "G4LossTableManager.hh"               << 
 32 #include "G4Material.hh"                       << 
 33 #include "G4MaterialCutsCouple.hh"             << 
 34 #include "G4ParticleChange.hh"                 << 
 35 #include "G4ParticleDefinition.hh"             << 
 36 #include "G4PhysicalConstants.hh"              << 
 37 #include "G4Step.hh"                               30 #include "G4Step.hh"
 38 #include "G4SystemOfUnits.hh"                  <<  31 #include "G4ParticleDefinition.hh"
 39 #include "G4VEmFluctuationModel.hh"            << 
 40 #include "G4VEmModel.hh"                           32 #include "G4VEmModel.hh"
 41 #include "G4VEnergyLossProcess.hh"             <<  33 #include "G4VEmFluctuationModel.hh"
 42 #include "G4VParticleChange.hh"                    34 #include "G4VParticleChange.hh"
                                                   >>  35 #include "G4UnitsTable.hh"
                                                   >>  36 #include "G4AdjointCSManager.hh"
                                                   >>  37 #include "G4LossTableManager.hh"
 43                                                    38 
 44 //////////////////////////////////////////////     39 ///////////////////////////////////////////////////////
 45 G4ContinuousGainOfEnergy::G4ContinuousGainOfEn <<  40 //
 46                                                <<  41 G4ContinuousGainOfEnergy::G4ContinuousGainOfEnergy(const G4String& name, 
 47   : G4VContinuousProcess(name, type)           <<  42   G4ProcessType type): G4VContinuousProcess(name, type)
 48 {}                                             <<  43 {
 49                                                    44 
 50 ////////////////////////////////////////////// <<  45 
 51 G4ContinuousGainOfEnergy::~G4ContinuousGainOfE <<  46   linLossLimit=0.05;
                                                   >>  47   lossFluctuationArePossible =true;
                                                   >>  48   lossFluctuationFlag=true;
                                                   >>  49   is_integral = false;
                                                   >>  50   
                                                   >>  51   //Will be properly set in SetDirectParticle()
                                                   >>  52   IsIon=false;
                                                   >>  53   massRatio =1.;
                                                   >>  54   chargeSqRatio=1.;
                                                   >>  55   preStepChargeSqRatio=1.;
                                                   >>  56   
                                                   >>  57   //Some initialization
                                                   >>  58   currentCoupleIndex=9999999;
                                                   >>  59   currentCutInRange=0.;
                                                   >>  60   currentMaterialIndex=9999999;
                                                   >>  61   currentTcut=0.;
                                                   >>  62   preStepKinEnergy=0.;
                                                   >>  63   preStepRange=0.;
                                                   >>  64   preStepScaledKinEnergy=0.;
                                                   >>  65   
                                                   >>  66   currentCouple=0;  
                                                   >>  67 }
 52                                                    68 
 53 //////////////////////////////////////////////     69 ///////////////////////////////////////////////////////
 54 void G4ContinuousGainOfEnergy::ProcessDescript <<  70 //
                                                   >>  71 G4ContinuousGainOfEnergy::~G4ContinuousGainOfEnergy()
 55 {                                                  72 {
 56   out << "Continuous process acting on adjoint <<  73  
 57          "continuous gain of energy of charged <<  74 }
 58          "tracked back.\n";                    <<  75 ///////////////////////////////////////////////////////
                                                   >>  76 //
                                                   >>  77 
                                                   >>  78 void G4ContinuousGainOfEnergy::PreparePhysicsTable(
                                                   >>  79      const G4ParticleDefinition& )
                                                   >>  80 {//theDirectEnergyLossProcess->PreparePhysicsTable(part);
                                                   >>  81 
                                                   >>  82 ; 
 59 }                                                  83 }
 60                                                    84 
 61 //////////////////////////////////////////////     85 ///////////////////////////////////////////////////////
 62 void G4ContinuousGainOfEnergy::SetDirectPartic <<  86 //
 63 {                                              <<  87 
 64   fDirectPartDef = p;                          <<  88 void G4ContinuousGainOfEnergy::BuildPhysicsTable(const G4ParticleDefinition&)
 65   if(fDirectPartDef->GetParticleType() == "nuc <<  89 {//theDirectEnergyLossProcess->BuildPhysicsTable(part);
 66   {                                            <<  90 ;
 67     fIsIon     = true;                         <<  91 }
 68     fMassRatio = proton_mass_c2 / fDirectPartD <<  92 
 69   }                                            <<  93 ///////////////////////////////////////////////////////
                                                   >>  94 //
                                                   >>  95 void  G4ContinuousGainOfEnergy::SetDirectParticle(G4ParticleDefinition* p)
                                                   >>  96 {theDirectPartDef=p;
                                                   >>  97  if (theDirectPartDef->GetParticleType()== "nucleus") {
                                                   >>  98    IsIon=true;
                                                   >>  99    massRatio = proton_mass_c2/theDirectPartDef->GetPDGMass();
                                                   >> 100    G4double q=theDirectPartDef->GetPDGCharge();
                                                   >> 101    chargeSqRatio=q*q;
                                                   >> 102   
                                                   >> 103    
                                                   >> 104  }
                                                   >> 105  
 70 }                                                 106 }
 71                                                   107 
 72 //////////////////////////////////////////////    108 ///////////////////////////////////////////////////////
                                                   >> 109 //
                                                   >> 110 // 
 73 G4VParticleChange* G4ContinuousGainOfEnergy::A    111 G4VParticleChange* G4ContinuousGainOfEnergy::AlongStepDoIt(const G4Track& track,
 74                                                << 112                                                        const G4Step& step)
 75 {                                                 113 {
 76   // Caution in this method the step length sh << 114    
 77   // A problem is that this is computed by the << 115   //Caution in this method the  step length should be the true step length
 78   // not know the energy at the end of the adj << 116   // A problem is that this is compute by the multiple scattering that does not know the energy at the end of the adjoint step. This energy is used during the 
 79   // during the forward sim. Nothing we can re << 117   //Forward sim. Nothing we can really do against that at this time. This is inherent to the MS method
 80   // time. This is inherent to the MS method   << 118   //
 81                                                << 119   
                                                   >> 120   
                                                   >> 121   
 82   aParticleChange.Initialize(track);              122   aParticleChange.Initialize(track);
 83                                                << 123   
 84   // Get the actual (true) Step length            124   // Get the actual (true) Step length
                                                   >> 125   //----------------------------------
 85   G4double length = step.GetStepLength();         126   G4double length = step.GetStepLength();
 86   G4double degain = 0.0;                       << 127   G4double degain  = 0.0;
 87                                                << 128   
                                                   >> 129   
                                                   >> 130  
 88   // Compute this for weight change after cont    131   // Compute this for weight change after continuous energy loss
 89   G4double DEDX_before =                       << 132   //-------------------------------------------------------------
 90     fDirectEnergyLossProcess->GetDEDX(fPreStep << 133   G4double DEDX_before = theDirectEnergyLossProcess->GetDEDX(preStepKinEnergy, currentCouple);
 91                                                << 134    
 92   // For the fluctuation we generate a new dyn << 135   
 93   // = preEnergy+egain and then compute the fl << 136   
 94   // case.                                     << 137   // For the fluctuation we generate a new dynamic particle with energy =preEnergy+egain
                                                   >> 138   // and then compute the fluctuation given in  the direct case.
                                                   >> 139   //-----------------------------------------------------------------------
 95   G4DynamicParticle* dynParticle = new G4Dynam    140   G4DynamicParticle* dynParticle = new G4DynamicParticle();
 96   *dynParticle                   = *(track.Get << 141   *dynParticle = *(track.GetDynamicParticle());
 97   dynParticle->SetDefinition(fDirectPartDef);  << 142   dynParticle->SetDefinition(theDirectPartDef);
 98   G4double Tkin = dynParticle->GetKineticEnerg << 143   G4double Tkin = dynParticle->GetKineticEnergy(); 
 99                                                << 144 
100   G4double dlength = length;                   << 145 
101   if(Tkin != fPreStepKinEnergy && fIsIon)      << 146   size_t n=1;
102   {                                            << 147   if (is_integral ) n=10;
103     G4double chargeSqRatio = fCurrentModel->Ge << 148   n=1;
104       fDirectPartDef, fCurrentMaterial, Tkin); << 149   G4double dlength= length/n; 
105     fDirectEnergyLossProcess->SetDynamicMassCh << 150   for (size_t i=0;i<n;i++) {
106   }                                            << 151     if (Tkin != preStepKinEnergy && IsIon) {
107                                                << 152       chargeSqRatio =  currentModel->GetChargeSquareRatio(theDirectPartDef,currentMaterial,Tkin);
108   G4double r = fDirectEnergyLossProcess->GetRa << 153     theDirectEnergyLossProcess->SetDynamicMassCharge(massRatio,chargeSqRatio); 
109   if(dlength <= fLinLossLimit * r)             << 154   
110   {                                            << 155     }
111     degain = DEDX_before * dlength;            << 156   
112   }                                            << 157     G4double r = theDirectEnergyLossProcess->GetRange(Tkin, currentCouple);
113   else                                         << 158     if( dlength <= linLossLimit * r ) {
114   {                                            << 159         degain = DEDX_before*dlength;   
115     G4double x = r + dlength;                  << 160   } 
116     G4double E = fDirectEnergyLossProcess->Get << 161     else {
117     if(fIsIon)                                 << 162         G4double x = r + dlength;
118     {                                          << 163         //degain = theDirectEnergyLossProcess->GetKineticEnergy(x,currentCouple) - theDirectEnergyLossProcess->GetKineticEnergy(r,currentCouple);
119       G4double chargeSqRatio = fCurrentModel-> << 164     G4double E = theDirectEnergyLossProcess->GetKineticEnergy(x,currentCouple);
120         fDirectPartDef, fCurrentMaterial, E);  << 165     if (IsIon){
121       fDirectEnergyLossProcess->SetDynamicMass << 166       chargeSqRatio =  currentModel->GetChargeSquareRatio(theDirectPartDef,currentMaterial,E);
122       G4double x1 = fDirectEnergyLossProcess-> << 167       theDirectEnergyLossProcess->SetDynamicMassCharge(massRatio,chargeSqRatio);
123                                                << 168       G4double x1= theDirectEnergyLossProcess->GetRange(E, currentCouple);
124       G4int ii              = 0;               << 169       while (std::abs(x-x1)>0.01*x) {
125       constexpr G4int iimax = 100;             << 170         E = theDirectEnergyLossProcess->GetKineticEnergy(x,currentCouple);
126       while(std::abs(x - x1) > 0.01 * x)       << 171         chargeSqRatio =  currentModel->GetChargeSquareRatio(theDirectPartDef,currentMaterial,E);
127       {                                        << 172         theDirectEnergyLossProcess->SetDynamicMassCharge(massRatio,chargeSqRatio);
128         E = fDirectEnergyLossProcess->GetKinet << 173         x1= theDirectEnergyLossProcess->GetRange(E, currentCouple);
129         chargeSqRatio = fCurrentModel->GetChar << 174       
130           fDirectPartDef, fCurrentMaterial, E) << 175       } 
131         fDirectEnergyLossProcess->SetDynamicMa << 176     }
132                                                << 177     
133         x1 = fDirectEnergyLossProcess->GetRang << 178     degain=E-Tkin;  
134         ++ii;                                  << 179     
135         if(ii >= iimax)                        << 180     
136         {                                      << 181     
137           break;                               << 182     }
138         }                                      << 183   //G4cout<<degain<<G4endl;
139       }                                        << 184     G4double tmax = currentModel->MaxSecondaryKinEnergy(dynParticle);
140     }                                          << 185   tmax = std::min(tmax,currentTcut);
141                                                << 186   
142     degain = E - Tkin;                         << 187   
143   }                                            << 188   dynParticle->SetKineticEnergy(Tkin+degain);
144   G4double tmax = fCurrentModel->MaxSecondaryK << 189 
145   fCurrentTcut = std::min(fCurrentTcut, tmax); << 190   // Corrections, which cannot be tabulated for ions
146                                                << 191   //----------------------------------------
147   dynParticle->SetKineticEnergy(Tkin + degain) << 192   G4double esecdep=0;//not used in most models
148                                                << 193   currentModel->CorrectionsAlongStep(currentCouple, dynParticle, degain,esecdep, dlength); 
149   // Corrections, which cannot be tabulated fo << 194 
150   fCurrentModel->CorrectionsAlongStep(fCurrent << 195     // Sample fluctuations
151                                                << 196     //-------------------
152   // Sample fluctuations                       << 197   
153   G4double deltaE = 0.;                        << 198     
154   if(fLossFluctuationFlag)                     << 199   G4double deltaE =0.;
155   {                                            << 200     if (lossFluctuationFlag ) {
156     deltaE = fCurrentModel->GetModelOfFluctuat << 201           deltaE = currentModel->GetModelOfFluctuations()->
157       fCurrentCouple, dynParticle, fCurrentTcu << 202                   SampleFluctuations(currentMaterial,dynParticle,tmax,dlength,degain)-degain;
158       - degain;                                << 203     }
159   }                                            << 204   
160                                                << 205   G4double egain=degain+deltaE;
161   G4double egain = degain + deltaE;            << 206   if (egain <=0) egain=degain;
162   if(egain <= 0.)                              << 207   Tkin+=egain;
163     egain = degain;                            << 208   dynParticle->SetKineticEnergy(Tkin);
164   Tkin += egain;                               << 209  }
165   dynParticle->SetKineticEnergy(Tkin);         << 210  
                                                   >> 211  
                                                   >> 212   
166                                                   213 
                                                   >> 214   
167   delete dynParticle;                             215   delete dynParticle;
168                                                << 216  
169   if(fIsIon)                                   << 217   if (IsIon){
170   {                                            << 218   chargeSqRatio =  currentModel->GetChargeSquareRatio(theDirectPartDef,currentMaterial,Tkin);
171     G4double chargeSqRatio = fCurrentModel->Ge << 219   theDirectEnergyLossProcess->SetDynamicMassCharge(massRatio,chargeSqRatio);
172       fDirectPartDef, fCurrentMaterial, Tkin); << 220     
173     fDirectEnergyLossProcess->SetDynamicMassCh << 
174   }                                               221   }
175                                                << 222   
176   G4double DEDX_after = fDirectEnergyLossProce << 223   G4double DEDX_after = theDirectEnergyLossProcess->GetDEDX(Tkin, currentCouple);
177   G4double weight_correction = DEDX_after / DE << 224   
178                                                << 225   
                                                   >> 226   G4double weight_correction=DEDX_after/DEDX_before;
                                                   >> 227  
                                                   >> 228   
179   aParticleChange.ProposeEnergy(Tkin);            229   aParticleChange.ProposeEnergy(Tkin);
180                                                << 230   
181   // Caution!!! It is important to select the  << 231   //we still need to register in the particleChange the modification of the weight of the particle 
182   // as the current weight and not the weight  << 232   G4double new_weight=weight_correction*track.GetWeight();
183   // the track is changed after having applied << 
184                                                << 
185   G4double new_weight =                        << 
186     weight_correction * step.GetPostStepPoint( << 
187   aParticleChange.SetParentWeightByProcess(fal    233   aParticleChange.SetParentWeightByProcess(false);
188   aParticleChange.ProposeParentWeight(new_weig    234   aParticleChange.ProposeParentWeight(new_weight);
                                                   >> 235   
189                                                   236 
190   return &aParticleChange;                        237   return &aParticleChange;
191 }                                              << 
192                                                   238 
                                                   >> 239 }
193 //////////////////////////////////////////////    240 ///////////////////////////////////////////////////////
                                                   >> 241 //
194 void G4ContinuousGainOfEnergy::SetLossFluctuat    242 void G4ContinuousGainOfEnergy::SetLossFluctuations(G4bool val)
195 {                                                 243 {
196   if(val && !fLossFluctuationArePossible)      << 244   if(val && !lossFluctuationArePossible) return;
197     return;                                    << 245   lossFluctuationFlag = val;
198   fLossFluctuationFlag = val;                  << 
199 }                                                 246 }
200                                                << 
201 //////////////////////////////////////////////    247 ///////////////////////////////////////////////////////
202 G4double G4ContinuousGainOfEnergy::GetContinuo << 248 //
203                                                << 
204                                                << 
205 {                                              << 
206   DefineMaterial(track.GetMaterialCutsCouple() << 
207                                                << 
208   fPreStepKinEnergy = track.GetKineticEnergy() << 
209   fCurrentModel     = fDirectEnergyLossProcess << 
210     track.GetKineticEnergy() * fMassRatio, fCu << 
211   G4double emax_model           = fCurrentMode << 
212   G4double preStepChargeSqRatio = 0.;          << 
213   if(fIsIon)                                   << 
214   {                                            << 
215     G4double chargeSqRatio = fCurrentModel->Ge << 
216       fDirectPartDef, fCurrentMaterial, fPreSt << 
217     preStepChargeSqRatio = chargeSqRatio;      << 
218     fDirectEnergyLossProcess->SetDynamicMassCh << 
219                                                << 
220   }                                            << 
221                                                << 
222   G4double maxE = 1.1 * fPreStepKinEnergy;     << 
223                                                << 
224   if(fPreStepKinEnergy < fCurrentTcut)         << 
225     maxE = std::min(fCurrentTcut, maxE);       << 
226                                                << 
227   maxE = std::min(emax_model * 1.001, maxE);   << 
228                                                << 
229   G4double preStepRange =                      << 
230     fDirectEnergyLossProcess->GetRange(fPreSte << 
231                                                << 
232   if(fIsIon)                                   << 
233   {                                            << 
234     G4double chargeSqRatioAtEmax = fCurrentMod << 
235       fDirectPartDef, fCurrentMaterial, maxE); << 
236     fDirectEnergyLossProcess->SetDynamicMassCh << 
237                                                << 
238   }                                            << 
239                                                   249 
240   G4double r1 = fDirectEnergyLossProcess->GetR << 
241                                                   250 
242   if(fIsIon)                                   << 
243     fDirectEnergyLossProcess->SetDynamicMassCh << 
244                                                << 
245                                                   251 
246   return std::max(r1 - preStepRange, 0.001 * m << 252 G4double G4ContinuousGainOfEnergy::GetContinuousStepLimit(const G4Track& track,
                                                   >> 253                 G4double , G4double , G4double& )
                                                   >> 254 { 
                                                   >> 255   G4double x = DBL_MAX;
                                                   >> 256   x=.1*mm;
                                                   >> 257  
                                                   >> 258  
                                                   >> 259   DefineMaterial(track.GetMaterialCutsCouple());
                                                   >> 260  
                                                   >> 261   preStepKinEnergy = track.GetKineticEnergy(); 
                                                   >> 262   preStepScaledKinEnergy = track.GetKineticEnergy()*massRatio;
                                                   >> 263   currentModel = theDirectEnergyLossProcess->SelectModelForMaterial(preStepScaledKinEnergy,currentCoupleIndex);
                                                   >> 264   G4double emax_model=currentModel->HighEnergyLimit();
                                                   >> 265   if (IsIon) {
                                                   >> 266     chargeSqRatio =  currentModel->GetChargeSquareRatio(theDirectPartDef,currentMaterial,preStepKinEnergy);
                                                   >> 267   preStepChargeSqRatio = chargeSqRatio;
                                                   >> 268   theDirectEnergyLossProcess->SetDynamicMassCharge(massRatio,preStepChargeSqRatio);
                                                   >> 269   } 
                                                   >> 270   
                                                   >> 271   
                                                   >> 272   G4double maxE =1.1*preStepKinEnergy;
                                                   >> 273   /*if (preStepKinEnergy< 0.05*MeV) maxE =2.*preStepKinEnergy;
                                                   >> 274   else if (preStepKinEnergy< 0.1*MeV) maxE =1.5*preStepKinEnergy;
                                                   >> 275   else if (preStepKinEnergy< 0.5*MeV) maxE =1.25*preStepKinEnergy;*/
                                                   >> 276    
                                                   >> 277   if (preStepKinEnergy < currentTcut) maxE = std::min(currentTcut,maxE);
                                                   >> 278  
                                                   >> 279   maxE=std::min(emax_model*1.001,maxE);
                                                   >> 280     
                                                   >> 281   preStepRange = theDirectEnergyLossProcess->GetRange(preStepKinEnergy, currentCouple);
                                                   >> 282   
                                                   >> 283   if (IsIon) {
                                                   >> 284     G4double chargeSqRatioAtEmax = currentModel->GetChargeSquareRatio(theDirectPartDef,currentMaterial,maxE);
                                                   >> 285   theDirectEnergyLossProcess->SetDynamicMassCharge(massRatio,chargeSqRatioAtEmax);
                                                   >> 286   } 
                                                   >> 287   
                                                   >> 288   G4double r1 = theDirectEnergyLossProcess->GetRange(maxE, currentCouple);
                                                   >> 289   
                                                   >> 290   if (IsIon) theDirectEnergyLossProcess->SetDynamicMassCharge(massRatio,preStepChargeSqRatio);
                                                   >> 291   
                                                   >> 292   
                                                   >> 293 
                                                   >> 294   x=r1-preStepRange;
                                                   >> 295   x=std::max(r1-preStepRange,0.001*mm);
                                                   >> 296  
                                                   >> 297   return x;
                                                   >> 298   
                                                   >> 299  
247 }                                                 300 }
248                                                << 301 #include "G4EmCorrections.hh"
249 //////////////////////////////////////////////    302 ///////////////////////////////////////////////////////
250 void G4ContinuousGainOfEnergy::SetDynamicMassC << 303 //
251                                                << 304 
252 {                                              << 305 void G4ContinuousGainOfEnergy::SetDynamicMassCharge(const G4Track& ,G4double energy)
253   G4double ChargeSqRatio =                     << 306 { 
254     G4LossTableManager::Instance()->EmCorrecti << 307 
255       fDirectPartDef, fCurrentMaterial, energy << 308   G4double ChargeSqRatio= G4LossTableManager::Instance()->EmCorrections()->EffectiveChargeSquareRatio(theDirectPartDef,currentMaterial,energy); 
256   if(fDirectEnergyLossProcess)                 << 309   if (theDirectEnergyLossProcess) theDirectEnergyLossProcess->SetDynamicMassCharge(massRatio,ChargeSqRatio);
257     fDirectEnergyLossProcess->SetDynamicMassCh << 
258 }                                                 310 }
259                                                   311