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

#include <G4ElNucleusSFcs.hh>

Inheritance diagram for G4ElNucleusSFcs:
G4VCrossSectionDataSet

Public Member Functions

virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
G4double ComputeCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=nullptr)
 
virtual void CrossSectionDescription (std::ostream &) const
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
bool ForAllAtomsAndEnergies () const
 
 G4ElNucleusSFcs ()
 
G4double GetCrossSection (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)
 
G4double GetMaxKinEnergy () const
 
G4double GetMinKinEnergy () const
 
const G4StringGetName () const
 
G4double GetRatio (G4int Z, G4int A)
 
virtual G4int GetVerboseLevel () const
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, const G4Material *)
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=nullptr, const G4Material *mat=nullptr)
 
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)
 
virtual void SetVerboseLevel (G4int value)
 
G4double ThresholdEnergy ()
 
virtual ~G4ElNucleusSFcs ()
 

Static Public Member Functions

static const char * Default_Name ()
 

Protected Attributes

G4String name
 
G4int verboseLevel
 

Private Attributes

G4ElectroNuclearCrossSectionfCHIPScs
 
G4bool isForAllAtomsAndEnergies
 
G4double maxKinEnergy
 
G4double minKinEnergy
 
G4CrossSectionDataSetRegistryregistry
 

Static Private Attributes

static const G4double fRR [19]
 
static const G4double fZZ [19]
 

Detailed Description

Definition at line 51 of file G4ElNucleusSFcs.hh.

Constructor & Destructor Documentation

◆ G4ElNucleusSFcs()

G4ElNucleusSFcs::G4ElNucleusSFcs ( )

Definition at line 42 of file G4ElNucleusSFcs.cc.

43{
45}
G4ElectroNuclearCrossSection * fCHIPScs
static const char * Default_Name()
G4VCrossSectionDataSet(const G4String &nam="")

References fCHIPScs.

◆ ~G4ElNucleusSFcs()

G4ElNucleusSFcs::~G4ElNucleusSFcs ( )
virtual

Definition at line 47 of file G4ElNucleusSFcs.cc.

48{
49 if( fCHIPScs != nullptr ) delete fCHIPScs;
50}

References fCHIPScs.

Member Function Documentation

◆ BuildPhysicsTable()

void G4VCrossSectionDataSet::BuildPhysicsTable ( const G4ParticleDefinition )
virtualinherited

◆ 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
int G4int
Definition: G4Types.hh:85
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 G4ElNucleusSFcs::CrossSectionDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 53 of file G4ElNucleusSFcs.cc.

54{
55 outFile << "G4ElNucleusSFcs provides the inelastic\n"
56 << "cross section for e- and e+ interactions with nuclei according to the structure function approach."
57 << "all energies.\n";
58}

◆ Default_Name()

static const char * G4ElNucleusSFcs::Default_Name ( )
inlinestatic

Definition at line 58 of file G4ElNucleusSFcs.hh.

58{return "ElectronNucleusSFcs";}

◆ DumpPhysicsTable()

void G4VCrossSectionDataSet::DumpPhysicsTable ( const G4ParticleDefinition )
virtualinherited

◆ ForAllAtomsAndEnergies()

bool G4VCrossSectionDataSet::ForAllAtomsAndEnergies ( ) const
inlineinherited

◆ GetCrossSection()

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().

◆ 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
static constexpr double MeV
Definition: G4SIunits.hh:200
#define G4endl
Definition: G4ios.hh:57
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
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 G4ElNucleusSFcs::GetIsoCrossSection ( const G4DynamicParticle aPart,
G4int  Z,
G4int  A,
const G4Isotope iso = nullptr,
const G4Element elm = nullptr,
const G4Material mat = nullptr 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 70 of file G4ElNucleusSFcs.cc.

72{
73 G4double xsc(0.), ratio(1.);
74
75 xsc = fCHIPScs->GetElementCrossSection( aPart, ZZ, nullptr ); //, const G4Material*);
76
77 ratio = GetRatio(ZZ,AA);
78
79 if( ratio > 0. ) xsc /= ratio;
80
81 return xsc;
82}
G4double GetRatio(G4int Z, G4int A)
virtual G4double GetElementCrossSection(const G4DynamicParticle *, G4int Z, const G4Material *mat)
static const G4double AA[5]
Definition: paraMaker.cc:44

References anonymous_namespace{paraMaker.cc}::AA, fCHIPScs, G4ElectroNuclearCrossSection::GetElementCrossSection(), and GetRatio().

◆ GetMaxKinEnergy()

G4double G4VCrossSectionDataSet::GetMaxKinEnergy ( ) const
inlineinherited

◆ GetMinKinEnergy()

G4double G4VCrossSectionDataSet::GetMinKinEnergy ( ) const
inlineinherited

◆ GetName()

const G4String & G4VCrossSectionDataSet::GetName ( ) const
inlineinherited

◆ GetRatio()

G4double G4ElNucleusSFcs::GetRatio ( G4int  Z,
G4int  A 
)

Definition at line 88 of file G4ElNucleusSFcs.cc.

89{
90 G4double ratio(1.);
91
92 if ( Z == 1 && A == 1 ) return 1.51;
93 else if ( Z == 1 && A == 2 ) return 0.33;
94 else if ( Z == 1 && A == 3 ) return 0.27;
95 else if ( Z == 2 && A == 4 ) return 1.81;
96 else if ( Z == 6 && A == 12 ) return 2.26;
97 else if ( Z == 7 && A == 14 ) return 2.47;
98 else if ( Z == 8 && A == 16 ) return 2.61;
99 else if ( Z == 13 && A == 27 ) return 2.57;
100 else if ( Z == 14 && A == 28 ) return 2.49;
101 else if ( Z == 18 && A == 40 ) return 2.72;
102 else if ( Z == 22 && A == 48 ) return 2.71;
103 else if ( Z == 26 && A == 56 ) return 2.79;
104 else if ( Z == 29 && A == 64 ) return 2.78;
105 else if ( Z == 32 && A == 73 ) return 2.87;
106 else if ( Z == 42 && A == 96 ) return 3.02;
107 else if ( Z == 46 && A == 106 ) return 3.02;
108 else if ( Z == 47 && A == 108 ) return 2.99;
109 else if ( Z == 48 && A == 112 ) return 3.00;
110 else if ( Z == 74 && A == 184 ) return 3.44;
111 else if ( Z == 79 && A == 200 ) return 3.49;
112 else if ( Z == 82 && A == 207 ) return 3.48;
113 else if ( Z == 92 && A == 238 ) return 3.88;
114 else
115 {
116 G4int it(0), iMax(19);
117 G4double zz = G4double(Z);
118
119 for ( it = 0; it < iMax; ++it ) if ( zz <= fZZ[it] ) break;
120
121 if ( it == 0 ) return fRR[0];
122 else if( it == iMax ) return fRR[iMax-1];
123 else
124 {
125 G4double x1 = fZZ[it-1];
126 G4double x2 = fZZ[it];
127 G4double y1 = fRR[it-1];
128 G4double y2 = fRR[it];
129
130 if( x1 >= x2 ) return fRR[it];
131 else
132 {
133 G4double angle = (y2-y1)/(x2-x1);
134 ratio = y1 + ( zz - x1 )*angle;
135 }
136 }
137 }
138 return ratio;
139}
static const G4double angle[DIMMOTT]
static const G4double fRR[19]
static const G4double fZZ[19]

References A, angle, fRR, fZZ, and Z.

Referenced by GetIsoCrossSection().

◆ GetVerboseLevel()

G4int G4VCrossSectionDataSet::GetVerboseLevel ( ) const
inlinevirtualinherited

◆ IsElementApplicable()

G4bool G4ElNucleusSFcs::IsElementApplicable ( const G4DynamicParticle aPart,
G4int  Z,
const G4Material  
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 60 of file G4ElNucleusSFcs.cc.

61{
62
63 G4double eKin = aPart->GetKineticEnergy();
64 G4double thTkin = ThresholdEnergy();
65
66 return ( Z > 0 && Z < 120 && eKin > thTkin );
67}
G4double ThresholdEnergy()

References G4DynamicParticle::GetKineticEnergy(), ThresholdEnergy(), and Z.

◆ IsIsoApplicable()

G4bool G4VCrossSectionDataSet::IsIsoApplicable ( const G4DynamicParticle ,
G4int  Z,
G4int  A,
const G4Element elm = nullptr,
const G4Material mat = nullptr 
)
virtualinherited

◆ 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 G4VCrossSectionDataSet::SetVerboseLevel ( G4int  value)
inlinevirtualinherited

◆ ThresholdEnergy()

G4double G4ElNucleusSFcs::ThresholdEnergy ( )

Definition at line 143 of file G4ElNucleusSFcs.cc.

144{
145 G4double thTkin = 134.9766*CLHEP::MeV; // Pi0 threshold for the nucleon
146
147 thTkin *= 1.3; // nucleon motion
148
149 return thTkin;
150}
static constexpr double MeV

References CLHEP::MeV.

Referenced by IsElementApplicable().

Field Documentation

◆ fCHIPScs

G4ElectroNuclearCrossSection* G4ElNucleusSFcs::fCHIPScs
private

Definition at line 78 of file G4ElNucleusSFcs.hh.

Referenced by G4ElNucleusSFcs(), GetIsoCrossSection(), and ~G4ElNucleusSFcs().

◆ fRR

const G4double G4ElNucleusSFcs::fRR
staticprivate
Initial value:
=
{
1.81, 2.26, 2.47, 2.61, 2.57, 2.49, 2.72, 2.71, 2.79, 2.78,
2.87, 3.02, 3.02, 2.99, 3., 3.44, 3.49, 3.48, 2.88
}

Definition at line 76 of file G4ElNucleusSFcs.hh.

Referenced by GetRatio().

◆ fZZ

const G4double G4ElNucleusSFcs::fZZ
staticprivate
Initial value:
=
{
2., 6., 7., 8., 13., 14., 18., 22., 26., 29.,
32., 42., 46., 47., 48., 74., 79., 82., 92.
}

Definition at line 75 of file G4ElNucleusSFcs.hh.

Referenced by GetRatio().

◆ isForAllAtomsAndEnergies

G4bool G4VCrossSectionDataSet::isForAllAtomsAndEnergies
privateinherited

◆ maxKinEnergy

G4double G4VCrossSectionDataSet::maxKinEnergy
privateinherited

◆ minKinEnergy

G4double G4VCrossSectionDataSet::minKinEnergy
privateinherited

◆ name

G4String G4VCrossSectionDataSet::name
protectedinherited

◆ registry

G4CrossSectionDataSetRegistry* G4VCrossSectionDataSet::registry
privateinherited

◆ 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().


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