G4PenelopeCrossSection Class Reference

#include <G4PenelopeCrossSection.hh>


Public Member Functions

 G4PenelopeCrossSection (size_t nOfEnergyPoints, size_t nOfShells=0)
 ~G4PenelopeCrossSection ()
G4double GetTotalCrossSection (G4double energy)
 Returns total cross section at the given energy.
G4double GetHardCrossSection (G4double energy)
 Returns hard cross section at the given energy.
G4double GetSoftStoppingPower (G4double energy)
 Returns the total stopping power due to soft collisions.
G4double GetShellCrossSection (size_t shellID, G4double energy)
 Returns the hard cross section for the given shell (per molecule).
G4double GetNormalizedShellCrossSection (size_t shellID, G4double energy)
 Returns the hard cross section for the given shell (normalized to 1).
size_t GetNumberOfShells ()
void AddCrossSectionPoint (size_t binNumber, G4double energy, G4double XH0, G4double XH1, G4double XH2, G4double XS0, G4double XS1, G4double XS2)
void AddShellCrossSectionPoint (size_t binNumber, size_t shellID, G4double energy, G4double xs)


Detailed Description

Definition at line 72 of file G4PenelopeCrossSection.hh.


Constructor & Destructor Documentation

G4PenelopeCrossSection::G4PenelopeCrossSection ( size_t  nOfEnergyPoints,
size_t  nOfShells = 0 
)

Definition at line 44 of file G4PenelopeCrossSection.cc.

References FatalException, G4endl, G4Exception(), and G4PhysicsTable::push_back().

00044                                                                              : 
00045   numberOfEnergyPoints(nPointsE),numberOfShells(nShells),softCrossSections(0),
00046   hardCrossSections(0),shellCrossSections(0),shellNormalizedCrossSections(0)
00047 {
00048   //check the number of points is not zero
00049   if (!numberOfEnergyPoints)
00050     {
00051       G4ExceptionDescription ed;
00052       ed << "G4PenelopeCrossSection: invalid number of energy points " << G4endl;
00053       G4Exception("G4PenelopeCrossSection::G4PenelopeCrossSection()",
00054                   "em2017",FatalException,ed);
00055     }
00056 
00057   isNormalized = false;
00058 
00059   // 1) soft XS table
00060   softCrossSections = new G4PhysicsTable();
00061   //the table contains 3 G4PhysicsFreeVectors, 
00062   //(softCrossSections)[0] -->  log XS0 vs. log E
00063   //(softCrossSections)[1] -->  log XS1 vs. log E
00064   //(softCrossSections)[2] -->  log XS2 vs. log E
00065   
00066   //everything is log-log
00067   for (size_t i=0;i<3;i++)
00068     softCrossSections->push_back(new G4PhysicsFreeVector(numberOfEnergyPoints));
00069   
00070   //2) hard XS table
00071   hardCrossSections = new G4PhysicsTable();
00072   //the table contains 3 G4PhysicsFreeVectors, 
00073   //(hardCrossSections)[0] -->  log XH0 vs. log E
00074   //(hardCrossSections)[1] -->  log XH1 vs. log E 
00075   //(hardCrossSections)[2] -->  log XH2 vs. log E
00076   
00077   //everything is log-log
00078   for (size_t i=0;i<3;i++)
00079     hardCrossSections->push_back(new G4PhysicsFreeVector(numberOfEnergyPoints));
00080   
00081   //3) shell XS table, if it is the case
00082   if (numberOfShells)
00083     {
00084       shellCrossSections = new G4PhysicsTable();
00085       shellNormalizedCrossSections = new G4PhysicsTable();
00086       //the table has to contain numberofShells G4PhysicsFreeVectors, 
00087       //(theTable)[ishell] --> cross section for shell #ishell
00088       for (size_t i=0;i<numberOfShells;i++)     
00089         {
00090           shellCrossSections->push_back(new G4PhysicsFreeVector(numberOfEnergyPoints));       
00091           shellNormalizedCrossSections->push_back(new G4PhysicsFreeVector(numberOfEnergyPoints));  
00092         }
00093     }
00094 }

G4PenelopeCrossSection::~G4PenelopeCrossSection (  ) 

Definition at line 97 of file G4PenelopeCrossSection.cc.

References G4PhysicsTable::clearAndDestroy().

00098 {
00099   //clean up tables
00100   if (shellCrossSections)
00101     {
00102       shellCrossSections->clearAndDestroy();
00103       delete shellCrossSections;          
00104     }
00105   if (shellNormalizedCrossSections)
00106     {
00107       shellNormalizedCrossSections->clearAndDestroy();
00108       delete shellNormalizedCrossSections;
00109     }
00110   if (softCrossSections)
00111     {
00112       softCrossSections->clearAndDestroy();
00113       delete softCrossSections;
00114     }
00115   if (hardCrossSections)
00116     {
00117       hardCrossSections->clearAndDestroy();
00118       delete hardCrossSections;
00119     }
00120 }


Member Function Documentation

void G4PenelopeCrossSection::AddCrossSectionPoint ( size_t  binNumber,
G4double  energy,
G4double  XH0,
G4double  XH1,
G4double  XH2,
G4double  XS0,
G4double  XS1,
G4double  XS2 
)

Definition at line 123 of file G4PenelopeCrossSection.cc.

References G4cout, G4endl, and G4PhysicsFreeVector::PutValue().

00128 {
00129   if (!softCrossSections || !hardCrossSections)
00130     {
00131       G4cout << "Something wrong in G4PenelopeCrossSection::AddCrossSectionPoint" <<
00132         G4endl;
00133       G4cout << "Trying to fill un-initialized tables" << G4endl;
00134       return;
00135     }
00136   
00137   //fill vectors
00138   G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*softCrossSections)[0];
00139  
00140   if (binNumber >= numberOfEnergyPoints)
00141     {
00142       G4cout << "Something wrong in G4PenelopeCrossSection::AddCrossSectionPoint" <<
00143         G4endl;
00144       G4cout << "Trying to register more points than originally declared" << G4endl;
00145       return;     
00146     }
00147    G4double logEne = std::log(energy);
00148 
00149    //XS0
00150    G4double val = std::log(std::max(XS0,1e-42*cm2)); //avoid log(0)
00151    theVector->PutValue(binNumber,logEne,val);
00152 
00153    //XS1
00154    theVector = (G4PhysicsFreeVector*) (*softCrossSections)[1];
00155    val =  std::log(std::max(XS1,1e-42*eV*cm2)); //avoid log(0)
00156    theVector->PutValue(binNumber,logEne,val);
00157 
00158    //XS2
00159    theVector = (G4PhysicsFreeVector*) (*softCrossSections)[2];
00160    val =  std::log(std::max(XS2,1e-42*eV*eV*cm2)); //avoid log(0)
00161    theVector->PutValue(binNumber,logEne,val);
00162 
00163    //XH0
00164    theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[0];
00165    val =  std::log(std::max(XH0,1e-42*cm2)); //avoid log(0)
00166    theVector->PutValue(binNumber,logEne,val);
00167 
00168    //XH1
00169    theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[1];
00170    val =  std::log(std::max(XH1,1e-42*eV*cm2)); //avoid log(0)
00171    theVector->PutValue(binNumber,logEne,val);
00172 
00173     //XH2
00174    theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[2];
00175    val =  std::log(std::max(XH2,1e-42*eV*eV*cm2)); //avoid log(0)
00176    theVector->PutValue(binNumber,logEne,val);
00177    
00178    return;
00179 }

void G4PenelopeCrossSection::AddShellCrossSectionPoint ( size_t  binNumber,
size_t  shellID,
G4double  energy,
G4double  xs 
)

Definition at line 183 of file G4PenelopeCrossSection.cc.

References G4cout, G4endl, and G4PhysicsFreeVector::PutValue().

00187 {
00188   if (!shellCrossSections)
00189     {
00190       G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
00191         G4endl;
00192       G4cout << "Trying to fill un-initialized table" << G4endl;
00193       return;
00194     }
00195   
00196   if (shellID >= numberOfShells)
00197     {
00198       G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
00199         G4endl;
00200       G4cout << "Trying to fill shell #" << shellID << " while the maximum is " 
00201              <<  numberOfShells-1 << G4endl;      
00202       return;    
00203     }
00204 
00205   //fill vector
00206   G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*shellCrossSections)[shellID];
00207  
00208   if (binNumber >= numberOfEnergyPoints)
00209     {
00210       G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
00211         G4endl;
00212       G4cout << "Trying to register more points than originally declared" << G4endl;
00213       return;     
00214     }
00215    G4double logEne = std::log(energy);
00216    G4double val = std::log(std::max(xs,1e-42*cm2)); //avoid log(0)
00217    theVector->PutValue(binNumber,logEne,val);
00218 
00219    return;
00220 }

G4double G4PenelopeCrossSection::GetHardCrossSection ( G4double  energy  ) 

Returns hard cross section at the given energy.

Definition at line 268 of file G4PenelopeCrossSection.cc.

References G4cout, G4endl, G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

Referenced by G4PenelopeIonisationModel::CrossSectionPerVolume(), and G4PenelopeBremsstrahlungModel::CrossSectionPerVolume().

00269 {
00270   G4double result = 0;
00271   //take here XH0
00272   if (!hardCrossSections)
00273     {
00274       G4cout << "Something wrong in G4PenelopeCrossSection::GetHardCrossSection" <<
00275         G4endl;
00276       G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
00277       return result;
00278     }
00279  
00280   G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[0];
00281   if (theVector->GetVectorLength() < numberOfEnergyPoints)
00282     {
00283       G4cout << "Something wrong in G4PenelopeCrossSection::GetHardCrossSection" <<
00284         G4endl;
00285       G4cout << "Hard cross section table looks not filled" << G4endl;
00286       return result;
00287     }
00288   G4double logene = std::log(energy);
00289   G4double logXS = theVector->Value(logene);
00290   result = std::exp(logXS);
00291 
00292   return result;
00293 }

G4double G4PenelopeCrossSection::GetNormalizedShellCrossSection ( size_t  shellID,
G4double  energy 
)

Returns the hard cross section for the given shell (normalized to 1).

Definition at line 363 of file G4PenelopeCrossSection.cc.

References G4cout, G4endl, G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

00364 {
00365   G4double result = 0;
00366   if (!shellNormalizedCrossSections)
00367     {
00368       G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
00369         G4endl;
00370       G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
00371       return result;
00372     }
00373 
00374   if (!isNormalized)
00375     NormalizeShellCrossSections();
00376 
00377   if (shellID >= numberOfShells)
00378     {
00379       G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
00380         G4endl;
00381       G4cout << "Trying to retrieve shell #" << shellID << " while the maximum is " 
00382              <<  numberOfShells-1 << G4endl;  
00383       return result;
00384     }
00385  
00386   G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*shellNormalizedCrossSections)[shellID];
00387 
00388   if (theVector->GetVectorLength() < numberOfEnergyPoints)
00389     {
00390       G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
00391         G4endl;
00392       G4cout << "Shell cross section table looks not filled" << G4endl;
00393       return result;
00394     }
00395   G4double logene = std::log(energy);
00396   G4double logXS = theVector->Value(logene);
00397   result = std::exp(logXS);
00398 
00399   return result;
00400 }

size_t G4PenelopeCrossSection::GetNumberOfShells (  )  [inline]

Definition at line 94 of file G4PenelopeCrossSection.hh.

00094 {return numberOfShells;};

G4double G4PenelopeCrossSection::GetShellCrossSection ( size_t  shellID,
G4double  energy 
)

Returns the hard cross section for the given shell (per molecule).

Definition at line 327 of file G4PenelopeCrossSection.cc.

References G4cout, G4endl, G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

Referenced by G4PenelopeIonisationCrossSection::CrossSection().

00328 {
00329   G4double result = 0;
00330   if (!shellCrossSections)
00331     {
00332       G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
00333         G4endl;
00334       G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
00335       return result;
00336     }
00337   if (shellID >= numberOfShells)
00338     {
00339       G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
00340         G4endl;
00341       G4cout << "Trying to retrieve shell #" << shellID << " while the maximum is " 
00342              <<  numberOfShells-1 << G4endl;  
00343       return result;
00344     }
00345 
00346   G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*shellCrossSections)[shellID];
00347 
00348   if (theVector->GetVectorLength() < numberOfEnergyPoints)
00349     {
00350       G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
00351         G4endl;
00352       G4cout << "Shell cross section table looks not filled" << G4endl;
00353       return result;
00354     }
00355   G4double logene = std::log(energy);
00356   G4double logXS = theVector->Value(logene);
00357   result = std::exp(logXS);
00358 
00359   return result;
00360 }

G4double G4PenelopeCrossSection::GetSoftStoppingPower ( G4double  energy  ) 

Returns the total stopping power due to soft collisions.

Definition at line 298 of file G4PenelopeCrossSection.cc.

References G4cout, G4endl, G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

Referenced by G4PenelopeIonisationModel::ComputeDEDXPerVolume(), and G4PenelopeBremsstrahlungModel::ComputeDEDXPerVolume().

00299 {
00300   G4double result = 0;
00301   //take here XH0
00302   if (!softCrossSections)
00303     {
00304       G4cout << "Something wrong in G4PenelopeCrossSection::GetSoftStoppingPower" <<
00305         G4endl;
00306       G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
00307       return result;
00308     }
00309  
00310   G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*softCrossSections)[1];
00311   if (theVector->GetVectorLength() < numberOfEnergyPoints)
00312     {
00313       G4cout << "Something wrong in G4PenelopeCrossSection::GetSoftStoppingPower" <<
00314         G4endl;
00315       G4cout << "Soft cross section table looks not filled" << G4endl;
00316       return result;
00317     }
00318   G4double logene = std::log(energy);
00319   G4double logXS = theVector->Value(logene);
00320   result = std::exp(logXS);
00321 
00322   return result;
00323 }

G4double G4PenelopeCrossSection::GetTotalCrossSection ( G4double  energy  ) 

Returns total cross section at the given energy.

Definition at line 224 of file G4PenelopeCrossSection.cc.

References G4cout, G4endl, G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

00225 {
00226   G4double result = 0;
00227   //take here XS0 + XH0
00228   if (!softCrossSections || !hardCrossSections)
00229     {
00230       G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
00231         G4endl;
00232       G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
00233       return result;
00234     }
00235   
00236   // 1) soft part
00237   G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*softCrossSections)[0];
00238   if (theVector->GetVectorLength() < numberOfEnergyPoints)
00239     {
00240       G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
00241         G4endl;
00242       G4cout << "Soft cross section table looks not filled" << G4endl;
00243       return result;
00244     }
00245   G4double logene = std::log(energy);
00246   G4double logXS = theVector->Value(logene);
00247   G4double softXS = std::exp(logXS);
00248 
00249    // 2) hard part
00250   theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[0];
00251   if (theVector->GetVectorLength() < numberOfEnergyPoints)
00252     {
00253       G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
00254         G4endl;
00255       G4cout << "Hard cross section table looks not filled" << G4endl;
00256       return result;
00257     }
00258   logXS = theVector->Value(logene);
00259   G4double hardXS = std::exp(logXS);
00260 
00261   result = hardXS + softXS;
00262   return result;  
00263 
00264 }


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