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

Last change on this file since 1036 was 1007, checked in by garnier, 17 years ago

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