Geant4 11.2.2
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 ComputeCrossSectionPerElement (G4double kinEnergy, G4double loge, const G4ParticleDefinition *, 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 G4double ComputeIsoCrossSection (G4double kinEnergy, G4double loge, const G4ParticleDefinition *, 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 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
 
void SetName (const G4String &nam)
 
G4VCrossSectionDataSetoperator= (const G4VCrossSectionDataSet &right)=delete
 
 G4VCrossSectionDataSet (const G4VCrossSectionDataSet &)=delete
 

Protected Attributes

G4int verboseLevel
 
G4String name
 

Detailed Description

Definition at line 69 of file G4VCrossSectionDataSet.hh.

Constructor & Destructor Documentation

◆ G4VCrossSectionDataSet() [1/2]

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

Definition at line 49 of file G4VCrossSectionDataSet.cc.

49 :
50 verboseLevel(0),name(nam),minKinEnergy(0.0),
51 maxKinEnergy(G4HadronicParameters::Instance()->GetMaxEnergy()),
52 isForAllAtomsAndEnergies(false)
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 *)

◆ G4VCrossSectionDataSet() [2/2]

G4VCrossSectionDataSet::G4VCrossSectionDataSet ( const G4VCrossSectionDataSet & )
delete

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 std::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 for (std::size_t j=0; j<nIso; ++j) {
102 const G4Isotope* iso = (*isoVector)[j];
103 G4int A = iso->GetN();
104 if(abundVector[j] > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) {
105 fact += abundVector[j];
106 xsec += abundVector[j]*GetIsoCrossSection(part, Z, A, iso, elm, mat);
107 }
108 }
109 return (fact > 0.0) ? xsec/fact : 0.0;
110}
std::vector< G4Isotope * > G4IsotopeVector
double G4double
Definition G4Types.hh:83
int G4int
Definition G4Types.hh:85
const G4double A[17]
G4double * GetRelativeAbundanceVector() const
Definition G4Element.hh:149
size_t GetNumberOfIsotopes() const
Definition G4Element.hh:143
G4int GetZasInt() const
Definition G4Element.hh:120
G4IsotopeVector * GetIsotopeVector() const
Definition G4Element.hh:146
G4int GetN() const
Definition G4Isotope.hh:83
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().

◆ ComputeCrossSectionPerElement()

G4double G4VCrossSectionDataSet::ComputeCrossSectionPerElement ( G4double kinEnergy,
G4double loge,
const G4ParticleDefinition * pd,
const G4Element * elm,
const G4Material * mat = nullptr )
virtual

Reimplemented in G4NeutronCaptureXS, G4NeutronElasticXS, G4NeutronInelasticXS, and G4ParticleInelasticXS.

Definition at line 113 of file G4VCrossSectionDataSet.cc.

117{
118 G4int Z = elm->GetZasInt();
119 std::size_t nIso = elm->GetNumberOfIsotopes();
120 G4double xsec = 0.0;
121 const G4IsotopeVector* isoVector = elm->GetIsotopeVector();
122 const G4double* abundVector = elm->GetRelativeAbundanceVector();
123 for (std::size_t j=0; j<nIso; ++j) {
124 const G4Isotope* iso = (*isoVector)[j];
125 G4int A = iso->GetN();
126 xsec += abundVector[j]*
127 ComputeIsoCrossSection(kinEnergy, loge, pd, Z, A, iso, elm, mat);
128 }
129 return xsec;
130}
virtual G4double ComputeIsoCrossSection(G4double kinEnergy, G4double loge, const G4ParticleDefinition *, G4int Z, G4int A, const G4Isotope *iso=nullptr, const G4Element *elm=nullptr, const G4Material *mat=nullptr)

◆ ComputeIsoCrossSection()

G4double G4VCrossSectionDataSet::ComputeIsoCrossSection ( G4double kinEnergy,
G4double loge,
const G4ParticleDefinition * pd,
G4int Z,
G4int A,
const G4Isotope * iso = nullptr,
const G4Element * elm = nullptr,
const G4Material * mat = nullptr )
virtual

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

Definition at line 168 of file G4VCrossSectionDataSet.cc.

174{
176 ed << "GetIsoCrossSection is not implemented in <" << name << ">\n"
177 << "Particle: " << pd->GetParticleName()
178 << " Ekin(MeV)= " << kinEnergy/CLHEP::MeV;
179 if(nullptr != mat) { ed << " material: " << mat->GetName(); }
180 if(nullptr != elm) { ed << " element: " << elm->GetName(); }
181 ed << " target Z= " << Z << " A= " << A << G4endl;
182 G4Exception("G4VCrossSectionDataSet::GetIsoCrossSection", "had001",
183 FatalException, ed);
184 return 0.0;
185}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
std::ostringstream G4ExceptionDescription
#define G4endl
Definition G4ios.hh:67
const G4String & GetName() const
Definition G4Element.hh:115
const G4String & GetName() const
const G4String & GetParticleName() const

Referenced by ComputeCrossSectionPerElement().

◆ CrossSectionDescription()

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

◆ DumpPhysicsTable()

◆ ForAllAtomsAndEnergies()

bool G4VCrossSectionDataSet::ForAllAtomsAndEnergies ( ) const
inline

Definition at line 234 of file G4VCrossSectionDataSet.hh.

235{
236 return isForAllAtomsAndEnergies;
237}

Referenced by G4CrossSectionDataStore::AddDataSet(), and G4CrossSectionDataStore::AddDataSet().

◆ GetCrossSection()

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

Definition at line 197 of file G4VCrossSectionDataSet.hh.

200{
201 return ComputeCrossSection(dp, elm, mat);
202}
G4double ComputeCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *mat=nullptr)

◆ GetElementCrossSection()

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

◆ 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 G4BGGNucleonElasticXS, G4BGGNucleonInelasticXS, G4BGGPionElasticXS, G4BGGPionInelasticXS, G4ChargeExchangeXS, G4ChipsAntiBaryonElasticXS, G4ChipsAntiBaryonInelasticXS, G4ChipsHyperonElasticXS, G4ChipsHyperonInelasticXS, G4ChipsKaonMinusElasticXS, G4ChipsKaonMinusInelasticXS, G4ChipsKaonPlusElasticXS, G4ChipsKaonPlusInelasticXS, G4ChipsKaonZeroElasticXS, G4ChipsKaonZeroInelasticXS, G4ChipsNeutronElasticXS, G4ChipsNeutronInelasticXS, G4ChipsPionMinusElasticXS, G4ChipsPionMinusInelasticXS, G4ChipsPionPlusElasticXS, G4ChipsPionPlusInelasticXS, G4ChipsProtonElasticXS, G4ChipsProtonInelasticXS, G4CrossSectionHP, G4ElNeutrinoNucleusTotXsc, G4ElNucleusSFcs, G4GammaNuclearXS, G4IonsShenCrossSection, G4LENDCombinedCrossSection, G4LENDCrossSection, G4LENDGammaCrossSection, G4MuNeutrinoNucleusTotXsc, G4NeutronCaptureXS, G4NeutronElasticXS, G4NeutronHPCaptureData, G4NeutronInelasticXS, G4ParticleHPElasticData, G4ParticleHPFissionData, G4ParticleHPInelasticData, G4ParticleHPThermalScatteringData, G4ParticleInelasticXS, G4PhotoNuclearCrossSection, and G4TauNeutrinoNucleusTotXsc.

Definition at line 149 of file G4VCrossSectionDataSet.cc.

154{
156 ed << "GetIsoCrossSection is not implemented in <" << name << ">\n"
157 << "Particle: " << dynPart->GetDefinition()->GetParticleName()
158 << " Ekin(MeV)= " << dynPart->GetKineticEnergy()/MeV;
159 if(nullptr != mat) { ed << " material: " << mat->GetName(); }
160 if(nullptr != elm) { ed << " element: " << elm->GetName(); }
161 ed << " target Z= " << Z << " A= " << A << G4endl;
162 G4Exception("G4VCrossSectionDataSet::GetIsoCrossSection", "had001",
163 FatalException, ed);
164 return 0.0;
165}

Referenced by ComputeCrossSection(), and G4GammaNuclearXS::GetIsoCrossSection().

◆ GetMaxKinEnergy()

◆ GetMinKinEnergy()

◆ GetName()

const G4String & G4VCrossSectionDataSet::GetName ( ) const
inline

◆ IsElementApplicable()

◆ IsIsoApplicable()

◆ operator=()

G4VCrossSectionDataSet & G4VCrossSectionDataSet::operator= ( const G4VCrossSectionDataSet & right)
delete

◆ SelectIsotope()

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

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

Definition at line 188 of file G4VCrossSectionDataSet.cc.

190{
191 G4int nIso = (G4int)anElement->GetNumberOfIsotopes();
192 const G4Isotope* iso = anElement->GetIsotope(0);
193
194 // more than 1 isotope
195 if(1 < nIso) {
196 const G4double* abundVector = anElement->GetRelativeAbundanceVector();
197 G4double sum = 0.0;
199 for (G4int j=0; j<nIso; ++j) {
200 sum += abundVector[j];
201 if(q <= sum) {
202 iso = anElement->GetIsotope(j);
203 break;
204 }
205 }
206 }
207 return iso;
208}
#define G4UniformRand()
Definition Randomize.hh:52
const G4Isotope * GetIsotope(G4int iso) const
Definition G4Element.hh:151

◆ SetForAllAtomsAndEnergies()

◆ SetMaxKinEnergy()

◆ SetMinKinEnergy()

◆ SetName()

void G4VCrossSectionDataSet::SetName ( const G4String & nam)
inline

Definition at line 244 of file G4VCrossSectionDataSet.hh.

245{
246 name = nam;
247}

Referenced by G4ParticleHPInelasticData::G4ParticleHPInelasticData().

◆ SetVerboseLevel()

void G4VCrossSectionDataSet::SetVerboseLevel ( G4int value)
inlinevirtual

Member Data Documentation

◆ name

G4String G4VCrossSectionDataSet::name
protected

◆ verboseLevel

G4int G4VCrossSectionDataSet::verboseLevel
protected

Definition at line 181 of file G4VCrossSectionDataSet.hh.

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


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