Geant4-11
Public Member Functions | Private Member Functions | Private Attributes
G4ElementSelector Class Reference

#include <G4ElementSelector.hh>

Public Member Functions

 G4ElementSelector ()
 
virtual const G4ElementSelectZandA (const G4Track &track, G4Nucleus *)
 
virtual ~G4ElementSelector ()
 

Private Member Functions

 G4ElementSelector (const G4ElementSelector &)
 
G4ElementSelectoroperator= (const G4ElementSelector &right)
 

Private Attributes

std::vector< G4doubleprob
 

Detailed Description

Definition at line 60 of file G4ElementSelector.hh.

Constructor & Destructor Documentation

◆ G4ElementSelector() [1/2]

G4ElementSelector::G4ElementSelector ( )

Definition at line 48 of file G4ElementSelector.cc.

49{}

◆ ~G4ElementSelector()

G4ElementSelector::~G4ElementSelector ( )
virtual

Definition at line 53 of file G4ElementSelector.cc.

54{}

◆ G4ElementSelector() [2/2]

G4ElementSelector::G4ElementSelector ( const G4ElementSelector )
private

Member Function Documentation

◆ operator=()

G4ElementSelector & G4ElementSelector::operator= ( const G4ElementSelector right)
private

◆ SelectZandA()

const G4Element * G4ElementSelector::SelectZandA ( const G4Track track,
G4Nucleus target 
)
virtual

Definition at line 59 of file G4ElementSelector.cc.

60{
61 // Fermi-Teller Z-low of mu- capture and exceptions
62 // for halogens and oxigen.
63 // N.C.Mukhopadhyay Phys. Rep. 30 (1977) 1.
64
65 size_t i = 0;
66 const G4Material* mat = track.GetMaterial();
67 size_t numberOfElements = mat->GetNumberOfElements();
68 const G4ElementVector* theElementVector = mat->GetElementVector();
69
70 if(1 < numberOfElements) {
71 if(numberOfElements > prob.size()) { prob.resize(numberOfElements, 0.0); }
72
73 const G4double* theAtomNumDensity = mat->GetAtomicNumDensityVector();
74
75 G4double sum = 0.0;
76 for (i=0; i < numberOfElements; ++i) {
77
78 G4int Z = (*theElementVector)[i]->GetZasInt();
79
80 // Halogens
81 if( (9 == Z) || (17 == Z) || (35 == Z) || (53 == Z) || (85 == Z) ) {
82 sum += 0.66 * Z * theAtomNumDensity[i];
83
84 // Oxigen
85 } else if( 8 == Z ) {
86 sum += 0.56 * Z * theAtomNumDensity[i];
87
88 // Others
89 } else {
90 sum += Z * theAtomNumDensity[i];
91 }
92 prob[i] = sum;
93 }
94
95 sum *= G4UniformRand();
96 for (i=0; i < numberOfElements; ++i) {
97 if(sum <= prob[i]) { break; }
98 }
99 }
100
101 const G4Element* elm = (*theElementVector)[i];
102 G4int Z = elm->GetZasInt();
103
104 // select isotope
105 const G4IsotopeVector* isv = elm->GetIsotopeVector();
106 size_t ni = isv->size();
107 i = 0;
108
109 if(1 < ni) {
110
111 const G4double* ab = elm->GetRelativeAbundanceVector();
113 for(i=0; i<ni; ++i) {
114 y -= ab[i];
115 if(y <= 0.0) { break; }
116 }
117 }
118
119 G4int A = elm->GetIsotope(i)->GetN();
120 target->SetParameters(A, Z);
121
122 return elm;
123}
static const G4double ab
std::vector< const G4Element * > G4ElementVector
std::vector< G4Isotope * > G4IsotopeVector
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
const G4int Z[17]
const G4double A[17]
#define G4UniformRand()
Definition: Randomize.hh:52
std::vector< G4double > prob
G4double * GetRelativeAbundanceVector() const
Definition: G4Element.hh:167
const G4Isotope * GetIsotope(G4int iso) const
Definition: G4Element.hh:170
G4int GetZasInt() const
Definition: G4Element.hh:132
G4IsotopeVector * GetIsotopeVector() const
Definition: G4Element.hh:163
G4int GetN() const
Definition: G4Isotope.hh:93
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:186
size_t GetNumberOfElements() const
Definition: G4Material.hh:182
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:212
void SetParameters(const G4double A, const G4double Z, const G4int numberOfLambdas=0)
Definition: G4Nucleus.cc:307
G4Material * GetMaterial() const

References A, ab, G4UniformRand, G4Material::GetAtomicNumDensityVector(), G4Material::GetElementVector(), G4Element::GetIsotope(), G4Element::GetIsotopeVector(), G4Track::GetMaterial(), G4Isotope::GetN(), G4Material::GetNumberOfElements(), G4Element::GetRelativeAbundanceVector(), G4Element::GetZasInt(), prob, G4Nucleus::SetParameters(), and Z.

Referenced by G4HadronStoppingProcess::AtRestDoIt(), and G4MuonMinusAtomicCapture::AtRestDoIt().

Field Documentation

◆ prob

std::vector<G4double> G4ElementSelector::prob
private

Definition at line 76 of file G4ElementSelector.hh.

Referenced by SelectZandA().


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