G4eIonisationCrossSectionHandler Class Reference

#include <G4eIonisationCrossSectionHandler.hh>

Inheritance diagram for G4eIonisationCrossSectionHandler:

G4VCrossSectionHandler

Public Member Functions

 G4eIonisationCrossSectionHandler (const G4VEnergySpectrum *spec, G4VDataSetAlgorithm *alg, G4double emin, G4double emax, G4int nbin)
 ~G4eIonisationCrossSectionHandler ()
G4double GetCrossSectionAboveThresholdForElement (G4double energy, G4double cutEnergy, G4int Z)

Protected Member Functions

std::vector< G4VEMDataSet * > * BuildCrossSectionsForMaterials (const G4DataVector &energyVector, const G4DataVector *energyCuts)

Detailed Description

Definition at line 63 of file G4eIonisationCrossSectionHandler.hh.


Constructor & Destructor Documentation

G4eIonisationCrossSectionHandler::G4eIonisationCrossSectionHandler ( const G4VEnergySpectrum spec,
G4VDataSetAlgorithm alg,
G4double  emin,
G4double  emax,
G4int  nbin 
)

Definition at line 73 of file G4eIonisationCrossSectionHandler.cc.

References G4VCrossSectionHandler::Initialise().

00076  :  G4VCrossSectionHandler(),
00077     theParam(spec),verbose(0)
00078 {
00079   G4VCrossSectionHandler::Initialise(alg, emin, emax, nbin);
00080   interp = new G4LinLogLogInterpolation();
00081 }

G4eIonisationCrossSectionHandler::~G4eIonisationCrossSectionHandler (  ) 

Definition at line 84 of file G4eIonisationCrossSectionHandler.cc.

00085 {
00086   delete interp;
00087 }


Member Function Documentation

std::vector< G4VEMDataSet * > * G4eIonisationCrossSectionHandler::BuildCrossSectionsForMaterials ( const G4DataVector energyVector,
const G4DataVector energyCuts 
) [protected, virtual]

Implements G4VCrossSectionHandler.

Definition at line 90 of file G4eIonisationCrossSectionHandler.cc.

References G4VEMDataSet::AddComponent(), G4VDataSetAlgorithm::Clone(), G4VCrossSectionHandler::FindValue(), G4cout, G4endl, G4Material::GetAtomicNumDensityVector(), G4Material::GetElementVector(), G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4Material::GetName(), G4Material::GetNumberOfElements(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), CLHEP::detail::n, G4VCrossSectionHandler::NumberOfComponents(), and G4VEnergySpectrum::Probability().

00093 {
00094   std::vector<G4VEMDataSet*>* set = new std::vector<G4VEMDataSet*>;
00095 
00096   G4DataVector* energies;
00097   G4DataVector* cs;
00098 
00099   G4DataVector* log_energies;
00100   G4DataVector* log_cs;
00101 
00102   G4int nOfBins = energyVector.size();
00103 
00104   const G4ProductionCutsTable* theCoupleTable=
00105         G4ProductionCutsTable::GetProductionCutsTable();
00106   size_t numOfCouples = theCoupleTable->GetTableSize();
00107 
00108   for (size_t mLocal=0; mLocal<numOfCouples; mLocal++) {
00109 
00110     const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(mLocal);
00111     const G4Material* material= couple->GetMaterial();
00112     const G4ElementVector* elementVector = material->GetElementVector();
00113     const G4double* nAtomsPerVolume = material->GetAtomicNumDensityVector();
00114     G4int nElements = material->GetNumberOfElements();
00115 
00116     if(verbose > 0) 
00117       {
00118         G4cout << "eIonisation CS for " << mLocal << "th material "
00119                << material->GetName()
00120                << "  eEl= " << nElements << G4endl;
00121       }
00122 
00123     G4double tcut  = (*energyCuts)[mLocal];
00124 
00125     G4VDataSetAlgorithm* algo = interp->Clone();
00126     G4VEMDataSet* setForMat = new G4CompositeEMDataSet(algo,1.,1.);
00127 
00128     for (G4int i=0; i<nElements; i++) {
00129 
00130       G4int Z = (G4int) (*elementVector)[i]->GetZ();
00131       G4int nShells = NumberOfComponents(Z);
00132 
00133       energies = new G4DataVector;
00134       cs       = new G4DataVector;
00135 
00136       log_energies = new G4DataVector;
00137       log_cs       = new G4DataVector;
00138 
00139       G4double density = nAtomsPerVolume[i];
00140 
00141       for (G4int bin=0; bin<nOfBins; bin++) {
00142 
00143         G4double e = energyVector[bin];
00144         energies->push_back(e);
00145         log_energies->push_back(std::log10(e));
00146         G4double value = 0.0;
00147         G4double log_value = -300;
00148 
00149         if(e > tcut) {
00150           for (G4int n=0; n<nShells; n++) {
00151             G4double cross = FindValue(Z, e, n);
00152             G4double p = theParam->Probability(Z, tcut, e, e, n);
00153             value += cross * p * density;
00154 
00155             if(verbose>0 && mLocal == 0 && e>=1. && e<=0.) 
00156             {
00157               G4cout << "G4eIonCrossSH: e(MeV)= " << e/MeV
00158                      << " n= " << n
00159                      << " cross= " << cross
00160                      << " p= " << p
00161                      << " value= " << value
00162                      << " tcut(MeV)= " << tcut/MeV
00163                      << " rho= " << density
00164                      << " Z= " << Z
00165                      << G4endl;
00166               }
00167 
00168           }
00169           if (value == 0.) value = 1e-300;
00170           log_value = std::log10(value);
00171         }
00172         cs->push_back(value);
00173         log_cs->push_back(log_value);
00174       }
00175       G4VDataSetAlgorithm* algoLocal = interp->Clone();
00176 
00177       //G4VEMDataSet* elSet = new G4EMDataSet(i,energies,cs,algoLocal,1.,1.);
00178 
00179       G4VEMDataSet* elSet = new G4EMDataSet(i,energies,cs,log_energies,log_cs,algoLocal,1.,1.);
00180 
00181       setForMat->AddComponent(elSet);
00182     }
00183     set->push_back(setForMat);
00184   }
00185 
00186   return set;
00187 }

G4double G4eIonisationCrossSectionHandler::GetCrossSectionAboveThresholdForElement ( G4double  energy,
G4double  cutEnergy,
G4int  Z 
)

Definition at line 189 of file G4eIonisationCrossSectionHandler.cc.

References G4VCrossSectionHandler::FindValue(), CLHEP::detail::n, G4VCrossSectionHandler::NumberOfComponents(), and G4VEnergySpectrum::Probability().

Referenced by G4LivermoreIonisationModel::ComputeCrossSectionPerAtom().

00192 {
00193   G4int nShells = NumberOfComponents(Z);
00194   G4double value = 0.;
00195   if(energy > cutEnergy) 
00196     {
00197       for (G4int n=0; n<nShells; n++) {
00198         G4double cross = FindValue(Z, energy, n);
00199         G4double p = theParam->Probability(Z, cutEnergy, energy, energy, n);
00200         value += cross * p;
00201       }
00202     }
00203   return value;
00204 }


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