G4PAIxSection Class Reference

#include <G4PAIxSection.hh>


Public Member Functions

 G4PAIxSection (G4MaterialCutsCouple *matCC)
 G4PAIxSection (G4int materialIndex, G4double maxEnergyTransfer)
 G4PAIxSection (G4int materialIndex, G4double maxEnergyTransfer, G4double betaGammaSq, G4double **photoAbsCof, G4int intNumber)
 G4PAIxSection (G4int materialIndex, G4double maxEnergyTransfer, G4double betaGammaSq)
 ~G4PAIxSection ()
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 GetPhotonRange (G4double energy)
G4double GetElectronRange (G4double energy)
G4double RePartDielectricConst (G4double energy)
G4double DifPAIxSection (G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxCerenkov (G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxMM (G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxPlasmon (G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxResonance (G4int intervalNumber, G4double betaGammaSq)
void IntegralPAIxSection ()
void IntegralCerenkov ()
void IntegralMM ()
void IntegralPlasmon ()
void IntegralResonance ()
G4double SumOverInterval (G4int intervalNumber)
G4double SumOverIntervaldEdx (G4int intervalNumber)
G4double SumOverInterCerenkov (G4int intervalNumber)
G4double SumOverInterMM (G4int intervalNumber)
G4double SumOverInterPlasmon (G4int intervalNumber)
G4double SumOverInterResonance (G4int intervalNumber)
G4double SumOverBorder (G4int intervalNumber, G4double energy)
G4double SumOverBorderdEdx (G4int intervalNumber, G4double energy)
G4double SumOverBordCerenkov (G4int intervalNumber, G4double energy)
G4double SumOverBordMM (G4int intervalNumber, G4double energy)
G4double SumOverBordPlasmon (G4int intervalNumber, G4double energy)
G4double SumOverBordResonance (G4int intervalNumber, G4double energy)
G4double GetStepEnergyLoss (G4double step)
G4double GetStepCerenkovLoss (G4double step)
G4double GetStepMMLoss (G4double step)
G4double GetStepPlasmonLoss (G4double step)
G4double GetStepResonanceLoss (G4double step)
G4double GetEnergyTransfer ()
G4double GetCerenkovEnergyTransfer ()
G4double GetMMEnergyTransfer ()
G4double GetPlasmonEnergyTransfer ()
G4double GetResonanceEnergyTransfer ()
G4double GetRutherfordEnergyTransfer ()
G4int GetNumberOfGammas () const
G4int GetSplineSize () const
G4int GetIntervalNumber () const
G4double GetEnergyInterval (G4int i)
G4double GetDifPAIxSection (G4int i)
G4double GetPAIdNdxCerenkov (G4int i)
G4double GetPAIdNdxMM (G4int i)
G4double GetPAIdNdxPlasmon (G4int i)
G4double GetPAIdNdxResonance (G4int i)
G4double GetMeanEnergyLoss () const
G4double GetMeanCerenkovLoss () const
G4double GetMeanMMLoss () const
G4double GetMeanPlasmonLoss () const
G4double GetMeanResonanceLoss () const
G4double GetNormalizationCof () const
G4double GetPAItable (G4int i, G4int j) const
G4double GetLorentzFactor (G4int i) const
G4double GetSplineEnergy (G4int i) const
G4double GetIntegralPAIxSection (G4int i) const
G4double GetIntegralPAIdEdx (G4int i) const
G4double GetIntegralCerenkov (G4int i) const
G4double GetIntegralMM (G4int i) const
G4double GetIntegralPlasmon (G4int i) const
G4double GetIntegralResonance (G4int i) const


Detailed Description

Definition at line 68 of file G4PAIxSection.hh.


Constructor & Destructor Documentation

G4PAIxSection::G4PAIxSection ( G4MaterialCutsCouple matCC  ) 

Definition at line 95 of file G4PAIxSection.cc.

References G4Material::GetDensity(), G4Material::GetIndex(), G4MaterialCutsCouple::GetMaterial(), G4SandiaTable::GetMaxInterval(), and G4SandiaTable::GetSandiaMatTable().

00096 {
00097   fDensity       = matCC->GetMaterial()->GetDensity();
00098   G4int matIndex = matCC->GetMaterial()->GetIndex();
00099   fMaterialIndex = matIndex;   
00100   fSandia        = new G4SandiaTable(matIndex);
00101 
00102   G4int i, j; 
00103   fMatSandiaMatrix = new G4OrderedTable();
00104  
00105   for (i = 0; i < fSandia->GetMaxInterval()-1; i++)
00106   {
00107      fMatSandiaMatrix->push_back(new G4DataVector(5,0.));
00108   }                     
00109   for (i = 0; i < fSandia->GetMaxInterval()-1; i++)
00110   {
00111     (*(*fMatSandiaMatrix)[i])[0] = fSandia->GetSandiaMatTable(i,0);
00112 
00113     for(j = 1; j < 5; j++)
00114     {
00115       (*(*fMatSandiaMatrix)[i])[j] = fSandia->GetSandiaMatTable(i,j)*fDensity;
00116     }     
00117   }                     
00118   fEnergyInterval = fA1 = fA2 = fA3 = fA4 = 0;
00119 }

G4PAIxSection::G4PAIxSection ( G4int  materialIndex,
G4double  maxEnergyTransfer 
)

Definition at line 123 of file G4PAIxSection.cc.

References G4Material::GetMaterialTable(), and InitPAI().

00125 {
00126    fSandia = 0;
00127    fMatSandiaMatrix = 0;
00128    const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable();
00129    G4int i, j;   
00130 
00131       fMaterialIndex   = materialIndex;   
00132       fDensity                = (*theMaterialTable)[materialIndex]->GetDensity();
00133       fElectronDensity        = (*theMaterialTable)[materialIndex]->
00134                              GetElectronDensity();
00135       fIntervalNumber         = (*theMaterialTable)[materialIndex]->
00136                              GetSandiaTable()->GetMatNbOfIntervals();
00137       fIntervalNumber--;      
00138       // G4cout<<fDensity<<"\t"<<fElectronDensity<<"\t"<<fIntervalNumber<<G4endl;
00139 
00140       fEnergyInterval = new G4double[fIntervalNumber+2];
00141       fA1             = new G4double[fIntervalNumber+2];
00142       fA2             = new G4double[fIntervalNumber+2];
00143       fA3             = new G4double[fIntervalNumber+2];
00144       fA4             = new G4double[fIntervalNumber+2];
00145 
00146       for(i = 1; i <= fIntervalNumber; i++ )
00147       {
00148          if(((*theMaterialTable)[materialIndex]->
00149     GetSandiaTable()->GetSandiaCofForMaterial(i-1,0) >= maxEnergyTransfer) ||
00150               i > fIntervalNumber               )
00151          {
00152             fEnergyInterval[i] = maxEnergyTransfer;
00153             fIntervalNumber = i;
00154             break;
00155          }
00156          fEnergyInterval[i] = (*theMaterialTable)[materialIndex]->
00157                               GetSandiaTable()->GetSandiaCofForMaterial(i-1,0);
00158          fA1[i]             = (*theMaterialTable)[materialIndex]->
00159                               GetSandiaTable()->GetSandiaCofForMaterial(i-1,1);
00160          fA2[i]             = (*theMaterialTable)[materialIndex]->
00161                               GetSandiaTable()->GetSandiaCofForMaterial(i-1,2);
00162          fA3[i]             = (*theMaterialTable)[materialIndex]->
00163                               GetSandiaTable()->GetSandiaCofForMaterial(i-1,3);
00164          fA4[i]             = (*theMaterialTable)[materialIndex]->
00165                               GetSandiaTable()->GetSandiaCofForMaterial(i-1,4);
00166          // G4cout<<i<<"\t"<<fEnergyInterval[i]<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
00167          //                               <<fA3[i]<<"\t"<<fA4[i]<<"\t"<<G4endl;
00168       }   
00169       if(fEnergyInterval[fIntervalNumber] != maxEnergyTransfer)
00170       {
00171          fIntervalNumber++;
00172          fEnergyInterval[fIntervalNumber] = maxEnergyTransfer;
00173       }
00174 
00175       // Now checking, if two borders are too close together
00176 
00177       for(i=1;i<fIntervalNumber;i++)
00178       {
00179         if(fEnergyInterval[i+1]-fEnergyInterval[i] >
00180            1.5*fDelta*(fEnergyInterval[i+1]+fEnergyInterval[i]))
00181         {
00182           continue;
00183         }
00184         else
00185         {
00186           for(j=i;j<fIntervalNumber;j++)
00187           {
00188             fEnergyInterval[j] = fEnergyInterval[j+1];
00189                         fA1[j] = fA1[j+1];
00190                         fA2[j] = fA2[j+1];
00191                         fA3[j] = fA3[j+1];
00192                         fA4[j] = fA4[j+1];
00193           }
00194           fIntervalNumber--;
00195           i--;
00196         }
00197       }
00198 
00199 
00200       /* *********************************
00201 
00202       fSplineEnergy          = new G4double[fMaxSplineSize];   
00203       fRePartDielectricConst = new G4double[fMaxSplineSize];   
00204       fImPartDielectricConst = new G4double[fMaxSplineSize];   
00205       fIntegralTerm          = new G4double[fMaxSplineSize];   
00206       fDifPAIxSection        = new G4double[fMaxSplineSize];   
00207       fIntegralPAIxSection   = new G4double[fMaxSplineSize];   
00208       
00209       for(i=0;i<fMaxSplineSize;i++)
00210       {
00211          fSplineEnergy[i]          = 0.0;   
00212          fRePartDielectricConst[i] = 0.0;   
00213          fImPartDielectricConst[i] = 0.0;   
00214          fIntegralTerm[i]          = 0.0;   
00215          fDifPAIxSection[i]        = 0.0;   
00216          fIntegralPAIxSection[i]   = 0.0;   
00217       }
00218       **************************************************  */   
00219 
00220       InitPAI();  // create arrays allocated above
00221       
00222       delete[] fEnergyInterval;
00223       delete[] fA1;
00224       delete[] fA2;
00225       delete[] fA3;
00226       delete[] fA4;    
00227 }

G4PAIxSection::G4PAIxSection ( G4int  materialIndex,
G4double  maxEnergyTransfer,
G4double  betaGammaSq,
G4double **  photoAbsCof,
G4int  intNumber 
)

Definition at line 233 of file G4PAIxSection.cc.

References DifPAIxSection(), G4Material::GetMaterialTable(), IntegralCerenkov(), IntegralMM(), IntegralPAIxSection(), IntegralPlasmon(), IntegralResonance(), NormShift(), PAIdNdxCerenkov(), PAIdNdxMM(), PAIdNdxPlasmon(), PAIdNdxResonance(), and SplainPAI().

00238 {
00239    fSandia = 0;
00240    fMatSandiaMatrix = 0;
00241    const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable();
00242    G4int i, j; 
00243   
00244       fMaterialIndex   = materialIndex;      
00245       fDensity                = (*theMaterialTable)[materialIndex]->GetDensity();
00246       fElectronDensity        = (*theMaterialTable)[materialIndex]->
00247                              GetElectronDensity();
00248 
00249       fIntervalNumber         = intNumber;
00250       fIntervalNumber--;
00251       //   G4cout<<fDensity<<"\t"<<fElectronDensity<<"\t"<<fIntervalNumber<<G4endl;
00252   
00253       fEnergyInterval = new G4double[fIntervalNumber+2];
00254       fA1             = new G4double[fIntervalNumber+2];
00255       fA2             = new G4double[fIntervalNumber+2];
00256       fA3             = new G4double[fIntervalNumber+2];
00257       fA4             = new G4double[fIntervalNumber+2];
00258 
00259       for( i = 1; i <= fIntervalNumber; i++ )
00260       {
00261          if( ( photoAbsCof[i-1][0] >= maxEnergyTransfer ) ||
00262              i > fIntervalNumber )
00263          {
00264             fEnergyInterval[i] = maxEnergyTransfer;
00265             fIntervalNumber = i;
00266             break;
00267          }
00268          fEnergyInterval[i] = photoAbsCof[i-1][0];
00269          fA1[i]             = photoAbsCof[i-1][1];
00270          fA2[i]             = photoAbsCof[i-1][2];
00271          fA3[i]             = photoAbsCof[i-1][3];
00272          fA4[i]             = photoAbsCof[i-1][4];
00273          // G4cout<<i<<"\t"<<fEnergyInterval[i]<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
00274          //        <<fA3[i]<<"\t"<<fA4[i]<<"\t"<<G4endl;
00275       }
00276   // G4cout<<"i last = "<<i<<"; "<<"fIntervalNumber = "<<fIntervalNumber<<G4endl;   
00277       if(fEnergyInterval[fIntervalNumber] != maxEnergyTransfer)
00278       {
00279          fIntervalNumber++;
00280          fEnergyInterval[fIntervalNumber] = maxEnergyTransfer;
00281       }
00282       for(i=1;i<=fIntervalNumber;i++)
00283       {
00284         //  G4cout<<i<<"\t"<<fEnergyInterval[i]<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
00285         //    <<fA3[i]<<"\t"<<fA4[i]<<"\t"<<G4endl;
00286       }
00287       // Now checking, if two borders are too close together
00288 
00289       for( i = 1; i < fIntervalNumber; i++ )
00290       {
00291         if(fEnergyInterval[i+1]-fEnergyInterval[i] >
00292            1.5*fDelta*(fEnergyInterval[i+1]+fEnergyInterval[i]))
00293         {
00294           continue;
00295         }
00296         else
00297         {
00298           for(j=i;j<fIntervalNumber;j++)
00299           {
00300             fEnergyInterval[j] = fEnergyInterval[j+1];
00301                         fA1[j] = fA1[j+1];
00302                         fA2[j] = fA2[j+1];
00303                         fA3[j] = fA3[j+1];
00304                         fA4[j] = fA4[j+1];
00305           }
00306           fIntervalNumber--;
00307           i--;
00308         }
00309       }
00310 
00311       // Preparation of fSplineEnergy array corresponding to min ionisation, G~4
00312       
00313       G4double   betaGammaSqRef = 
00314       fLorentzFactor[fRefGammaNumber]*fLorentzFactor[fRefGammaNumber] - 1;
00315 
00316       NormShift(betaGammaSqRef);             
00317       SplainPAI(betaGammaSqRef);
00318       
00319       // Preparation of integral PAI cross section for input betaGammaSq
00320    
00321       for(i = 1; i <= fSplineNumber; i++)
00322       {
00323          fdNdxCerenkov[i]   = PAIdNdxCerenkov(i,betaGammaSq);
00324          fdNdxMM[i]   = PAIdNdxMM(i,betaGammaSq);
00325          fdNdxPlasmon[i]    = PAIdNdxPlasmon(i,betaGammaSq);
00326          fdNdxResonance[i]  = PAIdNdxResonance(i,betaGammaSq);
00327          fDifPAIxSection[i] = DifPAIxSection(i,betaGammaSq);
00328 
00329          // G4cout<<i<<"; dNdxC = "<<fdNdxCerenkov[i]<<"; dNdxP = "<<fdNdxPlasmon[i]
00330          //    <<"; dNdxPAI = "<<fDifPAIxSection[i]<<G4endl;
00331       }
00332       IntegralCerenkov();
00333       IntegralMM();
00334       IntegralPlasmon();
00335       IntegralResonance();
00336       IntegralPAIxSection();
00337       
00338       delete[] fEnergyInterval;
00339       delete[] fA1;
00340       delete[] fA2;
00341       delete[] fA3;
00342       delete[] fA4;    
00343 }

G4PAIxSection::G4PAIxSection ( G4int  materialIndex,
G4double  maxEnergyTransfer,
G4double  betaGammaSq 
)

Definition at line 349 of file G4PAIxSection.cc.

References DifPAIxSection(), G4Material::GetMaterialTable(), G4SandiaTable::GetPhotoAbsorpCof(), IntegralCerenkov(), IntegralMM(), IntegralPAIxSection(), IntegralPlasmon(), IntegralResonance(), NormShift(), PAIdNdxCerenkov(), PAIdNdxMM(), PAIdNdxPlasmon(), PAIdNdxResonance(), G4SandiaTable::SandiaIntervals(), G4SandiaTable::SandiaMixing(), and SplainPAI().

00352 {
00353    fSandia = 0;
00354    fMatSandiaMatrix = 0;
00355    const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable();
00356 
00357    G4int i, j, numberOfElements;   
00358 
00359    fMaterialIndex   = materialIndex;   
00360    fDensity         = (*theMaterialTable)[materialIndex]->GetDensity();
00361    fElectronDensity = (*theMaterialTable)[materialIndex]->GetElectronDensity();
00362    numberOfElements = (*theMaterialTable)[materialIndex]->GetNumberOfElements();
00363 
00364    G4int* thisMaterialZ = new G4int[numberOfElements];
00365    
00366    for( i = 0; i < numberOfElements; i++ )
00367    {
00368          thisMaterialZ[i] = (G4int)(*theMaterialTable)[materialIndex]->
00369                                       GetElement(i)->GetZ();
00370    }
00371    // fSandia = new G4SandiaTable(materialIndex);
00372    fSandia = (*theMaterialTable)[materialIndex]->
00373      GetSandiaTable();
00374    G4SandiaTable     thisMaterialSandiaTable(materialIndex);
00375    fIntervalNumber = thisMaterialSandiaTable.SandiaIntervals
00376                            (thisMaterialZ,numberOfElements);   
00377    fIntervalNumber = thisMaterialSandiaTable.SandiaMixing
00378                            ( thisMaterialZ ,
00379                       (*theMaterialTable)[materialIndex]->GetFractionVector() ,
00380                              numberOfElements,fIntervalNumber);
00381 
00382    fIntervalNumber--;
00383 
00384       fEnergyInterval = new G4double[fIntervalNumber+2];
00385       fA1             = new G4double[fIntervalNumber+2];
00386       fA2             = new G4double[fIntervalNumber+2];
00387       fA3             = new G4double[fIntervalNumber+2];
00388       fA4             = new G4double[fIntervalNumber+2];
00389 
00390       for( i = 1; i <= fIntervalNumber; i++ )
00391       {
00392   if((thisMaterialSandiaTable.GetPhotoAbsorpCof(i,0) >= maxEnergyTransfer) ||
00393           i > fIntervalNumber)
00394          {
00395             fEnergyInterval[i] = maxEnergyTransfer;
00396             fIntervalNumber = i;
00397             break;
00398          }
00399    fEnergyInterval[i] = thisMaterialSandiaTable.GetPhotoAbsorpCof(i,0);
00400    fA1[i]             = thisMaterialSandiaTable.GetPhotoAbsorpCof(i,1)*fDensity;
00401    fA2[i]             = thisMaterialSandiaTable.GetPhotoAbsorpCof(i,2)*fDensity;
00402    fA3[i]             = thisMaterialSandiaTable.GetPhotoAbsorpCof(i,3)*fDensity;
00403    fA4[i]             = thisMaterialSandiaTable.GetPhotoAbsorpCof(i,4)*fDensity;
00404 
00405       }   
00406       if(fEnergyInterval[fIntervalNumber] != maxEnergyTransfer)
00407       {
00408          fIntervalNumber++;
00409          fEnergyInterval[fIntervalNumber] = maxEnergyTransfer;
00410          fA1[fIntervalNumber] = fA1[fIntervalNumber-1];
00411          fA2[fIntervalNumber] = fA2[fIntervalNumber-1];
00412          fA3[fIntervalNumber] = fA3[fIntervalNumber-1];
00413          fA4[fIntervalNumber] = fA4[fIntervalNumber-1];
00414       }
00415       for(i=1;i<=fIntervalNumber;i++)
00416       {
00417         // G4cout<<fEnergyInterval[i]<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
00418         //   <<fA3[i]<<"\t"<<fA4[i]<<"\t"<<G4endl;
00419       }
00420       // Now checking, if two borders are too close together
00421 
00422       for( i = 1; i < fIntervalNumber; i++ )
00423       {
00424         if(fEnergyInterval[i+1]-fEnergyInterval[i] >
00425            1.5*fDelta*(fEnergyInterval[i+1]+fEnergyInterval[i]))
00426         {
00427           continue;
00428         }
00429         else
00430         {
00431           for( j = i; j < fIntervalNumber; j++ )
00432           {
00433             fEnergyInterval[j] = fEnergyInterval[j+1];
00434                         fA1[j] = fA1[j+1];
00435                         fA2[j] = fA2[j+1];
00436                         fA3[j] = fA3[j+1];
00437                         fA4[j] = fA4[j+1];
00438           }
00439           fIntervalNumber--;
00440           i--;
00441         }
00442       }
00443 
00444       /* *********************************
00445       fSplineEnergy          = new G4double[fMaxSplineSize];   
00446       fRePartDielectricConst = new G4double[fMaxSplineSize];   
00447       fImPartDielectricConst = new G4double[fMaxSplineSize];   
00448       fIntegralTerm          = new G4double[fMaxSplineSize];   
00449       fDifPAIxSection        = new G4double[fMaxSplineSize];   
00450       fIntegralPAIxSection   = new G4double[fMaxSplineSize];   
00451       
00452       for(i=0;i<fMaxSplineSize;i++)
00453       {
00454          fSplineEnergy[i]          = 0.0;   
00455          fRePartDielectricConst[i] = 0.0;   
00456          fImPartDielectricConst[i] = 0.0;   
00457          fIntegralTerm[i]          = 0.0;   
00458          fDifPAIxSection[i]        = 0.0;   
00459          fIntegralPAIxSection[i]   = 0.0;   
00460       }
00461       */ 
00462 
00463       // Preparation of fSplineEnergy array corresponding to min ionisation, G~4
00464       
00465       G4double   betaGammaSqRef = 
00466       fLorentzFactor[fRefGammaNumber]*fLorentzFactor[fRefGammaNumber] - 1;
00467 
00468       NormShift(betaGammaSqRef);             
00469       SplainPAI(betaGammaSqRef);
00470       
00471       // Preparation of integral PAI cross section for input betaGammaSq
00472    
00473       for(i = 1; i <= fSplineNumber; i++)
00474       {
00475          fDifPAIxSection[i] = DifPAIxSection(i,betaGammaSq);
00476          fdNdxCerenkov[i]   = PAIdNdxCerenkov(i,betaGammaSq);
00477          fdNdxMM[i]   = PAIdNdxMM(i,betaGammaSq);
00478          fdNdxPlasmon[i]    = PAIdNdxPlasmon(i,betaGammaSq);
00479          fdNdxResonance[i]  = PAIdNdxResonance(i,betaGammaSq);
00480       }
00481       IntegralPAIxSection();
00482       IntegralCerenkov();
00483       IntegralMM();
00484       IntegralPlasmon();
00485       IntegralResonance();
00486       
00487       delete[] fEnergyInterval;
00488       delete[] fA1;
00489       delete[] fA2;
00490       delete[] fA3;
00491       delete[] fA4;    
00492       delete[] thisMaterialZ;
00493 }

G4PAIxSection::~G4PAIxSection (  ) 

Definition at line 500 of file G4PAIxSection.cc.

00501 {
00502    /* ************************
00503    delete[] fSplineEnergy         ;   
00504    delete[] fRePartDielectricConst;   
00505    delete[] fImPartDielectricConst;   
00506    delete[] fIntegralTerm         ;   
00507    delete[] fDifPAIxSection       ;   
00508    delete[] fIntegralPAIxSection  ;
00509    */ 
00510   delete fSandia;
00511   delete fMatSandiaMatrix;
00512 }


Member Function Documentation

G4double G4PAIxSection::DifPAIxSection ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 910 of file G4PAIxSection.cc.

References G4INCL::Math::pi.

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

00912 {        
00913    G4double cof,x1,x2,x3,x4,x5,x6,x7,x8,result;
00914 
00915    G4double betaBohr  = fine_structure_const;
00916    // G4double betaBohr2 = fine_structure_const*fine_structure_const;
00917    // G4double betaBohr3 = betaBohr*betaBohr2; // *4.0;
00918 
00919    G4double be2  = betaGammaSq/(1 + betaGammaSq);
00920    G4double beta = sqrt(be2);
00921    // G4double be3 = beta*be2;
00922 
00923    cof = 1.;
00924    x1  = log(2*electron_mass_c2/fSplineEnergy[i]);
00925 
00926    if( betaGammaSq < 0.01 ) x2 = log(be2);
00927    else
00928    {
00929      x2 = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
00930                 (1/betaGammaSq - fRePartDielectricConst[i]) + 
00931                 fImPartDielectricConst[i]*fImPartDielectricConst[i] )/2;
00932    }
00933    if( fImPartDielectricConst[i] == 0.0 ||betaGammaSq < 0.01 )
00934    {
00935      x6 = 0.;
00936    }
00937    else
00938    {
00939      x3 = -fRePartDielectricConst[i] + 1/betaGammaSq;
00940      x5 = -1 - fRePartDielectricConst[i] +
00941           be2*((1 +fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
00942           fImPartDielectricConst[i]*fImPartDielectricConst[i]);
00943 
00944      x7 = atan2(fImPartDielectricConst[i],x3);
00945      x6 = x5 * x7;
00946    }
00947     // if(fImPartDielectricConst[i] == 0) x6 = 0.;
00948    
00949    x4 = ((x1 + x2)*fImPartDielectricConst[i] + x6)/hbarc;
00950 
00951    //   if( x4 < 0.0 ) x4 = 0.0;
00952 
00953    x8 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) + 
00954         fImPartDielectricConst[i]*fImPartDielectricConst[i];
00955 
00956    result = (x4 + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i]);
00957 
00958    if( result < 1.0e-8 ) result = 1.0e-8;
00959 
00960    result *= fine_structure_const/be2/pi;
00961 
00962    // low energy correction
00963 
00964    G4double lowCof = 4.0 ; // Ar ~ 4.; -> fLowCof as f(Z1,Z2)? 
00965 
00966    result *= (1 - exp(-beta/betaBohr/lowCof));
00967 
00968    // result *= (1 - exp(-be2/betaBohr2/lowCof));
00969 
00970    // result *= (1 - exp(-be3/betaBohr3/lowCof)); // ~ be for be<<betaBohr
00971 
00972    // result *= (1 - exp(-be4/betaBohr4/lowCof));
00973 
00974    if(fDensity >= 0.1)
00975    { 
00976       result /= x8;
00977    }
00978    return result;
00979 
00980 } // end of DifPAIxSection 

G4double G4PAIxSection::GetCerenkovEnergyTransfer (  ) 

Definition at line 1976 of file G4PAIxSection.cc.

References G4UniformRand, and position.

Referenced by GetStepCerenkovLoss().

01977 {  
01978   G4int iTransfer ;
01979 
01980   G4double energyTransfer, position;
01981 
01982   position = fIntegralCerenkov[1]*G4UniformRand();
01983 
01984   for( iTransfer = 1; iTransfer <= fSplineNumber; iTransfer++ )
01985   {
01986         if( position >= fIntegralCerenkov[iTransfer] ) break;
01987   }
01988   if(iTransfer > fSplineNumber) iTransfer--;
01989  
01990   energyTransfer = fSplineEnergy[iTransfer];
01991 
01992   if(iTransfer > 1)
01993   {
01994     energyTransfer -= (fSplineEnergy[iTransfer]-fSplineEnergy[iTransfer-1])*G4UniformRand();
01995   }
01996   return energyTransfer;
01997 }

G4double G4PAIxSection::GetDifPAIxSection ( G4int  i  )  [inline]

Definition at line 185 of file G4PAIxSection.hh.

00185 { return fDifPAIxSection[i]; } 

G4double G4PAIxSection::GetElectronRange ( G4double  energy  ) 

Definition at line 814 of file G4PAIxSection.cc.

00815 {
00816   G4double range;
00817   /*
00818   const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable();
00819 
00820   G4double Z = (*theMaterialTable)[fMaterialIndex]->GetIonisation()->GetZeffective();
00821   G4double A = (*theMaterialTable)[fMaterialIndex]->GetA();
00822 
00823   energy /= keV; // energy in keV in parametrised formula
00824 
00825   if (energy < 10.)
00826   {
00827     range = 3.872e-3*A/Z;
00828     range *= pow(energy,1.492);
00829   }
00830   else
00831   {
00832     range = 6.97e-3*pow(energy,1.6);
00833   }
00834   */
00835   // Blum&Rolandi Particle Detection with Drift Chambers, p. 7
00836 
00837   G4double cofA = 5.37e-4*g/cm2/keV;
00838   G4double cofB = 0.9815;
00839   G4double cofC = 3.123e-3/keV;
00840   // energy /= keV;
00841 
00842   range = cofA*energy*( 1 - cofB/(1 + cofC*energy) ); 
00843 
00844   // range *= g/cm2;
00845   range /= fDensity;
00846 
00847   return range;
00848 }

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

Definition at line 183 of file G4PAIxSection.hh.

00183 { return fEnergyInterval[i]; } 

G4double G4PAIxSection::GetEnergyTransfer (  ) 

Definition at line 1897 of file G4PAIxSection.cc.

References G4UniformRand, and position.

Referenced by GetStepEnergyLoss().

01898 {  
01899   G4int iTransfer ;
01900 
01901   G4double energyTransfer, position;
01902 
01903   position = fIntegralPAIxSection[1]*G4UniformRand();
01904 
01905   for( iTransfer = 1; iTransfer <= fSplineNumber; iTransfer++ )
01906   {
01907         if( position >= fIntegralPAIxSection[iTransfer] ) break;
01908   }
01909   if(iTransfer > fSplineNumber) iTransfer--;
01910  
01911   energyTransfer = fSplineEnergy[iTransfer];
01912 
01913   if(iTransfer > 1)
01914   {
01915     energyTransfer -= (fSplineEnergy[iTransfer]-fSplineEnergy[iTransfer-1])*G4UniformRand();
01916   }
01917   return energyTransfer;
01918 }

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

Definition at line 316 of file G4PAIxSection.hh.

Referenced by G4PAIPhotonModel::BuildPAIonisationTable().

00317 {
00318   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralCerenkov"); }
00319   return fIntegralCerenkov[i];
00320 }

G4double G4PAIxSection::GetIntegralMM ( G4int  i  )  const [inline]

Definition at line 322 of file G4PAIxSection.hh.

00323 {
00324   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralMM"); }
00325   return fIntegralMM[i];
00326 }

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

Definition at line 310 of file G4PAIxSection.hh.

Referenced by G4PAIPhotonModel::BuildPAIonisationTable().

00311 {
00312   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIdEdx"); }
00313   return fIntegralPAIdEdx[i];
00314 }

G4double G4PAIxSection::GetIntegralPAIxSection ( G4int  i  )  const [inline]

Definition at line 304 of file G4PAIxSection.hh.

Referenced by G4PAIPhotonModel::BuildPAIonisationTable().

00305 {
00306   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIxSection"); }
00307   return fIntegralPAIxSection[i];
00308 }

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

Definition at line 328 of file G4PAIxSection.hh.

Referenced by G4PAIPhotonModel::BuildPAIonisationTable().

00329 {
00330   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPlasmon"); }
00331   return fIntegralPlasmon[i];
00332 }

G4double G4PAIxSection::GetIntegralResonance ( G4int  i  )  const [inline]

Definition at line 334 of file G4PAIxSection.hh.

00335 {
00336   if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralResonance"); }
00337   return fIntegralResonance[i];
00338 }

G4int G4PAIxSection::GetIntervalNumber (  )  const [inline]

Definition at line 181 of file G4PAIxSection.hh.

00181 { return fIntervalNumber; }

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

Definition at line 293 of file G4PAIxSection.hh.

00294 {
00295    return fLorentzFactor[j];
00296 }

G4double G4PAIxSection::GetMeanCerenkovLoss (  )  const [inline]

Definition at line 192 of file G4PAIxSection.hh.

00192 {return fIntegralCerenkov[0]; }

G4double G4PAIxSection::GetMeanEnergyLoss (  )  const [inline]

Definition at line 191 of file G4PAIxSection.hh.

Referenced by G4PAIPhotonModel::BuildPAIonisationTable().

00191 {return fIntegralPAIxSection[0]; }

G4double G4PAIxSection::GetMeanMMLoss (  )  const [inline]

Definition at line 193 of file G4PAIxSection.hh.

00193 {return fIntegralMM[0]; }

G4double G4PAIxSection::GetMeanPlasmonLoss (  )  const [inline]

Definition at line 194 of file G4PAIxSection.hh.

00194 {return fIntegralPlasmon[0]; }

G4double G4PAIxSection::GetMeanResonanceLoss (  )  const [inline]

Definition at line 195 of file G4PAIxSection.hh.

00195 {return fIntegralResonance[0]; }

G4double G4PAIxSection::GetMMEnergyTransfer (  ) 

Definition at line 2003 of file G4PAIxSection.cc.

References G4UniformRand, and position.

Referenced by GetStepMMLoss().

02004 {  
02005   G4int iTransfer ;
02006 
02007   G4double energyTransfer, position;
02008 
02009   position = fIntegralMM[1]*G4UniformRand();
02010 
02011   for( iTransfer = 1; iTransfer <= fSplineNumber; iTransfer++ )
02012   {
02013         if( position >= fIntegralMM[iTransfer] ) break;
02014   }
02015   if(iTransfer > fSplineNumber) iTransfer--;
02016  
02017   energyTransfer = fSplineEnergy[iTransfer];
02018 
02019   if(iTransfer > 1)
02020   {
02021     energyTransfer -= (fSplineEnergy[iTransfer]-fSplineEnergy[iTransfer-1])*G4UniformRand();
02022   }
02023   return energyTransfer;
02024 }

G4double G4PAIxSection::GetNormalizationCof (  )  const [inline]

Definition at line 197 of file G4PAIxSection.hh.

00197 { return fNormalizationCof; }

G4int G4PAIxSection::GetNumberOfGammas (  )  const [inline]

Definition at line 177 of file G4PAIxSection.hh.

00177 { return fNumberOfGammas; }

G4double G4PAIxSection::GetPAIdNdxCerenkov ( G4int  i  )  [inline]

Definition at line 186 of file G4PAIxSection.hh.

00186 { return fdNdxCerenkov[i]; } 

G4double G4PAIxSection::GetPAIdNdxMM ( G4int  i  )  [inline]

Definition at line 187 of file G4PAIxSection.hh.

00187 { return fdNdxMM[i]; } 

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

Definition at line 188 of file G4PAIxSection.hh.

00188 { return fdNdxPlasmon[i]; } 

G4double G4PAIxSection::GetPAIdNdxResonance ( G4int  i  )  [inline]

Definition at line 189 of file G4PAIxSection.hh.

00189 { return fdNdxResonance[i]; } 

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

Definition at line 288 of file G4PAIxSection.hh.

00289 {
00290    return fPAItable[i][j];
00291 }

G4double G4PAIxSection::GetPhotonRange ( G4double  energy  ) 

Definition at line 781 of file G4PAIxSection.cc.

References DBL_MAX, DBL_MIN, and G4InuclParticleNames::lambda.

00782 {
00783   G4int i;
00784   G4double energy2, energy3, energy4, result, lambda;
00785 
00786   energy2 = energy1*energy1;
00787   energy3 = energy2*energy1;
00788   energy4 = energy3*energy1;
00789 
00790   // G4double* SandiaCof = fSandia->GetSandiaCofForMaterialPAI(energy1);
00791   // result = SandiaCof[0]/energy1+SandiaCof[1]/energy2+SandiaCof[2]/energy3+SandiaCof[3]/energy4;
00792   // result *= fDensity;
00793 
00794   for( i = 1; i <= fIntervalNumber; i++ )
00795   {
00796      if( energy1 < fEnergyInterval[i]) break;
00797   }
00798   i--;
00799   if(i == 0) i = 1;
00800 
00801   result = fA1[i]/energy1+fA2[i]/energy2+fA3[i]/energy3+fA4[i]/energy4;  
00802 
00803   if( result > DBL_MIN ) lambda = 1./result;
00804   else                   lambda = DBL_MAX;
00805    
00806   return lambda;
00807 }  

G4double G4PAIxSection::GetPlasmonEnergyTransfer (  ) 

Definition at line 2056 of file G4PAIxSection.cc.

References G4UniformRand, and position.

Referenced by GetStepPlasmonLoss().

02057 {  
02058   G4int iTransfer ;
02059 
02060   G4double energyTransfer, position;
02061 
02062   position = fIntegralPlasmon[1]*G4UniformRand();
02063 
02064   for( iTransfer = 1; iTransfer <= fSplineNumber; iTransfer++ )
02065   {
02066         if( position >= fIntegralPlasmon[iTransfer] ) break;
02067   }
02068   if(iTransfer > fSplineNumber) iTransfer--;
02069  
02070   energyTransfer = fSplineEnergy[iTransfer];
02071 
02072   if(iTransfer > 1)
02073   {
02074     energyTransfer -= (fSplineEnergy[iTransfer]-fSplineEnergy[iTransfer-1])*G4UniformRand();
02075   }
02076   return energyTransfer;
02077 }

G4double G4PAIxSection::GetResonanceEnergyTransfer (  ) 

Definition at line 2110 of file G4PAIxSection.cc.

References G4UniformRand, and position.

Referenced by GetStepResonanceLoss().

02111 {  
02112   G4int iTransfer ;
02113 
02114   G4double energyTransfer, position;
02115 
02116   position = fIntegralResonance[1]*G4UniformRand();
02117 
02118   for( iTransfer = 1; iTransfer <= fSplineNumber; iTransfer++ )
02119   {
02120         if( position >= fIntegralResonance[iTransfer] ) break;
02121   }
02122   if(iTransfer > fSplineNumber) iTransfer--;
02123  
02124   energyTransfer = fSplineEnergy[iTransfer];
02125 
02126   if(iTransfer > 1)
02127   {
02128     energyTransfer -= (fSplineEnergy[iTransfer]-fSplineEnergy[iTransfer-1])*G4UniformRand();
02129   }
02130   return energyTransfer;
02131 }

G4double G4PAIxSection::GetRutherfordEnergyTransfer (  ) 

Definition at line 2138 of file G4PAIxSection.cc.

References G4UniformRand, and position.

02139 {  
02140   G4int iTransfer ;
02141 
02142   G4double energyTransfer, position;
02143 
02144   position = (fIntegralPlasmon[1]-fIntegralResonance[1])*G4UniformRand();
02145 
02146   for( iTransfer = 1; iTransfer <= fSplineNumber; iTransfer++ )
02147   {
02148         if( position >= (fIntegralPlasmon[iTransfer]-fIntegralResonance[iTransfer]) ) break;
02149   }
02150   if(iTransfer > fSplineNumber) iTransfer--;
02151  
02152   energyTransfer = fSplineEnergy[iTransfer];
02153 
02154   if(iTransfer > 1)
02155   {
02156     energyTransfer -= (fSplineEnergy[iTransfer]-fSplineEnergy[iTransfer-1])*G4UniformRand();
02157   }
02158   return energyTransfer;
02159 }

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

Definition at line 298 of file G4PAIxSection.hh.

Referenced by G4PAIPhotonModel::BuildPAIonisationTable().

00299 {
00300   if(i < 1 || i > fSplineNumber) { CallError(i, "GetSplineEnergy"); }
00301   return fSplineEnergy[i];
00302 }

G4int G4PAIxSection::GetSplineSize (  )  const [inline]

Definition at line 179 of file G4PAIxSection.hh.

Referenced by G4PAIPhotonModel::BuildPAIonisationTable().

00179 { return fSplineNumber; }

G4double G4PAIxSection::GetStepCerenkovLoss ( G4double  step  ) 

Definition at line 1924 of file G4PAIxSection.cc.

References G4Poisson(), and GetCerenkovEnergyTransfer().

01925 {  
01926   G4long numOfCollisions;
01927   G4double meanNumber, loss = 0.0;
01928 
01929   // G4cout<<" G4PAIxSection::GetStepCerenkovLoss "<<G4endl;
01930 
01931   meanNumber = fIntegralCerenkov[1]*step;
01932   numOfCollisions = G4Poisson(meanNumber);
01933 
01934   //   G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
01935 
01936   while(numOfCollisions)
01937   {
01938     loss += GetCerenkovEnergyTransfer();
01939     numOfCollisions--;
01940   }
01941   // G4cout<<"PAI Cerenkov loss = "<<loss/keV<<" keV"<<G4endl; 
01942 
01943   return loss;
01944 }

G4double G4PAIxSection::GetStepEnergyLoss ( G4double  step  ) 

Definition at line 1871 of file G4PAIxSection.cc.

References G4Poisson(), and GetEnergyTransfer().

01872 {  
01873   G4long numOfCollisions;
01874   G4double meanNumber, loss = 0.0;
01875 
01876   // G4cout<<" G4PAIxSection::GetStepEnergyLoss "<<G4endl;
01877 
01878   meanNumber = fIntegralPAIxSection[1]*step;
01879   numOfCollisions = G4Poisson(meanNumber);
01880 
01881   //   G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
01882 
01883   while(numOfCollisions)
01884   {
01885    loss += GetEnergyTransfer();
01886    numOfCollisions--;
01887   }
01888   // G4cout<<"PAI energy loss = "<<loss/keV<<" keV"<<G4endl; 
01889 
01890   return loss;
01891 }

G4double G4PAIxSection::GetStepMMLoss ( G4double  step  ) 

Definition at line 1950 of file G4PAIxSection.cc.

References G4Poisson(), and GetMMEnergyTransfer().

01951 {  
01952   G4long numOfCollisions;
01953   G4double meanNumber, loss = 0.0;
01954 
01955   // G4cout<<" G4PAIxSection::GetStepMMLoss "<<G4endl;
01956 
01957   meanNumber = fIntegralMM[1]*step;
01958   numOfCollisions = G4Poisson(meanNumber);
01959 
01960   //   G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
01961 
01962   while(numOfCollisions)
01963   {
01964     loss += GetMMEnergyTransfer();
01965     numOfCollisions--;
01966   }
01967   // G4cout<<"PAI MM-Cerenkov loss = "<<loss/keV<<" keV"<<G4endl; 
01968 
01969   return loss;
01970 }

G4double G4PAIxSection::GetStepPlasmonLoss ( G4double  step  ) 

Definition at line 2030 of file G4PAIxSection.cc.

References G4Poisson(), and GetPlasmonEnergyTransfer().

02031 {  
02032   G4long numOfCollisions;
02033   G4double  meanNumber, loss = 0.0;
02034 
02035   // G4cout<<" G4PAIxSection::GetStepPlasmonLoss "<<G4endl;
02036 
02037   meanNumber = fIntegralPlasmon[1]*step;
02038   numOfCollisions = G4Poisson(meanNumber);
02039 
02040   //   G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
02041 
02042   while(numOfCollisions)
02043   {
02044     loss += GetPlasmonEnergyTransfer();
02045     numOfCollisions--;
02046   }
02047   // G4cout<<"PAI Plasmon loss = "<<loss/keV<<" keV"<<G4endl; 
02048 
02049   return loss;
02050 }

G4double G4PAIxSection::GetStepResonanceLoss ( G4double  step  ) 

Definition at line 2083 of file G4PAIxSection.cc.

References G4Poisson(), and GetResonanceEnergyTransfer().

02084 {  
02085   G4long numOfCollisions;
02086   G4double meanNumber, loss = 0.0;
02087 
02088   // G4cout<<" G4PAIxSection::GetStepCreLosnkovs "<<G4endl;
02089 
02090   meanNumber = fIntegralResonance[1]*step;
02091   numOfCollisions = G4Poisson(meanNumber);
02092 
02093   //   G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
02094 
02095   while(numOfCollisions)
02096   {
02097     loss += GetResonanceEnergyTransfer();
02098     numOfCollisions--;
02099   }
02100   // G4cout<<"PAI resonance loss = "<<loss/keV<<" keV"<<G4endl; 
02101 
02102   return loss;
02103 }

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

Definition at line 761 of file G4PAIxSection.cc.

Referenced by NormShift(), and SplainPAI().

00763 {
00764    G4double energy2,energy3,energy4,result;
00765 
00766    energy2 = energy1*energy1;
00767    energy3 = energy2*energy1;
00768    energy4 = energy3*energy1;
00769    
00770    result = fA1[k]/energy1+fA2[k]/energy2+fA3[k]/energy3+fA4[k]/energy4;  
00771    result *=hbarc/energy1;
00772    
00773    return result;
00774 
00775 }  // end of ImPartDielectricConst 

void G4PAIxSection::InitPAI (  ) 

Definition at line 519 of file G4PAIxSection.cc.

References DifPAIxSection(), IntegralCerenkov(), IntegralMM(), IntegralPAIxSection(), IntegralPlasmon(), IntegralResonance(), NormShift(), PAIdNdxCerenkov(), PAIdNdxMM(), PAIdNdxPlasmon(), PAIdNdxResonance(), and SplainPAI().

Referenced by G4PAIxSection().

00520 {    
00521    G4int i;
00522    G4double betaGammaSq = fLorentzFactor[fRefGammaNumber]*
00523                           fLorentzFactor[fRefGammaNumber] - 1;
00524 
00525    // Preparation of integral PAI cross section for reference gamma
00526    
00527    NormShift(betaGammaSq);             
00528    SplainPAI(betaGammaSq);
00529 
00530    IntegralPAIxSection();
00531    IntegralCerenkov();
00532    IntegralMM();
00533    IntegralPlasmon();
00534    IntegralResonance();
00535 
00536    for(i = 0; i<= fSplineNumber; i++)
00537    {
00538       fPAItable[i][fRefGammaNumber] = fIntegralPAIxSection[i];
00539       if(i != 0) 
00540       {
00541          fPAItable[i][0] = fSplineEnergy[i];
00542       }
00543    }
00544    fPAItable[0][0] = fSplineNumber;
00545    
00546    for(G4int j = 1; j < 112; j++)       // for other gammas
00547    {
00548       if( j == fRefGammaNumber ) continue;
00549       
00550       betaGammaSq = fLorentzFactor[j]*fLorentzFactor[j] - 1;
00551       
00552       for(i = 1; i <= fSplineNumber; i++)
00553       {
00554          fDifPAIxSection[i] = DifPAIxSection(i,betaGammaSq);
00555          fdNdxCerenkov[i]   = PAIdNdxCerenkov(i,betaGammaSq);
00556          fdNdxMM[i]   = PAIdNdxMM(i,betaGammaSq);
00557          fdNdxPlasmon[i]    = PAIdNdxPlasmon(i,betaGammaSq);
00558          fdNdxResonance[i]  = PAIdNdxResonance(i,betaGammaSq);
00559       }
00560       IntegralPAIxSection();
00561       IntegralCerenkov();
00562       IntegralMM();
00563       IntegralPlasmon();
00564       IntegralResonance();
00565       
00566       for(i = 0; i <= fSplineNumber; i++)
00567       {
00568          fPAItable[i][j] = fIntegralPAIxSection[i];
00569       }
00570    } 
00571 
00572 }  

void G4PAIxSection::IntegralCerenkov (  ) 

Definition at line 1205 of file G4PAIxSection.cc.

References SumOverBordCerenkov(), and SumOverInterCerenkov().

Referenced by G4PAIxSection(), and InitPAI().

01206 {
01207   G4int i, k;
01208    fIntegralCerenkov[fSplineNumber] = 0;
01209    fIntegralCerenkov[0] = 0;
01210    k = fIntervalNumber -1;
01211 
01212    for( i = fSplineNumber-1; i >= 1; i-- )
01213    {
01214       if(fSplineEnergy[i] >= fEnergyInterval[k])
01215       {
01216         fIntegralCerenkov[i] = fIntegralCerenkov[i+1] + SumOverInterCerenkov(i);
01217         // G4cout<<"int: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
01218       }
01219       else
01220       {
01221         fIntegralCerenkov[i] = fIntegralCerenkov[i+1] + 
01222                                    SumOverBordCerenkov(i+1,fEnergyInterval[k]);
01223         k--;
01224         // G4cout<<"bord: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
01225       }
01226    }
01227 
01228 }   // end of IntegralCerenkov 

void G4PAIxSection::IntegralMM (  ) 

Definition at line 1236 of file G4PAIxSection.cc.

References SumOverBordMM(), and SumOverInterMM().

Referenced by G4PAIxSection(), and InitPAI().

01237 {
01238   G4int i, k;
01239    fIntegralMM[fSplineNumber] = 0;
01240    fIntegralMM[0] = 0;
01241    k = fIntervalNumber -1;
01242 
01243    for( i = fSplineNumber-1; i >= 1; i-- )
01244    {
01245       if(fSplineEnergy[i] >= fEnergyInterval[k])
01246       {
01247         fIntegralMM[i] = fIntegralMM[i+1] + SumOverInterMM(i);
01248         // G4cout<<"int: i = "<<i<<"; sumC = "<<fIntegralMM[i]<<G4endl;
01249       }
01250       else
01251       {
01252         fIntegralMM[i] = fIntegralMM[i+1] + 
01253                                    SumOverBordMM(i+1,fEnergyInterval[k]);
01254         k--;
01255         // G4cout<<"bord: i = "<<i<<"; sumC = "<<fIntegralMM[i]<<G4endl;
01256       }
01257    }
01258 
01259 }   // end of IntegralMM 

void G4PAIxSection::IntegralPAIxSection (  ) 

Definition at line 1174 of file G4PAIxSection.cc.

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

Referenced by G4PAIxSection(), and InitPAI().

01175 {
01176   fIntegralPAIxSection[fSplineNumber] = 0;
01177   fIntegralPAIdEdx[fSplineNumber]     = 0;
01178   fIntegralPAIxSection[0]             = 0;
01179   G4int k = fIntervalNumber -1;
01180 
01181   for(G4int i = fSplineNumber-1; i >= 1; i--)
01182   {
01183     if(fSplineEnergy[i] >= fEnergyInterval[k])
01184     {
01185       fIntegralPAIxSection[i] = fIntegralPAIxSection[i+1] + SumOverInterval(i);
01186       fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] + SumOverIntervaldEdx(i);
01187     }
01188     else
01189     {
01190       fIntegralPAIxSection[i] = fIntegralPAIxSection[i+1] + 
01191                                    SumOverBorder(i+1,fEnergyInterval[k]);
01192       fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] + 
01193                                    SumOverBorderdEdx(i+1,fEnergyInterval[k]);
01194       k--;
01195     }
01196   }
01197 }   // end of IntegralPAIxSection 

void G4PAIxSection::IntegralPlasmon (  ) 

Definition at line 1267 of file G4PAIxSection.cc.

References SumOverBordPlasmon(), and SumOverInterPlasmon().

Referenced by G4PAIxSection(), and InitPAI().

01268 {
01269    fIntegralPlasmon[fSplineNumber] = 0;
01270    fIntegralPlasmon[0] = 0;
01271    G4int k = fIntervalNumber -1;
01272    for(G4int i=fSplineNumber-1;i>=1;i--)
01273    {
01274       if(fSplineEnergy[i] >= fEnergyInterval[k])
01275       {
01276         fIntegralPlasmon[i] = fIntegralPlasmon[i+1] + SumOverInterPlasmon(i);
01277       }
01278       else
01279       {
01280         fIntegralPlasmon[i] = fIntegralPlasmon[i+1] + 
01281                                    SumOverBordPlasmon(i+1,fEnergyInterval[k]);
01282         k--;
01283       }
01284    }
01285 
01286 }   // end of IntegralPlasmon

void G4PAIxSection::IntegralResonance (  ) 

Definition at line 1294 of file G4PAIxSection.cc.

References SumOverBordResonance(), and SumOverInterResonance().

Referenced by G4PAIxSection(), and InitPAI().

01295 {
01296    fIntegralResonance[fSplineNumber] = 0;
01297    fIntegralResonance[0] = 0;
01298    G4int k = fIntervalNumber -1;
01299    for(G4int i=fSplineNumber-1;i>=1;i--)
01300    {
01301       if(fSplineEnergy[i] >= fEnergyInterval[k])
01302       {
01303         fIntegralResonance[i] = fIntegralResonance[i+1] + SumOverInterResonance(i);
01304       }
01305       else
01306       {
01307         fIntegralResonance[i] = fIntegralResonance[i+1] + 
01308                                    SumOverBordResonance(i+1,fEnergyInterval[k]);
01309         k--;
01310       }
01311    }
01312 
01313 }   // end of IntegralResonance

void G4PAIxSection::NormShift ( G4double  betaGammaSq  ) 

Definition at line 579 of file G4PAIxSection.cc.

References DifPAIxSection(), ImPartDielectricConst(), PAIdNdxCerenkov(), PAIdNdxMM(), PAIdNdxPlasmon(), PAIdNdxResonance(), G4INCL::Math::pi, RePartDielectricConst(), and RutherfordIntegral().

Referenced by G4PAIxSection(), and InitPAI().

00580 {
00581   G4int i, j;
00582 
00583   for( i = 1; i <= fIntervalNumber-1; i++ )
00584   {
00585     for( j = 1; j <= 2; j++ )
00586     {
00587       fSplineNumber = (i-1)*2 + j;
00588 
00589       if( j == 1 ) fSplineEnergy[fSplineNumber] = fEnergyInterval[i  ]*(1+fDelta);
00590       else         fSplineEnergy[fSplineNumber] = fEnergyInterval[i+1]*(1-fDelta); 
00591       //    G4cout<<"cn = "<<fSplineNumber<<"; "<<"energy = "
00592       //  <<fSplineEnergy[fSplineNumber]<<G4endl;
00593     }
00594   }
00595   fIntegralTerm[1]=RutherfordIntegral(1,fEnergyInterval[1],fSplineEnergy[1]);
00596 
00597   j = 1;
00598 
00599   for( i = 2; i <= fSplineNumber; i++ )
00600   {
00601     if(fSplineEnergy[i]<fEnergyInterval[j+1])
00602     {
00603          fIntegralTerm[i] = fIntegralTerm[i-1] + 
00604                             RutherfordIntegral(j,fSplineEnergy[i-1],
00605                                                  fSplineEnergy[i]   );
00606     }
00607     else
00608     {
00609        G4double x = RutherfordIntegral(j,fSplineEnergy[i-1],
00610                                            fEnergyInterval[j+1]   );
00611          j++;
00612          fIntegralTerm[i] = fIntegralTerm[i-1] + x + 
00613                             RutherfordIntegral(j,fEnergyInterval[j],
00614                                                  fSplineEnergy[i]    );
00615     }
00616     // G4cout<<i<<"\t"<<fSplineEnergy[i]<<"\t"<<fIntegralTerm[i]<<"\n"<<G4endl;
00617   } 
00618   fNormalizationCof = 2*pi*pi*hbarc*hbarc*fine_structure_const/electron_mass_c2;
00619   fNormalizationCof *= fElectronDensity/fIntegralTerm[fSplineNumber];
00620 
00621   // G4cout<<"fNormalizationCof = "<<fNormalizationCof<<G4endl;
00622 
00623           // Calculation of PAI differrential cross-section (1/(keV*cm))
00624           // in the energy points near borders of energy intervals
00625 
00626    for(G4int k = 1; k <= fIntervalNumber-1; k++ )
00627    {
00628       for( j = 1; j <= 2; j++ )
00629       {
00630          i = (k-1)*2 + j;
00631          fImPartDielectricConst[i] = fNormalizationCof*
00632                                      ImPartDielectricConst(k,fSplineEnergy[i]);
00633          fRePartDielectricConst[i] = fNormalizationCof*
00634                                      RePartDielectricConst(fSplineEnergy[i]);
00635          fIntegralTerm[i] *= fNormalizationCof;
00636 
00637          fDifPAIxSection[i] = DifPAIxSection(i,betaGammaSq);
00638          fdNdxCerenkov[i]   = PAIdNdxCerenkov(i,betaGammaSq);
00639          fdNdxMM[i]   = PAIdNdxMM(i,betaGammaSq);
00640          fdNdxPlasmon[i]    = PAIdNdxPlasmon(i,betaGammaSq);
00641          fdNdxResonance[i]    = PAIdNdxResonance(i,betaGammaSq);
00642       }
00643    }
00644 
00645 }  // end of NormShift 

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

Definition at line 986 of file G4PAIxSection.cc.

References G4INCL::Math::pi.

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

00988 {        
00989    G4double logarithm, x3, x5, argument, modul2, dNdxC; 
00990    G4double be2, be4, betaBohr2,betaBohr4,cofBetaBohr;
00991 
00992    cofBetaBohr = 4.0;
00993    betaBohr2   = fine_structure_const*fine_structure_const;
00994    betaBohr4   = betaBohr2*betaBohr2*cofBetaBohr;
00995 
00996    be2 = betaGammaSq/(1 + betaGammaSq);
00997    be4 = be2*be2;
00998 
00999    if( betaGammaSq < 0.01 ) logarithm = log(1.0+betaGammaSq); // 0.0;
01000    else
01001    {
01002      logarithm  = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
01003                         (1/betaGammaSq - fRePartDielectricConst[i]) + 
01004                         fImPartDielectricConst[i]*fImPartDielectricConst[i] )*0.5;
01005      logarithm += log(1+1.0/betaGammaSq);
01006    }
01007 
01008    if( fImPartDielectricConst[i] == 0.0 || betaGammaSq < 0.01 )
01009    {
01010      argument = 0.0;
01011    }
01012    else
01013    {
01014      x3 = -fRePartDielectricConst[i] + 1.0/betaGammaSq;
01015      x5 = -1.0 - fRePartDielectricConst[i] +
01016           be2*((1.0 +fRePartDielectricConst[i])*(1.0 + fRePartDielectricConst[i]) +
01017           fImPartDielectricConst[i]*fImPartDielectricConst[i]);
01018      if( x3 == 0.0 ) argument = 0.5*pi;
01019      else            argument = atan2(fImPartDielectricConst[i],x3);
01020      argument *= x5 ;
01021    }   
01022    dNdxC = ( logarithm*fImPartDielectricConst[i] + argument )/hbarc;
01023   
01024    if(dNdxC < 1.0e-8) dNdxC = 1.0e-8;
01025 
01026    dNdxC *= fine_structure_const/be2/pi;
01027 
01028    dNdxC *= (1-exp(-be4/betaBohr4));
01029 
01030    if(fDensity >= 0.1)
01031    { 
01032       modul2 = (1.0 + fRePartDielectricConst[i])*(1.0 + fRePartDielectricConst[i]) + 
01033                     fImPartDielectricConst[i]*fImPartDielectricConst[i];
01034       dNdxC /= modul2;
01035    }
01036    return dNdxC;
01037 
01038 } // end of PAIdNdxCerenkov 

G4double G4PAIxSection::PAIdNdxMM ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 1044 of file G4PAIxSection.cc.

References G4INCL::Math::pi.

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

01046 {        
01047    G4double logarithm, x3, x5, argument, dNdxC; 
01048    G4double be2, be4, betaBohr2,betaBohr4,cofBetaBohr;
01049 
01050    cofBetaBohr = 4.0;
01051    betaBohr2   = fine_structure_const*fine_structure_const;
01052    betaBohr4   = betaBohr2*betaBohr2*cofBetaBohr;
01053 
01054    be2 = betaGammaSq/(1 + betaGammaSq);
01055    be4 = be2*be2;
01056 
01057    if( betaGammaSq < 0.01 ) logarithm = log(1.0+betaGammaSq); // 0.0;
01058    else
01059    {
01060      logarithm  = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
01061                         (1/betaGammaSq - fRePartDielectricConst[i]) + 
01062                         fImPartDielectricConst[i]*fImPartDielectricConst[i] )*0.5;
01063      logarithm += log(1+1.0/betaGammaSq);
01064    }
01065 
01066    if( fImPartDielectricConst[i] == 0.0 || betaGammaSq < 0.01 )
01067    {
01068      argument = 0.0;
01069    }
01070    else
01071    {
01072      x3 = -fRePartDielectricConst[i] + 1.0/betaGammaSq;
01073      x5 = be2*( 1.0 + fRePartDielectricConst[i] ) - 1.0;
01074      if( x3 == 0.0 ) argument = 0.5*pi;
01075      else            argument = atan2(fImPartDielectricConst[i],x3);
01076      argument *= x5 ;
01077    }   
01078    dNdxC = ( logarithm*fImPartDielectricConst[i]*be2 + argument )/hbarc;
01079   
01080    if(dNdxC < 1.0e-8) dNdxC = 1.0e-8;
01081 
01082    dNdxC *= fine_structure_const/be2/pi;
01083 
01084    dNdxC *= (1-exp(-be4/betaBohr4));
01085    return dNdxC;
01086 
01087 } // end of PAIdNdxMM 

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

Definition at line 1094 of file G4PAIxSection.cc.

References G4INCL::Math::pi.

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

01096 {        
01097    G4double resonance, modul2, dNdxP, cof = 1.;
01098    G4double be2, be4, betaBohr2, betaBohr4, cofBetaBohr;
01099 
01100   
01101    cofBetaBohr = 4.0;
01102    betaBohr2   = fine_structure_const*fine_structure_const;
01103    betaBohr4   = betaBohr2*betaBohr2*cofBetaBohr;
01104 
01105    be2 = betaGammaSq/(1 + betaGammaSq);
01106    be4 = be2*be2;
01107  
01108    resonance = log(2*electron_mass_c2*be2/fSplineEnergy[i]);  
01109    resonance *= fImPartDielectricConst[i]/hbarc;
01110 
01111 
01112    dNdxP = ( resonance + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i] );
01113 
01114    if( dNdxP < 1.0e-8 ) dNdxP = 1.0e-8;
01115 
01116    dNdxP *= fine_structure_const/be2/pi;
01117    dNdxP *= (1-exp(-be4/betaBohr4));
01118 
01119    if( fDensity >= 0.1 )
01120    { 
01121      modul2 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) + 
01122         fImPartDielectricConst[i]*fImPartDielectricConst[i];
01123      dNdxP /= modul2;
01124    }
01125    return dNdxP;
01126 
01127 } // end of PAIdNdxPlasmon 

G4double G4PAIxSection::PAIdNdxResonance ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 1134 of file G4PAIxSection.cc.

References G4INCL::Math::pi.

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

01136 {        
01137    G4double resonance, modul2, dNdxP;
01138    G4double be2, be4, betaBohr2, betaBohr4, cofBetaBohr;
01139 
01140    cofBetaBohr = 4.0;
01141    betaBohr2   = fine_structure_const*fine_structure_const;
01142    betaBohr4   = betaBohr2*betaBohr2*cofBetaBohr;
01143 
01144    be2 = betaGammaSq/(1 + betaGammaSq);
01145    be4 = be2*be2;
01146  
01147    resonance = log(2*electron_mass_c2*be2/fSplineEnergy[i]);  
01148    resonance *= fImPartDielectricConst[i]/hbarc;
01149 
01150 
01151    dNdxP = resonance;
01152 
01153    if( dNdxP < 1.0e-8 ) dNdxP = 1.0e-8;
01154 
01155    dNdxP *= fine_structure_const/be2/pi;
01156    dNdxP *= (1-exp(-be4/betaBohr4));
01157 
01158    if( fDensity >= 0.1 )
01159    { 
01160      modul2 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) + 
01161         fImPartDielectricConst[i]*fImPartDielectricConst[i];
01162      dNdxP /= modul2;
01163    }
01164    return dNdxP;
01165 
01166 } // end of PAIdNdxResonance 

G4double G4PAIxSection::RePartDielectricConst ( G4double  energy  ) 

Definition at line 856 of file G4PAIxSection.cc.

References G4INCL::Math::pi.

Referenced by NormShift(), and SplainPAI().

00857 {       
00858    G4double x0, x02, x03, x04, x05, x1, x2, xx1 ,xx2 , xx12,
00859             c1, c2, c3, cof1, cof2, xln1, xln2, xln3, result;
00860 
00861    x0 = enb;
00862    result = 0;
00863    
00864    for(G4int i=1;i<=fIntervalNumber-1;i++)
00865    {
00866       x1 = fEnergyInterval[i];
00867       x2 = fEnergyInterval[i+1];
00868       xx1 = x1 - x0;
00869       xx2 = x2 - x0;
00870       xx12 = xx2/xx1;
00871       
00872       if(xx12<0)
00873       {
00874          xx12 = -xx12;
00875       }
00876       xln1 = log(x2/x1);
00877       xln2 = log(xx12);
00878       xln3 = log((x2 + x0)/(x1 + x0));
00879       x02 = x0*x0;
00880       x03 = x02*x0;
00881       x04 = x03*x0;
00882       x05 = x04*x0;
00883       c1  = (x2 - x1)/x1/x2;
00884       c2  = (x2 - x1)*(x2 +x1)/x1/x1/x2/x2;
00885       c3  = (x2 -x1)*(x1*x1 + x1*x2 + x2*x2)/x1/x1/x1/x2/x2/x2;
00886 
00887       result -= (fA1[i]/x02 + fA3[i]/x04)*xln1;
00888       result -= (fA2[i]/x02 + fA4[i]/x04)*c1;
00889       result -= fA3[i]*c2/2/x02;
00890       result -= fA4[i]*c3/3/x02;
00891 
00892       cof1 = fA1[i]/x02 + fA3[i]/x04;
00893       cof2 = fA2[i]/x03 + fA4[i]/x05;
00894 
00895       result += 0.5*(cof1 +cof2)*xln2;
00896       result += 0.5*(cof1 - cof2)*xln3;
00897    } 
00898    result *= 2*hbarc/pi;
00899    
00900    return result;
00901 
00902 }   // end of RePartDielectricConst 

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

Definition at line 740 of file G4PAIxSection.cc.

Referenced by NormShift(), and SplainPAI().

00743 {
00744    G4double  c1, c2, c3;
00745    // G4cout<<"RI: x1 = "<<x1<<"; "<<"x2 = "<<x2<<G4endl;   
00746    c1 = (x2 - x1)/x1/x2;
00747    c2 = (x2 - x1)*(x2 + x1)/x1/x1/x2/x2;
00748    c3 = (x2 - x1)*(x1*x1 + x1*x2 + x2*x2)/x1/x1/x1/x2/x2/x2;
00749    // G4cout<<" RI: c1 = "<<c1<<"; "<<"c2 = "<<c2<<"; "<<"c3 = "<<c3<<G4endl;   
00750    
00751    return  fA1[k]*log(x2/x1) + fA2[k]*c1 + fA3[k]*c2/2 + fA4[k]*c3/3;
00752 
00753 }   // end of RutherfordIntegral 

void G4PAIxSection::SplainPAI ( G4double  betaGammaSq  ) 

Definition at line 653 of file G4PAIxSection.cc.

References DifPAIxSection(), ImPartDielectricConst(), PAIdNdxCerenkov(), PAIdNdxMM(), PAIdNdxPlasmon(), PAIdNdxResonance(), RePartDielectricConst(), and RutherfordIntegral().

Referenced by G4PAIxSection(), and InitPAI().

00654 {
00655    G4int k = 1;
00656    G4int i = 1;
00657 
00658    while ( (i < fSplineNumber) && (fSplineNumber < fMaxSplineSize-1) )
00659    {
00660       if(fSplineEnergy[i+1] > fEnergyInterval[k+1])
00661       {
00662           k++;   // Here next energy point is in next energy interval
00663           i++;
00664           continue;
00665       }
00666                        // Shifting of arrayes for inserting the geometrical 
00667                        // average of 'i' and 'i+1' energy points to 'i+1' place
00668       fSplineNumber++;
00669 
00670       for(G4int j = fSplineNumber; j >= i+2; j-- )
00671       {
00672          fSplineEnergy[j]          = fSplineEnergy[j-1];
00673          fImPartDielectricConst[j] = fImPartDielectricConst[j-1];
00674          fRePartDielectricConst[j] = fRePartDielectricConst[j-1];
00675          fIntegralTerm[j]          = fIntegralTerm[j-1];
00676 
00677          fDifPAIxSection[j] = fDifPAIxSection[j-1];
00678          fdNdxCerenkov[j]   = fdNdxCerenkov[j-1];
00679          fdNdxMM[j]   = fdNdxMM[j-1];
00680          fdNdxPlasmon[j]    = fdNdxPlasmon[j-1];
00681          fdNdxResonance[j]  = fdNdxResonance[j-1];
00682       }
00683       G4double x1  = fSplineEnergy[i];
00684       G4double x2  = fSplineEnergy[i+1];
00685       G4double yy1 = fDifPAIxSection[i];
00686       G4double y2  = fDifPAIxSection[i+1];
00687 
00688       G4double en1 = sqrt(x1*x2);
00689       fSplineEnergy[i+1] = en1;
00690 
00691                  // Calculation of logarithmic linear approximation
00692                  // in this (enr) energy point, which number is 'i+1' now
00693 
00694       G4double a = log10(y2/yy1)/log10(x2/x1);
00695       G4double b = log10(yy1) - a*log10(x1);
00696       G4double y = a*log10(en1) + b;
00697       y = pow(10.,y);
00698 
00699                  // Calculation of the PAI dif. cross-section at this point
00700 
00701       fImPartDielectricConst[i+1] = fNormalizationCof*
00702                                     ImPartDielectricConst(k,fSplineEnergy[i+1]);
00703       fRePartDielectricConst[i+1] = fNormalizationCof*
00704                                     RePartDielectricConst(fSplineEnergy[i+1]);
00705       fIntegralTerm[i+1] = fIntegralTerm[i] + fNormalizationCof*
00706                            RutherfordIntegral(k,fSplineEnergy[i],
00707                                                 fSplineEnergy[i+1]);
00708 
00709       fDifPAIxSection[i+1] = DifPAIxSection(i+1,betaGammaSq);
00710       fdNdxCerenkov[i+1]   = PAIdNdxCerenkov(i+1,betaGammaSq);
00711       fdNdxMM[i+1]   = PAIdNdxMM(i+1,betaGammaSq);
00712       fdNdxPlasmon[i+1]    = PAIdNdxPlasmon(i+1,betaGammaSq);
00713       fdNdxResonance[i+1]  = PAIdNdxResonance(i+1,betaGammaSq);
00714 
00715                   // Condition for next division of this segment or to pass
00716                   // to higher energies
00717 
00718       G4double x = 2*(fDifPAIxSection[i+1] - y)/(fDifPAIxSection[i+1] + y);
00719 
00720       if( x < 0 ) 
00721       {
00722          x = -x;
00723       }
00724       if( x > fError && fSplineNumber < fMaxSplineSize-1 )
00725       {
00726          continue;  // next division
00727       }
00728       i += 2;  // pass to next segment
00729 
00730    }   // close 'while'
00731 
00732 }  // end of SplainPAI 

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

Definition at line 1638 of file G4PAIxSection.cc.

Referenced by IntegralCerenkov().

01640 {               
01641    G4double x0,x1,y0,yy1,a,b,e0,c,d,result;
01642 
01643    e0 = en0;
01644    x0 = fSplineEnergy[i];
01645    x1 = fSplineEnergy[i+1];
01646    y0 = fdNdxCerenkov[i];
01647    yy1 = fdNdxCerenkov[i+1];
01648 
01649    //  G4cout<<G4endl;
01650    //  G4cout<<"SumBordC, i = "<<i<<"; en0 = "<<en0<<"; x0 ="<<x0<<"; x1 = "<<x1
01651    //     <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
01652    c = x1/x0;
01653    d = e0/x0;
01654    a = log10(yy1/y0)/log10(c);
01655    // b0 = log10(y0) - a*log10(x0);
01656    b = y0/pow(x0,a); // pow(10.,b0);   
01657    
01658    a += 1.0;
01659    if( a == 0 ) result = b*log(x0/e0);
01660    else         result = y0*(x0 - e0*pow(d,a-1))/a;   
01661    a += 1.0;
01662 
01663    if( a == 0 ) fIntegralCerenkov[0] += b*log(x0/e0);
01664    else         fIntegralCerenkov[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01665 
01666 // G4cout<<"a = "<<a<<"; b0 = "<<b0<<"; b = "<<b<<"; result = "<<result<<G4endl;
01667    
01668    x0  = fSplineEnergy[i - 1];
01669    x1  = fSplineEnergy[i - 2];
01670    y0  = fdNdxCerenkov[i - 1];
01671    yy1 = fdNdxCerenkov[i - 2];
01672 
01673    // G4cout<<"x0 ="<<x0<<"; x1 = "<<x1
01674    //    <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
01675 
01676    c = x1/x0;
01677    d = e0/x0;
01678    a  = log10(yy1/y0)/log10(x1/x0);
01679    // b0 = log10(y0) - a*log10(x0);
01680    b  =  y0/pow(x0,a);  // pow(10.,b0);
01681 
01682    a += 1.0;
01683    if( a == 0 ) result += b*log(e0/x0);
01684    else         result += y0*(e0*pow(d,a-1) - x0 )/a;
01685    a += 1.0;
01686 
01687    if( a == 0 )   fIntegralCerenkov[0] += b*log(e0/x0);
01688    else           fIntegralCerenkov[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01689 
01690    // G4cout<<"a = "<<a<<"; b0 = "<<b0<<"; b = "
01691    // <<b<<"; result = "<<result<<G4endl;    
01692 
01693    return result;
01694 
01695 } 

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

Definition at line 1515 of file G4PAIxSection.cc.

Referenced by IntegralPAIxSection().

01517 {               
01518   G4double x0,x1,y0,yy1,a,b,/*c,*/d,e0,result;
01519 
01520    e0 = en0;
01521    x0 = fSplineEnergy[i];
01522    x1 = fSplineEnergy[i+1];
01523    y0 = fDifPAIxSection[i];
01524    yy1 = fDifPAIxSection[i+1];
01525 
01526    //c = x1/x0;
01527    d = e0/x0;   
01528    a = log10(yy1/y0)/log10(x1/x0);
01529    // b0 = log10(y0) - a*log10(x0);
01530    b = y0/pow(x0,a);  // pow(10.,b);
01531    
01532    a += 1.;
01533    if( std::fabs(a) < 1.e-6 )
01534    {
01535       result = b*log(x0/e0);
01536    }
01537    else
01538    {
01539       result = y0*(x0 - e0*pow(d,a-1))/a;
01540    }
01541    a += 1.;
01542    if( std::fabs(a) < 1.e-6 )
01543    {
01544       fIntegralPAIxSection[0] += b*log(x0/e0);
01545    }
01546    else 
01547    {
01548       fIntegralPAIxSection[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01549    }
01550    x0 = fSplineEnergy[i - 1];
01551    x1 = fSplineEnergy[i - 2];
01552    y0 = fDifPAIxSection[i - 1];
01553    yy1 = fDifPAIxSection[i - 2];
01554 
01555    //c = x1/x0;
01556    d = e0/x0;   
01557    a = log10(yy1/y0)/log10(x1/x0);
01558    //  b0 = log10(y0) - a*log10(x0);
01559    b = y0/pow(x0,a);
01560    a += 1.;
01561    if( std::fabs(a) < 1.e-6 )
01562    {
01563       result += b*log(e0/x0);
01564    }
01565    else
01566    {
01567       result += y0*(e0*pow(d,a-1) - x0)/a;
01568    }
01569    a += 1.;
01570    if( std::fabs(a) < 1.e-6 ) 
01571    {
01572       fIntegralPAIxSection[0] += b*log(e0/x0);
01573    }
01574    else
01575    {
01576       fIntegralPAIxSection[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01577    }
01578    return result;
01579 
01580 } 

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

Definition at line 1584 of file G4PAIxSection.cc.

Referenced by IntegralPAIxSection().

01586 {               
01587   G4double x0,x1,y0,yy1,a,b,/*c,*/d,e0,result;
01588 
01589    e0 = en0;
01590    x0 = fSplineEnergy[i];
01591    x1 = fSplineEnergy[i+1];
01592    y0 = fDifPAIxSection[i];
01593    yy1 = fDifPAIxSection[i+1];
01594 
01595    //c = x1/x0;
01596    d = e0/x0;   
01597    a = log10(yy1/y0)/log10(x1/x0);
01598    // b0 = log10(y0) - a*log10(x0);
01599    b = y0/pow(x0,a);  // pow(10.,b);
01600    
01601    a += 2;
01602    if(a == 0)
01603    {
01604       result = b*log(x0/e0);
01605    }
01606    else 
01607    {
01608       result = y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01609    }
01610    x0 = fSplineEnergy[i - 1];
01611    x1 = fSplineEnergy[i - 2];
01612    y0 = fDifPAIxSection[i - 1];
01613    yy1 = fDifPAIxSection[i - 2];
01614 
01615    // c = x1/x0;
01616    d = e0/x0;   
01617    a = log10(yy1/y0)/log10(x1/x0);
01618    //  b0 = log10(y0) - a*log10(x0);
01619    b = y0/pow(x0,a);
01620    a += 2;
01621    if(a == 0) 
01622    {
01623       result += b*log(e0/x0);
01624    }
01625    else
01626    {
01627       result += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01628    }
01629    return result;
01630 
01631 } 

G4double G4PAIxSection::SumOverBordMM ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1702 of file G4PAIxSection.cc.

Referenced by IntegralMM().

01704 {               
01705    G4double x0,x1,y0,yy1,a,b,e0,c,d,result;
01706 
01707    e0 = en0;
01708    x0 = fSplineEnergy[i];
01709    x1 = fSplineEnergy[i+1];
01710    y0 = fdNdxMM[i];
01711    yy1 = fdNdxMM[i+1];
01712 
01713    //  G4cout<<G4endl;
01714    //  G4cout<<"SumBordC, i = "<<i<<"; en0 = "<<en0<<"; x0 ="<<x0<<"; x1 = "<<x1
01715    //     <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
01716    c = x1/x0;
01717    d = e0/x0;
01718    a = log10(yy1/y0)/log10(c);
01719    // b0 = log10(y0) - a*log10(x0);
01720    b = y0/pow(x0,a); // pow(10.,b0);   
01721    
01722    a += 1.0;
01723    if( a == 0 ) result = b*log(x0/e0);
01724    else         result = y0*(x0 - e0*pow(d,a-1))/a;   
01725    a += 1.0;
01726 
01727    if( a == 0 ) fIntegralMM[0] += b*log(x0/e0);
01728    else         fIntegralMM[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01729 
01730 // G4cout<<"a = "<<a<<"; b0 = "<<b0<<"; b = "<<b<<"; result = "<<result<<G4endl;
01731    
01732    x0  = fSplineEnergy[i - 1];
01733    x1  = fSplineEnergy[i - 2];
01734    y0  = fdNdxMM[i - 1];
01735    yy1 = fdNdxMM[i - 2];
01736 
01737    // G4cout<<"x0 ="<<x0<<"; x1 = "<<x1
01738    //    <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
01739 
01740    c = x1/x0;
01741    d = e0/x0;
01742    a  = log10(yy1/y0)/log10(x1/x0);
01743    // b0 = log10(y0) - a*log10(x0);
01744    b  =  y0/pow(x0,a);  // pow(10.,b0);
01745 
01746    a += 1.0;
01747    if( a == 0 ) result += b*log(e0/x0);
01748    else         result += y0*(e0*pow(d,a-1) - x0 )/a;
01749    a += 1.0;
01750 
01751    if( a == 0 )   fIntegralMM[0] += b*log(e0/x0);
01752    else           fIntegralMM[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01753 
01754    // G4cout<<"a = "<<a<<"; b0 = "<<b0<<"; b = "
01755    // <<b<<"; result = "<<result<<G4endl;    
01756 
01757    return result;
01758 
01759 } 

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

Definition at line 1766 of file G4PAIxSection.cc.

Referenced by IntegralPlasmon().

01768 {               
01769    G4double x0,x1,y0,yy1,a,b,c,d,e0,result;
01770 
01771    e0 = en0;
01772    x0 = fSplineEnergy[i];
01773    x1 = fSplineEnergy[i+1];
01774    y0 = fdNdxPlasmon[i];
01775    yy1 = fdNdxPlasmon[i+1];
01776 
01777    c = x1/x0;
01778    d = e0/x0;   
01779    a = log10(yy1/y0)/log10(c);
01780    //  b0 = log10(y0) - a*log10(x0);
01781    b = y0/pow(x0,a); //pow(10.,b);
01782    
01783    a += 1.0;
01784    if( a == 0 ) result = b*log(x0/e0);
01785    else         result = y0*(x0 - e0*pow(d,a-1))/a;   
01786    a += 1.0;
01787 
01788    if( a == 0 ) fIntegralPlasmon[0] += b*log(x0/e0);
01789    else         fIntegralPlasmon[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01790    
01791    x0 = fSplineEnergy[i - 1];
01792    x1 = fSplineEnergy[i - 2];
01793    y0 = fdNdxPlasmon[i - 1];
01794    yy1 = fdNdxPlasmon[i - 2];
01795 
01796    c = x1/x0;
01797    d = e0/x0;
01798    a = log10(yy1/y0)/log10(c);
01799    // b0 = log10(y0) - a*log10(x0);
01800    b = y0/pow(x0,a);// pow(10.,b0);
01801 
01802    a += 1.0;
01803    if( a == 0 ) result += b*log(e0/x0);
01804    else         result += y0*(e0*pow(d,a-1) - x0)/a;
01805    a += 1.0;
01806 
01807    if( a == 0 )   fIntegralPlasmon[0] += b*log(e0/x0);
01808    else           fIntegralPlasmon[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01809    
01810    return result;
01811 
01812 } 

G4double G4PAIxSection::SumOverBordResonance ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1819 of file G4PAIxSection.cc.

Referenced by IntegralResonance().

01821 {               
01822    G4double x0,x1,y0,yy1,a,b,c,d,e0,result;
01823 
01824    e0 = en0;
01825    x0 = fSplineEnergy[i];
01826    x1 = fSplineEnergy[i+1];
01827    y0 = fdNdxResonance[i];
01828    yy1 = fdNdxResonance[i+1];
01829 
01830    c = x1/x0;
01831    d = e0/x0;   
01832    a = log10(yy1/y0)/log10(c);
01833    //  b0 = log10(y0) - a*log10(x0);
01834    b = y0/pow(x0,a); //pow(10.,b);
01835    
01836    a += 1.0;
01837    if( a == 0 ) result = b*log(x0/e0);
01838    else         result = y0*(x0 - e0*pow(d,a-1))/a;   
01839    a += 1.0;
01840 
01841    if( a == 0 ) fIntegralResonance[0] += b*log(x0/e0);
01842    else         fIntegralResonance[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
01843    
01844    x0 = fSplineEnergy[i - 1];
01845    x1 = fSplineEnergy[i - 2];
01846    y0 = fdNdxResonance[i - 1];
01847    yy1 = fdNdxResonance[i - 2];
01848 
01849    c = x1/x0;
01850    d = e0/x0;
01851    a = log10(yy1/y0)/log10(c);
01852    // b0 = log10(y0) - a*log10(x0);
01853    b = y0/pow(x0,a);// pow(10.,b0);
01854 
01855    a += 1.0;
01856    if( a == 0 ) result += b*log(e0/x0);
01857    else         result += y0*(e0*pow(d,a-1) - x0)/a;
01858    a += 1.0;
01859 
01860    if( a == 0 )   fIntegralResonance[0] += b*log(e0/x0);
01861    else           fIntegralResonance[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
01862    
01863    return result;
01864 
01865 } 

G4double G4PAIxSection::SumOverInterCerenkov ( G4int  intervalNumber  ) 

Definition at line 1388 of file G4PAIxSection.cc.

Referenced by IntegralCerenkov().

01389 {         
01390    G4double x0,x1,y0,yy1,a,b,c,result;
01391 
01392    x0  = fSplineEnergy[i];
01393    x1  = fSplineEnergy[i+1];
01394    y0  = fdNdxCerenkov[i];
01395    yy1 = fdNdxCerenkov[i+1];
01396    // G4cout<<"SumC, i = "<<i<<"; x0 ="<<x0<<"; x1 = "<<x1
01397    //   <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
01398 
01399    c = x1/x0;
01400    a = log10(yy1/y0)/log10(c);
01401    b = y0/pow(x0,a);
01402 
01403    a += 1.0;
01404    if(a == 0) result = b*log(c);
01405    else       result = y0*(x1*pow(c,a-1) - x0)/a;   
01406    a += 1.0;
01407 
01408    if( a == 0 ) fIntegralCerenkov[0] += b*log(x1/x0);
01409    else         fIntegralCerenkov[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
01410    //  G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;   
01411    return result;
01412 
01413 } //  end of SumOverInterCerenkov

G4double G4PAIxSection::SumOverInterMM ( G4int  intervalNumber  ) 

Definition at line 1421 of file G4PAIxSection.cc.

Referenced by IntegralMM().

01422 {         
01423    G4double x0,x1,y0,yy1,a,b,c,result;
01424 
01425    x0  = fSplineEnergy[i];
01426    x1  = fSplineEnergy[i+1];
01427    y0  = fdNdxMM[i];
01428    yy1 = fdNdxMM[i+1];
01429    // G4cout<<"SumC, i = "<<i<<"; x0 ="<<x0<<"; x1 = "<<x1
01430    //   <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
01431 
01432    c = x1/x0;
01433    a = log10(yy1/y0)/log10(c);
01434    b = y0/pow(x0,a);
01435 
01436    a += 1.0;
01437    if(a == 0) result = b*log(c);
01438    else       result = y0*(x1*pow(c,a-1) - x0)/a;   
01439    a += 1.0;
01440 
01441    if( a == 0 ) fIntegralMM[0] += b*log(x1/x0);
01442    else         fIntegralMM[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
01443    //  G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;   
01444    return result;
01445 
01446 } //  end of SumOverInterMM

G4double G4PAIxSection::SumOverInterPlasmon ( G4int  intervalNumber  ) 

Definition at line 1454 of file G4PAIxSection.cc.

Referenced by IntegralPlasmon().

01455 {         
01456    G4double x0,x1,y0,yy1,a,b,c,result;
01457 
01458    x0  = fSplineEnergy[i];
01459    x1  = fSplineEnergy[i+1];
01460    y0  = fdNdxPlasmon[i];
01461    yy1 = fdNdxPlasmon[i+1];
01462    c =x1/x0;
01463    a = log10(yy1/y0)/log10(c);
01464    // b = log10(y0) - a*log10(x0);
01465    b = y0/pow(x0,a);
01466 
01467    a += 1.0;
01468    if(a == 0) result = b*log(x1/x0);
01469    else       result = y0*(x1*pow(c,a-1) - x0)/a;   
01470    a += 1.0;
01471 
01472    if( a == 0 ) fIntegralPlasmon[0] += b*log(x1/x0);
01473    else         fIntegralPlasmon[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
01474    
01475    return result;
01476 
01477 } //  end of SumOverInterPlasmon

G4double G4PAIxSection::SumOverInterResonance ( G4int  intervalNumber  ) 

Definition at line 1485 of file G4PAIxSection.cc.

Referenced by IntegralResonance().

01486 {         
01487    G4double x0,x1,y0,yy1,a,b,c,result;
01488 
01489    x0  = fSplineEnergy[i];
01490    x1  = fSplineEnergy[i+1];
01491    y0  = fdNdxResonance[i];
01492    yy1 = fdNdxResonance[i+1];
01493    c =x1/x0;
01494    a = log10(yy1/y0)/log10(c);
01495    // b = log10(y0) - a*log10(x0);
01496    b = y0/pow(x0,a);
01497 
01498    a += 1.0;
01499    if(a == 0) result = b*log(x1/x0);
01500    else       result = y0*(x1*pow(c,a-1) - x0)/a;   
01501    a += 1.0;
01502 
01503    if( a == 0 ) fIntegralResonance[0] += b*log(x1/x0);
01504    else         fIntegralResonance[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
01505    
01506    return result;
01507 
01508 } //  end of SumOverInterResonance

G4double G4PAIxSection::SumOverInterval ( G4int  intervalNumber  ) 

Definition at line 1321 of file G4PAIxSection.cc.

Referenced by IntegralPAIxSection().

01322 {         
01323    G4double x0,x1,y0,yy1,a,b,c,result;
01324 
01325    x0 = fSplineEnergy[i];
01326    x1 = fSplineEnergy[i+1];
01327    y0 = fDifPAIxSection[i];
01328    yy1 = fDifPAIxSection[i+1];
01329    c = x1/x0;
01330    a = log10(yy1/y0)/log10(c);
01331    // b = log10(y0) - a*log10(x0);
01332    b = y0/pow(x0,a);
01333    a += 1.;
01334    if( std::fabs(a) < 1.e-6 ) 
01335    {
01336       result = b*log(x1/x0);
01337    }
01338    else
01339    {
01340       result = y0*(x1*pow(c,a-1) - x0)/a;
01341    }
01342    a += 1.;
01343    if( std::fabs(a) < 1.e-6 ) 
01344    {
01345       fIntegralPAIxSection[0] += b*log(x1/x0);
01346    }
01347    else
01348    {
01349       fIntegralPAIxSection[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
01350    }
01351    return result;
01352 
01353 } //  end of SumOverInterval

G4double G4PAIxSection::SumOverIntervaldEdx ( G4int  intervalNumber  ) 

Definition at line 1357 of file G4PAIxSection.cc.

Referenced by IntegralPAIxSection().

01358 {         
01359    G4double x0,x1,y0,yy1,a,b,c,result;
01360 
01361    x0 = fSplineEnergy[i];
01362    x1 = fSplineEnergy[i+1];
01363    y0 = fDifPAIxSection[i];
01364    yy1 = fDifPAIxSection[i+1];
01365    c = x1/x0;
01366    a = log10(yy1/y0)/log10(c);
01367    // b = log10(y0) - a*log10(x0);
01368    b = y0/pow(x0,a);
01369    a += 2;
01370    if(a == 0) 
01371    {
01372      result = b*log(x1/x0);
01373    }
01374    else
01375    {
01376      result = y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
01377    }
01378    return result;
01379 
01380 } //  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