Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/lowenergy/src/G4eIonisationSpectrum.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/lowenergy/src/G4eIonisationSpectrum.cc (Version 11.3.0) and /processes/electromagnetic/lowenergy/src/G4eIonisationSpectrum.cc (Version 3.1)


  1 //                                                  1 
  2 // *******************************************    
  3 // * License and Disclaimer                       
  4 // *                                              
  5 // * The  Geant4 software  is  copyright of th    
  6 // * the Geant4 Collaboration.  It is provided    
  7 // * conditions of the Geant4 Software License    
  8 // * LICENSE and available at  http://cern.ch/    
  9 // * include a list of copyright holders.         
 10 // *                                              
 11 // * Neither the authors of this software syst    
 12 // * institutes,nor the agencies providing fin    
 13 // * work  make  any representation or  warran    
 14 // * regarding  this  software system or assum    
 15 // * use.  Please see the license in the file     
 16 // * for the full disclaimer and the limitatio    
 17 // *                                              
 18 // * This  code  implementation is the result     
 19 // * technical work of the GEANT4 collaboratio    
 20 // * By using,  copying,  modifying or  distri    
 21 // * any work based  on the software)  you  ag    
 22 // * use  in  resulting  scientific  publicati    
 23 // * acceptance of all terms of the Geant4 Sof    
 24 // *******************************************    
 25 //                                                
 26 //                                                
 27 // -------------------------------------------    
 28 //                                                
 29 // GEANT4 Class file                              
 30 //                                                
 31 //                                                
 32 // File name:     G4eIonisationSpectrum           
 33 //                                                
 34 // Author:        V.Ivanchenko (Vladimir.Ivanc    
 35 //                                                
 36 // Creation date: 29 September 2001               
 37 //                                                
 38 // Modifications:                                 
 39 // 10.10.2001 MGP           Revision to improv    
 40 //                          consistency with d    
 41 // 02.11.2001 VI            Optimize sampling     
 42 // 29.11.2001 VI            New parametrisatio    
 43 // 19.04.2002 VI            Add protection in     
 44 // 30.05.2002 VI            Update to 24-param    
 45 // 11.07.2002 VI            Fix in integration    
 46 // 23.03.2009 LP            Added protection a    
 47 //                          faulty database fi    
 48 //                                                
 49 // -------------------------------------------    
 50 //                                                
 51                                                   
 52 #include "G4eIonisationSpectrum.hh"               
 53 #include "G4AtomicTransitionManager.hh"           
 54 #include "G4AtomicShell.hh"                       
 55 #include "G4DataVector.hh"                        
 56 #include "Randomize.hh"                           
 57 #include "G4PhysicalConstants.hh"                 
 58 #include "G4SystemOfUnits.hh"                     
 59 #include "G4Exp.hh"                               
 60                                                   
 61 G4eIonisationSpectrum::G4eIonisationSpectrum()    
 62   lowestE(0.1*eV),                                
 63   factor(1.3),                                    
 64   iMax(24),                                       
 65   verbose(0)                                      
 66 {                                                 
 67   theParam = new G4eIonisationParameters();       
 68 }                                                 
 69                                                   
 70                                                   
 71 G4eIonisationSpectrum::~G4eIonisationSpectrum(    
 72 {                                                 
 73   delete theParam;                                
 74 }                                                 
 75                                                   
 76                                                   
 77 G4double G4eIonisationSpectrum::Probability(G4    
 78                     G4double tMin,                
 79               G4double tMax,                      
 80               G4double e,                         
 81               G4int shell,                        
 82               const G4ParticleDefinition* ) co    
 83 {                                                 
 84   // Please comment what Probability does and     
 85   // functions mentioned below                    
 86   // Describe the algorithms used                 
 87                                                   
 88   G4double eMax = MaxEnergyOfSecondaries(e);      
 89   G4double t0 = std::max(tMin, lowestE);          
 90   G4double tm = std::min(tMax, eMax);             
 91   if(t0 >= tm) return 0.0;                        
 92                                                   
 93   G4double bindingEnergy = (G4AtomicTransition    
 94                            Shell(Z, shell)->Bi    
 95                                                   
 96   if(e <= bindingEnergy) return 0.0;              
 97                                                   
 98   G4double energy = e + bindingEnergy;            
 99                                                   
100   G4double x1 = std::min(0.5,(t0 + bindingEner    
101   G4double x2 = std::min(0.5,(tm + bindingEner    
102                                                   
103   if(verbose > 1 || (Z==4 && e>= 1.0 && e<= 0.    
104     G4cout << "G4eIonisationSpectrum::Probabil    
105            << "; shell= " << shell                
106            << "; E(keV)= " << e/keV               
107            << "; Eb(keV)= " << bindingEnergy/k    
108            << "; x1= " << x1                      
109            << "; x2= " << x2                      
110            << G4endl;                             
111                                                   
112   }                                               
113                                                   
114   G4DataVector p;                                 
115                                                   
116   // Access parameters                            
117   for (G4int i=0; i<iMax; i++)                    
118   {                                               
119     G4double x = theParam->Parameter(Z, shell,    
120     if(i<4) x /= energy;                          
121     p.push_back(x);                               
122   }                                               
123                                                   
124   if(p[3] > 0.5) p[3] = 0.5;                      
125                                                   
126   G4double gLocal = energy/electron_mass_c2 +     
127   p.push_back((2.0*gLocal - 1.0)/(gLocal*gLoca    
128                                                   
129   //Add protection against division by zero: a    
130   //parameter p[3] appears in the denominator.    
131   if (p[3] > 0)                                   
132     p[iMax-1] = Function(p[3], p);                
133   else                                            
134     {                                             
135       G4cout << "WARNING: G4eIonisationSpectru    
136        << "parameter p[3] <= 0. G4LEDATA dabat    
137        << Z << ". Please check and/or update i    
138     }                                             
139                                                   
140   if(e >= 1. && e <= 0. && Z == 4) p.push_back    
141                                                   
142                                                   
143   G4double val = IntSpectrum(x1, x2, p);          
144   G4double x0  = (lowestE + bindingEnergy)/ene    
145   G4double nor = IntSpectrum(x0, 0.5, p);         
146                                                   
147   if(verbose > 1 || (Z==4 && e>= 1.0 && e<= 0.    
148     G4cout << "tcut= " << tMin                    
149            << "; tMax= " << tMax                  
150            << "; x0= " << x0                      
151            << "; x1= " << x1                      
152            << "; x2= " << x2                      
153            << "; val= " << val                    
154            << "; nor= " << nor                    
155            << "; sum= " << p[0]                   
156            << "; a= " << p[1]                     
157            << "; b= " << p[2]                     
158            << "; c= " << p[3]                     
159            << G4endl;                             
160     if(shell == 1) G4cout << "============" <<    
161   }                                               
162                                                   
163   p.clear();                                      
164                                                   
165   if(nor > 0.0) val /= nor;                       
166   else          val  = 0.0;                       
167                                                   
168   return val;                                     
169 }                                                 
170                                                   
171                                                   
172 G4double G4eIonisationSpectrum::AverageEnergy(    
173                       G4double tMin,              
174                 G4double tMax,                    
175                 G4double e,                       
176                 G4int shell,                      
177                 const G4ParticleDefinition* )     
178 {                                                 
179   // Please comment what AverageEnergy does an    
180   // functions mentioned below                    
181   // Describe the algorithms used                 
182                                                   
183   G4double eMax = MaxEnergyOfSecondaries(e);      
184   G4double t0 = std::max(tMin, lowestE);          
185   G4double tm = std::min(tMax, eMax);             
186   if(t0 >= tm) return 0.0;                        
187                                                   
188   G4double bindingEnergy = (G4AtomicTransition    
189                            Shell(Z, shell)->Bi    
190                                                   
191   if(e <= bindingEnergy) return 0.0;              
192                                                   
193   G4double energy = e + bindingEnergy;            
194                                                   
195   G4double x1 = std::min(0.5,(t0 + bindingEner    
196   G4double x2 = std::min(0.5,(tm + bindingEner    
197                                                   
198   if(verbose > 1) {                               
199     G4cout << "G4eIonisationSpectrum::AverageE    
200            << "; shell= " << shell                
201            << "; E(keV)= " << e/keV               
202            << "; bindingE(keV)= " << bindingEn    
203            << "; x1= " << x1                      
204            << "; x2= " << x2                      
205            << G4endl;                             
206   }                                               
207                                                   
208   G4DataVector p;                                 
209                                                   
210   // Access parameters                            
211   for (G4int i=0; i<iMax; i++)                    
212   {                                               
213     G4double x = theParam->Parameter(Z, shell,    
214     if(i<4) x /= energy;                          
215     p.push_back(x);                               
216   }                                               
217                                                   
218   if(p[3] > 0.5) p[3] = 0.5;                      
219                                                   
220   G4double gLocal2 = energy/electron_mass_c2 +    
221   p.push_back((2.0*gLocal2 - 1.0)/(gLocal2*gLo    
222                                                   
223                                                   
224   //Add protection against division by zero: a    
225   //parameter p[3] appears in the denominator.    
226   if (p[3] > 0)                                   
227     p[iMax-1] = Function(p[3], p);                
228   else                                            
229     {                                             
230       G4cout << "WARNING: G4eIonisationSpectru    
231        << "parameter p[3] <= 0. G4LEDATA dabat    
232        << Z << ". Please check and/or update i    
233     }                                             
234                                                   
235   G4double val = AverageValue(x1, x2, p);         
236   G4double x0  = (lowestE + bindingEnergy)/ene    
237   G4double nor = IntSpectrum(x0, 0.5, p);         
238   val *= energy;                                  
239                                                   
240   if(verbose > 1) {                               
241     G4cout << "tcut(MeV)= " << tMin/MeV           
242            << "; tMax(MeV)= " << tMax/MeV         
243            << "; x0= " << x0                      
244            << "; x1= " << x1                      
245            << "; x2= " << x2                      
246            << "; val= " << val                    
247            << "; nor= " << nor                    
248            << "; sum= " << p[0]                   
249            << "; a= " << p[1]                     
250            << "; b= " << p[2]                     
251            << "; c= " << p[3]                     
252            << G4endl;                             
253   }                                               
254                                                   
255   p.clear();                                      
256                                                   
257   if(nor > 0.0) val /= nor;                       
258   else          val  = 0.0;                       
259                                                   
260   return val;                                     
261 }                                                 
262                                                   
263                                                   
264 G4double G4eIonisationSpectrum::SampleEnergy(G    
265                G4double tMin,                     
266                G4double tMax,                     
267                      G4double e,                  
268                      G4int shell,                 
269                const G4ParticleDefinition* ) c    
270 {                                                 
271   // Please comment what SampleEnergy does        
272   G4double tDelta = 0.0;                          
273   G4double t0 = std::max(tMin, lowestE);          
274   G4double tm = std::min(tMax, MaxEnergyOfSeco    
275   if(t0 > tm) return tDelta;                      
276                                                   
277   G4double bindingEnergy = (G4AtomicTransition    
278                            Shell(Z, shell)->Bi    
279                                                   
280   if(e <= bindingEnergy) return 0.0;              
281                                                   
282   G4double energy = e + bindingEnergy;            
283                                                   
284   G4double x1 = std::min(0.5,(t0 + bindingEner    
285   G4double x2 = std::min(0.5,(tm + bindingEner    
286   if(x1 >= x2) return tDelta;                     
287                                                   
288   if(verbose > 1) {                               
289     G4cout << "G4eIonisationSpectrum::SampleEn    
290            << "; shell= " << shell                
291            << "; E(keV)= " << e/keV               
292            << G4endl;                             
293   }                                               
294                                                   
295   // Access parameters                            
296   G4DataVector p;                                 
297                                                   
298   // Access parameters                            
299   for (G4int i=0; i<iMax; i++)                    
300   {                                               
301     G4double x = theParam->Parameter(Z, shell,    
302     if(i<4) x /= energy;                          
303     p.push_back(x);                               
304   }                                               
305                                                   
306   if(p[3] > 0.5) p[3] = 0.5;                      
307                                                   
308   G4double gLocal3 = energy/electron_mass_c2 +    
309   p.push_back((2.0*gLocal3 - 1.0)/(gLocal3*gLo    
310                                                   
311                                                   
312   //Add protection against division by zero: a    
313   //parameter p[3] appears in the denominator.    
314   if (p[3] > 0)                                   
315     p[iMax-1] = Function(p[3], p);                
316   else                                            
317     {                                             
318       G4cout << "WARNING: G4eIonisationSpectru    
319        << "parameter p[3] <= 0. G4LEDATA dabat    
320        << Z << ". Please check and/or update i    
321     }                                             
322                                                   
323   G4double aria1 = 0.0;                           
324   G4double a1 = std::max(x1,p[1]);                
325   G4double a2 = std::min(x2,p[3]);                
326   if(a1 < a2) aria1 = IntSpectrum(a1, a2, p);     
327   G4double aria2 = 0.0;                           
328   G4double a3 = std::max(x1,p[3]);                
329   G4double a4 = x2;                               
330   if(a3 < a4) aria2 = IntSpectrum(a3, a4, p);     
331                                                   
332   G4double aria = (aria1 + aria2)*G4UniformRan    
333   G4double amaj, fun, q, x, z1, z2, dx, dx1;      
334                                                   
335   //======= First aria to sample =====            
336                                                   
337   if(aria <= aria1) {                             
338                                                   
339     amaj = p[4];                                  
340     for (G4int j=5; j<iMax; j++) {                
341       if(p[j] > amaj) amaj = p[j];                
342     }                                             
343                                                   
344     a1 = 1./a1;                                   
345     a2 = 1./a2;                                   
346                                                   
347     G4int i;                                      
348     do {                                          
349                                                   
350       x = 1./(a2 + G4UniformRand()*(a1 - a2));    
351       z1 = p[1];                                  
352       z2 = p[3];                                  
353       dx = (p[2] - p[1]) / 3.0;                   
354       dx1= G4Exp(std::log(p[3]/p[2]) / 16.0);     
355       for (i=4; i<iMax-1; i++) {                  
356                                                   
357         if (i < 7) {                              
358           z2 = z1 + dx;                           
359         } else if(iMax-2 == i) {                  
360           z2 = p[3];                              
361           break;                                  
362         } else {                                  
363           z2 = z1*dx1;                            
364         }                                         
365         if(x >= z1 && x <= z2) break;             
366         z1 = z2;                                  
367       }                                           
368       fun = p[i] + (x - z1) * (p[i+1] - p[i])/    
369                                                   
370       if(fun > amaj) {                            
371           G4cout << "WARNING in G4eIonisationS    
372                  << " Majoranta " << amaj         
373                  << " < " << fun                  
374                  << " in the first aria at x=     
375                  << G4endl;                       
376       }                                           
377                                                   
378       q = amaj*G4UniformRand();                   
379                                                   
380     } while (q >= fun);                           
381                                                   
382   //======= Second aria to sample =====           
383                                                   
384   } else {                                        
385                                                   
386     amaj = std::max(p[iMax-1], Function(0.5, p    
387     a1 = 1./a3;                                   
388     a2 = 1./a4;                                   
389                                                   
390     do {                                          
391                                                   
392       x = 1./(a2 + G4UniformRand()*(a1 - a2));    
393       fun = Function(x, p);                       
394                                                   
395       if(fun > amaj) {                            
396           G4cout << "WARNING in G4eIonisationS    
397                  << " Majoranta " << amaj         
398                  << " < " << fun                  
399                  << " in the second aria at x=    
400                  << G4endl;                       
401       }                                           
402                                                   
403       q = amaj*G4UniformRand();                   
404                                                   
405     } while (q >= fun);                           
406                                                   
407   }                                               
408                                                   
409   p.clear();                                      
410                                                   
411   tDelta = x*energy - bindingEnergy;              
412                                                   
413   if(verbose > 1) {                               
414     G4cout << "tcut(MeV)= " << tMin/MeV           
415            << "; tMax(MeV)= " << tMax/MeV         
416            << "; x1= " << x1                      
417            << "; x2= " << x2                      
418            << "; a1= " << a1                      
419            << "; a2= " << a2                      
420            << "; x= " << x                        
421            << "; be= " << bindingEnergy           
422            << "; e= " << e                        
423            << "; tDelta= " << tDelta              
424            << G4endl;                             
425   }                                               
426                                                   
427                                                   
428   return tDelta;                                  
429 }                                                 
430                                                   
431                                                   
432 G4double G4eIonisationSpectrum::IntSpectrum(G4    
433               G4double xMax,                      
434               const G4DataVector& p) const        
435 {                                                 
436   // Please comment what IntSpectrum does         
437   G4double sum = 0.0;                             
438   if(xMin >= xMax) return sum;                    
439                                                   
440   G4double x1, x2, xs1, xs2, y1, y2, ys1, ys2,    
441                                                   
442   // Integral over interpolation aria             
443   if(xMin < p[3]) {                               
444                                                   
445     x1 = p[1];                                    
446     y1 = p[4];                                    
447                                                   
448     G4double dx = (p[2] - p[1]) / 3.0;            
449     G4double dx1= G4Exp(std::log(p[3]/p[2]) /     
450                                                   
451     for (size_t i=0; i<19; i++) {                 
452                                                   
453       q = 0.0;                                    
454       if (i < 3) {                                
455         x2 = x1 + dx;                             
456       } else if(18 == i) {                        
457         x2 = p[3];                                
458       } else {                                    
459         x2 = x1*dx1;                              
460       }                                           
461                                                   
462       y2 = p[5 + i];                              
463                                                   
464       if (xMax <= x1) {                           
465         break;                                    
466       } else if (xMin < x2) {                     
467                                                   
468         xs1 = x1;                                 
469         xs2 = x2;                                 
470         ys1 = y1;                                 
471         ys2 = y2;                                 
472                                                   
473         if (x2 > x1) {                            
474           if (xMin > x1) {                        
475             xs1 = xMin;                           
476             ys1 += (xs1 - x1)*(y2 - y1)/(x2 -     
477     }                                             
478           if (xMax < x2) {                        
479             xs2 = xMax;                           
480             ys2 += (xs2 - x2)*(y1 - y2)/(x1 -     
481     }                                             
482           if (xs2 > xs1) {                        
483             q = (ys1*xs2 - ys2*xs1)/(xs1*xs2)     
484               +  std::log(xs2/xs1)*(ys2 - ys1)    
485             sum += q;                             
486             if(p.size() == 26) G4cout << "i= "    
487     }                                             
488   }                                               
489       }                                           
490       x1 = x2;                                    
491       y1 = y2;                                    
492     }                                             
493   }                                               
494                                                   
495   // Integral over aria with parametrised form    
496                                                   
497   x1 = std::max(xMin, p[3]);                      
498   if(x1 >= xMax) return sum;                      
499   x2 = xMax;                                      
500                                                   
501   xs1 = 1./x1;                                    
502   xs2 = 1./x2;                                    
503   q = (xs1 - xs2)*(1.0 - p[0])                    
504        - p[iMax]*std::log(x2/x1)                  
505        + (1. - p[iMax])*(x2 - x1)                 
506        + 1./(1. - x2) - 1./(1. - x1)              
507        + p[iMax]*std::log((1. - x2)/(1. - x1))    
508        + 0.25*p[0]*(xs1*xs1 - xs2*xs2);           
509   sum += q;                                       
510   if(p.size() == 26) G4cout << "param...  q= "    
511                                                   
512   return sum;                                     
513 }                                                 
514                                                   
515                                                   
516 G4double G4eIonisationSpectrum::AverageValue(G    
517                      G4double xMax,               
518                const G4DataVector& p) const       
519 {                                                 
520   G4double sum = 0.0;                             
521   if(xMin >= xMax) return sum;                    
522                                                   
523   G4double x1, x2, xs1, xs2, y1, y2, ys1, ys2;    
524                                                   
525   // Integral over interpolation aria             
526   if(xMin < p[3]) {                               
527                                                   
528     x1 = p[1];                                    
529     y1 = p[4];                                    
530                                                   
531     G4double dx = (p[2] - p[1]) / 3.0;            
532     G4double dx1= G4Exp(std::log(p[3]/p[2]) /     
533                                                   
534     for (size_t i=0; i<19; i++) {                 
535                                                   
536       if (i < 3) {                                
537         x2 = x1 + dx;                             
538       } else if(18 == i) {                        
539         x2 = p[3];                                
540       } else {                                    
541         x2 = x1*dx1;                              
542       }                                           
543                                                   
544       y2 = p[5 + i];                              
545                                                   
546       if (xMax <= x1) {                           
547         break;                                    
548       } else if (xMin < x2) {                     
549                                                   
550         xs1 = x1;                                 
551         xs2 = x2;                                 
552         ys1 = y1;                                 
553         ys2 = y2;                                 
554                                                   
555         if (x2 > x1) {                            
556           if (xMin > x1) {                        
557             xs1 = xMin;                           
558             ys1 += (xs1 - x1)*(y2 - y1)/(x2 -     
559     }                                             
560           if (xMax < x2) {                        
561             xs2 = xMax;                           
562             ys2 += (xs2 - x2)*(y1 - y2)/(x1 -     
563     }                                             
564           if (xs2 > xs1) {                        
565             sum += std::log(xs2/xs1)*(ys1*xs2     
566                 +  ys2 - ys1;                     
567     }                                             
568   }                                               
569       }                                           
570       x1 = x2;                                    
571       y1 = y2;                                    
572                                                   
573     }                                             
574   }                                               
575                                                   
576   // Integral over aria with parametrised form    
577                                                   
578   x1 = std::max(xMin, p[3]);                      
579   if(x1 >= xMax) return sum;                      
580   x2 = xMax;                                      
581                                                   
582   xs1 = 1./x1;                                    
583   xs2 = 1./x2;                                    
584                                                   
585   sum  += std::log(x2/x1)*(1.0 - p[0])            
586         + 0.5*(1. - p[iMax])*(x2*x2 - x1*x1)      
587         + 1./(1. - x2) - 1./(1. - x1)             
588         + (1. + p[iMax])*std::log((1. - x2)/(1    
589         + 0.5*p[0]*(xs1 - xs2);                   
590                                                   
591   return sum;                                     
592 }                                                 
593                                                   
594                                                   
595 void G4eIonisationSpectrum::PrintData() const     
596 {                                                 
597   theParam->PrintData();                          
598 }                                                 
599                                                   
600 G4double G4eIonisationSpectrum::MaxEnergyOfSec    
601                    G4int, // Z = 0,               
602                    const G4ParticleDefinition*    
603 {                                                 
604   return 0.5 * kineticEnergy;                     
605 }                                                 
606