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

#include <G4RDBremsstrahlungCrossSectionHandler.hh>

Inheritance diagram for G4RDBremsstrahlungCrossSectionHandler:
Collaboration diagram for G4RDBremsstrahlungCrossSectionHandler:

Public Member Functions

 G4RDBremsstrahlungCrossSectionHandler (const G4RDVEnergySpectrum *spectrum, G4RDVDataSetAlgorithm *interpolation)
 
 ~G4RDBremsstrahlungCrossSectionHandler ()
 
- 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 62 of file G4RDBremsstrahlungCrossSectionHandler.hh.

Constructor & Destructor Documentation

G4RDBremsstrahlungCrossSectionHandler::G4RDBremsstrahlungCrossSectionHandler ( const G4RDVEnergySpectrum spectrum,
G4RDVDataSetAlgorithm interpolation 
)

Definition at line 57 of file G4RDBremsstrahlungCrossSectionHandler.cc.

59  : theBR(spec)
60 {
61  interp = new G4RDSemiLogInterpolation();
62 }
G4RDBremsstrahlungCrossSectionHandler::~G4RDBremsstrahlungCrossSectionHandler ( )

Definition at line 65 of file G4RDBremsstrahlungCrossSectionHandler.cc.

66 {
67  delete interp;
68 }

Member Function Documentation

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

Implements G4RDVCrossSectionHandler.

Definition at line 72 of file G4RDBremsstrahlungCrossSectionHandler.cc.

74 {
75  std::vector<G4RDVEMDataSet*>* set = new std::vector<G4RDVEMDataSet*>;
76 
77  G4DataVector* energies;
78  G4DataVector* cs;
79  G4int nOfBins = energyVector.size();
80 
81  const G4ProductionCutsTable* theCoupleTable=
83  size_t numOfCouples = theCoupleTable->GetTableSize();
84 
85  for (size_t m=0; m<numOfCouples; m++) {
86 
87  const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(m);
88  const G4Material* material= couple->GetMaterial();
89  const G4ElementVector* elementVector = material->GetElementVector();
90  const G4double* nAtomsPerVolume = material->GetVecNbOfAtomsPerVolume();
91  G4int nElements = material->GetNumberOfElements();
92 
93  G4double tcut = (*energyCuts)[m];
94 
95  G4RDVDataSetAlgorithm* algo = interp->Clone();
96  G4RDVEMDataSet* setForMat = new G4RDCompositeEMDataSet(algo,1.,1.);
97 
98  for (G4int i=0; i<nElements; i++) {
99 
100  G4int Z = (G4int) ((*elementVector)[i]->GetZ());
101  energies = new G4DataVector;
102  cs = new G4DataVector;
103  G4double density = nAtomsPerVolume[i];
104 
105  for (G4int bin=0; bin<nOfBins; bin++) {
106 
107  G4double e = energyVector[bin];
108  energies->push_back(e);
109  G4double value = 0.0;
110 
111  if(e > tcut) {
112  G4double elemCs = FindValue(Z, e);
113 
114  value = theBR->Probability(Z, tcut, e, e);
115 
116  value *= elemCs*density;
117  }
118  cs->push_back(value);
119  }
120  G4RDVDataSetAlgorithm* algol = interp->Clone();
121  G4RDVEMDataSet* elSet = new G4RDEMDataSet(i,energies,cs,algol,1.,1.);
122  setForMat->AddComponent(elSet);
123  }
124  set->push_back(setForMat);
125  }
126 
127  return set;
128 }
virtual G4RDVDataSetAlgorithm * Clone() const =0
std::vector< G4Element * > G4ElementVector
tuple bin
Definition: plottest35.py:22
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
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:206
static constexpr double m
Definition: G4SIunits.hh:129
const XML_Char int const XML_Char * value
Definition: expat.h:331
static G4ProductionCutsTable * GetProductionCutsTable()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
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: