Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4VCrossSectionDataSet Class Reference

#include <G4VCrossSectionDataSet.hh>

+ Inheritance diagram for G4VCrossSectionDataSet:

Public Member Functions

 G4VCrossSectionDataSet (const G4String &nam="")
 
virtual ~G4VCrossSectionDataSet ()
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=0, const G4Material *mat=0)
 
G4double GetCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
G4double ComputeCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
virtual G4double GetIsoCrossSection (const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
 
virtual G4IsotopeSelectIsotope (const G4Element *, G4double kinEnergy)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
virtual void CrossSectionDescription (std::ostream &) const
 
void SetVerboseLevel (G4int value)
 
G4double GetMinKinEnergy () const
 
void SetMinKinEnergy (G4double value)
 
G4double GetMaxKinEnergy () const
 
void SetMaxKinEnergy (G4double value)
 
const G4StringGetName () const
 

Protected Member Functions

void SetName (const G4String &)
 

Protected Attributes

G4int verboseLevel
 

Detailed Description

Definition at line 71 of file G4VCrossSectionDataSet.hh.

Constructor & Destructor Documentation

◆ G4VCrossSectionDataSet()

G4VCrossSectionDataSet::G4VCrossSectionDataSet ( const G4String nam = "")

Definition at line 54 of file G4VCrossSectionDataSet.cc.

54 :
55 verboseLevel(0),minKinEnergy(0.0),maxKinEnergy(100*TeV),name(nam)
56{
58}
void Register(G4VCrossSectionDataSet *)
static G4CrossSectionDataSetRegistry * Instance()

◆ ~G4VCrossSectionDataSet()

G4VCrossSectionDataSet::~G4VCrossSectionDataSet ( )
virtual

Definition at line 60 of file G4VCrossSectionDataSet.cc.

Member Function Documentation

◆ BuildPhysicsTable()

◆ ComputeCrossSection()

G4double G4VCrossSectionDataSet::ComputeCrossSection ( const G4DynamicParticle part,
const G4Element elm,
const G4Material mat = 0 
)

Definition at line 83 of file G4VCrossSectionDataSet.cc.

86{
87 G4int Z = G4lrint(elm->GetZ());
88
89 if (IsElementApplicable(part, Z, mat)) {
90 return GetElementCrossSection(part, Z, mat);
91 }
92
93 // isotope-wise cross section making sum over available
94 // isotope cross sections, which may be incomplete, so
95 // the result is corrected
96 G4int nIso = elm->GetNumberOfIsotopes();
97 G4double fact = 0.0;
98 G4double xsec = 0.0;
99 G4Isotope* iso = 0;
100
101 if (0 < nIso) {
102
103 // user-defined isotope abundances
104 G4IsotopeVector* isoVector = elm->GetIsotopeVector();
105 G4double* abundVector = elm->GetRelativeAbundanceVector();
106
107 for (G4int j = 0; j<nIso; ++j) {
108 iso = (*isoVector)[j];
109 G4int A = iso->GetN();
110 if(abundVector[j] > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) {
111 fact += abundVector[j];
112 xsec += abundVector[j]*GetIsoCrossSection(part, Z, A, iso, elm, mat);
113 }
114 }
115
116 } else {
117
118 // natural isotope abundances
120 G4int n0 = nist->GetNistFirstIsotopeN(Z);
122 for (G4int A = n0; A < n0+nn; ++A) {
123 G4double abund = nist->GetIsotopeAbundance(Z, A);
124 if(abund > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) {
125 fact += abund;
126 xsec += abund*GetIsoCrossSection(part, Z, A, iso, elm, mat);
127 }
128 }
129 }
130 if(fact > 0.0) { xsec /= fact; }
131 return xsec;
132}
std::vector< G4Isotope * > G4IsotopeVector
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
G4double * GetRelativeAbundanceVector() const
Definition: G4Element.hh:166
G4double GetZ() const
Definition: G4Element.hh:131
size_t GetNumberOfIsotopes() const
Definition: G4Element.hh:158
G4IsotopeVector * GetIsotopeVector() const
Definition: G4Element.hh:162
G4int GetN() const
Definition: G4Isotope.hh:94
G4int GetNumberOfNistIsotopes(G4int Z) const
G4int GetNistFirstIsotopeN(G4int Z) const
static G4NistManager * Instance()
G4double GetIsotopeAbundance(G4int Z, G4int N) const
virtual G4bool IsIsoApplicable(const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=0, const G4Material *mat=0)
virtual G4double GetElementCrossSection(const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
virtual G4bool IsElementApplicable(const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
virtual G4double GetIsoCrossSection(const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
int G4lrint(double ad)
Definition: templates.hh:163

Referenced by GetCrossSection().

◆ CrossSectionDescription()

◆ DumpPhysicsTable()

◆ GetCrossSection()

G4double G4VCrossSectionDataSet::GetCrossSection ( const G4DynamicParticle dp,
const G4Element elm,
const G4Material mat = 0 
)
inline

Definition at line 179 of file G4VCrossSectionDataSet.hh.

182{
183 return ComputeCrossSection(dp, elm, mat);
184}
G4double ComputeCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)

Referenced by G4ElectroNuclearReaction::ApplyYourself().

◆ GetElementCrossSection()

G4double G4VCrossSectionDataSet::GetElementCrossSection ( const G4DynamicParticle dynPart,
G4int  Z,
const G4Material mat = 0 
)
virtual

Reimplemented in G4EMDissociationCrossSection, G4GGNuclNuclCrossSection, G4IonsKoxCrossSection, G4IonsShenCrossSection, G4IonsSihverCrossSection, G4NeutronInelasticCrossSection, G4ProtonInelasticCrossSection, G4TripathiCrossSection, G4BGGNucleonElasticXS, G4BGGNucleonInelasticXS, G4BGGPionElasticXS, G4BGGPionInelasticXS, G4CrossSectionElastic, G4CrossSectionInelastic, G4CrossSectionPairGG, G4NeutronCaptureXS, G4NeutronElasticXS, G4NeutronInelasticXS, G4IonProtonCrossSection, G4HadronCaptureDataSet, G4HadronElasticDataSet, G4HadronFissionDataSet, G4HadronInelasticDataSet, G4NucleonNuclearCrossSection, G4KokoulinMuonNuclearXS, G4PiNuclearCrossSection, G4GeneralSpaceNNCrossSection, and G4TripathiLightCrossSection.

Definition at line 135 of file G4VCrossSectionDataSet.cc.

138{
139 G4cout << "G4VCrossSectionDataSet::GetCrossSection per element ERROR: "
140 << " there is no cross section for "
141 << dynPart->GetDefinition()->GetParticleName()
142 << " E(MeV)= " << dynPart->GetKineticEnergy()/MeV;
143 if(mat) { G4cout << " inside " << mat->GetName(); }
144 G4cout << " for Z= " << Z << G4endl;
145 throw G4HadronicException(__FILE__, __LINE__,
146 "G4VCrossSectionDataSet::GetElementCrossSection is absent");
147 return 0.0;
148}
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4String & GetName() const
Definition: G4Material.hh:177
const G4String & GetParticleName() const

Referenced by G4CrossSectionPairGG::BuildPhysicsTable(), ComputeCrossSection(), and G4CrossSectionPairGG::GetElementCrossSection().

◆ GetIsoCrossSection()

G4double G4VCrossSectionDataSet::GetIsoCrossSection ( const G4DynamicParticle dynPart,
G4int  Z,
G4int  A,
const G4Isotope iso = 0,
const G4Element elm = 0,
const G4Material mat = 0 
)
virtual

Reimplemented in G4ChipsAntiBaryonElasticXS, G4ChipsAntiBaryonInelasticXS, G4ChipsHyperonElasticXS, G4ChipsHyperonInelasticXS, G4ChipsKaonMinusElasticXS, G4ChipsKaonMinusInelasticXS, G4ChipsKaonPlusElasticXS, G4ChipsKaonPlusInelasticXS, G4ChipsKaonZeroElasticXS, G4ChipsKaonZeroInelasticXS, G4ChipsNeutronElasticXS, G4ChipsNeutronInelasticXS, G4ChipsPionMinusElasticXS, G4ChipsPionMinusInelasticXS, G4ChipsPionPlusElasticXS, G4ChipsPionPlusInelasticXS, G4ChipsProtonElasticXS, G4ChipsProtonInelasticXS, G4NeutronCaptureXS, G4BGGNucleonElasticXS, G4BGGNucleonInelasticXS, G4BGGPionElasticXS, G4BGGPionInelasticXS, G4GlauberGribovCrossSection, G4IonsShenCrossSection, G4PhotoNuclearCrossSection, G4CHIPSElasticXS, G4LENDCrossSection, G4NeutronHPCaptureData, G4NeutronHPElasticData, G4NeutronHPFissionData, G4NeutronHPInelasticData, G4NeutronHPorLCaptureData, G4NeutronHPorLEInelasticData, G4NeutronHPorLElasticData, G4NeutronHPorLFissionData, G4NeutronHPThermalScatteringData, G4ElectroNuclearCrossSection, G4QHadronElasticDataSet, and G4QHadronInelasticDataSet.

Definition at line 151 of file G4VCrossSectionDataSet.cc.

156{
157 G4cout << "G4VCrossSectionDataSet::GetCrossSection per isotope ERROR: "
158 << " there is no cross section for "
159 << dynPart->GetDefinition()->GetParticleName()
160 << " E(MeV)= " << dynPart->GetKineticEnergy()/MeV;
161 if(mat) { G4cout << " inside " << mat->GetName(); }
162 if(elm) { G4cout << " for " << elm->GetName(); }
163 G4cout << " Z= " << Z << " A= " << A << G4endl;
164 throw G4HadronicException(__FILE__, __LINE__,
165 "G4VCrossSectionDataSet::GetIsoCrossSection is absent");
166 return 0.0;
167}
const G4String & GetName() const
Definition: G4Element.hh:127

Referenced by G4QMDReaction::ApplyYourself(), and ComputeCrossSection().

◆ GetMaxKinEnergy()

◆ GetMinKinEnergy()

◆ GetName()

◆ IsElementApplicable()

◆ IsIsoApplicable()

◆ SelectIsotope()

G4Isotope * G4VCrossSectionDataSet::SelectIsotope ( const G4Element anElement,
G4double  kinEnergy 
)
virtual

Reimplemented in G4NeutronCaptureXS.

Definition at line 170 of file G4VCrossSectionDataSet.cc.

171{
172 G4int nIso = anElement->GetNumberOfIsotopes();
173 G4IsotopeVector* isoVector = anElement->GetIsotopeVector();
174 G4Isotope* iso = (*isoVector)[0];
175
176 // more than 1 isotope
177 if(1 < nIso) {
178 G4double* abundVector = anElement->GetRelativeAbundanceVector();
179 G4double sum = 0.0;
181 for (G4int j = 0; j<nIso; ++j) {
182 sum += abundVector[j];
183 if(q <= sum) {
184 iso = (*isoVector)[j];
185 break;
186 }
187 }
188 }
189 return iso;
190}
#define G4UniformRand()
Definition: Randomize.hh:53

◆ SetMaxKinEnergy()

◆ SetMinKinEnergy()

◆ SetName()

void G4VCrossSectionDataSet::SetName ( const G4String nam)
inlineprotected

Definition at line 216 of file G4VCrossSectionDataSet.hh.

217{
218 name = nam;
219}

◆ SetVerboseLevel()

void G4VCrossSectionDataSet::SetVerboseLevel ( G4int  value)
inline

Definition at line 186 of file G4VCrossSectionDataSet.hh.

187{
188 verboseLevel = value;
189}

Member Data Documentation

◆ verboseLevel

G4int G4VCrossSectionDataSet::verboseLevel
protected

Definition at line 165 of file G4VCrossSectionDataSet.hh.

Referenced by G4BGGNucleonElasticXS::BuildPhysicsTable(), G4BGGNucleonInelasticXS::BuildPhysicsTable(), G4BGGPionElasticXS::BuildPhysicsTable(), G4BGGPionInelasticXS::BuildPhysicsTable(), G4CrossSectionPairGG::BuildPhysicsTable(), G4NeutronCaptureXS::BuildPhysicsTable(), G4NeutronElasticXS::BuildPhysicsTable(), G4NeutronInelasticXS::BuildPhysicsTable(), G4LENDCrossSection::create_used_target_map(), G4BGGNucleonElasticXS::G4BGGNucleonElasticXS(), G4BGGNucleonInelasticXS::G4BGGNucleonInelasticXS(), G4BGGPionElasticXS::G4BGGPionElasticXS(), G4BGGPionInelasticXS::G4BGGPionInelasticXS(), G4CrossSectionPairGG::G4CrossSectionPairGG(), G4NeutronCaptureXS::G4NeutronCaptureXS(), G4NeutronElasticXS::G4NeutronElasticXS(), G4NeutronInelasticXS::G4NeutronInelasticXS(), G4BGGNucleonElasticXS::GetElementCrossSection(), G4BGGNucleonInelasticXS::GetElementCrossSection(), G4BGGPionElasticXS::GetElementCrossSection(), G4BGGPionInelasticXS::GetElementCrossSection(), G4CrossSectionPairGG::GetElementCrossSection(), G4NeutronCaptureXS::GetElementCrossSection(), G4NeutronElasticXS::GetElementCrossSection(), G4NeutronInelasticXS::GetElementCrossSection(), G4GeneralSpaceNNCrossSection::GetElementCrossSection(), G4NeutronCaptureXS::GetIsoCrossSection(), G4BGGNucleonElasticXS::GetIsoCrossSection(), G4BGGNucleonInelasticXS::GetIsoCrossSection(), G4BGGPionElasticXS::GetIsoCrossSection(), G4BGGPionInelasticXS::GetIsoCrossSection(), and SetVerboseLevel().


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