Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4RDeIonisationCrossSectionHandler Class Reference

#include <G4RDeIonisationCrossSectionHandler.hh>

Inheritance diagram for G4RDeIonisationCrossSectionHandler:
Collaboration diagram for G4RDeIonisationCrossSectionHandler:

Public Member Functions

 G4RDeIonisationCrossSectionHandler (const G4RDVEnergySpectrum *spec, G4RDVDataSetAlgorithm *alg, G4double emin, G4double emax, G4int nbin)
 
 ~G4RDeIonisationCrossSectionHandler ()
 
- Public Member Functions inherited from G4RDVCrossSectionHandler
 G4RDVCrossSectionHandler ()
 
 G4RDVCrossSectionHandler (G4RDVDataSetAlgorithm *interpolation, G4double minE=250 *CLHEP::eV, G4double maxE=100 *CLHEP::GeV, G4int nBins=200, G4double unitE=CLHEP::MeV, G4double unitData=CLHEP::barn, G4int minZ=1, G4int maxZ=99)
 
virtual ~G4RDVCrossSectionHandler ()
 
void Initialise (G4RDVDataSetAlgorithm *interpolation=0, G4double minE=250 *CLHEP::eV, G4double maxE=100 *CLHEP::GeV, G4int numberOfBins=200, G4double unitE=CLHEP::MeV, G4double unitData=CLHEP::barn, G4int minZ=1, G4int maxZ=99)
 
G4int SelectRandomAtom (const G4MaterialCutsCouple *couple, G4double e) const
 
const G4ElementSelectRandomElement (const G4MaterialCutsCouple *material, G4double e) const
 
G4int SelectRandomShell (G4int Z, G4double e) const
 
G4RDVEMDataSetBuildMeanFreePathForMaterials (const G4DataVector *energyCuts=0)
 
G4double FindValue (G4int Z, G4double e) const
 
G4double FindValue (G4int Z, G4double e, G4int shellIndex) const
 
G4double ValueForMaterial (const G4Material *material, G4double e) const
 
void LoadData (const G4String &dataFile)
 
void LoadShellData (const G4String &dataFile)
 
void PrintData () const
 
void Clear ()
 

Protected Member Functions

std::vector< G4RDVEMDataSet * > * BuildCrossSectionsForMaterials (const G4DataVector &energyVector, const G4DataVector *energyCuts)
 
- Protected Member Functions inherited from G4RDVCrossSectionHandler
G4int NumberOfComponents (G4int Z) const
 
void ActiveElements ()
 
virtual G4RDVDataSetAlgorithmCreateInterpolation ()
 
const G4RDVDataSetAlgorithmGetInterpolation () const
 

Detailed Description

Definition at line 63 of file G4RDeIonisationCrossSectionHandler.hh.

Constructor & Destructor Documentation

G4RDeIonisationCrossSectionHandler::G4RDeIonisationCrossSectionHandler ( const G4RDVEnergySpectrum spec,
G4RDVDataSetAlgorithm alg,
G4double  emin,
G4double  emax,
G4int  nbin 
)

Definition at line 62 of file G4RDeIonisationCrossSectionHandler.cc.

66  theParam(spec)
67 {
69  interp = new G4RDLinLogLogInterpolation();
70 }
void Initialise(G4RDVDataSetAlgorithm *interpolation=0, G4double minE=250 *CLHEP::eV, G4double maxE=100 *CLHEP::GeV, G4int numberOfBins=200, G4double unitE=CLHEP::MeV, G4double unitData=CLHEP::barn, G4int minZ=1, G4int maxZ=99)
static const G4double emax

Here is the call graph for this function:

G4RDeIonisationCrossSectionHandler::~G4RDeIonisationCrossSectionHandler ( )

Definition at line 73 of file G4RDeIonisationCrossSectionHandler.cc.

74 {
75  delete interp;
76 }

Member Function Documentation

std::vector< G4RDVEMDataSet * > * G4RDeIonisationCrossSectionHandler::BuildCrossSectionsForMaterials ( const G4DataVector energyVector,
const G4DataVector energyCuts 
)
protectedvirtual

Implements G4RDVCrossSectionHandler.

Definition at line 79 of file G4RDeIonisationCrossSectionHandler.cc.

82 {
83  G4int verbose = 0;
84  std::vector<G4RDVEMDataSet*>* set = new std::vector<G4RDVEMDataSet*>;
85 
86  G4DataVector* energies;
87  G4DataVector* cs;
88  G4int nOfBins = energyVector.size();
89 
90  const G4ProductionCutsTable* theCoupleTable=
92  size_t numOfCouples = theCoupleTable->GetTableSize();
93 
94  for (size_t m=0; m<numOfCouples; m++) {
95 
96  const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(m);
97  const G4Material* material= couple->GetMaterial();
98  const G4ElementVector* elementVector = material->GetElementVector();
99  const G4double* nAtomsPerVolume = material->GetAtomicNumDensityVector();
100  G4int nElements = material->GetNumberOfElements();
101 
102  if(verbose > 0) {
103  G4cout << "eIonisation CS for " << m << "th material "
104  << material->GetName()
105  << " eEl= " << nElements << G4endl;
106  }
107 
108  G4double tcut = (*energyCuts)[m];
109 
110  G4RDVDataSetAlgorithm* algo = interp->Clone();
111  G4RDVEMDataSet* setForMat = new G4RDCompositeEMDataSet(algo,1.,1.);
112 
113  for (G4int i=0; i<nElements; i++) {
114 
115  G4int Z = (G4int) (*elementVector)[i]->GetZ();
116  G4int nShells = NumberOfComponents(Z);
117  energies = new G4DataVector;
118  cs = new G4DataVector;
119  G4double density = nAtomsPerVolume[i];
120 
121  for (G4int bin=0; bin<nOfBins; bin++) {
122 
123  G4double e = energyVector[bin];
124  energies->push_back(e);
125  G4double value = 0.0;
126 
127  if(e > tcut) {
128  for (G4int n=0; n<nShells; n++) {
129  G4double cross = FindValue(Z, e, n);
130  G4double p = theParam->Probability(Z, tcut, e, e, n);
131  value += cross * p * density;
132 
133  if(verbose>0 && m == 0 && e>=1. && e<=0.) {
134  G4cout << "G4eIonCrossSH: e(MeV)= " << e/MeV
135  << " n= " << n
136  << " cross= " << cross
137  << " p= " << p
138  << " value= " << value
139  << " tcut(MeV)= " << tcut/MeV
140  << " rho= " << density
141  << " Z= " << Z
142  << G4endl;
143  }
144 
145  }
146  }
147  cs->push_back(value);
148  }
149  G4RDVDataSetAlgorithm* algo = interp->Clone();
150  G4RDVEMDataSet* elSet = new G4RDEMDataSet(i,energies,cs,algo,1.,1.);
151  setForMat->AddComponent(elSet);
152  }
153  set->push_back(setForMat);
154  }
155 
156  return set;
157 }
virtual G4RDVDataSetAlgorithm * Clone() const =0
std::vector< G4Element * > G4ElementVector
tuple bin
Definition: plottest35.py:22
const char * p
Definition: xmltok.h:285
const G4String & GetName() const
Definition: G4Material.hh:178
G4double FindValue(G4int Z, G4double e) const
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
string material
Definition: eplot.py:19
G4GLOB_DLL std::ostream G4cout
static constexpr double m
Definition: G4SIunits.hh:129
const XML_Char int const XML_Char * value
Definition: expat.h:331
const G4int n
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:216
static G4ProductionCutsTable * GetProductionCutsTable()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
G4int NumberOfComponents(G4int Z) const
double G4double
Definition: G4Types.hh:76
virtual void AddComponent(G4RDVEMDataSet *dataSet)=0
const G4Material * GetMaterial() const
virtual G4double Probability(G4int Z, G4double minKineticEnergy, G4double maxKineticEnergy, G4double kineticEnergy, G4int shell=0, const G4ParticleDefinition *pd=0) const =0

Here is the call graph for this function:


The documentation for this class was generated from the following files: