G4PAIySection Class Reference

#include <G4PAIySection.hh>


Public Member Functions

 G4PAIySection ()
 ~G4PAIySection ()
void Initialize (const G4Material *material, G4double maxEnergyTransfer, G4double betaGammaSq)
void ComputeLowEnergyCof (const G4Material *material)
void InitPAI ()
void NormShift (G4double betaGammaSq)
void SplainPAI (G4double betaGammaSq)
G4double RutherfordIntegral (G4int intervalNumber, G4double limitLow, G4double limitHigh)
G4double ImPartDielectricConst (G4int intervalNumber, G4double energy)
G4double RePartDielectricConst (G4double energy)
G4double DifPAIySection (G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxCerenkov (G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxPlasmon (G4int intervalNumber, G4double betaGammaSq)
void IntegralPAIySection ()
void IntegralCerenkov ()
void IntegralPlasmon ()
G4double SumOverInterval (G4int intervalNumber)
G4double SumOverIntervaldEdx (G4int intervalNumber)
G4double SumOverInterCerenkov (G4int intervalNumber)
G4double SumOverInterPlasmon (G4int intervalNumber)
G4double SumOverBorder (G4int intervalNumber, G4double energy)
G4double SumOverBorderdEdx (G4int intervalNumber, G4double energy)
G4double SumOverBordCerenkov (G4int intervalNumber, G4double energy)
G4double SumOverBordPlasmon (G4int intervalNumber, G4double energy)
G4double GetStepEnergyLoss (G4double step)
G4double GetStepCerenkovLoss (G4double step)
G4double GetStepPlasmonLoss (G4double step)
G4int GetNumberOfGammas () const
G4int GetSplineSize () const
G4int GetIntervalNumber () const
G4double GetEnergyInterval (G4int i)
G4double GetDifPAIySection (G4int i)
G4double GetPAIdNdxCrenkov (G4int i)
G4double GetPAIdNdxPlasmon (G4int i)
G4double GetMeanEnergyLoss () const
G4double GetMeanCerenkovLoss () const
G4double GetMeanPlasmonLoss () const
G4double GetNormalizationCof () const
G4double GetPAItable (G4int i, G4int j) const
G4double GetLorentzFactor (G4int i) const
G4double GetSplineEnergy (G4int i) const
G4double GetIntegralPAIySection (G4int i) const
G4double GetIntegralPAIdEdx (G4int i) const
G4double GetIntegralCerenkov (G4int i) const
G4double GetIntegralPlasmon (G4int i) const
void SetVerbose (G4int v)


Detailed Description

Definition at line 52 of file G4PAIySection.hh.


Constructor & Destructor Documentation

G4PAIySection::G4PAIySection (  ) 

Definition at line 73 of file G4PAIySection.cc.

00074 {
00075   fSandia = 0;
00076   fDensity = fElectronDensity = fNormalizationCof = fLowEnergyCof = 0.0;
00077   fIntervalNumber = fSplineNumber = 0;
00078   fVerbose = 0;
00079   for(G4int i=0; i<500; ++i) {
00080     fSplineEnergy[i] = 0.0;
00081     fRePartDielectricConst[i] = 0.0;
00082     fImPartDielectricConst[i] = 0.0;
00083     fIntegralTerm[i] = 0.0;
00084     fDifPAIySection[i] = 0.0;
00085     fdNdxCerenkov[i] = 0.0;
00086     fdNdxPlasmon[i] = 0.0;
00087     fIntegralPAIySection[i] = 0.0;
00088     fIntegralPAIdEdx[i] = 0.0;
00089     fIntegralCerenkov[i] = 0.0;
00090     fIntegralPlasmon[i] = 0.0;
00091     for(G4int j=0; j<112; ++j) { fPAItable[i][j] = 0.0; }
00092   }
00093 }

G4PAIySection::~G4PAIySection (  ) 

Definition at line 99 of file G4PAIySection.cc.

00100 {}


Member Function Documentation

void G4PAIySection::ComputeLowEnergyCof ( const G4Material material  ) 

Definition at line 215 of file G4PAIySection.cc.

References G4Material::GetElement(), G4Material::GetNumberOfElements(), and G4Element::GetZ().

Referenced by Initialize().

00216 {    
00217   G4int i, numberOfElements = material->GetNumberOfElements();
00218   G4double sumZ = 0., sumCof = 0.; 
00219 
00220   const G4double p0 =  1.20923e+00; 
00221   const G4double p1 =  3.53256e-01; 
00222   const G4double p2 = -1.45052e-03; 
00223   
00224   G4double* thisMaterialZ   = new G4double[numberOfElements];
00225   G4double* thisMaterialCof = new G4double[numberOfElements];
00226    
00227   for( i = 0; i < numberOfElements; i++ )
00228   {
00229     thisMaterialZ[i] = material->GetElement(i)->GetZ();
00230     sumZ += thisMaterialZ[i];
00231     thisMaterialCof[i] = p0+p1*thisMaterialZ[i]+p2*thisMaterialZ[i]*thisMaterialZ[i];   
00232   }
00233   for( i = 0; i < numberOfElements; i++ )
00234   {
00235     sumCof += thisMaterialCof[i]*thisMaterialZ[i]/sumZ;
00236   }
00237   fLowEnergyCof = sumCof;
00238   delete [] thisMaterialZ;
00239   delete [] thisMaterialCof;
00240   // G4cout<<"fLowEnergyCof = "<<fLowEnergyCof<<G4endl;
00241 }

G4double G4PAIySection::DifPAIySection ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 558 of file G4PAIySection.cc.

References G4INCL::Math::pi.

Referenced by Initialize(), InitPAI(), NormShift(), and SplainPAI().

00560 {        
00561   G4double beta, be2,cof,x1,x2,x3,x4,x5,x6,x7,x8,result;
00562    //G4double beta, be4;
00563    //G4double be4;
00564    G4double betaBohr = fine_structure_const;
00565    // G4double betaBohr2 = fine_structure_const*fine_structure_const;
00566    // G4double betaBohr4 = betaBohr2*betaBohr2*4.0;
00567    be2 = betaGammaSq/(1 + betaGammaSq);
00568    //be4 = be2*be2;
00569    beta = sqrt(be2);
00570    cof = 1;
00571    x1 = log(2*electron_mass_c2/fSplineEnergy[i]);
00572 
00573    if( betaGammaSq < 0.01 ) x2 = log(be2);
00574    else
00575    {
00576      x2 = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
00577                 (1/betaGammaSq - fRePartDielectricConst[i]) + 
00578                 fImPartDielectricConst[i]*fImPartDielectricConst[i] )/2;
00579    }
00580    if( fImPartDielectricConst[i] == 0.0 ||betaGammaSq < 0.01 )
00581    {
00582      x6=0;
00583    }
00584    else
00585    {
00586      x3 = -fRePartDielectricConst[i] + 1/betaGammaSq;
00587      x5 = -1 - fRePartDielectricConst[i] +
00588           be2*((1 +fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
00589           fImPartDielectricConst[i]*fImPartDielectricConst[i]);
00590 
00591      x7 = atan2(fImPartDielectricConst[i],x3);
00592      x6 = x5 * x7;
00593    }
00594     // if(fImPartDielectricConst[i] == 0) x6 = 0;
00595    
00596    x4 = ((x1 + x2)*fImPartDielectricConst[i] + x6)/hbarc;
00597    //   if( x4 < 0.0 ) x4 = 0.0;
00598    x8 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) + 
00599         fImPartDielectricConst[i]*fImPartDielectricConst[i];
00600 
00601    result = (x4 + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i]);
00602    if(result < 1.0e-8) result = 1.0e-8;
00603    result *= fine_structure_const/be2/pi;
00604    // low energy correction
00605 
00606    G4double lowCof = fLowEnergyCof; // 6.0 ; // Ar ~ 4.; -> fLowCof as f(Z1,Z2)? 
00607 
00608    result *= (1 - exp(-beta/betaBohr/lowCof));
00609 
00610    //   result *= (1-exp(-beta/betaBohr))*(1-exp(-beta/betaBohr));
00611    //  result *= (1-exp(-be2/betaBohr2));
00612    // result *= (1-exp(-be4/betaBohr4));
00613    //   if(fDensity >= 0.1)
00614    if(x8 > 0.)
00615    { 
00616       result /= x8;
00617    }
00618    return result;
00619 
00620 } // end of DifPAIySection 

G4double G4PAIySection::GetDifPAIySection ( G4int  i  )  [inline]

Definition at line 123 of file G4PAIySection.hh.

00123 { return fDifPAIySection[i]; } 

G4double G4PAIySection::GetEnergyInterval ( G4int  i  )  [inline]

Definition at line 121 of file G4PAIySection.hh.

00121 { return fEnergyInterval[i]; } 

G4double G4PAIySection::GetIntegralCerenkov ( G4int  i  )  const [inline]

Definition at line 230 of file G4PAIySection.hh.

00231 {
00232   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralCerenkov"); }
00233   return fIntegralCerenkov[i];
00234 }

G4double G4PAIySection::GetIntegralPAIdEdx ( G4int  i  )  const [inline]

Definition at line 224 of file G4PAIySection.hh.

Referenced by G4PAIModel::BuildPAIonisationTable().

00225 {
00226   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIdEdx"); }
00227   return fIntegralPAIdEdx[i];
00228 }

G4double G4PAIySection::GetIntegralPAIySection ( G4int  i  )  const [inline]

Definition at line 218 of file G4PAIySection.hh.

Referenced by G4PAIModel::BuildPAIonisationTable().

00219 {
00220   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIySection"); }
00221   return fIntegralPAIySection[i];
00222 }

G4double G4PAIySection::GetIntegralPlasmon ( G4int  i  )  const [inline]

Definition at line 236 of file G4PAIySection.hh.

00237 {
00238   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPlasmon"); }
00239   return fIntegralPlasmon[i];
00240 }

G4int G4PAIySection::GetIntervalNumber (  )  const [inline]

Definition at line 119 of file G4PAIySection.hh.

00119 { return fIntervalNumber; }

G4double G4PAIySection::GetLorentzFactor ( G4int  i  )  const [inline]

Definition at line 207 of file G4PAIySection.hh.

00208 {
00209    return fLorentzFactor[j];
00210 }

G4double G4PAIySection::GetMeanCerenkovLoss (  )  const [inline]

Definition at line 128 of file G4PAIySection.hh.

00128 {return fIntegralCerenkov[0]; }

G4double G4PAIySection::GetMeanEnergyLoss (  )  const [inline]

Definition at line 127 of file G4PAIySection.hh.

Referenced by G4PAIModel::BuildPAIonisationTable().

00127 {return fIntegralPAIySection[0]; }

G4double G4PAIySection::GetMeanPlasmonLoss (  )  const [inline]

Definition at line 129 of file G4PAIySection.hh.

00129 {return fIntegralPlasmon[0]; }

G4double G4PAIySection::GetNormalizationCof (  )  const [inline]

Definition at line 131 of file G4PAIySection.hh.

00131 { return fNormalizationCof; }

G4int G4PAIySection::GetNumberOfGammas (  )  const [inline]

Definition at line 115 of file G4PAIySection.hh.

00115 { return fNumberOfGammas; }

G4double G4PAIySection::GetPAIdNdxCrenkov ( G4int  i  )  [inline]

Definition at line 124 of file G4PAIySection.hh.

00124 { return fdNdxCerenkov[i]; } 

G4double G4PAIySection::GetPAIdNdxPlasmon ( G4int  i  )  [inline]

Definition at line 125 of file G4PAIySection.hh.

00125 { return fdNdxPlasmon[i]; } 

G4double G4PAIySection::GetPAItable ( G4int  i,
G4int  j 
) const [inline]

Definition at line 202 of file G4PAIySection.hh.

00203 {
00204    return fPAItable[i][j];
00205 }

G4double G4PAIySection::GetSplineEnergy ( G4int  i  )  const [inline]

Definition at line 212 of file G4PAIySection.hh.

Referenced by G4PAIModel::BuildPAIonisationTable().

00213 {
00214   if(i < 1 || i > fSplineNumber) { CallError(i, "GetSplineEnergy"); }
00215   return fSplineEnergy[i];
00216 }

G4int G4PAIySection::GetSplineSize (  )  const [inline]

Definition at line 117 of file G4PAIySection.hh.

Referenced by G4PAIModel::BuildPAIonisationTable().

00117 { return fSplineNumber; }

G4double G4PAIySection::GetStepCerenkovLoss ( G4double  step  ) 

Definition at line 1238 of file G4PAIySection.cc.

References G4Poisson(), G4UniformRand, and position.

01239 {  
01240   G4int iTransfer ;
01241   G4long numOfCollisions;
01242   G4double loss = 0.0;
01243   G4double meanNumber, position;
01244 
01245   // G4cout<<" G4PAIySection::GetStepCreLosnkovs "<<G4endl;
01246 
01247 
01248 
01249   meanNumber = fIntegralCerenkov[1]*step;
01250   numOfCollisions = G4Poisson(meanNumber);
01251 
01252   //   G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
01253 
01254   while(numOfCollisions)
01255   {
01256     position = fIntegralCerenkov[1]*G4UniformRand();
01257 
01258     for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
01259     {
01260         if( position >= fIntegralCerenkov[iTransfer] ) break;
01261     }
01262     loss += fSplineEnergy[iTransfer] ;
01263     numOfCollisions--;
01264   }
01265   // G4cout<<"PAI Cerenkov loss = "<<loss/keV<<" keV"<<G4endl; 
01266 
01267   return loss;
01268 }

G4double G4PAIySection::GetStepEnergyLoss ( G4double  step  ) 

Definition at line 1202 of file G4PAIySection.cc.

References G4Poisson(), G4UniformRand, and position.

01203 {  
01204   G4int iTransfer ;
01205   G4long numOfCollisions;
01206   G4double loss = 0.0;
01207   G4double meanNumber, position;
01208 
01209   // G4cout<<" G4PAIySection::GetStepEnergyLoss "<<G4endl;
01210 
01211 
01212 
01213   meanNumber = fIntegralPAIySection[1]*step;
01214   numOfCollisions = G4Poisson(meanNumber);
01215 
01216   //   G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
01217 
01218   while(numOfCollisions)
01219   {
01220     position = fIntegralPAIySection[1]*G4UniformRand();
01221 
01222     for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
01223     {
01224         if( position >= fIntegralPAIySection[iTransfer] ) break;
01225     }
01226     loss += fSplineEnergy[iTransfer] ;
01227     numOfCollisions--;
01228   }
01229   // G4cout<<"PAI energy loss = "<<loss/keV<<" keV"<<G4endl; 
01230 
01231   return loss;
01232 }

G4double G4PAIySection::GetStepPlasmonLoss ( G4double  step  ) 

Definition at line 1274 of file G4PAIySection.cc.

References G4Poisson(), G4UniformRand, and position.

01275 {  
01276   G4int iTransfer ;
01277   G4long numOfCollisions;
01278   G4double loss = 0.0;
01279   G4double meanNumber, position;
01280 
01281   // G4cout<<" G4PAIySection::GetStepCreLosnkovs "<<G4endl;
01282 
01283 
01284 
01285   meanNumber = fIntegralPlasmon[1]*step;
01286   numOfCollisions = G4Poisson(meanNumber);
01287 
01288   //   G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
01289 
01290   while(numOfCollisions)
01291   {
01292     position = fIntegralPlasmon[1]*G4UniformRand();
01293 
01294     for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
01295     {
01296         if( position >= fIntegralPlasmon[iTransfer] ) break;
01297     }
01298     loss += fSplineEnergy[iTransfer] ;
01299     numOfCollisions--;
01300   }
01301   // G4cout<<"PAI Plasmon loss = "<<loss/keV<<" keV"<<G4endl; 
01302 
01303   return loss;
01304 }

G4double G4PAIySection::ImPartDielectricConst ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 481 of file G4PAIySection.cc.

Referenced by NormShift(), and SplainPAI().

00483 {
00484    G4double energy2,energy3,energy4,result;
00485 
00486    energy2 = energy1*energy1;
00487    energy3 = energy2*energy1;
00488    energy4 = energy3*energy1;
00489    
00490    result = fA1[k]/energy1+fA2[k]/energy2+fA3[k]/energy3+fA4[k]/energy4;  
00491    result *=hbarc/energy1;
00492    
00493    return result;
00494 
00495 }  // end of ImPartDielectricConst 

void G4PAIySection::Initialize ( const G4Material material,
G4double  maxEnergyTransfer,
G4double  betaGammaSq 
)

Definition at line 106 of file G4PAIySection.cc.

References ComputeLowEnergyCof(), DifPAIySection(), G4cout, G4endl, G4Material::GetDensity(), G4Material::GetElectronDensity(), G4SandiaTable::GetMaxInterval(), G4Material::GetName(), G4SandiaTable::GetSandiaMatTablePAI(), G4Material::GetSandiaTable(), IntegralCerenkov(), IntegralPAIySection(), IntegralPlasmon(), NormShift(), PAIdNdxCerenkov(), PAIdNdxPlasmon(), and SplainPAI().

Referenced by G4PAIModel::BuildPAIonisationTable().

00109 {
00110   G4int i, j;
00111   G4double energy;   
00112   // fVerbose = 1;   
00113   fDensity         = material->GetDensity();
00114   fElectronDensity = material->GetElectronDensity();
00115   //G4int numberOfElements = material->GetNumberOfElements();
00116 
00117   fSandia = material->GetSandiaTable();
00118 
00119   fIntervalNumber = fSandia->GetMaxInterval();
00120 
00121   fIntervalNumber--;
00122 
00123   for( i = 1; i <= fIntervalNumber; i++ )
00124   {
00125     energy = fSandia->GetSandiaMatTablePAI(i-1,0); //vmg 20.11.10
00126 
00127     if( energy >= maxEnergyTransfer || i > fIntervalNumber )
00128     {
00129       fEnergyInterval[i] = maxEnergyTransfer;
00130       fIntervalNumber = i;
00131       break;
00132     }
00133     fEnergyInterval[i] = energy;
00134     fA1[i]             = fSandia->GetSandiaMatTablePAI(i-1,1);
00135     fA2[i]             = fSandia->GetSandiaMatTablePAI(i-1,2);
00136     fA3[i]             = fSandia->GetSandiaMatTablePAI(i-1,3);
00137     fA4[i]             = fSandia->GetSandiaMatTablePAI(i-1,4);
00138 
00139     if( fVerbose > 0 && std::fabs( betaGammaSq - 8. ) < 0.4 )
00140     {
00141         G4cout<<i<<"\t"<<fEnergyInterval[i]/keV<<" keV \t"<<fA1[i]<<"\t"<<fA2[i] <<"\t"<<fA3[i] <<"\t"<<fA4[i]<<G4endl;
00142     }
00143   } 
00144 
00145   
00146   if( fEnergyInterval[fIntervalNumber] != maxEnergyTransfer )
00147   {
00148     fIntervalNumber++;
00149     fEnergyInterval[fIntervalNumber] = maxEnergyTransfer;
00150     fA1[fIntervalNumber] = fA1[fIntervalNumber-1];
00151     fA2[fIntervalNumber] = fA2[fIntervalNumber-1];
00152     fA3[fIntervalNumber] = fA3[fIntervalNumber-1];
00153     fA4[fIntervalNumber] = fA4[fIntervalNumber-1];
00154   }
00155 
00156    // Now checking, if two borders are too close together
00157   for( i = 1; i < fIntervalNumber; i++ )
00158   {
00159         // G4cout<<fEnergyInterval[i]<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
00160         //   <<fA3[i]<<"\t"<<fA4[i]<<G4endl;
00161     if(fEnergyInterval[i+1]-fEnergyInterval[i] <
00162            1.5*fDelta*(fEnergyInterval[i+1]+fEnergyInterval[i]))
00163     {
00164       for( j = i; j < fIntervalNumber; j++ )
00165       {
00166         fEnergyInterval[j] = fEnergyInterval[j+1];
00167         fA1[j] = fA1[j+1];
00168         fA2[j] = fA2[j+1];
00169         fA3[j] = fA3[j+1];
00170         fA4[j] = fA4[j+1];
00171       }
00172       fIntervalNumber--;
00173       i--;
00174     }
00175   }
00176   if( fVerbose > 0 && std::fabs( betaGammaSq - 8. ) < 0.4 )
00177   {
00178     G4cout<<"Sandia cofs in G4PAIySection::Initialize(), mat = "<<material->GetName()<<G4endl;
00179     G4cout<<"for bg2 = "<<betaGammaSq<<", Tmax = "<< maxEnergyTransfer/keV<<" keV"<<G4endl;
00180     G4cout<<"energy \t"<<"a1 \t"<<"a2 \t"<<"a3 \t"<<"a4 \t"<<G4endl;
00181 
00182       for( j = 1; j < fIntervalNumber; j++ )
00183       {
00184         G4cout<<j<<"\t"<<fEnergyInterval[j]/keV<<" keV \t"<<fA1[j]<<"\t"<<fA2[j] <<"\t"<<fA3[j] <<"\t"<<fA4[j]<<G4endl;
00185       }
00186   }
00187 
00188       // Preparation of fSplineEnergy array corresponding to min ionisation, G~4
00189       
00190   G4double   betaGammaSqRef = 
00191      fLorentzFactor[fRefGammaNumber]*fLorentzFactor[fRefGammaNumber] - 1;
00192       
00193   ComputeLowEnergyCof(material);            
00194   NormShift(betaGammaSqRef);             
00195   SplainPAI(betaGammaSqRef);
00196       
00197    // Preparation of integral PAI cross section for input betaGammaSq
00198    
00199   for( i = 1; i <= fSplineNumber; i++ )
00200   {
00201     fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
00202     fdNdxCerenkov[i]   = PAIdNdxCerenkov(i,betaGammaSq);
00203     fdNdxPlasmon[i]    = PAIdNdxPlasmon(i,betaGammaSq);
00204   }
00205   IntegralPAIySection();
00206   IntegralCerenkov();
00207   IntegralPlasmon();
00208 }

void G4PAIySection::InitPAI (  ) 

Definition at line 251 of file G4PAIySection.cc.

References DifPAIySection(), IntegralCerenkov(), IntegralPAIySection(), IntegralPlasmon(), NormShift(), PAIdNdxCerenkov(), PAIdNdxPlasmon(), and SplainPAI().

00252 {    
00253    G4int i;
00254    G4double betaGammaSq = fLorentzFactor[fRefGammaNumber]*
00255                           fLorentzFactor[fRefGammaNumber] - 1;
00256 
00257    // Preparation of integral PAI cross section for reference gamma
00258    
00259    NormShift(betaGammaSq);             
00260    SplainPAI(betaGammaSq);
00261 
00262    IntegralPAIySection();
00263    IntegralCerenkov();
00264    IntegralPlasmon();
00265 
00266    for(i = 0; i<=fSplineNumber; i++)
00267    {
00268       fPAItable[i][fRefGammaNumber] = fIntegralPAIySection[i];
00269       if(i != 0) 
00270       {
00271          fPAItable[i][0] = fSplineEnergy[i];
00272       }
00273    }
00274    fPAItable[0][0] = fSplineNumber;
00275    
00276    for(G4int j = 1; j < 112; j++)       // for other gammas
00277    {
00278       if( j == fRefGammaNumber ) continue;
00279       
00280       betaGammaSq = fLorentzFactor[j]*fLorentzFactor[j] - 1;
00281       
00282       for(i = 1; i <= fSplineNumber; i++)
00283       {
00284          fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
00285          fdNdxCerenkov[i]   = PAIdNdxCerenkov(i,betaGammaSq);
00286          fdNdxPlasmon[i]    = PAIdNdxPlasmon(i,betaGammaSq);
00287       }
00288       IntegralPAIySection();
00289       IntegralCerenkov();
00290       IntegralPlasmon();
00291       
00292       for(i = 0; i <= fSplineNumber; i++)
00293       {
00294          fPAItable[i][j] = fIntegralPAIySection[i];
00295       }
00296    } 
00297 
00298 }  

void G4PAIySection::IntegralCerenkov (  ) 

Definition at line 762 of file G4PAIySection.cc.

References SumOverBordCerenkov(), and SumOverInterCerenkov().

Referenced by Initialize(), and InitPAI().

00763 {
00764   G4int i, k;
00765    fIntegralCerenkov[fSplineNumber] = 0;
00766    fIntegralCerenkov[0] = 0;
00767    k = fIntervalNumber -1;
00768 
00769    for( i = fSplineNumber-1; i >= 1; i-- )
00770    {
00771       if(fSplineEnergy[i] >= fEnergyInterval[k])
00772       {
00773         fIntegralCerenkov[i] = fIntegralCerenkov[i+1] + SumOverInterCerenkov(i);
00774         // G4cout<<"int: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
00775       }
00776       else
00777       {
00778         fIntegralCerenkov[i] = fIntegralCerenkov[i+1] + 
00779                                    SumOverBordCerenkov(i+1,fEnergyInterval[k]);
00780         k--;
00781         // G4cout<<"bord: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
00782       }
00783    }
00784 
00785 }   // end of IntegralCerenkov 

void G4PAIySection::IntegralPAIySection (  ) 

Definition at line 731 of file G4PAIySection.cc.

References SumOverBorder(), SumOverBorderdEdx(), SumOverInterval(), and SumOverIntervaldEdx().

Referenced by Initialize(), and InitPAI().

00732 {
00733   fIntegralPAIySection[fSplineNumber] = 0;
00734   fIntegralPAIdEdx[fSplineNumber]     = 0;
00735   fIntegralPAIySection[0]             = 0;
00736   G4int k = fIntervalNumber -1;
00737 
00738   for(G4int i = fSplineNumber-1; i >= 1; i--)
00739   {
00740     if(fSplineEnergy[i] >= fEnergyInterval[k])
00741     {
00742       fIntegralPAIySection[i] = fIntegralPAIySection[i+1] + SumOverInterval(i);
00743       fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] + SumOverIntervaldEdx(i);
00744     }
00745     else
00746     {
00747       fIntegralPAIySection[i] = fIntegralPAIySection[i+1] + 
00748                                    SumOverBorder(i+1,fEnergyInterval[k]);
00749       fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] + 
00750                                    SumOverBorderdEdx(i+1,fEnergyInterval[k]);
00751       k--;
00752     }
00753   }
00754 }   // end of IntegralPAIySection 

void G4PAIySection::IntegralPlasmon (  ) 

Definition at line 793 of file G4PAIySection.cc.

References SumOverBordPlasmon(), and SumOverInterPlasmon().

Referenced by Initialize(), and InitPAI().

00794 {
00795    fIntegralPlasmon[fSplineNumber] = 0;
00796    fIntegralPlasmon[0] = 0;
00797    G4int k = fIntervalNumber -1;
00798    for(G4int i=fSplineNumber-1;i>=1;i--)
00799    {
00800       if(fSplineEnergy[i] >= fEnergyInterval[k])
00801       {
00802         fIntegralPlasmon[i] = fIntegralPlasmon[i+1] + SumOverInterPlasmon(i);
00803       }
00804       else
00805       {
00806         fIntegralPlasmon[i] = fIntegralPlasmon[i+1] + 
00807                                    SumOverBordPlasmon(i+1,fEnergyInterval[k]);
00808         k--;
00809       }
00810    }
00811 
00812 }   // end of IntegralPlasmon

void G4PAIySection::NormShift ( G4double  betaGammaSq  ) 

Definition at line 305 of file G4PAIySection.cc.

References DifPAIySection(), ImPartDielectricConst(), PAIdNdxCerenkov(), PAIdNdxPlasmon(), G4INCL::Math::pi, RePartDielectricConst(), and RutherfordIntegral().

Referenced by Initialize(), and InitPAI().

00306 {
00307   G4int i, j;
00308 
00309   for( i = 1; i <= fIntervalNumber-1; i++ )
00310   {
00311     for( j = 1; j <= 2; j++ )
00312     {
00313       fSplineNumber = (i-1)*2 + j;
00314 
00315       if( j == 1 ) fSplineEnergy[fSplineNumber] = fEnergyInterval[i  ]*(1+fDelta);
00316       else         fSplineEnergy[fSplineNumber] = fEnergyInterval[i+1]*(1-fDelta); 
00317       //    G4cout<<"cn = "<<fSplineNumber<<"; "<<"energy = "
00318       //  <<fSplineEnergy[fSplineNumber]<<G4endl;
00319     }
00320   }
00321   fIntegralTerm[1]=RutherfordIntegral(1,fEnergyInterval[1],fSplineEnergy[1]);
00322 
00323   j = 1;
00324 
00325   for(i=2;i<=fSplineNumber;i++)
00326   {
00327     if(fSplineEnergy[i]<fEnergyInterval[j+1])
00328     {
00329          fIntegralTerm[i] = fIntegralTerm[i-1] + 
00330                             RutherfordIntegral(j,fSplineEnergy[i-1],
00331                                                  fSplineEnergy[i]   );
00332     }
00333     else
00334     {
00335        G4double x = RutherfordIntegral(j,fSplineEnergy[i-1],
00336                                            fEnergyInterval[j+1]   );
00337          j++;
00338          fIntegralTerm[i] = fIntegralTerm[i-1] + x + 
00339                             RutherfordIntegral(j,fEnergyInterval[j],
00340                                                  fSplineEnergy[i]    );
00341     }
00342     // G4cout<<i<<"\t"<<fSplineEnergy[i]<<"\t"<<fIntegralTerm[i]<<"\n"<<G4endl;
00343   } 
00344   fNormalizationCof = 2*pi*pi*hbarc*hbarc*fine_structure_const/electron_mass_c2;
00345   fNormalizationCof *= fElectronDensity/fIntegralTerm[fSplineNumber];
00346 
00347   // G4cout<<"fNormalizationCof = "<<fNormalizationCof<<G4endl;
00348 
00349           // Calculation of PAI differrential cross-section (1/(keV*cm))
00350           // in the energy points near borders of energy intervals
00351 
00352    for(G4int k=1;k<=fIntervalNumber-1;k++)
00353    {
00354       for(j=1;j<=2;j++)
00355       {
00356          i = (k-1)*2 + j;
00357          fImPartDielectricConst[i] = fNormalizationCof*
00358                                      ImPartDielectricConst(k,fSplineEnergy[i]);
00359          fRePartDielectricConst[i] = fNormalizationCof*
00360                                      RePartDielectricConst(fSplineEnergy[i]);
00361          fIntegralTerm[i] *= fNormalizationCof;
00362 
00363          fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
00364          fdNdxCerenkov[i]   = PAIdNdxCerenkov(i,betaGammaSq);
00365          fdNdxPlasmon[i]    = PAIdNdxPlasmon(i,betaGammaSq);
00366       }
00367    }
00368 
00369 }  // end of NormShift 

G4double G4PAIySection::PAIdNdxCerenkov ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 626 of file G4PAIySection.cc.

References G4INCL::Math::pi.

Referenced by Initialize(), InitPAI(), NormShift(), and SplainPAI().

00628 {        
00629    G4double logarithm, x3, x5, argument, modul2, dNdxC; 
00630    G4double be2, be4, betaBohr2,betaBohr4,cofBetaBohr;
00631 
00632    //G4double cof         = 1.0;
00633    cofBetaBohr = 4.0;
00634    betaBohr2   = fine_structure_const*fine_structure_const;
00635    betaBohr4   = betaBohr2*betaBohr2*cofBetaBohr;
00636 
00637    be2 = betaGammaSq/(1 + betaGammaSq);
00638    be4 = be2*be2;
00639 
00640    if( betaGammaSq < 0.01 ) logarithm = log(1.0+betaGammaSq); // 0.0;
00641    else
00642    {
00643      logarithm  = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
00644                         (1/betaGammaSq - fRePartDielectricConst[i]) + 
00645                         fImPartDielectricConst[i]*fImPartDielectricConst[i] )*0.5;
00646      logarithm += log(1+1.0/betaGammaSq);
00647    }
00648 
00649    if( fImPartDielectricConst[i] == 0.0 || betaGammaSq < 0.01 )
00650    {
00651      argument = 0.0;
00652    }
00653    else
00654    {
00655      x3 = -fRePartDielectricConst[i] + 1.0/betaGammaSq;
00656      x5 = -1.0 - fRePartDielectricConst[i] +
00657           be2*((1.0 +fRePartDielectricConst[i])*(1.0 + fRePartDielectricConst[i]) +
00658           fImPartDielectricConst[i]*fImPartDielectricConst[i]);
00659      if( x3 == 0.0 ) argument = 0.5*pi;
00660      else            argument = atan2(fImPartDielectricConst[i],x3);
00661      argument *= x5 ;
00662    }   
00663    dNdxC = ( logarithm*fImPartDielectricConst[i] + argument )/hbarc;
00664   
00665    if(dNdxC < 1.0e-8) dNdxC = 1.0e-8;
00666 
00667    dNdxC *= fine_structure_const/be2/pi;
00668 
00669    dNdxC *= (1-exp(-be4/betaBohr4));
00670 
00671    //   if(fDensity >= 0.1)
00672    // { 
00673    modul2 = (1.0 + fRePartDielectricConst[i])*(1.0 + fRePartDielectricConst[i]) + 
00674                     fImPartDielectricConst[i]*fImPartDielectricConst[i];
00675    if(modul2 > 0.)
00676      {
00677        dNdxC /= modul2;
00678      }
00679    return dNdxC;
00680 
00681 } // end of PAIdNdxCerenkov 

G4double G4PAIySection::PAIdNdxPlasmon ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 688 of file G4PAIySection.cc.

References G4INCL::Math::pi.

Referenced by Initialize(), InitPAI(), NormShift(), and SplainPAI().

00690 {        
00691    G4double cof, resonance, modul2, dNdxP;
00692    G4double be2, be4, betaBohr2, betaBohr4, cofBetaBohr;
00693 
00694    cof = 1;
00695    cofBetaBohr = 4.0;
00696    betaBohr2   = fine_structure_const*fine_structure_const;
00697    betaBohr4   = betaBohr2*betaBohr2*cofBetaBohr;
00698 
00699    be2 = betaGammaSq/(1 + betaGammaSq);
00700    be4 = be2*be2;
00701  
00702    resonance = log(2*electron_mass_c2*be2/fSplineEnergy[i]);  
00703    resonance *= fImPartDielectricConst[i]/hbarc;
00704 
00705 
00706    dNdxP = ( resonance + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i] );
00707 
00708    if( dNdxP < 1.0e-8 ) dNdxP = 1.0e-8;
00709 
00710    dNdxP *= fine_structure_const/be2/pi;
00711    dNdxP *= (1-exp(-be4/betaBohr4));
00712 
00713 //   if( fDensity >= 0.1 )
00714 //   { 
00715    modul2 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) + 
00716      fImPartDielectricConst[i]*fImPartDielectricConst[i];
00717    if(modul2 > 0.)
00718      { 
00719        dNdxP /= modul2;
00720      }
00721    return dNdxP;
00722 
00723 } // end of PAIdNdxPlasmon 

G4double G4PAIySection::RePartDielectricConst ( G4double  energy  ) 

Definition at line 504 of file G4PAIySection.cc.

References G4INCL::Math::pi.

Referenced by NormShift(), and SplainPAI().

00505 {       
00506    G4double x0, x02, x03, x04, x05, x1, x2, xx1 ,xx2 , xx12,
00507             c1, c2, c3, cof1, cof2, xln1, xln2, xln3, result;
00508 
00509    x0 = enb;
00510    result = 0;
00511    
00512    for(G4int i=1;i<=fIntervalNumber-1;i++)
00513    {
00514       x1 = fEnergyInterval[i];
00515       x2 = fEnergyInterval[i+1];
00516       xx1 = x1 - x0;
00517       xx2 = x2 - x0;
00518       xx12 = xx2/xx1;
00519       
00520       if(xx12<0)
00521       {
00522          xx12 = -xx12;
00523       }
00524       xln1 = log(x2/x1);
00525       xln2 = log(xx12);
00526       xln3 = log((x2 + x0)/(x1 + x0));
00527       x02 = x0*x0;
00528       x03 = x02*x0;
00529       x04 = x03*x0;
00530       x05 = x04*x0;
00531       c1  = (x2 - x1)/x1/x2;
00532       c2  = (x2 - x1)*(x2 +x1)/x1/x1/x2/x2;
00533       c3  = (x2 -x1)*(x1*x1 + x1*x2 + x2*x2)/x1/x1/x1/x2/x2/x2;
00534 
00535       result -= (fA1[i]/x02 + fA3[i]/x04)*xln1;
00536       result -= (fA2[i]/x02 + fA4[i]/x04)*c1;
00537       result -= fA3[i]*c2/2/x02;
00538       result -= fA4[i]*c3/3/x02;
00539 
00540       cof1 = fA1[i]/x02 + fA3[i]/x04;
00541       cof2 = fA2[i]/x03 + fA4[i]/x05;
00542 
00543       result += 0.5*(cof1 +cof2)*xln2;
00544       result += 0.5*(cof1 - cof2)*xln3;
00545    } 
00546    result *= 2*hbarc/pi;
00547    
00548    return result;
00549 
00550 }   // end of RePartDielectricConst 

G4double G4PAIySection::RutherfordIntegral ( G4int  intervalNumber,
G4double  limitLow,
G4double  limitHigh 
)

Definition at line 460 of file G4PAIySection.cc.

Referenced by NormShift(), and SplainPAI().

00463 {
00464    G4double  c1, c2, c3;
00465    // G4cout<<"RI: x1 = "<<x1<<"; "<<"x2 = "<<x2<<G4endl;   
00466    c1 = (x2 - x1)/x1/x2;
00467    c2 = (x2 - x1)*(x2 + x1)/x1/x1/x2/x2;
00468    c3 = (x2 - x1)*(x1*x1 + x1*x2 + x2*x2)/x1/x1/x1/x2/x2/x2;
00469    // G4cout<<" RI: c1 = "<<c1<<"; "<<"c2 = "<<c2<<"; "<<"c3 = "<<c3<<G4endl;   
00470    
00471    return  fA1[k]*log(x2/x1) + fA2[k]*c1 + fA3[k]*c2/2 + fA4[k]*c3/3;
00472 
00473 }   // end of RutherfordIntegral 

void G4PAIySection::SetVerbose ( G4int  v  )  [inline]

Definition at line 144 of file G4PAIySection.hh.

Referenced by G4PAIModel::Initialise().

00144 { fVerbose = v; };

void G4PAIySection::SplainPAI ( G4double  betaGammaSq  ) 

Definition at line 377 of file G4PAIySection.cc.

References DifPAIySection(), ImPartDielectricConst(), PAIdNdxCerenkov(), PAIdNdxPlasmon(), RePartDielectricConst(), and RutherfordIntegral().

Referenced by Initialize(), and InitPAI().

00378 {
00379    G4int k = 1;
00380    G4int i = 1;
00381 
00382    while ( (i < fSplineNumber) && (fSplineNumber < fMaxSplineSize-1) )
00383    {
00384       if(fSplineEnergy[i+1] > fEnergyInterval[k+1])
00385       {
00386           k++;   // Here next energy point is in next energy interval
00387           i++;
00388           continue;
00389       }
00390                        // Shifting of arrayes for inserting the geometrical 
00391                        // average of 'i' and 'i+1' energy points to 'i+1' place
00392       fSplineNumber++;
00393 
00394       for(G4int j = fSplineNumber; j >= i+2; j-- )
00395       {
00396          fSplineEnergy[j]          = fSplineEnergy[j-1];
00397          fImPartDielectricConst[j] = fImPartDielectricConst[j-1];
00398          fRePartDielectricConst[j] = fRePartDielectricConst[j-1];
00399          fIntegralTerm[j]          = fIntegralTerm[j-1];
00400 
00401          fDifPAIySection[j] = fDifPAIySection[j-1];
00402          fdNdxCerenkov[j]   = fdNdxCerenkov[j-1];
00403          fdNdxPlasmon[j]    = fdNdxPlasmon[j-1];
00404       }
00405       G4double x1  = fSplineEnergy[i];
00406       G4double x2  = fSplineEnergy[i+1];
00407       G4double yy1 = fDifPAIySection[i];
00408       G4double y2  = fDifPAIySection[i+1];
00409 
00410       G4double en1 = sqrt(x1*x2);
00411       fSplineEnergy[i+1] = en1;
00412 
00413                  // Calculation of logarithmic linear approximation
00414                  // in this (enr) energy point, which number is 'i+1' now
00415 
00416       G4double a = log10(y2/yy1)/log10(x2/x1);
00417       G4double b = log10(yy1) - a*log10(x1);
00418       G4double y = a*log10(en1) + b;
00419       y = pow(10.,y);
00420 
00421                  // Calculation of the PAI dif. cross-section at this point
00422 
00423       fImPartDielectricConst[i+1] = fNormalizationCof*
00424                                     ImPartDielectricConst(k,fSplineEnergy[i+1]);
00425       fRePartDielectricConst[i+1] = fNormalizationCof*
00426                                     RePartDielectricConst(fSplineEnergy[i+1]);
00427       fIntegralTerm[i+1] = fIntegralTerm[i] + fNormalizationCof*
00428                            RutherfordIntegral(k,fSplineEnergy[i],
00429                                                 fSplineEnergy[i+1]);
00430 
00431       fDifPAIySection[i+1] = DifPAIySection(i+1,betaGammaSq);
00432       fdNdxCerenkov[i+1]   = PAIdNdxCerenkov(i+1,betaGammaSq);
00433       fdNdxPlasmon[i+1]    = PAIdNdxPlasmon(i+1,betaGammaSq);
00434 
00435                   // Condition for next division of this segment or to pass
00436                   // to higher energies
00437 
00438       G4double x = 2*(fDifPAIySection[i+1] - y)/(fDifPAIySection[i+1] + y);
00439 
00440       if( x < 0 ) 
00441       {
00442          x = -x;
00443       }
00444       if( x > fError && fSplineNumber < fMaxSplineSize-1 )
00445       {
00446          continue;  // next division
00447       }
00448       i += 2;  // pass to next segment
00449 
00450    }   // close 'while'
00451 
00452 }  // end of SplainPAI 

G4double G4PAIySection::SumOverBordCerenkov ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1078 of file G4PAIySection.cc.

Referenced by IntegralCerenkov().

01080 {               
01081    G4double x0,x1,y0,yy1,a,e0,c,d,result;
01082 
01083    e0 = en0;
01084    x0 = fSplineEnergy[i];
01085    x1 = fSplineEnergy[i+1];
01086    y0 = fdNdxCerenkov[i];
01087    yy1 = fdNdxCerenkov[i+1];
01088 
01089    //  G4cout<<G4endl;
01090    //G4cout<<"SumBordC, i = "<<i<<"; en0 = "<<en0<<"; x0 ="<<x0<<"; x1 = "<<x1
01091    //     <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
01092    c = x1/x0;
01093    d = e0/x0;
01094    a = log10(yy1/y0)/log10(c);
01095  
01096    G4double b = 0.0;
01097    if(a < 20.) b = y0/pow(x0,a);
01098    
01099    a += 1.0;
01100    if( a == 0 ) result = b*log(x0/e0);
01101    else         result = y0*(x0 - e0*pow(d,a-1))/a;   
01102    a += 1.0;
01103 
01104    if( a == 0 ) fIntegralCerenkov[0] += b*log(x0/e0);
01105    else         fIntegralCerenkov[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01106 
01107    //G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
01108    
01109    x0  = fSplineEnergy[i - 1];
01110    x1  = fSplineEnergy[i - 2];
01111    y0  = fdNdxCerenkov[i - 1];
01112    yy1 = fdNdxCerenkov[i - 2];
01113 
01114    //G4cout<<"x0 ="<<x0<<"; x1 = "<<x1
01115    //    <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
01116 
01117    c = x1/x0;
01118    d = e0/x0;
01119    a  = log10(yy1/y0)/log10(x1/x0);
01120   
01121    //   G4cout << "a= " << a << G4endl;
01122    if(a < 20.) b = y0/pow(x0,a);
01123 
01124    if(a > 20.0) b = 0.0;
01125    else         b = y0/pow(x0,a);  // pow(10.,b0);
01126 
01127    //G4cout << "b= " << b << G4endl;
01128 
01129    a += 1.0;
01130    if( a == 0 ) result += b*log(e0/x0);
01131    else         result += y0*(e0*pow(d,a-1) - x0 )/a;
01132    a += 1.0;
01133    //G4cout << "result= " << result << G4endl;
01134 
01135    if( a == 0 )   fIntegralCerenkov[0] += b*log(e0/x0);
01136    else           fIntegralCerenkov[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01137 
01138    //G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;    
01139 
01140    return result;
01141 
01142 } 

G4double G4PAIySection::SumOverBorder ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 952 of file G4PAIySection.cc.

Referenced by IntegralPAIySection().

00954 {               
00955   G4double x0,x1,y0,yy1,a,/*c,*/d,e0,result;
00956 
00957    e0 = en0;
00958    x0 = fSplineEnergy[i];
00959    x1 = fSplineEnergy[i+1];
00960    y0 = fDifPAIySection[i];
00961    yy1 = fDifPAIySection[i+1];
00962 
00963    //c = x1/x0;
00964    d = e0/x0;   
00965    a = log10(yy1/y0)/log10(x1/x0);
00966 
00967    G4double b = 0.0;
00968    if(a < 20.) b = y0/pow(x0,a);
00969    
00970    a += 1;
00971    if(a == 0)
00972    {
00973       result = b*log(x0/e0);
00974    }
00975    else
00976    {
00977       result = y0*(x0 - e0*pow(d,a-1))/a;
00978    }
00979    a++;
00980    if(a == 0)
00981    {
00982       fIntegralPAIySection[0] += b*log(x0/e0);
00983    }
00984    else 
00985    {
00986       fIntegralPAIySection[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
00987    }
00988    x0 = fSplineEnergy[i - 1];
00989    x1 = fSplineEnergy[i - 2];
00990    y0 = fDifPAIySection[i - 1];
00991    yy1 = fDifPAIySection[i - 2];
00992 
00993    //c = x1/x0;
00994    d = e0/x0;   
00995    a = log10(yy1/y0)/log10(x1/x0);
00996    //  b0 = log10(y0) - a*log10(x0);
00997    b = y0/pow(x0,a);
00998    a += 1;
00999    if(a == 0)
01000    {
01001       result += b*log(e0/x0);
01002    }
01003    else
01004    {
01005       result += y0*(e0*pow(d,a-1) - x0)/a;
01006    }
01007    a++;
01008    if(a == 0) 
01009    {
01010       fIntegralPAIySection[0] += b*log(e0/x0);
01011    }
01012    else
01013    {
01014       fIntegralPAIySection[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01015    }
01016    return result;
01017 
01018 } 

G4double G4PAIySection::SumOverBorderdEdx ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1022 of file G4PAIySection.cc.

Referenced by IntegralPAIySection().

01024 {               
01025   G4double x0,x1,y0,yy1,a,/*c,*/d,e0,result;
01026 
01027    e0 = en0;
01028    x0 = fSplineEnergy[i];
01029    x1 = fSplineEnergy[i+1];
01030    y0 = fDifPAIySection[i];
01031    yy1 = fDifPAIySection[i+1];
01032 
01033    //c = x1/x0;
01034    d = e0/x0;   
01035    a = log10(yy1/y0)/log10(x1/x0);
01036    
01037    G4double b = 0.0;
01038    if(a < 20.) b = y0/pow(x0,a);
01039    
01040    a += 2;
01041    if(a == 0)
01042    {
01043       result = b*log(x0/e0);
01044    }
01045    else 
01046    {
01047       result = y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01048    }
01049    x0 = fSplineEnergy[i - 1];
01050    x1 = fSplineEnergy[i - 2];
01051    y0 = fDifPAIySection[i - 1];
01052    yy1 = fDifPAIySection[i - 2];
01053 
01054    //c = x1/x0;
01055    d = e0/x0;   
01056    a = log10(yy1/y0)/log10(x1/x0);
01057 
01058    if(a < 20.) b = y0/pow(x0,a);
01059 
01060    a += 2;
01061    if(a == 0) 
01062    {
01063       result += b*log(e0/x0);
01064    }
01065    else
01066    {
01067       result += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01068    }
01069    return result;
01070 
01071 } 

G4double G4PAIySection::SumOverBordPlasmon ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1149 of file G4PAIySection.cc.

Referenced by IntegralPlasmon().

01151 {               
01152    G4double x0,x1,y0,yy1,a,c,d,e0,result;
01153 
01154    e0 = en0;
01155    x0 = fSplineEnergy[i];
01156    x1 = fSplineEnergy[i+1];
01157    y0 = fdNdxPlasmon[i];
01158    yy1 = fdNdxPlasmon[i+1];
01159 
01160    c = x1/x0;
01161    d = e0/x0;   
01162    a = log10(yy1/y0)/log10(c);
01163 
01164    G4double b = 0.0;
01165    if(a < 20.) b = y0/pow(x0,a);
01166    
01167    a += 1.0;
01168    if( a == 0 ) result = b*log(x0/e0);
01169    else         result = y0*(x0 - e0*pow(d,a-1))/a;   
01170    a += 1.0;
01171 
01172    if( a == 0 ) fIntegralPlasmon[0] += b*log(x0/e0);
01173    else         fIntegralPlasmon[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01174    
01175    x0 = fSplineEnergy[i - 1];
01176    x1 = fSplineEnergy[i - 2];
01177    y0 = fdNdxPlasmon[i - 1];
01178    yy1 = fdNdxPlasmon[i - 2];
01179 
01180    c = x1/x0;
01181    d = e0/x0;
01182    a = log10(yy1/y0)/log10(c);
01183  
01184    if(a < 20.) b = y0/pow(x0,a);
01185 
01186    a += 1.0;
01187    if( a == 0 ) result += b*log(e0/x0);
01188    else         result += y0*(e0*pow(d,a-1) - x0)/a;
01189    a += 1.0;
01190 
01191    if( a == 0 )   fIntegralPlasmon[0] += b*log(e0/x0);
01192    else           fIntegralPlasmon[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01193    
01194    return result;
01195 
01196 } 

G4double G4PAIySection::SumOverInterCerenkov ( G4int  intervalNumber  ) 

Definition at line 887 of file G4PAIySection.cc.

Referenced by IntegralCerenkov().

00888 {         
00889    G4double x0,x1,y0,yy1,a,c,result;
00890 
00891    x0  = fSplineEnergy[i];
00892    x1  = fSplineEnergy[i+1];
00893    y0  = fdNdxCerenkov[i];
00894    yy1 = fdNdxCerenkov[i+1];
00895    // G4cout<<"SumC, i = "<<i<<"; x0 ="<<x0<<"; x1 = "<<x1
00896    //   <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
00897 
00898    c = x1/x0;
00899    a = log10(yy1/y0)/log10(c);
00900    G4double b = 0.0;
00901    if(a < 20.) b = y0/pow(x0,a);
00902 
00903    a += 1.0;
00904    if(a == 0) result = b*log(c);
00905    else       result = y0*(x1*pow(c,a-1) - x0)/a;   
00906    a += 1.0;
00907 
00908    if( a == 0 ) fIntegralCerenkov[0] += b*log(x1/x0);
00909    else         fIntegralCerenkov[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
00910    //  G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;   
00911    return result;
00912 
00913 } //  end of SumOverInterCerenkov

G4double G4PAIySection::SumOverInterPlasmon ( G4int  intervalNumber  ) 

Definition at line 921 of file G4PAIySection.cc.

Referenced by IntegralPlasmon().

00922 {         
00923   G4double x0,x1,y0,yy1,a,c,result;
00924 
00925    x0  = fSplineEnergy[i];
00926    x1  = fSplineEnergy[i+1];
00927    y0  = fdNdxPlasmon[i];
00928    yy1 = fdNdxPlasmon[i+1];
00929    c = x1/x0;
00930    a = log10(yy1/y0)/log10(c);
00931 
00932    G4double b = 0.0;
00933    if(a < 20.) b = y0/pow(x0,a);
00934 
00935    a += 1.0;
00936    if(a == 0) result = b*log(x1/x0);
00937    else       result = y0*(x1*pow(c,a-1) - x0)/a;   
00938    a += 1.0;
00939 
00940    if( a == 0 ) fIntegralPlasmon[0] += b*log(x1/x0);
00941    else         fIntegralPlasmon[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
00942    
00943    return result;
00944 
00945 } //  end of SumOverInterPlasmon

G4double G4PAIySection::SumOverInterval ( G4int  intervalNumber  ) 

Definition at line 820 of file G4PAIySection.cc.

Referenced by IntegralPAIySection().

00821 {         
00822    G4double x0,x1,y0,yy1,a,b,c,result;
00823 
00824    x0 = fSplineEnergy[i];
00825    x1 = fSplineEnergy[i+1];
00826    y0 = fDifPAIySection[i];
00827    yy1 = fDifPAIySection[i+1];
00828    c = x1/x0;
00829    a = log10(yy1/y0)/log10(c);
00830    // b = log10(y0) - a*log10(x0);
00831    b = y0/pow(x0,a);
00832    a += 1;
00833    if(a == 0) 
00834    {
00835       result = b*log(x1/x0);
00836    }
00837    else
00838    {
00839       result = y0*(x1*pow(c,a-1) - x0)/a;
00840    }
00841    a++;
00842    if(a == 0) 
00843    {
00844       fIntegralPAIySection[0] += b*log(x1/x0);
00845    }
00846    else
00847    {
00848       fIntegralPAIySection[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
00849    }
00850    return result;
00851 
00852 } //  end of SumOverInterval

G4double G4PAIySection::SumOverIntervaldEdx ( G4int  intervalNumber  ) 

Definition at line 856 of file G4PAIySection.cc.

Referenced by IntegralPAIySection().

00857 {         
00858    G4double x0,x1,y0,yy1,a,b,c,result;
00859 
00860    x0 = fSplineEnergy[i];
00861    x1 = fSplineEnergy[i+1];
00862    y0 = fDifPAIySection[i];
00863    yy1 = fDifPAIySection[i+1];
00864    c = x1/x0;
00865    a = log10(yy1/y0)/log10(c);
00866    // b = log10(y0) - a*log10(x0);
00867    b = y0/pow(x0,a);
00868    a += 2;
00869    if(a == 0) 
00870    {
00871      result = b*log(x1/x0);
00872    }
00873    else
00874    {
00875      result = y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
00876    }
00877    return result;
00878 
00879 } //  end of SumOverInterval


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:48 2013 for Geant4 by  doxygen 1.4.7