Geant4 10.7.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=nullptr)
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=nullptr, const G4Material *mat=nullptr)
 
G4double GetCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=nullptr)
 
G4double ComputeCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=nullptr)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *mat=nullptr)
 
virtual G4double GetIsoCrossSection (const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=nullptr, const G4Element *elm=nullptr, const G4Material *mat=nullptr)
 
virtual const G4IsotopeSelectIsotope (const G4Element *, G4double kinEnergy, G4double logE)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
virtual void CrossSectionDescription (std::ostream &) const
 
virtual G4int GetVerboseLevel () const
 
virtual void SetVerboseLevel (G4int value)
 
G4double GetMinKinEnergy () const
 
void SetMinKinEnergy (G4double value)
 
G4double GetMaxKinEnergy () const
 
void SetMaxKinEnergy (G4double value)
 
bool ForAllAtomsAndEnergies () const
 
void SetForAllAtomsAndEnergies (G4bool val)
 
const G4StringGetName () const
 

Protected Member Functions

void SetName (const G4String &)
 

Protected Attributes

G4int verboseLevel
 

Detailed Description

Definition at line 69 of file G4VCrossSectionDataSet.hh.

Constructor & Destructor Documentation

◆ G4VCrossSectionDataSet()

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

Definition at line 49 of file G4VCrossSectionDataSet.cc.

49 :
50 verboseLevel(0),minKinEnergy(0.0),
51 maxKinEnergy(G4HadronicParameters::Instance()->GetMaxEnergy()),
52 isForAllAtomsAndEnergies(false),name(nam)
53{
55 registry->Register(this);
56}
void Register(G4VCrossSectionDataSet *)
static G4CrossSectionDataSetRegistry * Instance()
static G4HadronicParameters * Instance()

◆ ~G4VCrossSectionDataSet()

G4VCrossSectionDataSet::~G4VCrossSectionDataSet ( )
virtual

Definition at line 58 of file G4VCrossSectionDataSet.cc.

59{
60 registry->DeRegister(this);
61}
void DeRegister(G4VCrossSectionDataSet *)

Member Function Documentation

◆ BuildPhysicsTable()

◆ ComputeCrossSection()

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

Definition at line 81 of file G4VCrossSectionDataSet.cc.

84{
85 G4int Z = elm->GetZasInt();
86
87 if (IsElementApplicable(part, Z, mat)) {
88 return GetElementCrossSection(part, Z, mat);
89 }
90
91 // isotope-wise cross section making sum over available
92 // isotope cross sections, which may be incomplete, so
93 // the result is corrected
94 size_t nIso = elm->GetNumberOfIsotopes();
95 G4double fact = 0.0;
96 G4double xsec = 0.0;
97
98 // user-defined isotope abundances
99 const G4IsotopeVector* isoVector = elm->GetIsotopeVector();
100 const G4double* abundVector = elm->GetRelativeAbundanceVector();
101
102 for (size_t j=0; j<nIso; ++j) {
103 const G4Isotope* iso = (*isoVector)[j];
104 G4int A = iso->GetN();
105 if(abundVector[j] > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) {
106 fact += abundVector[j];
107 xsec += abundVector[j]*GetIsoCrossSection(part, Z, A, iso, elm, mat);
108 }
109 }
110 return (fact > 0.0) ? xsec/fact : 0.0;
111}
double A(double temperature)
std::vector< G4Isotope * > G4IsotopeVector
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
G4double * GetRelativeAbundanceVector() const
Definition: G4Element.hh:166
size_t GetNumberOfIsotopes() const
Definition: G4Element.hh:158
G4int GetZasInt() const
Definition: G4Element.hh:131
G4IsotopeVector * GetIsotopeVector() const
Definition: G4Element.hh:162
G4int GetN() const
Definition: G4Isotope.hh:93
virtual G4double GetElementCrossSection(const G4DynamicParticle *, G4int Z, const G4Material *mat=nullptr)
virtual G4double GetIsoCrossSection(const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=nullptr, const G4Element *elm=nullptr, const G4Material *mat=nullptr)
virtual G4bool IsIsoApplicable(const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=nullptr, const G4Material *mat=nullptr)
virtual G4bool IsElementApplicable(const G4DynamicParticle *, G4int Z, const G4Material *mat=nullptr)

Referenced by GetCrossSection().

◆ CrossSectionDescription()

void G4VCrossSectionDataSet::CrossSectionDescription ( std::ostream &  outFile) const
virtual

Reimplemented in G4ChipsAntiBaryonElasticXS, G4ChipsAntiBaryonInelasticXS, G4ChipsHyperonElasticXS, G4ChipsHyperonInelasticXS, G4ChipsKaonMinusElasticXS, G4ChipsKaonMinusInelasticXS, G4ChipsKaonPlusElasticXS, G4ChipsKaonPlusInelasticXS, G4ChipsKaonZeroElasticXS, G4ChipsKaonZeroInelasticXS, G4ChipsNeutronElasticXS, G4ChipsNeutronInelasticXS, G4ChipsPionMinusElasticXS, G4ChipsPionMinusInelasticXS, G4ChipsPionPlusElasticXS, G4ChipsPionPlusInelasticXS, G4ChipsProtonElasticXS, G4ChipsProtonInelasticXS, G4CrossSectionElastic, G4CrossSectionInelastic, G4CrossSectionPairGG, G4ElectroNuclearCrossSection, G4HadronCaptureDataSet, G4HadronElasticDataSet, G4HadronFissionDataSet, G4HadronInelasticDataSet, G4IonsKoxCrossSection, G4IonsShenCrossSection, G4IonsSihverCrossSection, G4KokoulinMuonNuclearXS, G4NeutronInelasticCrossSection, G4PhotoNuclearCrossSection, G4PiNuclearCrossSection, G4ZeroXS, G4ParticleHPCaptureData, G4ParticleHPElasticData, G4ParticleHPFissionData, G4ParticleHPInelasticData, G4ParticleHPThermalScatteringData, G4BGGNucleonElasticXS, G4BGGPionElasticXS, G4BGGPionInelasticXS, G4GammaNuclearXS, G4NeutronCaptureXS, G4NeutronElasticXS, G4NeutronInelasticXS, G4NucleonNuclearCrossSection, G4ParticleInelasticXS, G4UPiNuclearCrossSection, G4BGGNucleonInelasticXS, G4IonProtonCrossSection, and G4GeneralSpaceNNCrossSection.

Definition at line 177 of file G4VCrossSectionDataSet.cc.

178{
179 outFile << "The description for this cross section data set has not been written yet.\n";
180}

Referenced by G4CrossSectionDataStore::PrintCrossSectionHtml().

◆ DumpPhysicsTable()

◆ ForAllAtomsAndEnergies()

bool G4VCrossSectionDataSet::ForAllAtomsAndEnergies ( ) const
inline

Definition at line 231 of file G4VCrossSectionDataSet.hh.

232{
233 return isForAllAtomsAndEnergies;
234}

Referenced by G4CrossSectionDataStore::AddDataSet().

◆ GetCrossSection()

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

Definition at line 188 of file G4VCrossSectionDataSet.hh.

191{
192 return ComputeCrossSection(dp, elm, mat);
193}
G4double ComputeCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *mat=nullptr)

◆ GetElementCrossSection()

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

Reimplemented in G4EMDissociationCrossSection, G4IonsKoxCrossSection, G4IonsShenCrossSection, G4IonsSihverCrossSection, G4NeutrinoElectronCcXsc, G4NeutrinoElectronNcXsc, G4NeutrinoElectronTotXsc, G4NeutronElectronElXsc, G4NeutronInelasticCrossSection, G4PhotoNuclearCrossSection, G4ProtonInelasticCrossSection, G4TripathiCrossSection, G4NeutronCaptureXS, G4NeutronElasticXS, G4NeutronInelasticXS, G4ElectroNuclearCrossSection, G4BGGNucleonElasticXS, G4BGGPionElasticXS, G4BGGPionInelasticXS, G4GammaNuclearXS, G4BGGNucleonInelasticXS, G4CrossSectionElastic, G4CrossSectionInelastic, G4CrossSectionPairGG, G4ParticleInelasticXS, G4ZeroXS, G4IonProtonCrossSection, G4HadronCaptureDataSet, G4HadronElasticDataSet, G4HadronFissionDataSet, G4HadronInelasticDataSet, G4NucleonNuclearCrossSection, G4MuNeutrinoNucleusTotXsc, G4KokoulinMuonNuclearXS, G4PiNuclearCrossSection, G4GeneralSpaceNNCrossSection, and G4TripathiLightCrossSection.

Definition at line 114 of file G4VCrossSectionDataSet.cc.

117{
119 ed << "GetElementCrossSection is not implemented in <" << name << ">\n"
120 << "Particle: " << dynPart->GetDefinition()->GetParticleName()
121 << " Ekin(MeV)= " << dynPart->GetKineticEnergy()/MeV;
122 if(mat) { ed << " material: " << mat->GetName(); }
123 ed << " target Z= " << Z << G4endl;
124 G4Exception("G4VCrossSectionDataSet::GetElementCrossSection", "had001",
125 FatalException, ed);
126 return 0.0;
127}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
#define G4endl
Definition: G4ios.hh:57
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4String & GetName() const
Definition: G4Material.hh:175
const G4String & GetParticleName() const

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

◆ GetIsoCrossSection()

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

Reimplemented in G4ChipsAntiBaryonElasticXS, G4ChipsAntiBaryonInelasticXS, G4ChipsHyperonElasticXS, G4ChipsHyperonInelasticXS, G4ChipsKaonMinusElasticXS, G4ChipsKaonMinusInelasticXS, G4ChipsKaonPlusElasticXS, G4ChipsKaonPlusInelasticXS, G4ChipsKaonZeroElasticXS, G4ChipsKaonZeroInelasticXS, G4ChipsNeutronElasticXS, G4ChipsNeutronInelasticXS, G4ChipsPionMinusElasticXS, G4ChipsPionMinusInelasticXS, G4ChipsPionPlusElasticXS, G4ChipsPionPlusInelasticXS, G4ChipsProtonElasticXS, G4ChipsProtonInelasticXS, G4GammaNuclearXS, G4NeutronCaptureXS, G4NeutronElasticXS, G4NeutronInelasticXS, G4IonsShenCrossSection, G4BGGNucleonElasticXS, G4BGGPionElasticXS, G4BGGPionInelasticXS, G4ParticleInelasticXS, G4BGGNucleonInelasticXS, G4IonProtonCrossSection, G4LENDCombinedCrossSection, G4LENDCrossSection, G4LENDGammaCrossSection, G4ParticleHPCaptureData, G4ParticleHPElasticData, G4ParticleHPFissionData, G4ParticleHPInelasticData, G4ParticleHPThermalScatteringData, G4MuNeutrinoNucleusTotXsc, G4ElNeutrinoNucleusTotXsc, and G4PhotoNuclearCrossSection.

Definition at line 130 of file G4VCrossSectionDataSet.cc.

135{
137 ed << "GetIsoCrossSection is not implemented in <" << name << ">\n"
138 << "Particle: " << dynPart->GetDefinition()->GetParticleName()
139 << " Ekin(MeV)= " << dynPart->GetKineticEnergy()/MeV;
140 if(mat) { ed << " material: " << mat->GetName(); }
141 if(elm) { ed << " element: " << elm->GetName(); }
142 ed << " target Z= " << Z << " A= " << A << G4endl;
143 G4Exception("G4VCrossSectionDataSet::GetIsoCrossSection", "had001",
144 FatalException, ed);
145 return 0.0;
146}
const G4String & GetName() const
Definition: G4Element.hh:126

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

◆ GetMaxKinEnergy()

◆ GetMinKinEnergy()

◆ GetName()

◆ GetVerboseLevel()

G4int G4VCrossSectionDataSet::GetVerboseLevel ( ) const
inlinevirtual

◆ IsElementApplicable()

◆ IsIsoApplicable()

◆ SelectIsotope()

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

Reimplemented in G4GammaNuclearXS, G4NeutronCaptureXS, G4NeutronElasticXS, G4NeutronInelasticXS, and G4ParticleInelasticXS.

Definition at line 149 of file G4VCrossSectionDataSet.cc.

151{
152 size_t nIso = anElement->GetNumberOfIsotopes();
153 const G4Isotope* iso = anElement->GetIsotope(0);
154
155 // more than 1 isotope
156 if(1 < nIso) {
157 const G4double* abundVector = anElement->GetRelativeAbundanceVector();
158 G4double sum = 0.0;
160 for (size_t j=0; j<nIso; ++j) {
161 sum += abundVector[j];
162 if(q <= sum) {
163 iso = anElement->GetIsotope(j);
164 break;
165 }
166 }
167 }
168 return iso;
169}
#define G4UniformRand()
Definition: Randomize.hh:52
const G4Isotope * GetIsotope(G4int iso) const
Definition: G4Element.hh:169

◆ SetForAllAtomsAndEnergies()

◆ SetMaxKinEnergy()

◆ SetMinKinEnergy()

◆ SetName()

void G4VCrossSectionDataSet::SetName ( const G4String nam)
inlineprotected

Definition at line 241 of file G4VCrossSectionDataSet.hh.

242{
243 name = nam;
244}

Referenced by G4ParticleHPInelasticData::G4ParticleHPInelasticData().

◆ SetVerboseLevel()

void G4VCrossSectionDataSet::SetVerboseLevel ( G4int  value)
inlinevirtual

Member Data Documentation

◆ verboseLevel

G4int G4VCrossSectionDataSet::verboseLevel
protected

Definition at line 170 of file G4VCrossSectionDataSet.hh.

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


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