Geant4-11
Public Member Functions | Protected Attributes | Private Attributes
G4ParticleHPCaptureData Class Reference

#include <G4ParticleHPCaptureData.hh>

Inheritance diagram for G4ParticleHPCaptureData:
G4VCrossSectionDataSet

Public Member Functions

void BuildPhysicsTable (const G4ParticleDefinition &)
 
G4double ComputeCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=nullptr)
 
virtual void CrossSectionDescription (std::ostream &) const
 
void DumpPhysicsTable (const G4ParticleDefinition &)
 
bool ForAllAtomsAndEnergies () const
 
 G4ParticleHPCaptureData ()
 
G4double GetCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=nullptr)
 
G4double GetCrossSection (const G4DynamicParticle *, const G4Element *, G4double aT)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *mat=nullptr)
 
G4double GetIsoCrossSection (const G4DynamicParticle *, G4int, G4int, const G4Isotope *, const G4Element *, const G4Material *)
 
G4double GetMaxKinEnergy () const
 
G4double GetMinKinEnergy () const
 
const G4StringGetName () const
 
G4int GetVerboseLevel () const
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, const G4Material *mat=nullptr)
 
G4bool IsIsoApplicable (const G4DynamicParticle *, G4int, G4int, const G4Element *, const G4Material *)
 
virtual const G4IsotopeSelectIsotope (const G4Element *, G4double kinEnergy, G4double logE)
 
void SetForAllAtomsAndEnergies (G4bool val)
 
void SetMaxKinEnergy (G4double value)
 
void SetMinKinEnergy (G4double value)
 
void SetName (const G4String &nam)
 
void SetVerboseLevel (G4int)
 
 ~G4ParticleHPCaptureData ()
 

Protected Attributes

G4String name
 
G4int verboseLevel
 

Private Attributes

const G4Elementelement_cache
 
G4bool instanceOfWorker
 
G4bool isForAllAtomsAndEnergies
 
G4double ke_cache
 
const G4Materialmaterial_cache
 
G4double maxKinEnergy
 
G4double minKinEnergy
 
G4CrossSectionDataSetRegistryregistry
 
G4PhysicsTabletheCrossSections
 
G4double xs_cache
 

Detailed Description

Definition at line 51 of file G4ParticleHPCaptureData.hh.

Constructor & Destructor Documentation

◆ G4ParticleHPCaptureData()

G4ParticleHPCaptureData::G4ParticleHPCaptureData ( )

Definition at line 51 of file G4ParticleHPCaptureData.cc.

52:G4VCrossSectionDataSet("NeutronHPCaptureXS")
53{
54 SetMinKinEnergy( 0*MeV );
55 SetMaxKinEnergy( 20*MeV );
56
58
59 instanceOfWorker = false;
61 instanceOfWorker = true;
62 }
63
64 element_cache = NULL;
65 material_cache = NULL;
66 ke_cache = 0.0;
67 xs_cache = 0.0;
68
69 //BuildPhysicsTable(*G4Neutron::Neutron());
70}
static constexpr double MeV
Definition: G4SIunits.hh:200
G4VCrossSectionDataSet(const G4String &nam="")
void SetMaxKinEnergy(G4double value)
void SetMinKinEnergy(G4double value)
G4bool IsWorkerThread()
Definition: G4Threading.cc:123

References element_cache, instanceOfWorker, G4Threading::IsWorkerThread(), ke_cache, material_cache, MeV, G4VCrossSectionDataSet::SetMaxKinEnergy(), G4VCrossSectionDataSet::SetMinKinEnergy(), theCrossSections, and xs_cache.

◆ ~G4ParticleHPCaptureData()

G4ParticleHPCaptureData::~G4ParticleHPCaptureData ( )

Definition at line 72 of file G4ParticleHPCaptureData.cc.

73{
74 if ( theCrossSections != NULL && instanceOfWorker != true ) {
76 delete theCrossSections;
77 theCrossSections = NULL;
78 }
79}
void clearAndDestroy()

References G4PhysicsTable::clearAndDestroy(), instanceOfWorker, and theCrossSections.

Member Function Documentation

◆ BuildPhysicsTable()

void G4ParticleHPCaptureData::BuildPhysicsTable ( const G4ParticleDefinition aP)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 120 of file G4ParticleHPCaptureData.cc.

121{
122 if(&aP!=G4Neutron::Neutron())
123 throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!");
124
127 return;
128 }
129
130 size_t numberOfElements = G4Element::GetNumberOfElements();
131 // G4cout << "CALLED G4ParticleHPCaptureData::BuildPhysicsTable "<<numberOfElements<<G4endl;
132 // TKDB
133 //if ( theCrossSections == 0 ) theCrossSections = new G4PhysicsTable( numberOfElements );
134 if ( theCrossSections == NULL )
135 theCrossSections = new G4PhysicsTable( numberOfElements );
136 else
138
139 // make a PhysicsVector for each element
140
141 static G4ThreadLocal G4ElementTable *theElementTable = 0 ; if (!theElementTable) theElementTable= G4Element::GetElementTable();
142 for( size_t i=0; i<numberOfElements; ++i )
143 {
144 #ifdef G4VERBOSE
145 if(std::getenv("CaptureDataIndexDebug"))
146 {
147 G4int index_debug = ((*theElementTable)[i])->GetIndex();
148 if ( G4HadronicParameters::Instance()->GetVerboseLevel() > 0 ) G4cout << "IndexDebug "<< i <<" "<<index_debug<<G4endl;
149 }
150 #endif
152 Instance(G4Neutron::Neutron())->MakePhysicsVector((*theElementTable)[i], this);
153 theCrossSections->push_back(physVec);
154 }
155
157}
std::vector< G4Element * > G4ElementTable
int G4int
Definition: G4Types.hh:85
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:397
static size_t GetNumberOfElements()
Definition: G4Element.cc:404
static G4HadronicParameters * Instance()
static G4Neutron * Neutron()
Definition: G4Neutron.cc:103
static G4ParticleHPData * Instance(G4ParticleDefinition *projectile)
G4PhysicsVector * MakePhysicsVector(G4Element *thE, G4ParticleHPFissionData *theP)
G4PhysicsTable * GetCaptureCrossSections()
void RegisterCaptureCrossSections(G4PhysicsTable *val)
static G4ParticleHPManager * GetInstance()
void push_back(G4PhysicsVector *)
#define G4ThreadLocal
Definition: tls.hh:77

References G4PhysicsTable::clearAndDestroy(), G4cout, G4endl, G4ThreadLocal, G4ParticleHPManager::GetCaptureCrossSections(), G4Element::GetElementTable(), G4ParticleHPManager::GetInstance(), G4Element::GetNumberOfElements(), G4HadronicParameters::GetVerboseLevel(), G4HadronicParameters::Instance(), G4ParticleHPData::Instance(), G4Threading::IsWorkerThread(), G4ParticleHPData::MakePhysicsVector(), G4Neutron::Neutron(), G4PhysicsTable::push_back(), G4ParticleHPManager::RegisterCaptureCrossSections(), and theCrossSections.

◆ ComputeCrossSection()

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

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}
std::vector< G4Isotope * > G4IsotopeVector
double G4double
Definition: G4Types.hh:83
const G4int Z[17]
const G4double A[17]
G4double * GetRelativeAbundanceVector() const
Definition: G4Element.hh:167
size_t GetNumberOfIsotopes() const
Definition: G4Element.hh:159
G4int GetZasInt() const
Definition: G4Element.hh:132
G4IsotopeVector * GetIsotopeVector() const
Definition: G4Element.hh:163
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)

References A, G4VCrossSectionDataSet::GetElementCrossSection(), G4VCrossSectionDataSet::GetIsoCrossSection(), G4Element::GetIsotopeVector(), G4Isotope::GetN(), G4Element::GetNumberOfIsotopes(), G4Element::GetRelativeAbundanceVector(), G4Element::GetZasInt(), G4VCrossSectionDataSet::IsElementApplicable(), G4VCrossSectionDataSet::IsIsoApplicable(), and Z.

Referenced by G4VCrossSectionDataSet::GetCrossSection().

◆ CrossSectionDescription()

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

Reimplemented from G4VCrossSectionDataSet.

Definition at line 295 of file G4ParticleHPCaptureData.cc.

296{
297 outFile << "High Precision cross data based on Evaluated Nuclear Data Files (ENDF) for radiative capture reaction of neutrons below 20MeV\n" ;
298}

◆ DumpPhysicsTable()

void G4ParticleHPCaptureData::DumpPhysicsTable ( const G4ParticleDefinition aP)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 159 of file G4ParticleHPCaptureData.cc.

160{
161 if(&aP!=G4Neutron::Neutron())
162 throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!");
163
164 #ifdef G4VERBOSE
165 if ( G4HadronicParameters::Instance()->GetVerboseLevel() == 0 ) return;
166
167//
168// Dump element based cross section
169// range 10e-5 eV to 20 MeV
170// 10 point per decade
171// in barn
172//
173
174 G4cout << G4endl;
175 G4cout << G4endl;
176 G4cout << "Capture Cross Section of Neutron HP"<< G4endl;
177 G4cout << "(Pointwise cross-section at 0 Kelvin.)" << G4endl;
178 G4cout << G4endl;
179 G4cout << "Name of Element" << G4endl;
180 G4cout << "Energy[eV] XS[barn]" << G4endl;
181 G4cout << G4endl;
182
183 size_t numberOfElements = G4Element::GetNumberOfElements();
184 static G4ThreadLocal G4ElementTable *theElementTable = 0 ; if (!theElementTable) theElementTable= G4Element::GetElementTable();
185
186 for ( size_t i = 0 ; i < numberOfElements ; ++i )
187 {
188
189 G4cout << (*theElementTable)[i]->GetName() << G4endl;
190
191 G4int ie = 0;
192
193 for ( ie = 0 ; ie < 130 ; ie++ )
194 {
195 G4double eKinetic = 1.0e-5 * G4Pow::GetInstance()->powA ( 10.0 , ie/10.0 ) *eV;
196 G4bool outOfRange = false;
197
198 if ( eKinetic < 20*MeV )
199 {
200 G4cout << eKinetic/eV << " " << (*((*theCrossSections)(i))).GetValue(eKinetic, outOfRange)/barn << G4endl;
201 }
202
203 }
204
205 G4cout << G4endl;
206 }
207
208 //G4cout << "G4ParticleHPCaptureData::DumpPhysicsTable still to be implemented"<<G4endl;
209 #endif
210}
static constexpr double barn
Definition: G4SIunits.hh:85
static constexpr double eV
Definition: G4SIunits.hh:201
bool G4bool
Definition: G4Types.hh:86
static G4Pow * GetInstance()
Definition: G4Pow.cc:41
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:230

References barn, eV, G4cout, G4endl, G4ThreadLocal, G4Element::GetElementTable(), G4Pow::GetInstance(), G4Element::GetNumberOfElements(), GetVerboseLevel(), G4HadronicParameters::Instance(), MeV, G4Neutron::Neutron(), and G4Pow::powA().

◆ ForAllAtomsAndEnergies()

bool G4VCrossSectionDataSet::ForAllAtomsAndEnergies ( ) const
inlineinherited

◆ GetCrossSection() [1/2]

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

Definition at line 187 of file G4VCrossSectionDataSet.hh.

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

References G4VCrossSectionDataSet::ComputeCrossSection().

◆ GetCrossSection() [2/2]

G4double G4ParticleHPCaptureData::GetCrossSection ( const G4DynamicParticle aP,
const G4Element anE,
G4double  aT 
)

Definition at line 214 of file G4ParticleHPCaptureData.cc.

216{
217 G4double result = 0;
218 G4bool outOfRange;
219 G4int index = anE->GetIndex();
220
221 // prepare neutron
222 G4double eKinetic = aP->GetKineticEnergy();
223
224 if ( G4ParticleHPManager::GetInstance()->GetNeglectDoppler() )
225 {
226 //NEGLECT_DOPPLER
227 G4double factor = 1.0;
228 if ( eKinetic < aT * k_Boltzmann )
229 {
230 // below 0.1 eV neutrons
231 // Have to do some, but now just igonre.
232 // Will take care after performance check.
233 // factor = factor * targetV;
234 }
235 return ( (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) )* factor;
236 }
237
238 G4ReactionProduct theNeutron( aP->GetDefinition() );
239 theNeutron.SetMomentum( aP->GetMomentum() );
240 theNeutron.SetKineticEnergy( eKinetic );
241
242 // prepare thermal nucleus
243 G4Nucleus aNuc;
244 G4double eps = 0.0001;
245 G4double theA = anE->GetN();
246 G4double theZ = anE->GetZ();
247 G4double eleMass;
248 eleMass = G4NucleiProperties::GetNuclearMass( static_cast<G4int>(theA+eps) , static_cast<G4int>(theZ+eps) ) / G4Neutron::Neutron()->GetPDGMass();
249
250 G4ReactionProduct boosted;
251 G4double aXsection;
252
253 // MC integration loop
254 G4int counter = 0;
255 G4double buffer = 0;
256 G4int size = G4int(std::max(10., aT/60*kelvin));
257 G4ThreeVector neutronVelocity = 1./G4Neutron::Neutron()->GetPDGMass()*theNeutron.GetMomentum();
258 G4double neutronVMag = neutronVelocity.mag();
259
260 while(counter == 0 || std::abs(buffer-result/std::max(1,counter)) > 0.03*buffer) // Loop checking, 11.05.2015, T. Koi
261 {
262 if(counter) buffer = result/counter;
263 while (counter<size) // Loop checking, 11.05.2015, T. Koi
264 {
265 counter ++;
266 G4ReactionProduct aThermalNuc = aNuc.GetThermalNucleus(eleMass, aT);
267 boosted.Lorentz(theNeutron, aThermalNuc);
268 G4double theEkin = boosted.GetKineticEnergy();
269 aXsection = (*((*theCrossSections)(index))).GetValue(theEkin, outOfRange);
270 // velocity correction, or luminosity factor...
271 G4ThreeVector targetVelocity = 1./aThermalNuc.GetMass()*aThermalNuc.GetMomentum();
272 aXsection *= (targetVelocity-neutronVelocity).mag()/neutronVMag;
273 result += aXsection;
274 }
275 size += size;
276 }
277 result /= counter;
278/*
279 // Checking impact of G4NEUTRONHP_NEGLECT_DOPPLER
280 G4cout << " result " << result << " "
281 << (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) << " "
282 << (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) /result << G4endl;
283*/
284 return result;
285}
static const G4double eps
static constexpr double kelvin
Definition: G4SIunits.hh:274
double mag() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
G4ThreeVector GetMomentum() const
G4double GetZ() const
Definition: G4Element.hh:131
size_t GetIndex() const
Definition: G4Element.hh:182
G4double GetN() const
Definition: G4Element.hh:135
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4ReactionProduct GetThermalNucleus(G4double aMass, G4double temp=-1) const
Definition: G4Nucleus.cc:236
G4double GetKineticEnergy() const
G4ThreeVector GetMomentum() const
void Lorentz(const G4ReactionProduct &p1, const G4ReactionProduct &p2)
G4double GetMass() const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
float k_Boltzmann
Definition: hepunit.py:298
#define buffer
Definition: xmlparse.cc:628

References buffer, eps, G4DynamicParticle::GetDefinition(), G4Element::GetIndex(), G4ParticleHPManager::GetInstance(), G4DynamicParticle::GetKineticEnergy(), G4ReactionProduct::GetKineticEnergy(), G4ReactionProduct::GetMass(), G4DynamicParticle::GetMomentum(), G4ReactionProduct::GetMomentum(), G4Element::GetN(), G4NucleiProperties::GetNuclearMass(), G4ParticleDefinition::GetPDGMass(), G4Nucleus::GetThermalNucleus(), G4Element::GetZ(), source.hepunit::k_Boltzmann, kelvin, G4ReactionProduct::Lorentz(), CLHEP::Hep3Vector::mag(), G4INCL::Math::max(), G4Neutron::Neutron(), G4ReactionProduct::SetKineticEnergy(), G4ReactionProduct::SetMomentum(), and theCrossSections.

Referenced by GetIsoCrossSection().

◆ GetElementCrossSection()

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

Reimplemented in G4EMDissociationCrossSection, G4IonsShenCrossSection, G4NeutrinoElectronCcXsc, G4NeutrinoElectronNcXsc, G4NeutrinoElectronTotXsc, G4NeutronElectronElXsc, G4PhotoNuclearCrossSection, G4NeutronCaptureXS, G4NeutronElasticXS, G4NeutronInelasticXS, G4ElectroNuclearCrossSection, G4BGGNucleonElasticXS, G4BGGPionElasticXS, G4BGGPionInelasticXS, G4BGGNucleonInelasticXS, G4CrossSectionElastic, G4CrossSectionInelastic, G4GammaNuclearXS, G4ParticleInelasticXS, G4ZeroXS, G4NucleonNuclearCrossSection, G4MuNeutrinoNucleusTotXsc, and G4KokoulinMuonNuclearXS.

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
const G4String & GetName() const
Definition: G4Material.hh:173
const G4String & GetParticleName() const

References FatalException, G4endl, G4Exception(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetKineticEnergy(), G4Material::GetName(), G4ParticleDefinition::GetParticleName(), MeV, G4VCrossSectionDataSet::name, and Z.

Referenced by G4QMDReaction::ApplyYourself(), G4VCrossSectionDataSet::ComputeCrossSection(), G4GammaNuclearXS::GetElementCrossSection(), G4GammaNuclearXS::GetIsoCrossSection(), and G4GammaNuclearXS::Initialise().

◆ GetIsoCrossSection()

G4double G4ParticleHPCaptureData::GetIsoCrossSection ( const G4DynamicParticle dp,
G4int  ,
G4int  ,
const G4Isotope ,
const G4Element element,
const G4Material material 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 94 of file G4ParticleHPCaptureData.cc.

99{
100 if ( dp->GetKineticEnergy() == ke_cache && element == element_cache && material == material_cache ) return xs_cache;
101
103 element_cache = element;
105 G4double xs = GetCrossSection( dp , element , material->GetTemperature() );
106 xs_cache = xs;
107 return xs;
108}
G4double GetCrossSection(const G4DynamicParticle *, const G4Element *, G4double aT)
string material
Definition: eplot.py:19

References element_cache, GetCrossSection(), G4DynamicParticle::GetKineticEnergy(), ke_cache, eplot::material, material_cache, and xs_cache.

◆ GetMaxKinEnergy()

G4double G4VCrossSectionDataSet::GetMaxKinEnergy ( ) const
inlineinherited

◆ GetMinKinEnergy()

G4double G4VCrossSectionDataSet::GetMinKinEnergy ( ) const
inlineinherited

◆ GetName()

const G4String & G4VCrossSectionDataSet::GetName ( ) const
inlineinherited

◆ GetVerboseLevel()

G4int G4ParticleHPCaptureData::GetVerboseLevel ( ) const
virtual

◆ IsElementApplicable()

G4bool G4VCrossSectionDataSet::IsElementApplicable ( const G4DynamicParticle ,
G4int  Z,
const G4Material mat = nullptr 
)
virtualinherited

◆ IsIsoApplicable()

G4bool G4ParticleHPCaptureData::IsIsoApplicable ( const G4DynamicParticle dp,
G4int  ,
G4int  ,
const G4Element ,
const G4Material  
)
virtual

◆ SelectIsotope()

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

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:170

References G4UniformRand, G4Element::GetIsotope(), G4Element::GetNumberOfIsotopes(), and G4Element::GetRelativeAbundanceVector().

◆ SetForAllAtomsAndEnergies()

void G4VCrossSectionDataSet::SetForAllAtomsAndEnergies ( G4bool  val)
inlineinherited

◆ SetMaxKinEnergy()

void G4VCrossSectionDataSet::SetMaxKinEnergy ( G4double  value)
inlineinherited

◆ SetMinKinEnergy()

void G4VCrossSectionDataSet::SetMinKinEnergy ( G4double  value)
inlineinherited

◆ SetName()

void G4VCrossSectionDataSet::SetName ( const G4String nam)
inlineinherited

Definition at line 240 of file G4VCrossSectionDataSet.hh.

241{
242 name = nam;
243}

References G4VCrossSectionDataSet::name.

Referenced by G4ParticleHPInelasticData::G4ParticleHPInelasticData().

◆ SetVerboseLevel()

void G4ParticleHPCaptureData::SetVerboseLevel ( G4int  newValue)
virtual

Field Documentation

◆ element_cache

const G4Element* G4ParticleHPCaptureData::element_cache
private

Definition at line 95 of file G4ParticleHPCaptureData.hh.

Referenced by G4ParticleHPCaptureData(), and GetIsoCrossSection().

◆ instanceOfWorker

G4bool G4ParticleHPCaptureData::instanceOfWorker
private

Definition at line 91 of file G4ParticleHPCaptureData.hh.

Referenced by G4ParticleHPCaptureData(), and ~G4ParticleHPCaptureData().

◆ isForAllAtomsAndEnergies

G4bool G4VCrossSectionDataSet::isForAllAtomsAndEnergies
privateinherited

◆ ke_cache

G4double G4ParticleHPCaptureData::ke_cache
private

Definition at line 93 of file G4ParticleHPCaptureData.hh.

Referenced by G4ParticleHPCaptureData(), and GetIsoCrossSection().

◆ material_cache

const G4Material* G4ParticleHPCaptureData::material_cache
private

Definition at line 96 of file G4ParticleHPCaptureData.hh.

Referenced by G4ParticleHPCaptureData(), and GetIsoCrossSection().

◆ maxKinEnergy

G4double G4VCrossSectionDataSet::maxKinEnergy
privateinherited

◆ minKinEnergy

G4double G4VCrossSectionDataSet::minKinEnergy
privateinherited

◆ name

G4String G4VCrossSectionDataSet::name
protectedinherited

◆ registry

G4CrossSectionDataSetRegistry* G4VCrossSectionDataSet::registry
privateinherited

◆ theCrossSections

G4PhysicsTable* G4ParticleHPCaptureData::theCrossSections
private

◆ verboseLevel

G4int G4VCrossSectionDataSet::verboseLevel
protectedinherited

Definition at line 168 of file G4VCrossSectionDataSet.hh.

Referenced by 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(), G4GammaNuclearXS::G4GammaNuclearXS(), G4NeutronCaptureXS::G4NeutronCaptureXS(), G4NeutronElasticXS::G4NeutronElasticXS(), G4NeutronInelasticXS::G4NeutronInelasticXS(), G4ParticleInelasticXS::G4ParticleInelasticXS(), G4NeutronCaptureXS::GetElementCrossSection(), G4NeutronElasticXS::GetElementCrossSection(), G4NeutronInelasticXS::GetElementCrossSection(), G4BGGNucleonElasticXS::GetElementCrossSection(), G4BGGPionElasticXS::GetElementCrossSection(), G4BGGPionInelasticXS::GetElementCrossSection(), G4BGGNucleonInelasticXS::GetElementCrossSection(), G4GammaNuclearXS::GetElementCrossSection(), G4ParticleInelasticXS::GetElementCrossSection(), G4BGGNucleonElasticXS::GetIsoCrossSection(), G4BGGPionElasticXS::GetIsoCrossSection(), G4BGGPionInelasticXS::GetIsoCrossSection(), G4GammaNuclearXS::GetIsoCrossSection(), G4BGGNucleonInelasticXS::GetIsoCrossSection(), G4VCrossSectionDataSet::GetVerboseLevel(), G4NeutronElasticXS::Initialise(), G4ParticleInelasticXS::IsoCrossSection(), G4NeutronCaptureXS::IsoCrossSection(), G4NeutronInelasticXS::IsoCrossSection(), G4GammaNuclearXS::RetrieveVector(), G4NeutronCaptureXS::RetrieveVector(), G4NeutronInelasticXS::RetrieveVector(), G4ParticleInelasticXS::RetrieveVector(), and G4VCrossSectionDataSet::SetVerboseLevel().

◆ xs_cache

G4double G4ParticleHPCaptureData::xs_cache
private

Definition at line 94 of file G4ParticleHPCaptureData.hh.

Referenced by G4ParticleHPCaptureData(), and GetIsoCrossSection().


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