source: trunk/source/processes/electromagnetic/standard/src/G4PAIySection.cc@ 1201

Last change on this file since 1201 was 1196, checked in by garnier, 16 years ago

update CVS release candidate geant4.9.3.01
File size: 36.2 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26//
27// $Id: G4PAIySection.cc,v 1.4 2009/07/26 15:51:01 vnivanch Exp $
28// GEANT4 tag $Name: geant4-09-03-cand-01 $
29//
30//
31// G4PAIySection.cc -- class implementation file
32//
33// GEANT 4 class implementation file
34//
35// For information related to this code, please, contact
36// the Geant4 Collaboration.
37//
38// R&D: Vladimir.Grichine@cern.ch
39//
40// History:
41//
42// 01.10.07, V.Ivanchenko create using V.Grichine G4PAIxSection class
43// 26.07.09, V.Ivanchenko added protection for mumerical exceptions for
44// low-density materials
45//
46
47#include "G4PAIySection.hh"
48
49#include "globals.hh"
50#include "G4ios.hh"
51#include "G4Poisson.hh"
52#include "G4Material.hh"
53#include "G4MaterialCutsCouple.hh"
54#include "G4SandiaTable.hh"
55
56using namespace std;
57
58// Local class constants
59
60const G4double G4PAIySection::fDelta = 0.005 ; // energy shift from interval border
61const G4double G4PAIySection::fError = 0.005 ; // error in lin-log approximation
62
63const G4int G4PAIySection::fMaxSplineSize = 500 ; // Max size of output spline
64 // arrays
65
66//////////////////////////////////////////////////////////////////
67//
68// Constructor
69//
70
71G4PAIySection::G4PAIySection()
72{}
73
74////////////////////////////////////////////////////////////////////////////
75//
76// Destructor
77
78G4PAIySection::~G4PAIySection()
79{}
80
81////////////////////////////////////////////////////////////////////////
82//
83// Test Constructor with beta*gamma square value
84
85void G4PAIySection::Initialize( const G4Material* material,
86 G4double maxEnergyTransfer,
87 G4double betaGammaSq)
88{
89 G4int i, j, numberOfElements ;
90
91 fDensity = material->GetDensity();
92 fElectronDensity = material->GetElectronDensity() ;
93 numberOfElements = material->GetNumberOfElements() ;
94
95 fSandia = material->GetSandiaTable();
96
97 fIntervalNumber = fSandia->GetMaxInterval();
98
99 fIntervalNumber--;
100
101 for(i=1;i<=fIntervalNumber;i++)
102 {
103 G4double e = fSandia->GetSandiaMatTablePAI(i,0);
104 if(e >= maxEnergyTransfer || i > fIntervalNumber)
105 {
106 fEnergyInterval[i] = maxEnergyTransfer ;
107 fIntervalNumber = i ;
108 break;
109 }
110 fEnergyInterval[i] = e;
111 fA1[i] = fSandia->GetSandiaMatTablePAI(i,1);
112 fA2[i] = fSandia->GetSandiaMatTablePAI(i,2);
113 fA3[i] = fSandia->GetSandiaMatTablePAI(i,3);
114 fA4[i] = fSandia->GetSandiaMatTablePAI(i,4);
115
116 }
117 if(fEnergyInterval[fIntervalNumber] != maxEnergyTransfer)
118 {
119 fIntervalNumber++;
120 fEnergyInterval[fIntervalNumber] = maxEnergyTransfer ;
121 fA1[fIntervalNumber] = fA1[fIntervalNumber-1] ;
122 fA2[fIntervalNumber] = fA2[fIntervalNumber-1] ;
123 fA3[fIntervalNumber] = fA3[fIntervalNumber-1] ;
124 fA4[fIntervalNumber] = fA4[fIntervalNumber-1] ;
125 }
126
127 // Now checking, if two borders are too close together
128 for(i=1;i<fIntervalNumber;i++)
129 {
130 // G4cout<<fEnergyInterval[i]<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
131 // <<fA3[i]<<"\t"<<fA4[i]<<"\t"<<G4endl ;
132 if(fEnergyInterval[i+1]-fEnergyInterval[i] <
133 1.5*fDelta*(fEnergyInterval[i+1]+fEnergyInterval[i]))
134 {
135 for(j=i;j<fIntervalNumber;j++)
136 {
137 fEnergyInterval[j] = fEnergyInterval[j+1] ;
138 fA1[j] = fA1[j+1] ;
139 fA2[j] = fA2[j+1] ;
140 fA3[j] = fA3[j+1] ;
141 fA4[j] = fA4[j+1] ;
142 }
143 fIntervalNumber-- ;
144 i-- ;
145 }
146 }
147
148 // Preparation of fSplineEnergy array corresponding to min ionisation, G~4
149
150 G4double betaGammaSqRef =
151 fLorentzFactor[fRefGammaNumber]*fLorentzFactor[fRefGammaNumber] - 1;
152
153 NormShift(betaGammaSqRef) ;
154 SplainPAI(betaGammaSqRef) ;
155
156 // Preparation of integral PAI cross section for input betaGammaSq
157
158 for(i = 1 ; i <= fSplineNumber ; i++)
159 {
160 fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
161 fdNdxCerenkov[i] = PAIdNdxCerenkov(i,betaGammaSq);
162 fdNdxPlasmon[i] = PAIdNdxPlasmon(i,betaGammaSq);
163 }
164 IntegralPAIySection() ;
165 IntegralCerenkov() ;
166 IntegralPlasmon() ;
167}
168
169/////////////////////////////////////////////////////////////////////////
170//
171// General control function for class G4PAIySection
172//
173
174void G4PAIySection::InitPAI()
175{
176 G4int i ;
177 G4double betaGammaSq = fLorentzFactor[fRefGammaNumber]*
178 fLorentzFactor[fRefGammaNumber] - 1;
179
180 // Preparation of integral PAI cross section for reference gamma
181
182 NormShift(betaGammaSq) ;
183 SplainPAI(betaGammaSq) ;
184
185 IntegralPAIySection() ;
186 IntegralCerenkov() ;
187 IntegralPlasmon() ;
188
189 for(i = 0 ; i<=fSplineNumber ; i++)
190 {
191 fPAItable[i][fRefGammaNumber] = fIntegralPAIySection[i] ;
192 if(i != 0)
193 {
194 fPAItable[i][0] = fSplineEnergy[i] ;
195 }
196 }
197 fPAItable[0][0] = fSplineNumber ;
198
199 for(G4int j = 1 ; j < 112 ; j++) // for other gammas
200 {
201 if( j == fRefGammaNumber ) continue ;
202
203 betaGammaSq = fLorentzFactor[j]*fLorentzFactor[j] - 1 ;
204
205 for(i = 1 ; i <= fSplineNumber ; i++)
206 {
207 fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
208 fdNdxCerenkov[i] = PAIdNdxCerenkov(i,betaGammaSq);
209 fdNdxPlasmon[i] = PAIdNdxPlasmon(i,betaGammaSq);
210 }
211 IntegralPAIySection() ;
212 IntegralCerenkov() ;
213 IntegralPlasmon() ;
214
215 for(i = 0 ; i <= fSplineNumber ; i++)
216 {
217 fPAItable[i][j] = fIntegralPAIySection[i] ;
218 }
219 }
220
221}
222
223///////////////////////////////////////////////////////////////////////
224//
225// Shifting from borders to intervals Creation of first energy points
226//
227
228void G4PAIySection::NormShift(G4double betaGammaSq)
229{
230 G4int i, j ;
231
232 for( i = 1 ; i <= fIntervalNumber-1 ; i++ )
233 {
234 for( j = 1 ; j <= 2 ; j++ )
235 {
236 fSplineNumber = (i-1)*2 + j ;
237
238 if( j == 1 ) fSplineEnergy[fSplineNumber] = fEnergyInterval[i ]*(1+fDelta);
239 else fSplineEnergy[fSplineNumber] = fEnergyInterval[i+1]*(1-fDelta);
240 // G4cout<<"cn = "<<fSplineNumber<<"; "<<"energy = "
241 // <<fSplineEnergy[fSplineNumber]<<G4endl;
242 }
243 }
244 fIntegralTerm[1]=RutherfordIntegral(1,fEnergyInterval[1],fSplineEnergy[1]);
245
246 j = 1 ;
247
248 for(i=2;i<=fSplineNumber;i++)
249 {
250 if(fSplineEnergy[i]<fEnergyInterval[j+1])
251 {
252 fIntegralTerm[i] = fIntegralTerm[i-1] +
253 RutherfordIntegral(j,fSplineEnergy[i-1],
254 fSplineEnergy[i] ) ;
255 }
256 else
257 {
258 G4double x = RutherfordIntegral(j,fSplineEnergy[i-1],
259 fEnergyInterval[j+1] ) ;
260 j++;
261 fIntegralTerm[i] = fIntegralTerm[i-1] + x +
262 RutherfordIntegral(j,fEnergyInterval[j],
263 fSplineEnergy[i] ) ;
264 }
265 // G4cout<<i<<"\t"<<fSplineEnergy[i]<<"\t"<<fIntegralTerm[i]<<"\n"<<G4endl;
266 }
267 fNormalizationCof = 2*pi*pi*hbarc*hbarc*fine_structure_const/electron_mass_c2 ;
268 fNormalizationCof *= fElectronDensity/fIntegralTerm[fSplineNumber] ;
269
270 // G4cout<<"fNormalizationCof = "<<fNormalizationCof<<G4endl ;
271
272 // Calculation of PAI differrential cross-section (1/(keV*cm))
273 // in the energy points near borders of energy intervals
274
275 for(G4int k=1;k<=fIntervalNumber-1;k++)
276 {
277 for(j=1;j<=2;j++)
278 {
279 i = (k-1)*2 + j ;
280 fImPartDielectricConst[i] = fNormalizationCof*
281 ImPartDielectricConst(k,fSplineEnergy[i]);
282 fRePartDielectricConst[i] = fNormalizationCof*
283 RePartDielectricConst(fSplineEnergy[i]);
284 fIntegralTerm[i] *= fNormalizationCof;
285
286 fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
287 fdNdxCerenkov[i] = PAIdNdxCerenkov(i,betaGammaSq);
288 fdNdxPlasmon[i] = PAIdNdxPlasmon(i,betaGammaSq);
289 }
290 }
291
292} // end of NormShift
293
294/////////////////////////////////////////////////////////////////////////
295//
296// Creation of new energy points as geometrical mean of existing
297// one, calculation PAI_cs for them, while the error of logarithmic
298// linear approximation would be smaller than 'fError'
299
300void G4PAIySection::SplainPAI(G4double betaGammaSq)
301{
302 G4int k = 1 ;
303 G4int i = 1 ;
304
305 while ( (i < fSplineNumber) && (fSplineNumber < fMaxSplineSize-1) )
306 {
307 if(fSplineEnergy[i+1] > fEnergyInterval[k+1])
308 {
309 k++ ; // Here next energy point is in next energy interval
310 i++;
311 continue;
312 }
313 // Shifting of arrayes for inserting the geometrical
314 // average of 'i' and 'i+1' energy points to 'i+1' place
315 fSplineNumber++;
316
317 for(G4int j = fSplineNumber; j >= i+2 ; j-- )
318 {
319 fSplineEnergy[j] = fSplineEnergy[j-1];
320 fImPartDielectricConst[j] = fImPartDielectricConst[j-1];
321 fRePartDielectricConst[j] = fRePartDielectricConst[j-1];
322 fIntegralTerm[j] = fIntegralTerm[j-1];
323
324 fDifPAIySection[j] = fDifPAIySection[j-1];
325 fdNdxCerenkov[j] = fdNdxCerenkov[j-1];
326 fdNdxPlasmon[j] = fdNdxPlasmon[j-1];
327 }
328 G4double x1 = fSplineEnergy[i];
329 G4double x2 = fSplineEnergy[i+1];
330 G4double yy1 = fDifPAIySection[i];
331 G4double y2 = fDifPAIySection[i+1];
332
333 G4double en1 = sqrt(x1*x2);
334 fSplineEnergy[i+1] = en1;
335
336 // Calculation of logarithmic linear approximation
337 // in this (enr) energy point, which number is 'i+1' now
338
339 G4double a = log10(y2/yy1)/log10(x2/x1);
340 G4double b = log10(yy1) - a*log10(x1);
341 G4double y = a*log10(en1) + b ;
342 y = pow(10.,y);
343
344 // Calculation of the PAI dif. cross-section at this point
345
346 fImPartDielectricConst[i+1] = fNormalizationCof*
347 ImPartDielectricConst(k,fSplineEnergy[i+1]);
348 fRePartDielectricConst[i+1] = fNormalizationCof*
349 RePartDielectricConst(fSplineEnergy[i+1]);
350 fIntegralTerm[i+1] = fIntegralTerm[i] + fNormalizationCof*
351 RutherfordIntegral(k,fSplineEnergy[i],
352 fSplineEnergy[i+1]);
353
354 fDifPAIySection[i+1] = DifPAIySection(i+1,betaGammaSq);
355 fdNdxCerenkov[i+1] = PAIdNdxCerenkov(i+1,betaGammaSq);
356 fdNdxPlasmon[i+1] = PAIdNdxPlasmon(i+1,betaGammaSq);
357
358 // Condition for next division of this segment or to pass
359 // to higher energies
360
361 G4double x = 2*(fDifPAIySection[i+1] - y)/(fDifPAIySection[i+1] + y);
362
363 if( x < 0 )
364 {
365 x = -x ;
366 }
367 if( x > fError && fSplineNumber < fMaxSplineSize-1 )
368 {
369 continue; // next division
370 }
371 i += 2; // pass to next segment
372
373 } // close 'while'
374
375} // end of SplainPAI
376
377
378////////////////////////////////////////////////////////////////////
379//
380// Integration over electrons that could be considered
381// quasi-free at energy transfer of interest
382
383G4double G4PAIySection::RutherfordIntegral( G4int k,
384 G4double x1,
385 G4double x2 )
386{
387 G4double c1, c2, c3 ;
388 // G4cout<<"RI: x1 = "<<x1<<"; "<<"x2 = "<<x2<<G4endl;
389 c1 = (x2 - x1)/x1/x2 ;
390 c2 = (x2 - x1)*(x2 + x1)/x1/x1/x2/x2 ;
391 c3 = (x2 - x1)*(x1*x1 + x1*x2 + x2*x2)/x1/x1/x1/x2/x2/x2 ;
392 // G4cout<<" RI: c1 = "<<c1<<"; "<<"c2 = "<<c2<<"; "<<"c3 = "<<c3<<G4endl;
393
394 return fA1[k]*log(x2/x1) + fA2[k]*c1 + fA3[k]*c2/2 + fA4[k]*c3/3 ;
395
396} // end of RutherfordIntegral
397
398
399/////////////////////////////////////////////////////////////////
400//
401// Imaginary part of dielectric constant
402// (G4int k - interval number, G4double en1 - energy point)
403
404G4double G4PAIySection::ImPartDielectricConst( G4int k ,
405 G4double energy1 )
406{
407 G4double energy2,energy3,energy4,result;
408
409 energy2 = energy1*energy1;
410 energy3 = energy2*energy1;
411 energy4 = energy3*energy1;
412
413 result = fA1[k]/energy1+fA2[k]/energy2+fA3[k]/energy3+fA4[k]/energy4 ;
414 result *=hbarc/energy1 ;
415
416 return result ;
417
418} // end of ImPartDielectricConst
419
420
421//////////////////////////////////////////////////////////////////////////////
422//
423// Real part of dielectric constant minus unit: epsilon_1 - 1
424// (G4double enb - energy point)
425//
426
427G4double G4PAIySection::RePartDielectricConst(G4double enb)
428{
429 G4double x0, x02, x03, x04, x05, x1, x2, xx1 ,xx2 , xx12,
430 c1, c2, c3, cof1, cof2, xln1, xln2, xln3, result ;
431
432 x0 = enb ;
433 result = 0 ;
434
435 for(G4int i=1;i<=fIntervalNumber-1;i++)
436 {
437 x1 = fEnergyInterval[i] ;
438 x2 = fEnergyInterval[i+1] ;
439 xx1 = x1 - x0 ;
440 xx2 = x2 - x0 ;
441 xx12 = xx2/xx1 ;
442
443 if(xx12<0)
444 {
445 xx12 = -xx12;
446 }
447 xln1 = log(x2/x1) ;
448 xln2 = log(xx12) ;
449 xln3 = log((x2 + x0)/(x1 + x0)) ;
450 x02 = x0*x0 ;
451 x03 = x02*x0 ;
452 x04 = x03*x0 ;
453 x05 = x04*x0;
454 c1 = (x2 - x1)/x1/x2 ;
455 c2 = (x2 - x1)*(x2 +x1)/x1/x1/x2/x2 ;
456 c3 = (x2 -x1)*(x1*x1 + x1*x2 + x2*x2)/x1/x1/x1/x2/x2/x2 ;
457
458 result -= (fA1[i]/x02 + fA3[i]/x04)*xln1 ;
459 result -= (fA2[i]/x02 + fA4[i]/x04)*c1 ;
460 result -= fA3[i]*c2/2/x02 ;
461 result -= fA4[i]*c3/3/x02 ;
462
463 cof1 = fA1[i]/x02 + fA3[i]/x04 ;
464 cof2 = fA2[i]/x03 + fA4[i]/x05 ;
465
466 result += 0.5*(cof1 +cof2)*xln2 ;
467 result += 0.5*(cof1 - cof2)*xln3 ;
468 }
469 result *= 2*hbarc/pi ;
470
471 return result ;
472
473} // end of RePartDielectricConst
474
475//////////////////////////////////////////////////////////////////////
476//
477// PAI differential cross-section in terms of
478// simplified Allison's equation
479//
480
481G4double G4PAIySection::DifPAIySection( G4int i ,
482 G4double betaGammaSq )
483{
484 G4double be2,cof,x1,x2,x3,x4,x5,x6,x7,x8,result ;
485 //G4double beta, be4 ;
486 G4double be4 ;
487 G4double betaBohr2 = fine_structure_const*fine_structure_const ;
488 G4double betaBohr4 = betaBohr2*betaBohr2*4.0 ;
489 be2 = betaGammaSq/(1 + betaGammaSq) ;
490 be4 = be2*be2 ;
491 // beta = sqrt(be2) ;
492 cof = 1 ;
493 x1 = log(2*electron_mass_c2/fSplineEnergy[i]) ;
494
495 if( betaGammaSq < 0.01 ) x2 = log(be2) ;
496 else
497 {
498 x2 = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
499 (1/betaGammaSq - fRePartDielectricConst[i]) +
500 fImPartDielectricConst[i]*fImPartDielectricConst[i] )/2 ;
501 }
502 if( fImPartDielectricConst[i] == 0.0 ||betaGammaSq < 0.01 )
503 {
504 x6=0 ;
505 }
506 else
507 {
508 x3 = -fRePartDielectricConst[i] + 1/betaGammaSq ;
509 x5 = -1 - fRePartDielectricConst[i] +
510 be2*((1 +fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
511 fImPartDielectricConst[i]*fImPartDielectricConst[i]) ;
512
513 x7 = atan2(fImPartDielectricConst[i],x3) ;
514 x6 = x5 * x7 ;
515 }
516 // if(fImPartDielectricConst[i] == 0) x6 = 0 ;
517
518 x4 = ((x1 + x2)*fImPartDielectricConst[i] + x6)/hbarc ;
519 // if( x4 < 0.0 ) x4 = 0.0 ;
520 x8 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
521 fImPartDielectricConst[i]*fImPartDielectricConst[i] ;
522
523 result = (x4 + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i]) ;
524 if(result < 1.0e-8) result = 1.0e-8 ;
525 result *= fine_structure_const/be2/pi ;
526 // result *= (1-exp(-beta/betaBohr))*(1-exp(-beta/betaBohr)) ;
527 // result *= (1-exp(-be2/betaBohr2)) ;
528 result *= (1-exp(-be4/betaBohr4)) ;
529 // if(fDensity >= 0.1)
530 if(x8 > 0.)
531 {
532 result /= x8 ;
533 }
534 return result ;
535
536} // end of DifPAIySection
537
538//////////////////////////////////////////////////////////////////////////
539//
540// Calculation od dN/dx of collisions with creation of Cerenkov pseudo-photons
541
542G4double G4PAIySection::PAIdNdxCerenkov( G4int i ,
543 G4double betaGammaSq )
544{
545 G4double cof, logarithm, x3, x5, argument, modul2, dNdxC ;
546 G4double be2, be4, betaBohr2,betaBohr4,cofBetaBohr ;
547
548 cof = 1.0 ;
549 cofBetaBohr = 4.0 ;
550 betaBohr2 = fine_structure_const*fine_structure_const ;
551 betaBohr4 = betaBohr2*betaBohr2*cofBetaBohr ;
552
553 be2 = betaGammaSq/(1 + betaGammaSq) ;
554 be4 = be2*be2 ;
555
556 if( betaGammaSq < 0.01 ) logarithm = log(1.0+betaGammaSq) ; // 0.0 ;
557 else
558 {
559 logarithm = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
560 (1/betaGammaSq - fRePartDielectricConst[i]) +
561 fImPartDielectricConst[i]*fImPartDielectricConst[i] )*0.5 ;
562 logarithm += log(1+1.0/betaGammaSq) ;
563 }
564
565 if( fImPartDielectricConst[i] == 0.0 || betaGammaSq < 0.01 )
566 {
567 argument = 0.0 ;
568 }
569 else
570 {
571 x3 = -fRePartDielectricConst[i] + 1.0/betaGammaSq ;
572 x5 = -1.0 - fRePartDielectricConst[i] +
573 be2*((1.0 +fRePartDielectricConst[i])*(1.0 + fRePartDielectricConst[i]) +
574 fImPartDielectricConst[i]*fImPartDielectricConst[i]) ;
575 if( x3 == 0.0 ) argument = 0.5*pi;
576 else argument = atan2(fImPartDielectricConst[i],x3) ;
577 argument *= x5 ;
578 }
579 dNdxC = ( logarithm*fImPartDielectricConst[i] + argument )/hbarc ;
580
581 if(dNdxC < 1.0e-8) dNdxC = 1.0e-8 ;
582
583 dNdxC *= fine_structure_const/be2/pi ;
584
585 dNdxC *= (1-exp(-be4/betaBohr4)) ;
586
587 // if(fDensity >= 0.1)
588 // {
589 modul2 = (1.0 + fRePartDielectricConst[i])*(1.0 + fRePartDielectricConst[i]) +
590 fImPartDielectricConst[i]*fImPartDielectricConst[i] ;
591 if(modul2 > 0.)
592 {
593 dNdxC /= modul2 ;
594 }
595 return dNdxC ;
596
597} // end of PAIdNdxCerenkov
598
599//////////////////////////////////////////////////////////////////////////
600//
601// Calculation od dN/dx of collisions with creation of longitudinal EM
602// excitations (plasmons, delta-electrons)
603
604G4double G4PAIySection::PAIdNdxPlasmon( G4int i ,
605 G4double betaGammaSq )
606{
607 G4double cof, resonance, modul2, dNdxP ;
608 G4double be2, be4, betaBohr2, betaBohr4, cofBetaBohr ;
609
610 cof = 1 ;
611 cofBetaBohr = 4.0 ;
612 betaBohr2 = fine_structure_const*fine_structure_const ;
613 betaBohr4 = betaBohr2*betaBohr2*cofBetaBohr ;
614
615 be2 = betaGammaSq/(1 + betaGammaSq) ;
616 be4 = be2*be2 ;
617
618 resonance = log(2*electron_mass_c2*be2/fSplineEnergy[i]) ;
619 resonance *= fImPartDielectricConst[i]/hbarc ;
620
621
622 dNdxP = ( resonance + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i] ) ;
623
624 if( dNdxP < 1.0e-8 ) dNdxP = 1.0e-8 ;
625
626 dNdxP *= fine_structure_const/be2/pi ;
627 dNdxP *= (1-exp(-be4/betaBohr4)) ;
628
629// if( fDensity >= 0.1 )
630// {
631 modul2 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
632 fImPartDielectricConst[i]*fImPartDielectricConst[i] ;
633 if(modul2 > 0.)
634 {
635 dNdxP /= modul2 ;
636 }
637 return dNdxP ;
638
639} // end of PAIdNdxPlasmon
640
641////////////////////////////////////////////////////////////////////////
642//
643// Calculation of the PAI integral cross-section
644// fIntegralPAIySection[1] = specific primary ionisation, 1/cm
645// and fIntegralPAIySection[0] = mean energy loss per cm in keV/cm
646
647void G4PAIySection::IntegralPAIySection()
648{
649 fIntegralPAIySection[fSplineNumber] = 0 ;
650 fIntegralPAIdEdx[fSplineNumber] = 0 ;
651 fIntegralPAIySection[0] = 0 ;
652 G4int k = fIntervalNumber -1 ;
653
654 for(G4int i = fSplineNumber-1 ; i >= 1 ; i--)
655 {
656 if(fSplineEnergy[i] >= fEnergyInterval[k])
657 {
658 fIntegralPAIySection[i] = fIntegralPAIySection[i+1] + SumOverInterval(i) ;
659 fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] + SumOverIntervaldEdx(i) ;
660 }
661 else
662 {
663 fIntegralPAIySection[i] = fIntegralPAIySection[i+1] +
664 SumOverBorder(i+1,fEnergyInterval[k]) ;
665 fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] +
666 SumOverBorderdEdx(i+1,fEnergyInterval[k]) ;
667 k-- ;
668 }
669 }
670} // end of IntegralPAIySection
671
672////////////////////////////////////////////////////////////////////////
673//
674// Calculation of the PAI Cerenkov integral cross-section
675// fIntegralCrenkov[1] = specific Crenkov ionisation, 1/cm
676// and fIntegralCerenkov[0] = mean Cerenkov loss per cm in keV/cm
677
678void G4PAIySection::IntegralCerenkov()
679{
680 G4int i, k ;
681 fIntegralCerenkov[fSplineNumber] = 0 ;
682 fIntegralCerenkov[0] = 0 ;
683 k = fIntervalNumber -1 ;
684
685 for( i = fSplineNumber-1 ; i >= 1 ; i-- )
686 {
687 if(fSplineEnergy[i] >= fEnergyInterval[k])
688 {
689 fIntegralCerenkov[i] = fIntegralCerenkov[i+1] + SumOverInterCerenkov(i) ;
690 // G4cout<<"int: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
691 }
692 else
693 {
694 fIntegralCerenkov[i] = fIntegralCerenkov[i+1] +
695 SumOverBordCerenkov(i+1,fEnergyInterval[k]) ;
696 k-- ;
697 // G4cout<<"bord: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
698 }
699 }
700
701} // end of IntegralCerenkov
702
703////////////////////////////////////////////////////////////////////////
704//
705// Calculation of the PAI Plasmon integral cross-section
706// fIntegralPlasmon[1] = splasmon primary ionisation, 1/cm
707// and fIntegralPlasmon[0] = mean plasmon loss per cm in keV/cm
708
709void G4PAIySection::IntegralPlasmon()
710{
711 fIntegralPlasmon[fSplineNumber] = 0 ;
712 fIntegralPlasmon[0] = 0 ;
713 G4int k = fIntervalNumber -1 ;
714 for(G4int i=fSplineNumber-1;i>=1;i--)
715 {
716 if(fSplineEnergy[i] >= fEnergyInterval[k])
717 {
718 fIntegralPlasmon[i] = fIntegralPlasmon[i+1] + SumOverInterPlasmon(i) ;
719 }
720 else
721 {
722 fIntegralPlasmon[i] = fIntegralPlasmon[i+1] +
723 SumOverBordPlasmon(i+1,fEnergyInterval[k]) ;
724 k-- ;
725 }
726 }
727
728} // end of IntegralPlasmon
729
730//////////////////////////////////////////////////////////////////////
731//
732// Calculation the PAI integral cross-section inside
733// of interval of continuous values of photo-ionisation
734// cross-section. Parameter 'i' is the number of interval.
735
736G4double G4PAIySection::SumOverInterval( G4int i )
737{
738 G4double x0,x1,y0,yy1,a,b,c,result ;
739
740 x0 = fSplineEnergy[i] ;
741 x1 = fSplineEnergy[i+1] ;
742 y0 = fDifPAIySection[i] ;
743 yy1 = fDifPAIySection[i+1];
744 c = x1/x0;
745 a = log10(yy1/y0)/log10(c) ;
746 // b = log10(y0) - a*log10(x0) ;
747 b = y0/pow(x0,a) ;
748 a += 1 ;
749 if(a == 0)
750 {
751 result = b*log(x1/x0) ;
752 }
753 else
754 {
755 result = y0*(x1*pow(c,a-1) - x0)/a ;
756 }
757 a++;
758 if(a == 0)
759 {
760 fIntegralPAIySection[0] += b*log(x1/x0) ;
761 }
762 else
763 {
764 fIntegralPAIySection[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a ;
765 }
766 return result ;
767
768} // end of SumOverInterval
769
770/////////////////////////////////
771
772G4double G4PAIySection::SumOverIntervaldEdx( G4int i )
773{
774 G4double x0,x1,y0,yy1,a,b,c,result ;
775
776 x0 = fSplineEnergy[i] ;
777 x1 = fSplineEnergy[i+1] ;
778 y0 = fDifPAIySection[i] ;
779 yy1 = fDifPAIySection[i+1];
780 c = x1/x0;
781 a = log10(yy1/y0)/log10(c) ;
782 // b = log10(y0) - a*log10(x0) ;
783 b = y0/pow(x0,a) ;
784 a += 2 ;
785 if(a == 0)
786 {
787 result = b*log(x1/x0) ;
788 }
789 else
790 {
791 result = y0*(x1*x1*pow(c,a-2) - x0*x0)/a ;
792 }
793 return result ;
794
795} // end of SumOverInterval
796
797//////////////////////////////////////////////////////////////////////
798//
799// Calculation the PAI Cerenkov integral cross-section inside
800// of interval of continuous values of photo-ionisation Cerenkov
801// cross-section. Parameter 'i' is the number of interval.
802
803G4double G4PAIySection::SumOverInterCerenkov( G4int i )
804{
805 G4double x0,x1,y0,yy1,a,c,result ;
806
807 x0 = fSplineEnergy[i] ;
808 x1 = fSplineEnergy[i+1] ;
809 y0 = fdNdxCerenkov[i] ;
810 yy1 = fdNdxCerenkov[i+1];
811 // G4cout<<"SumC, i = "<<i<<"; x0 ="<<x0<<"; x1 = "<<x1
812 // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
813
814 c = x1/x0;
815 a = log10(yy1/y0)/log10(c) ;
816 G4double b = 0.0;
817 if(a < 20.) b = y0/pow(x0,a) ;
818
819 a += 1.0 ;
820 if(a == 0) result = b*log(c) ;
821 else result = y0*(x1*pow(c,a-1) - x0)/a ;
822 a += 1.0 ;
823
824 if( a == 0 ) fIntegralCerenkov[0] += b*log(x1/x0) ;
825 else fIntegralCerenkov[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a ;
826 // G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
827 return result ;
828
829} // end of SumOverInterCerenkov
830
831//////////////////////////////////////////////////////////////////////
832//
833// Calculation the PAI Plasmon integral cross-section inside
834// of interval of continuous values of photo-ionisation Plasmon
835// cross-section. Parameter 'i' is the number of interval.
836
837G4double G4PAIySection::SumOverInterPlasmon( G4int i )
838{
839 G4double x0,x1,y0,yy1,a,c,result ;
840
841 x0 = fSplineEnergy[i] ;
842 x1 = fSplineEnergy[i+1] ;
843 y0 = fdNdxPlasmon[i] ;
844 yy1 = fdNdxPlasmon[i+1];
845 c =x1/x0;
846 a = log10(yy1/y0)/log10(c) ;
847
848 G4double b = 0.0;
849 if(a < 20.) b = y0/pow(x0,a) ;
850
851 a += 1.0 ;
852 if(a == 0) result = b*log(x1/x0) ;
853 else result = y0*(x1*pow(c,a-1) - x0)/a ;
854 a += 1.0 ;
855
856 if( a == 0 ) fIntegralPlasmon[0] += b*log(x1/x0) ;
857 else fIntegralPlasmon[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a ;
858
859 return result ;
860
861} // end of SumOverInterPlasmon
862
863///////////////////////////////////////////////////////////////////////////////
864//
865// Integration of PAI cross-section for the case of
866// passing across border between intervals
867
868G4double G4PAIySection::SumOverBorder( G4int i ,
869 G4double en0 )
870{
871 G4double x0,x1,y0,yy1,a,c,d,e0,result ;
872
873 e0 = en0 ;
874 x0 = fSplineEnergy[i] ;
875 x1 = fSplineEnergy[i+1] ;
876 y0 = fDifPAIySection[i] ;
877 yy1 = fDifPAIySection[i+1] ;
878
879 c = x1/x0;
880 d = e0/x0;
881 a = log10(yy1/y0)/log10(x1/x0) ;
882
883 G4double b = 0.0;
884 if(a < 20.) b = y0/pow(x0,a) ;
885
886 a += 1 ;
887 if(a == 0)
888 {
889 result = b*log(x0/e0) ;
890 }
891 else
892 {
893 result = y0*(x0 - e0*pow(d,a-1))/a ;
894 }
895 a++ ;
896 if(a == 0)
897 {
898 fIntegralPAIySection[0] += b*log(x0/e0) ;
899 }
900 else
901 {
902 fIntegralPAIySection[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a ;
903 }
904 x0 = fSplineEnergy[i - 1] ;
905 x1 = fSplineEnergy[i - 2] ;
906 y0 = fDifPAIySection[i - 1] ;
907 yy1 = fDifPAIySection[i - 2] ;
908
909 c = x1/x0;
910 d = e0/x0;
911 a = log10(yy1/y0)/log10(x1/x0) ;
912 // b0 = log10(y0) - a*log10(x0) ;
913 b = y0/pow(x0,a) ;
914 a += 1 ;
915 if(a == 0)
916 {
917 result += b*log(e0/x0) ;
918 }
919 else
920 {
921 result += y0*(e0*pow(d,a-1) - x0)/a ;
922 }
923 a++ ;
924 if(a == 0)
925 {
926 fIntegralPAIySection[0] += b*log(e0/x0) ;
927 }
928 else
929 {
930 fIntegralPAIySection[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a ;
931 }
932 return result ;
933
934}
935
936///////////////////////////////////////////////////////////////////////
937
938G4double G4PAIySection::SumOverBorderdEdx( G4int i ,
939 G4double en0 )
940{
941 G4double x0,x1,y0,yy1,a,c,d,e0,result ;
942
943 e0 = en0 ;
944 x0 = fSplineEnergy[i] ;
945 x1 = fSplineEnergy[i+1] ;
946 y0 = fDifPAIySection[i] ;
947 yy1 = fDifPAIySection[i+1] ;
948
949 c = x1/x0;
950 d = e0/x0;
951 a = log10(yy1/y0)/log10(x1/x0) ;
952
953 G4double b = 0.0;
954 if(a < 20.) b = y0/pow(x0,a) ;
955
956 a += 2 ;
957 if(a == 0)
958 {
959 result = b*log(x0/e0) ;
960 }
961 else
962 {
963 result = y0*(x0*x0 - e0*e0*pow(d,a-2))/a ;
964 }
965 x0 = fSplineEnergy[i - 1] ;
966 x1 = fSplineEnergy[i - 2] ;
967 y0 = fDifPAIySection[i - 1] ;
968 yy1 = fDifPAIySection[i - 2] ;
969
970 c = x1/x0;
971 d = e0/x0;
972 a = log10(yy1/y0)/log10(x1/x0) ;
973
974 if(a < 20.) b = y0/pow(x0,a) ;
975
976 a += 2 ;
977 if(a == 0)
978 {
979 result += b*log(e0/x0) ;
980 }
981 else
982 {
983 result += y0*(e0*e0*pow(d,a-2) - x0*x0)/a ;
984 }
985 return result ;
986
987}
988
989///////////////////////////////////////////////////////////////////////////////
990//
991// Integration of Cerenkov cross-section for the case of
992// passing across border between intervals
993
994G4double G4PAIySection::SumOverBordCerenkov( G4int i ,
995 G4double en0 )
996{
997 G4double x0,x1,y0,yy1,a,e0,c,d,result ;
998
999 e0 = en0 ;
1000 x0 = fSplineEnergy[i] ;
1001 x1 = fSplineEnergy[i+1] ;
1002 y0 = fdNdxCerenkov[i] ;
1003 yy1 = fdNdxCerenkov[i+1] ;
1004
1005 // G4cout<<G4endl;
1006 //G4cout<<"SumBordC, i = "<<i<<"; en0 = "<<en0<<"; x0 ="<<x0<<"; x1 = "<<x1
1007 // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
1008 c = x1/x0 ;
1009 d = e0/x0 ;
1010 a = log10(yy1/y0)/log10(c) ;
1011
1012 G4double b = 0.0;
1013 if(a < 20.) b = y0/pow(x0,a) ;
1014
1015 a += 1.0 ;
1016 if( a == 0 ) result = b*log(x0/e0) ;
1017 else result = y0*(x0 - e0*pow(d,a-1))/a ;
1018 a += 1.0 ;
1019
1020 if( a == 0 ) fIntegralCerenkov[0] += b*log(x0/e0) ;
1021 else fIntegralCerenkov[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a ;
1022
1023 //G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
1024
1025 x0 = fSplineEnergy[i - 1] ;
1026 x1 = fSplineEnergy[i - 2] ;
1027 y0 = fdNdxCerenkov[i - 1] ;
1028 yy1 = fdNdxCerenkov[i - 2] ;
1029
1030 //G4cout<<"x0 ="<<x0<<"; x1 = "<<x1
1031 // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
1032
1033 c = x1/x0 ;
1034 d = e0/x0 ;
1035 a = log10(yy1/y0)/log10(x1/x0) ;
1036
1037 // G4cout << "a= " << a << G4endl;
1038 if(a < 20.) b = y0/pow(x0,a) ;
1039
1040 if(a > 20.0) b = 0.0;
1041 else b = y0/pow(x0,a); // pow(10.,b0) ;
1042
1043 //G4cout << "b= " << b << G4endl;
1044
1045 a += 1.0 ;
1046 if( a == 0 ) result += b*log(e0/x0) ;
1047 else result += y0*(e0*pow(d,a-1) - x0 )/a ;
1048 a += 1.0 ;
1049 //G4cout << "result= " << result << G4endl;
1050
1051 if( a == 0 ) fIntegralCerenkov[0] += b*log(e0/x0) ;
1052 else fIntegralCerenkov[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a ;
1053
1054 //G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
1055
1056 return result ;
1057
1058}
1059
1060///////////////////////////////////////////////////////////////////////////////
1061//
1062// Integration of Plasmon cross-section for the case of
1063// passing across border between intervals
1064
1065G4double G4PAIySection::SumOverBordPlasmon( G4int i ,
1066 G4double en0 )
1067{
1068 G4double x0,x1,y0,yy1,a,c,d,e0,result ;
1069
1070 e0 = en0 ;
1071 x0 = fSplineEnergy[i] ;
1072 x1 = fSplineEnergy[i+1] ;
1073 y0 = fdNdxPlasmon[i] ;
1074 yy1 = fdNdxPlasmon[i+1] ;
1075
1076 c = x1/x0 ;
1077 d = e0/x0 ;
1078 a = log10(yy1/y0)/log10(c) ;
1079
1080 G4double b = 0.0;
1081 if(a < 20.) b = y0/pow(x0,a) ;
1082
1083 a += 1.0 ;
1084 if( a == 0 ) result = b*log(x0/e0) ;
1085 else result = y0*(x0 - e0*pow(d,a-1))/a ;
1086 a += 1.0 ;
1087
1088 if( a == 0 ) fIntegralPlasmon[0] += b*log(x0/e0) ;
1089 else fIntegralPlasmon[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a ;
1090
1091 x0 = fSplineEnergy[i - 1] ;
1092 x1 = fSplineEnergy[i - 2] ;
1093 y0 = fdNdxPlasmon[i - 1] ;
1094 yy1 = fdNdxPlasmon[i - 2] ;
1095
1096 c = x1/x0 ;
1097 d = e0/x0 ;
1098 a = log10(yy1/y0)/log10(c) ;
1099
1100 if(a < 20.) b = y0/pow(x0,a) ;
1101
1102 a += 1.0 ;
1103 if( a == 0 ) result += b*log(e0/x0) ;
1104 else result += y0*(e0*pow(d,a-1) - x0)/a ;
1105 a += 1.0 ;
1106
1107 if( a == 0 ) fIntegralPlasmon[0] += b*log(e0/x0) ;
1108 else fIntegralPlasmon[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a ;
1109
1110 return result ;
1111
1112}
1113
1114/////////////////////////////////////////////////////////////////////////
1115//
1116//
1117
1118G4double G4PAIySection::GetStepEnergyLoss( G4double step )
1119{
1120 G4int iTransfer ;
1121 G4long numOfCollisions ;
1122 G4double loss = 0.0 ;
1123 G4double meanNumber, position ;
1124
1125 // G4cout<<" G4PAIySection::GetStepEnergyLoss "<<G4endl ;
1126
1127
1128
1129 meanNumber = fIntegralPAIySection[1]*step ;
1130 numOfCollisions = G4Poisson(meanNumber) ;
1131
1132 // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl ;
1133
1134 while(numOfCollisions)
1135 {
1136 position = fIntegralPAIySection[1]*G4UniformRand() ;
1137
1138 for( iTransfer=1 ; iTransfer<=fSplineNumber ; iTransfer++ )
1139 {
1140 if( position >= fIntegralPAIySection[iTransfer] ) break ;
1141 }
1142 loss += fSplineEnergy[iTransfer] ;
1143 numOfCollisions-- ;
1144 }
1145 // G4cout<<"PAI energy loss = "<<loss/keV<<" keV"<<G4endl ;
1146
1147 return loss ;
1148}
1149
1150/////////////////////////////////////////////////////////////////////////
1151//
1152//
1153
1154G4double G4PAIySection::GetStepCerenkovLoss( G4double step )
1155{
1156 G4int iTransfer ;
1157 G4long numOfCollisions ;
1158 G4double loss = 0.0 ;
1159 G4double meanNumber, position ;
1160
1161 // G4cout<<" G4PAIySection::GetStepCreLosnkovs "<<G4endl ;
1162
1163
1164
1165 meanNumber = fIntegralCerenkov[1]*step ;
1166 numOfCollisions = G4Poisson(meanNumber) ;
1167
1168 // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl ;
1169
1170 while(numOfCollisions)
1171 {
1172 position = fIntegralCerenkov[1]*G4UniformRand() ;
1173
1174 for( iTransfer=1 ; iTransfer<=fSplineNumber ; iTransfer++ )
1175 {
1176 if( position >= fIntegralCerenkov[iTransfer] ) break ;
1177 }
1178 loss += fSplineEnergy[iTransfer] ;
1179 numOfCollisions-- ;
1180 }
1181 // G4cout<<"PAI Cerenkov loss = "<<loss/keV<<" keV"<<G4endl ;
1182
1183 return loss ;
1184}
1185
1186/////////////////////////////////////////////////////////////////////////
1187//
1188//
1189
1190G4double G4PAIySection::GetStepPlasmonLoss( G4double step )
1191{
1192 G4int iTransfer ;
1193 G4long numOfCollisions ;
1194 G4double loss = 0.0 ;
1195 G4double meanNumber, position ;
1196
1197 // G4cout<<" G4PAIySection::GetStepCreLosnkovs "<<G4endl ;
1198
1199
1200
1201 meanNumber = fIntegralPlasmon[1]*step ;
1202 numOfCollisions = G4Poisson(meanNumber) ;
1203
1204 // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl ;
1205
1206 while(numOfCollisions)
1207 {
1208 position = fIntegralPlasmon[1]*G4UniformRand() ;
1209
1210 for( iTransfer=1 ; iTransfer<=fSplineNumber ; iTransfer++ )
1211 {
1212 if( position >= fIntegralPlasmon[iTransfer] ) break ;
1213 }
1214 loss += fSplineEnergy[iTransfer] ;
1215 numOfCollisions-- ;
1216 }
1217 // G4cout<<"PAI Plasmon loss = "<<loss/keV<<" keV"<<G4endl ;
1218
1219 return loss ;
1220}
1221
1222
1223
1224/////////////////////////////////////////////////////////////////////////////
1225//
1226// Init array of Lorentz factors
1227//
1228
1229G4int G4PAIySection::fNumberOfGammas = 111 ;
1230
1231const G4double G4PAIySection::fLorentzFactor[112] = // fNumberOfGammas+1
1232{
12330.0,
12341.094989e+00, 1.107813e+00, 1.122369e+00, 1.138890e+00, 1.157642e+00,
12351.178925e+00, 1.203082e+00, 1.230500e+00, 1.261620e+00, 1.296942e+00, // 10
12361.337032e+00, 1.382535e+00, 1.434181e+00, 1.492800e+00, 1.559334e+00,
12371.634850e+00, 1.720562e+00, 1.817845e+00, 1.928263e+00, 2.053589e+00, // 20
12382.195835e+00, 2.357285e+00, 2.540533e+00, 2.748522e+00, 2.984591e+00,
12393.252533e+00, 3.556649e+00, 3.901824e+00, 4.293602e+00, 4.738274e+00, // 30
12405.242981e+00, 5.815829e+00, 6.466019e+00, 7.203990e+00, 8.041596e+00,
12418.992288e+00, 1.007133e+01, 1.129606e+01, 1.268614e+01, 1.426390e+01, // 40
12421.605467e+01, 1.808721e+01, 2.039417e+01, 2.301259e+01, 2.598453e+01,
12432.935771e+01, 3.318630e+01, 3.753180e+01, 4.246399e+01, 4.806208e+01, // 50
12445.441597e+01, 6.162770e+01, 6.981310e+01, 7.910361e+01, 8.964844e+01,
12451.016169e+02, 1.152013e+02, 1.306197e+02, 1.481198e+02, 1.679826e+02, // 60
12461.905270e+02, 2.161152e+02, 2.451581e+02, 2.781221e+02, 3.155365e+02,
12473.580024e+02, 4.062016e+02, 4.609081e+02, 5.230007e+02, 5.934765e+02, // 70
12486.734672e+02, 7.642575e+02, 8.673056e+02, 9.842662e+02, 1.117018e+03,
12491.267692e+03, 1.438709e+03, 1.632816e+03, 1.853128e+03, 2.103186e+03, // 80
12502.387004e+03, 2.709140e+03, 3.074768e+03, 3.489760e+03, 3.960780e+03,
12514.495394e+03, 5.102185e+03, 5.790900e+03, 6.572600e+03, 7.459837e+03, // 90
12528.466860e+03, 9.609843e+03, 1.090714e+04, 1.237959e+04, 1.405083e+04,
12531.594771e+04, 1.810069e+04, 2.054434e+04, 2.331792e+04, 2.646595e+04, // 100
12543.003901e+04, 3.409446e+04, 3.869745e+04, 4.392189e+04, 4.985168e+04,
12555.658206e+04, 6.422112e+04, 7.289153e+04, 8.273254e+04, 9.390219e+04, // 110
12561.065799e+05
1257} ;
1258
1259///////////////////////////////////////////////////////////////////////
1260//
1261// The number of gamma for creation of spline (near ion-min , G ~ 4 )
1262//
1263
1264const
1265G4int G4PAIySection::fRefGammaNumber = 29 ;
1266
1267
1268//
1269// end of G4PAIySection implementation file
1270//
1271////////////////////////////////////////////////////////////////////////////
1272
Note: See TracBrowser for help on using the repository browser.