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

#include <G4ECDecay.hh>

Inheritance diagram for G4ECDecay:
G4NuclearDecay G4VDecayChannel

Public Member Functions

virtual G4DecayProductsDecayIt (G4double)
 
void DumpInfo ()
 
virtual void DumpNuclearInfo ()
 
 G4ECDecay (const G4ParticleDefinition *theParentNucleus, const G4double &theBR, const G4double &Qvalue, const G4double &excitation, const G4Ions::G4FloatLevelBase &flb, const G4RadioactiveDecayMode &mode)
 
G4int GetAngularMomentum ()
 
G4double GetBR () const
 
G4ParticleDefinitionGetDaughter (G4int anIndex)
 
G4double GetDaughterExcitation ()
 
G4double GetDaughterMass (G4int anIndex) const
 
const G4StringGetDaughterName (G4int anIndex) const
 
G4ParticleDefinitionGetDaughterNucleus ()
 
G4RadioactiveDecayMode GetDecayMode ()
 
G4Ions::G4FloatLevelBase GetFloatingLevel ()
 
G4double GetHLThreshold ()
 
const G4StringGetKinematicsName () const
 
G4int GetNumberOfDaughters () const
 
G4ParticleDefinitionGetParent ()
 
G4double GetParentMass () const
 
const G4StringGetParentName () const
 
const G4ThreeVectorGetPolarization () const
 
G4double GetRangeMass () const
 
G4int GetVerboseLevel () const
 
virtual G4bool IsOKWithParentMass (G4double parentMass)
 
G4bool operator!= (const G4VDecayChannel &r) const
 
G4bool operator< (const G4VDecayChannel &right) const
 
G4bool operator== (const G4VDecayChannel &r) const
 
void SetARM (G4bool onoff)
 
void SetBR (G4double value)
 
void SetDaughter (G4int anIndex, const G4ParticleDefinition *particle_type)
 
void SetDaughter (G4int anIndex, const G4String &particle_name)
 
void SetHLThreshold (G4double HLT)
 
void SetNumberOfDaughters (G4int value)
 
void SetParent (const G4ParticleDefinition *particle_type)
 
void SetParent (const G4String &particle_name)
 
void SetPolarization (const G4ThreeVector &)
 
void SetRangeMass (G4double val)
 
void SetVerboseLevel (G4int value)
 
virtual ~G4ECDecay ()
 

Protected Member Functions

void CheckAndFillDaughters ()
 
void CheckAndFillParent ()
 
void ClearDaughtersName ()
 
G4double DynamicalMass (G4double massPDG, G4double width, G4double maxDev=1.0) const
 

Protected Attributes

G4String ** daughters_name = nullptr
 
G4Mutex daughtersMutex
 
G4ParticleDefinition ** G4MT_daughters = nullptr
 
G4doubleG4MT_daughters_mass = nullptr
 
G4doubleG4MT_daughters_width = nullptr
 
G4ParticleDefinitionG4MT_parent = nullptr
 
G4double G4MT_parent_mass = 0.0
 
G4String kinematics_name = ""
 
G4int numberOfDaughters = 0
 
G4Stringparent_name = nullptr
 
G4ThreeVector parent_polarization
 
G4Mutex parentMutex
 
G4ParticleTableparticletable = nullptr
 
G4double rangeMass = 2.5
 
G4double rbranch = 0.0
 
const G4RadioactiveDecayMode theMode
 
G4int verboseLevel = 1
 

Static Protected Attributes

static const G4String noName = " "
 

Private Member Functions

void DefineSubshellProbabilities (G4int Z, G4int A)
 
void FillDaughters ()
 
void FillParent ()
 
const G4StringGetNoName () const
 

Private Attributes

G4bool applyARM
 
const G4double daughterEx
 
const G4Ions::G4FloatLevelBase floatingLevel
 
G4double halflifeThreshold
 
G4double PL1
 
G4double PL2
 
G4double PM1
 
G4double PM2
 
G4double PN1
 
G4double PN2
 
const G4double transitionQ
 

Static Private Attributes

static const G4double PL2overPL1 [100]
 
static const G4double PM2overPM1 [100]
 
static const G4double PN2overPN1 [100]
 

Detailed Description

Definition at line 42 of file G4ECDecay.hh.

Constructor & Destructor Documentation

◆ G4ECDecay()

G4ECDecay::G4ECDecay ( const G4ParticleDefinition theParentNucleus,
const G4double theBR,
const G4double Qvalue,
const G4double excitation,
const G4Ions::G4FloatLevelBase flb,
const G4RadioactiveDecayMode mode 
)

Definition at line 47 of file G4ECDecay.cc.

52 : G4NuclearDecay("electron capture", mode, excitationE, flb), transitionQ(Qvalue),
53 applyARM(true)
54{
55 SetParent(theParentNucleus); // Store name of parent nucleus, delete G4MT_parent
56 SetBR(branch);
57
59 G4IonTable* theIonTable =
61 G4int daughterZ = theParentNucleus->GetAtomicNumber() - 1;
62 G4int daughterA = theParentNucleus->GetAtomicMass();
63 SetDaughter(0, theIonTable->GetIon(daughterZ, daughterA, excitationE, flb) );
64 SetDaughter(1, "nu_e");
65 DefineSubshellProbabilities(daughterZ, daughterZ);
66}
int G4int
Definition: G4Types.hh:85
G4bool applyARM
Definition: G4ECDecay.hh:64
void DefineSubshellProbabilities(G4int Z, G4int A)
Definition: G4ECDecay.cc:230
const G4double transitionQ
Definition: G4ECDecay.hh:63
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:522
G4NuclearDecay(const G4String &channelName, const G4RadioactiveDecayMode &mode, const G4double &excitation, const G4Ions::G4FloatLevelBase &floatingLevel)
G4int GetAtomicNumber() const
G4int GetAtomicMass() const
G4IonTable * GetIonTable() const
static G4ParticleTable * GetParticleTable()
void SetBR(G4double value)
void SetNumberOfDaughters(G4int value)
void SetDaughter(G4int anIndex, const G4ParticleDefinition *particle_type)
void SetParent(const G4ParticleDefinition *particle_type)

References DefineSubshellProbabilities(), G4ParticleDefinition::GetAtomicMass(), G4ParticleDefinition::GetAtomicNumber(), G4IonTable::GetIon(), G4ParticleTable::GetIonTable(), G4ParticleTable::GetParticleTable(), G4VDecayChannel::SetBR(), G4VDecayChannel::SetDaughter(), G4VDecayChannel::SetNumberOfDaughters(), and G4VDecayChannel::SetParent().

◆ ~G4ECDecay()

G4ECDecay::~G4ECDecay ( )
virtual

Definition at line 69 of file G4ECDecay.cc.

70{}

Member Function Documentation

◆ CheckAndFillDaughters()

void G4VDecayChannel::CheckAndFillDaughters ( )
inlineprotectedinherited

◆ CheckAndFillParent()

void G4VDecayChannel::CheckAndFillParent ( )
inlineprotectedinherited

◆ ClearDaughtersName()

void G4VDecayChannel::ClearDaughtersName ( )
protectedinherited

Definition at line 183 of file G4VDecayChannel.cc.

184{
186 if ( daughters_name != nullptr)
187 {
188 if (numberOfDaughters>0)
189 {
190#ifdef G4VERBOSE
191 if (verboseLevel>1)
192 {
193 G4cout << "G4VDecayChannel::ClearDaughtersName() "
194 << " for " << *parent_name << G4endl;
195 }
196#endif
197 for (G4int index=0; index<numberOfDaughters; ++index)
198 {
199 delete daughters_name[index];
200 }
201 }
202 delete [] daughters_name;
203 daughters_name = nullptr;
204 }
205
206 delete [] G4MT_daughters;
207 delete [] G4MT_daughters_mass;
208 delete [] G4MT_daughters_width;
209 G4MT_daughters_width = nullptr;
210 G4MT_daughters = nullptr;
211 G4MT_daughters_mass = nullptr;
212
214}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
G4String * parent_name
G4String ** daughters_name
G4double * G4MT_daughters_mass
G4double * G4MT_daughters_width

References G4VDecayChannel::daughters_name, G4VDecayChannel::daughtersMutex, G4cout, G4endl, G4VDecayChannel::G4MT_daughters, G4VDecayChannel::G4MT_daughters_mass, G4VDecayChannel::G4MT_daughters_width, G4VDecayChannel::numberOfDaughters, G4VDecayChannel::parent_name, and G4VDecayChannel::verboseLevel.

Referenced by G4DalitzDecayChannel::operator=(), G4KL3DecayChannel::operator=(), G4MuonDecayChannel::operator=(), G4MuonDecayChannelWithSpin::operator=(), G4MuonRadiativeDecayChannelWithSpin::operator=(), G4NeutronBetaDecayChannel::operator=(), G4PionRadiativeDecayChannel::operator=(), G4TauLeptonicDecayChannel::operator=(), G4VDecayChannel::operator=(), G4VDecayChannel::SetNumberOfDaughters(), and G4VDecayChannel::~G4VDecayChannel().

◆ DecayIt()

G4DecayProducts * G4ECDecay::DecayIt ( G4double  )
virtual

Implements G4VDecayChannel.

Definition at line 73 of file G4ECDecay.cc.

74{
75 // Fill G4MT_parent with theParentNucleus (stored by SetParent in ctor)
77
78 // Fill G4MT_daughters with alpha and residual nucleus (stored by SetDaughter)
80
81 // Get shell number of captured electron
82 G4int shellIndex = -1;
83 G4double ran;
84 switch (theMode)
85 {
86 case KshellEC:
87 shellIndex = 0;
88 break;
89 case LshellEC: // PL1+PL2+PL3=1
90 ran=G4UniformRand();
91 if (ran <= PL1) shellIndex =1;
92 else if (ran<= (PL1+PL2)) shellIndex =2;
93 else shellIndex =3;
94 break;
95 case MshellEC: // PM1+PM2+PM3=1
96 ran=G4UniformRand();
97 if (ran < PM1) shellIndex =4;
98 else if (ran< (PM1+PM2)) shellIndex =5;
99 else shellIndex = 6;
100 break;
101 case NshellEC: // PN1+PN2+PN3=1
102 ran=G4UniformRand();
103 if (ran < PN1) shellIndex = 9;
104 else if (ran<= (PN1+PN2)) shellIndex =2;
105 else shellIndex =10;
106 break;
107 default:
108 G4Exception("G4ECDecay::DecayIt()", "HAD_RDM_009",
109 FatalException, "Invalid electron shell selected");
110 }
111
112 // Initialize decay products with parent nucleus at rest
113 G4DynamicParticle parentParticle(G4MT_parent, G4ThreeVector(0,0,0), 0.0);
114 G4DecayProducts* products = new G4DecayProducts(parentParticle);
115 G4double eBind = 0.0;
116
117 // G4LossTableManager must already be initialized with G4UAtomicDeexcitation
118 // This is currently done in G4RadioactiveDecay::BuildPhysicsTable
119 G4VAtomDeexcitation* atomDeex =
121 std::vector<G4DynamicParticle*> armProducts;
122
123 if (applyARM) {
124 if (atomDeex) {
127 if (shellIndex >= nShells) shellIndex = nShells;
129 const G4AtomicShell* shell = atomDeex->GetAtomicShell(aZ, as);
130 eBind = shell->BindingEnergy();
131 if (atomDeex->IsFluoActive() && aZ > 5 && aZ < 100) {
132 // Do atomic relaxation
133 // VI, SI
134 // Allows fixing of Bugzilla 1727
135 //const G4double deexLimit = 0.1*keV;
136 G4double deexLimit = 0.1*keV;
137 if (G4EmParameters::Instance()->DeexcitationIgnoreCut()) deexLimit =0.;
138 //
139 atomDeex->GenerateParticles(&armProducts, shell, aZ, deexLimit, deexLimit);
140 }
141
142 G4double productEnergy = 0.;
143 for (G4int i = 0; i < G4int(armProducts.size()); i++)
144 productEnergy += armProducts[i]->GetKineticEnergy();
145
146 G4double deficit = shell->BindingEnergy() - productEnergy;
147 if (deficit > 0.0) {
148 // Add a dummy electron to make up extra energy
149 G4double cosTh = 1.-2.*G4UniformRand();
150 G4double sinTh = std::sqrt(1.- cosTh*cosTh);
152
153 G4ThreeVector electronDirection(sinTh*std::sin(phi),
154 sinTh*std::cos(phi), cosTh);
155 G4DynamicParticle* extra =
156 new G4DynamicParticle(G4Electron::Electron(), electronDirection,
157 deficit);
158 armProducts.push_back(extra);
159 }
160 } // atomDeex
161 } // applyARM
162
163 G4double daughterMass = G4MT_daughters[0]->GetPDGMass();
164
165 // CM momentum using Q value corrected for binding energy of captured electron
166 G4double Q = transitionQ - eBind;
167 G4double cmMomentum = Q*(Q + 2.*daughterMass)/(Q + daughterMass)/2.;
168
169 G4double costheta = 2.*G4UniformRand() - 1.0;
170 G4double sintheta = std::sqrt(1.0 - costheta*costheta);
172 G4ThreeVector direction(sintheta*std::cos(phi),sintheta*std::sin(phi),
173 costheta);
174 G4double KE = cmMomentum;
175 G4DynamicParticle* daughterParticle =
176 new G4DynamicParticle(G4MT_daughters[1], direction, KE, 0.0);
177 products->PushProducts(daughterParticle);
178
179 KE = std::sqrt(cmMomentum*cmMomentum + daughterMass*daughterMass) - daughterMass;
180 daughterParticle =
181 new G4DynamicParticle(G4MT_daughters[0], -1.0*direction, KE, daughterMass);
182 products->PushProducts(daughterParticle);
183
184 G4int nArm = armProducts.size();
185 if (nArm > 0) {
186 G4ThreeVector bst = daughterParticle->Get4Momentum().boostVector();
187 for (G4int i = 0; i < nArm; ++i) {
188 G4DynamicParticle* dp = armProducts[i];
189 G4LorentzVector lv = dp->Get4Momentum().boost(bst);
190 dp->Set4Momentum(lv);
191 products->PushProducts(dp);
192 }
193 }
194
195 // Energy conservation check
196 /*
197 G4int newSize = products->entries();
198 G4DynamicParticle* temp = 0;
199 G4double KEsum = 0.0;
200 for (G4int i = 0; i < newSize; i++) {
201 temp = products->operator[](i);
202 KEsum += temp->GetKineticEnergy();
203 }
204
205 G4double eCons = (transitionQ - KEsum)/keV;
206 G4cout << " EC check: Ediff (keV) = " << eCons << G4endl;
207 */
208 return products;
209}
G4AtomicShellEnumerator
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
static constexpr double twopi
Definition: G4SIunits.hh:56
static constexpr double rad
Definition: G4SIunits.hh:129
static constexpr double keV
Definition: G4SIunits.hh:202
CLHEP::Hep3Vector G4ThreeVector
double G4double
Definition: G4Types.hh:83
#define G4UniformRand()
Definition: Randomize.hh:52
Hep3Vector boostVector() const
HepLorentzVector & boost(double, double, double)
G4double BindingEnergy() const
static G4int GetNumberOfShells(G4int Z)
G4int PushProducts(G4DynamicParticle *aParticle)
G4LorentzVector Get4Momentum() const
void Set4Momentum(const G4LorentzVector &momentum)
G4double PM1
Definition: G4ECDecay.hh:65
G4double PM2
Definition: G4ECDecay.hh:65
G4double PN2
Definition: G4ECDecay.hh:65
G4double PL2
Definition: G4ECDecay.hh:65
G4double PL1
Definition: G4ECDecay.hh:65
G4double PN1
Definition: G4ECDecay.hh:65
static G4Electron * Electron()
Definition: G4Electron.cc:93
static G4EmParameters * Instance()
static G4LossTableManager * Instance()
G4VAtomDeexcitation * AtomDeexcitation()
const G4RadioactiveDecayMode theMode
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
void GenerateParticles(std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
void CheckAndFillDaughters()
static double Q[]

References applyARM, G4LossTableManager::AtomDeexcitation(), G4AtomicShell::BindingEnergy(), CLHEP::HepLorentzVector::boost(), CLHEP::HepLorentzVector::boostVector(), G4VDecayChannel::CheckAndFillDaughters(), G4VDecayChannel::CheckAndFillParent(), G4Electron::Electron(), FatalException, G4Exception(), G4VDecayChannel::G4MT_daughters, G4VDecayChannel::G4MT_parent, G4UniformRand, G4VAtomDeexcitation::GenerateParticles(), G4DynamicParticle::Get4Momentum(), G4ParticleDefinition::GetAtomicNumber(), G4VAtomDeexcitation::GetAtomicShell(), G4AtomicShells::GetNumberOfShells(), G4ParticleDefinition::GetPDGMass(), G4EmParameters::Instance(), G4LossTableManager::Instance(), G4VAtomDeexcitation::IsFluoActive(), keV, KshellEC, LshellEC, MshellEC, NshellEC, PL1, PL2, PM1, PM2, PN1, PN2, G4DecayProducts::PushProducts(), Q, rad, G4DynamicParticle::Set4Momentum(), G4NuclearDecay::theMode, transitionQ, and twopi.

◆ DefineSubshellProbabilities()

void G4ECDecay::DefineSubshellProbabilities ( G4int  Z,
G4int  A 
)
private

Definition at line 230 of file G4ECDecay.cc.

231{ //Implementation for the case of allowed transitions
232 //PL1+PL2=1. , PM1+PM2=1., PN1+PN2=1.
233 PL1 = 1./(1+PL2overPL1[Z-1]);
234 PL2 = PL1*PL2overPL1[Z-1];
235 PM1 = 1./(1+PM2overPM1[Z-1]);
236 PM2 = PM1*PM2overPM1[Z-1];
237 PN1 = 1./(1+PN2overPN1[Z-1]);
238 PN2 = PN1*PN2overPN1[Z-1];
239}
const G4int Z[17]
static const G4double PM2overPM1[100]
Definition: G4ECDecay.hh:69
static const G4double PN2overPN1[100]
Definition: G4ECDecay.hh:70
static const G4double PL2overPL1[100]
Definition: G4ECDecay.hh:68

References PL1, PL2, PL2overPL1, PM1, PM2, PM2overPM1, PN1, PN2, PN2overPN1, and Z.

Referenced by G4ECDecay().

◆ DumpInfo()

void G4VDecayChannel::DumpInfo ( )
inherited

Definition at line 560 of file G4VDecayChannel.cc.

561{
562 G4cout << " BR: " << rbranch << " [" << kinematics_name << "]";
563 G4cout << " : " ;
564 for (G4int index=0; index<numberOfDaughters; ++index)
565 {
566 if(daughters_name[index] != nullptr)
567 {
568 G4cout << " " << *(daughters_name[index]);
569 }
570 else
571 {
572 G4cout << " not defined ";
573 }
574 }
575 G4cout << G4endl;
576}
G4String kinematics_name

References G4VDecayChannel::daughters_name, G4cout, G4endl, G4VDecayChannel::kinematics_name, G4VDecayChannel::numberOfDaughters, and G4VDecayChannel::rbranch.

Referenced by G4GeneralPhaseSpaceDecay::DecayIt(), G4PhaseSpaceDecayChannel::DecayIt(), and G4KL3DecayChannel::G4KL3DecayChannel().

◆ DumpNuclearInfo()

void G4ECDecay::DumpNuclearInfo ( )
virtual

Implements G4NuclearDecay.

Definition at line 212 of file G4ECDecay.cc.

213{
214 G4cout << " G4ECDecay of parent nucleus " << GetParentName() << " from ";
215 if (theMode == 3) {
216 G4cout << "K shell";
217 } else if (theMode == 4) {
218 G4cout << "L shell";
219 } else if (theMode == 5) {
220 G4cout << "M shell";
221 }
222 else if (theMode == 6) {
223 G4cout << "N shell";
224 }
225 G4cout << G4endl;
226 G4cout << " to " << GetDaughterName(0) << " + " << GetDaughterName(1)
227 << " with branching ratio " << GetBR() << "% and Q value "
228 << transitionQ << G4endl;
229}
G4double GetBR() const
const G4String & GetParentName() const
const G4String & GetDaughterName(G4int anIndex) const

References G4cout, G4endl, G4VDecayChannel::GetBR(), G4VDecayChannel::GetDaughterName(), G4VDecayChannel::GetParentName(), G4NuclearDecay::theMode, and transitionQ.

◆ DynamicalMass()

G4double G4VDecayChannel::DynamicalMass ( G4double  massPDG,
G4double  width,
G4double  maxDev = 1.0 
) const
protectedinherited

Definition at line 585 of file G4VDecayChannel.cc.

587{
588 if (width<=0.0) return massPDG;
589 if (maxDev >rangeMass) maxDev = rangeMass;
590 if (maxDev <=-1.*rangeMass) return massPDG; // cannot calculate
591
592 G4double x = G4UniformRand()*(maxDev+rangeMass) - rangeMass;
594 const std::size_t MAX_LOOP=10000;
595 for (std::size_t loop_counter=0; loop_counter <MAX_LOOP; ++loop_counter)
596 {
597 if ( y * (width*width*x*x + massPDG*massPDG*width*width)
598 <= massPDG*massPDG*width*width ) break;
599 x = G4UniformRand()*(maxDev+rangeMass) - rangeMass;
600 y = G4UniformRand();
601 }
602 G4double mass = massPDG + x*width;
603 return mass;
604}

References G4UniformRand, and G4VDecayChannel::rangeMass.

Referenced by G4PhaseSpaceDecayChannel::ThreeBodyDecayIt(), and G4PhaseSpaceDecayChannel::TwoBodyDecayIt().

◆ FillDaughters()

void G4VDecayChannel::FillDaughters ( )
privateinherited

Definition at line 308 of file G4VDecayChannel.cc.

309{
311
312 // Double check, check again if another thread has already filled this, in
313 // case do not need to do anything
314 if ( G4MT_daughters != nullptr ) return;
315
316 G4int index;
317
318#ifdef G4VERBOSE
319 if (verboseLevel>1) G4cout << "G4VDecayChannel::FillDaughters()" << G4endl;
320#endif
321 if (G4MT_daughters != nullptr)
322 {
323 delete [] G4MT_daughters;
324 G4MT_daughters = nullptr;
325 }
326
327 // parent mass
329 G4double parentmass = G4MT_parent->GetPDGMass();
330
331 //
332 G4double sumofdaughtermass = 0.0;
333 G4double sumofdaughterwidthsq = 0.0;
334
335 if ((numberOfDaughters <=0) || (daughters_name == nullptr) )
336 {
337#ifdef G4VERBOSE
338 if (verboseLevel>0)
339 {
340 G4cout << "G4VDecayChannel::FillDaughters() - "
341 << "[ " << G4MT_parent->GetParticleName() << " ]"
342 << "numberOfDaughters is not defined yet";
343 }
344#endif
345 G4MT_daughters = nullptr;
346 G4Exception("G4VDecayChannel::FillDaughters()",
347 "PART011", FatalException,
348 "Cannot fill daughters: numberOfDaughters is not defined yet");
349 }
350
351 // create and set the array of pointers to daughter particles
353 if (G4MT_daughters_mass != nullptr) delete [] G4MT_daughters_mass;
354 if (G4MT_daughters_width != nullptr) delete [] G4MT_daughters_width;
357 // loop over all daughters
358 for (index=0; index<numberOfDaughters; ++index)
359 {
360 if (daughters_name[index] == nullptr)
361 {
362 // daughter name is not defined
363#ifdef G4VERBOSE
364 if (verboseLevel>0)
365 {
366 G4cout << "G4VDecayChannel::FillDaughters() - "
367 << "[ " << G4MT_parent->GetParticleName() << " ]"
368 << index << "-th daughter is not defined yet" << G4endl;
369 }
370#endif
371 G4MT_daughters[index] = nullptr;
372 G4Exception("G4VDecayChannel::FillDaughters()",
373 "PART011", FatalException,
374 "Cannot fill daughters: name of daughter is not defined yet");
375 }
376 // search daughter particles in the particle table
378 if (G4MT_daughters[index] == nullptr )
379 {
380 // cannot find the daughter particle
381#ifdef G4VERBOSE
382 if (verboseLevel>0)
383 {
384 G4cout << "G4VDecayChannel::FillDaughters() - "
385 << "[ " << G4MT_parent->GetParticleName() << " ]"
386 << index << ":" << *daughters_name[index]
387 << " is not defined !!" << G4endl;
388 G4cout << " The BR of this decay mode is set to zero." << G4endl;
389 }
390#endif
391 SetBR(0.0);
392 return;
393 }
394#ifdef G4VERBOSE
395 if (verboseLevel>1)
396 {
397 G4cout << index << ":" << *daughters_name[index];
398 G4cout << ":" << G4MT_daughters[index] << G4endl;
399 }
400#endif
402 G4double d_width = G4MT_daughters[index]->GetPDGWidth();
403 G4MT_daughters_width[index] = d_width;
404 sumofdaughtermass += G4MT_daughters[index]->GetPDGMass();
405 sumofdaughterwidthsq += d_width*d_width;
406 } // end loop over all daughters
407
408 // check sum of daghter mass
409 G4double widthMass = std::sqrt(G4MT_parent->GetPDGWidth()
411 +sumofdaughterwidthsq);
412 if ( (G4MT_parent->GetParticleType() != "nucleus")
413 && (numberOfDaughters !=1)
414 && (sumofdaughtermass > parentmass + rangeMass*widthMass) )
415 {
416 // !!! illegal mass !!!
417#ifdef G4VERBOSE
418 if (GetVerboseLevel()>0)
419 {
420 G4cout << "G4VDecayChannel::FillDaughters() - "
421 << "[ " << G4MT_parent->GetParticleName() << " ]"
422 << " Energy/Momentum conserevation breaks " << G4endl;
423 if (GetVerboseLevel()>1)
424 {
425 G4cout << " parent:" << *parent_name
426 << " mass:" << parentmass/GeV << "[GeV/c/c]" << G4endl;
427 for (index=0; index<numberOfDaughters; ++index)
428 {
429 G4cout << " daughter " << index << ":" << *daughters_name[index]
430 << " mass:" << G4MT_daughters[index]->GetPDGMass()/GeV
431 << "[GeV/c/c]" << G4endl;
432 }
433 }
434 }
435#endif
436 }
437}
static constexpr double GeV
Definition: G4SIunits.hh:203
const G4String & GetParticleType() const
G4double GetPDGWidth() const
const G4String & GetParticleName() const
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
G4int GetVerboseLevel() const
G4ParticleTable * particletable

References G4VDecayChannel::CheckAndFillParent(), G4VDecayChannel::daughters_name, G4VDecayChannel::daughtersMutex, FatalException, G4ParticleTable::FindParticle(), G4cout, G4endl, G4Exception(), G4VDecayChannel::G4MT_daughters, G4VDecayChannel::G4MT_daughters_mass, G4VDecayChannel::G4MT_daughters_width, G4VDecayChannel::G4MT_parent, G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetParticleType(), G4ParticleDefinition::GetPDGMass(), G4ParticleDefinition::GetPDGWidth(), G4VDecayChannel::GetVerboseLevel(), GeV, G4VDecayChannel::numberOfDaughters, G4VDecayChannel::parent_name, G4VDecayChannel::particletable, G4VDecayChannel::rangeMass, G4VDecayChannel::SetBR(), and G4VDecayChannel::verboseLevel.

Referenced by G4VDecayChannel::CheckAndFillDaughters().

◆ FillParent()

void G4VDecayChannel::FillParent ( )
privateinherited

Definition at line 440 of file G4VDecayChannel.cc.

441{
442 G4AutoLock lock(&parentMutex);
443
444 // Double check, check again if another thread has already filled this, in
445 // case do not need to do anything
446 if ( G4MT_parent != nullptr ) return;
447
448 if (parent_name == nullptr)
449 {
450 // parent name is not defined
451#ifdef G4VERBOSE
452 if (verboseLevel>0)
453 {
454 G4cout << "G4VDecayChannel::FillParent() - "
455 << "parent name is not defined !!" << G4endl;
456 }
457#endif
458 G4MT_parent = nullptr;
459 G4Exception("G4VDecayChannel::FillParent()",
460 "PART012", FatalException,
461 "Cannot fill parent: parent name is not defined yet");
462 return;
463 }
464
465 // search parent particle in the particle table
467 if (G4MT_parent == nullptr)
468 {
469 // parent particle does not exist
470#ifdef G4VERBOSE
471 if (verboseLevel>0)
472 {
473 G4cout << "G4VDecayChannel::FillParent() - "
474 << *parent_name << " does not exist !!" << G4endl;
475 }
476#endif
477 G4Exception("G4VDecayChannel::FillParent()",
478 "PART012", FatalException,
479 "Cannot fill parent: parent does not exist");
480 return;
481 }
483}
G4double G4MT_parent_mass

References FatalException, G4ParticleTable::FindParticle(), G4cout, G4endl, G4Exception(), G4VDecayChannel::G4MT_parent, G4VDecayChannel::G4MT_parent_mass, G4ParticleDefinition::GetPDGMass(), G4VDecayChannel::parent_name, G4VDecayChannel::parentMutex, G4VDecayChannel::particletable, and G4VDecayChannel::verboseLevel.

Referenced by G4VDecayChannel::CheckAndFillParent().

◆ GetAngularMomentum()

G4int G4VDecayChannel::GetAngularMomentum ( )
inherited

Definition at line 492 of file G4VDecayChannel.cc.

493{
494 // determine angular momentum
495
496 // fill pointers to daughter particles if not yet set
498
499 const G4int PiSpin = G4MT_parent->GetPDGiSpin();
500 const G4int PParity = G4MT_parent->GetPDGiParity();
501 if (2==numberOfDaughters) // up to now we can only handle two particle decays
502 {
503 const G4int D1iSpin = G4MT_daughters[0]->GetPDGiSpin();
504 const G4int D1Parity = G4MT_daughters[0]->GetPDGiParity();
505 const G4int D2iSpin = G4MT_daughters[1]->GetPDGiSpin();
506 const G4int D2Parity = G4MT_daughters[1]->GetPDGiParity();
507 const G4int MiniSpin = std::abs (D1iSpin - D2iSpin);
508 const G4int MaxiSpin = D1iSpin + D2iSpin;
509 const G4int lMax = (PiSpin+D1iSpin+D2iSpin)/2; // l is always int
510 G4int lMin;
511#ifdef G4VERBOSE
512 if (verboseLevel>1)
513 {
514 G4cout << "iSpin: " << PiSpin << " -> " << D1iSpin
515 << " + " << D2iSpin << G4endl;
516 G4cout << "2*jmin, 2*jmax, lmax " << MiniSpin << " " << MaxiSpin
517 << " " << lMax << G4endl;
518 }
519#endif
520 for (G4int j=MiniSpin; j<=MaxiSpin; j+=2) // loop over all possible
521 { // spin couplings
522 lMin = std::abs(PiSpin-j)/2;
523#ifdef G4VERBOSE
524 if (verboseLevel>1)
525 G4cout << "-> checking 2*j=" << j << G4endl;
526#endif
527 for (G4int l=lMin; l<=lMax; ++l)
528 {
529#ifdef G4VERBOSE
530 if (verboseLevel>1)
531 G4cout << " checking l=" << l << G4endl;
532#endif
533 if (l%2==0)
534 {
535 if (PParity == D1Parity*D2Parity) // check parity for this l
536 return l;
537 }
538 else
539 {
540 if (PParity == -1*D1Parity*D2Parity) // check parity for this l
541 return l;
542 }
543 }
544 }
545 }
546 else
547 {
548 G4Exception("G4VDecayChannel::GetAngularMomentum()",
549 "PART111", JustWarning,
550 "Sorry, can't handle 3 particle decays (up to now)");
551 return 0;
552 }
553 G4Exception ("G4VDecayChannel::GetAngularMomentum()",
554 "PART111", JustWarning,
555 "Can't find angular momentum for this decay");
556 return 0;
557}
@ JustWarning

References G4VDecayChannel::CheckAndFillDaughters(), G4cout, G4endl, G4Exception(), G4VDecayChannel::G4MT_daughters, G4VDecayChannel::G4MT_parent, G4ParticleDefinition::GetPDGiParity(), G4ParticleDefinition::GetPDGiSpin(), JustWarning, G4VDecayChannel::numberOfDaughters, and G4VDecayChannel::verboseLevel.

◆ GetBR()

G4double G4VDecayChannel::GetBR ( ) const
inlineinherited

◆ GetDaughter()

G4ParticleDefinition * G4VDecayChannel::GetDaughter ( G4int  anIndex)
inlineinherited

Definition at line 201 of file G4VDecayChannel.hh.

202{
203 // pointers to daughter particles are filled, if they are not set yet
205
206 // get the pointer to a daughter particle
207 if ( (anIndex>=0) && (anIndex<numberOfDaughters) )
208 {
209 return G4MT_daughters[anIndex];
210 }
211 else
212 {
213 if (verboseLevel>0)
214 G4cout << "G4VDecayChannel::GetDaughter index out of range "
215 << anIndex << G4endl;
216 return nullptr;
217 }
218}

References G4VDecayChannel::CheckAndFillDaughters(), G4cout, G4endl, G4VDecayChannel::G4MT_daughters, G4VDecayChannel::numberOfDaughters, and G4VDecayChannel::verboseLevel.

Referenced by G4IonTable::CreateIon(), G4KineticTrack::Decay(), G4KineticTrack::G4KineticTrack(), G4HtmlPPReporter::GeneratePropertyTable(), G4TextPPReporter::GeneratePropertyTable(), G4NuclearDecay::GetDaughterNucleus(), G4MTRunManagerKernel::SetUpDecayChannels(), and G4TaskRunManagerKernel::SetUpDecayChannels().

◆ GetDaughterExcitation()

G4double G4NuclearDecay::GetDaughterExcitation ( )
inlineinherited

Definition at line 55 of file G4NuclearDecay.hh.

55{return daughterEx;}
const G4double daughterEx

References G4NuclearDecay::daughterEx.

Referenced by G4Radioactivation::CalculateChainsFromParent().

◆ GetDaughterMass()

G4double G4VDecayChannel::GetDaughterMass ( G4int  anIndex) const
inlineinherited

Definition at line 239 of file G4VDecayChannel.hh.

240{
241 if ( (anIndex>=0) && (anIndex<numberOfDaughters) )
242 {
243 return G4MT_daughters_mass[anIndex];
244 }
245 else
246 {
247 if (verboseLevel>0)
248 {
249 G4cout << "G4VDecayChannel::GetDaughterMass ";
250 G4cout << "index out of range " << anIndex << G4endl;
251 }
252 return 0.0;
253 }
254}

References G4cout, G4endl, G4VDecayChannel::G4MT_daughters_mass, G4VDecayChannel::numberOfDaughters, and G4VDecayChannel::verboseLevel.

◆ GetDaughterName()

const G4String & G4VDecayChannel::GetDaughterName ( G4int  anIndex) const
inlineinherited

◆ GetDaughterNucleus()

G4ParticleDefinition * G4NuclearDecay::GetDaughterNucleus ( )
inlineinherited

Definition at line 59 of file G4NuclearDecay.hh.

59{return GetDaughter(0);}
G4ParticleDefinition * GetDaughter(G4int anIndex)

References G4VDecayChannel::GetDaughter().

Referenced by G4Radioactivation::CalculateChainsFromParent().

◆ GetDecayMode()

G4RadioactiveDecayMode G4NuclearDecay::GetDecayMode ( )
inlineinherited

◆ GetFloatingLevel()

G4Ions::G4FloatLevelBase G4NuclearDecay::GetFloatingLevel ( )
inlineinherited

Definition at line 57 of file G4NuclearDecay.hh.

57{return floatingLevel;}
const G4Ions::G4FloatLevelBase floatingLevel

References G4NuclearDecay::floatingLevel.

◆ GetHLThreshold()

G4double G4NuclearDecay::GetHLThreshold ( )
inlineinherited

Definition at line 62 of file G4NuclearDecay.hh.

62{return halflifeThreshold;}
G4double halflifeThreshold

References G4NuclearDecay::halflifeThreshold.

◆ GetKinematicsName()

const G4String & G4VDecayChannel::GetKinematicsName ( ) const
inlineinherited

◆ GetNoName()

const G4String & G4VDecayChannel::GetNoName ( ) const
privateinherited

Definition at line 579 of file G4VDecayChannel.cc.

580{
581 return noName;
582}
static const G4String noName

References G4VDecayChannel::noName.

Referenced by G4VDecayChannel::GetDaughterName().

◆ GetNumberOfDaughters()

G4int G4VDecayChannel::GetNumberOfDaughters ( ) const
inlineinherited

◆ GetParent()

G4ParticleDefinition * G4VDecayChannel::GetParent ( )
inlineinherited

Definition at line 257 of file G4VDecayChannel.hh.

258{
259 // the pointer to the parent particle is filled, if it is not set yet
261 // get the pointer to the parent particle
262 return G4MT_parent;
263}

References G4VDecayChannel::CheckAndFillParent(), and G4VDecayChannel::G4MT_parent.

Referenced by G4DecayTable::Insert().

◆ GetParentMass()

G4double G4VDecayChannel::GetParentMass ( ) const
inlineinherited

Definition at line 272 of file G4VDecayChannel.hh.

273{
274 return G4MT_parent_mass;
275}

References G4VDecayChannel::G4MT_parent_mass.

◆ GetParentName()

const G4String & G4VDecayChannel::GetParentName ( ) const
inlineinherited

◆ GetPolarization()

const G4ThreeVector & G4VDecayChannel::GetPolarization ( ) const
inlineinherited

Definition at line 331 of file G4VDecayChannel.hh.

332{
333 return parent_polarization;
334}
G4ThreeVector parent_polarization

References G4VDecayChannel::parent_polarization.

◆ GetRangeMass()

G4double G4VDecayChannel::GetRangeMass ( ) const
inlineinherited

Definition at line 316 of file G4VDecayChannel.hh.

317{
318 return rangeMass;
319}

References G4VDecayChannel::rangeMass.

◆ GetVerboseLevel()

G4int G4VDecayChannel::GetVerboseLevel ( ) const
inlineinherited

◆ IsOKWithParentMass()

G4bool G4VDecayChannel::IsOKWithParentMass ( G4double  parentMass)
virtualinherited

Reimplemented in G4PhaseSpaceDecayChannel.

Definition at line 607 of file G4VDecayChannel.cc.

608{
609 G4double sumOfDaughterMassMin = 0.0;
612 // skip one body decay
613 if (numberOfDaughters==1) return true;
614
615 for (G4int index=0; index<numberOfDaughters; ++index)
616 {
617 sumOfDaughterMassMin +=
619 }
620 return (parentMass >= sumOfDaughterMassMin);
621}

References G4VDecayChannel::CheckAndFillDaughters(), G4VDecayChannel::CheckAndFillParent(), G4VDecayChannel::G4MT_daughters_mass, G4VDecayChannel::G4MT_daughters_width, G4VDecayChannel::numberOfDaughters, and G4VDecayChannel::rangeMass.

Referenced by G4Decay::DecayIt(), G4MuonicAtomDecay::DecayIt(), and G4PhaseSpaceDecayChannel::IsOKWithParentMass().

◆ operator!=()

G4bool G4VDecayChannel::operator!= ( const G4VDecayChannel r) const
inlineinherited

Definition at line 69 of file G4VDecayChannel.hh.

69{ return (this != &r); }

◆ operator<()

G4bool G4VDecayChannel::operator< ( const G4VDecayChannel right) const
inlineinherited

Definition at line 195 of file G4VDecayChannel.hh.

196{
197 return (this->rbranch < right.rbranch);
198}

References G4VDecayChannel::rbranch.

◆ operator==()

G4bool G4VDecayChannel::operator== ( const G4VDecayChannel r) const
inlineinherited

Definition at line 68 of file G4VDecayChannel.hh.

68{ return (this == &r); }

◆ SetARM()

void G4ECDecay::SetARM ( G4bool  onoff)
inline

Definition at line 56 of file G4ECDecay.hh.

56{applyARM = onoff;}

References applyARM.

Referenced by G4RadioactiveDecay::LoadDecayTable().

◆ SetBR()

void G4VDecayChannel::SetBR ( G4double  value)
inherited

◆ SetDaughter() [1/2]

void G4VDecayChannel::SetDaughter ( G4int  anIndex,
const G4ParticleDefinition particle_type 
)
inherited

◆ SetDaughter() [2/2]

void G4VDecayChannel::SetDaughter ( G4int  anIndex,
const G4String particle_name 
)
inherited

Definition at line 234 of file G4VDecayChannel.cc.

235{
236 // check numberOfDaughters is positive
237 if (numberOfDaughters<=0)
238 {
239#ifdef G4VERBOSE
240 if (verboseLevel>0)
241 {
242 G4cout << "G4VDecayChannel::SetDaughter() - "
243 << "Number of daughters is not defined" << G4endl;
244 }
245#endif
246 return;
247 }
248
249 // An analysis of this code, shows that this method is called
250 // only in the constructor of derived classes.
251 // The general idea of this method is probably to support
252 // the possibility to re-define daughters on the fly, however
253 // this design is extremely problematic for MT mode, we thus
254 // require (as practically happens) that the method is called only
255 // at construction, i.e. when G4MT_daughters == 0
256 // moreover this method can be called only after SetNumberOfDaugthers()
257 // has been called (see previous if), in such a case daughters_name != nullptr
258 //
259 if ( daughters_name == nullptr )
260 {
261 G4Exception("G4VDecayChannel::SetDaughter()", "PART112", FatalException,
262 "Trying to add a daughter without specifying number of secondaries!");
263 return;
264 }
265 if ( G4MT_daughters != nullptr )
266 {
267 G4Exception("G4VDecayChannel::SetDaughter()", "PART111", FatalException,
268 "Trying to modify a daughter of a decay channel, \
269 but decay channel already has daughters.");
270 return;
271 }
272
273 // check an index
274 if ( (anIndex<0) || (anIndex>=numberOfDaughters) )
275 {
276#ifdef G4VERBOSE
277 if (verboseLevel>0)
278 {
279 G4cout << "G4VDecayChannel::SetDaughter() - "
280 << "index out of range " << anIndex << G4endl;
281 }
282#endif
283 }
284 else
285 {
286 // fill the name
287 daughters_name[anIndex] = new G4String(particle_name);
288#ifdef G4VERBOSE
289 if (verboseLevel>1)
290 {
291 G4cout << "G4VDecayChannel::SetDaughter[" << anIndex <<"] :";
292 G4cout << daughters_name[anIndex] << ":"
293 << *daughters_name[anIndex]<<G4endl;
294 }
295#endif
296 }
297}

References G4VDecayChannel::daughters_name, FatalException, G4cout, G4endl, G4Exception(), G4VDecayChannel::G4MT_daughters, G4VDecayChannel::numberOfDaughters, and G4VDecayChannel::verboseLevel.

◆ SetHLThreshold()

void G4NuclearDecay::SetHLThreshold ( G4double  HLT)
inlineinherited

Definition at line 61 of file G4NuclearDecay.hh.

61{halflifeThreshold = HLT;}

References G4NuclearDecay::halflifeThreshold.

◆ SetNumberOfDaughters()

void G4VDecayChannel::SetNumberOfDaughters ( G4int  value)
inherited

◆ SetParent() [1/2]

void G4VDecayChannel::SetParent ( const G4ParticleDefinition particle_type)
inherited

◆ SetParent() [2/2]

void G4VDecayChannel::SetParent ( const G4String particle_name)
inlineinherited

Definition at line 278 of file G4VDecayChannel.hh.

279{
280 if (parent_name != nullptr) delete parent_name;
281 parent_name = new G4String(particle_name);
282 G4MT_parent = nullptr;
283}

References G4VDecayChannel::G4MT_parent, and G4VDecayChannel::parent_name.

◆ SetPolarization()

void G4VDecayChannel::SetPolarization ( const G4ThreeVector polar)
inlineinherited

◆ SetRangeMass()

void G4VDecayChannel::SetRangeMass ( G4double  val)
inlineinherited

Definition at line 322 of file G4VDecayChannel.hh.

322{ if(val>=0.) rangeMass=val; }

References G4VDecayChannel::rangeMass.

◆ SetVerboseLevel()

void G4VDecayChannel::SetVerboseLevel ( G4int  value)
inlineinherited

Definition at line 304 of file G4VDecayChannel.hh.

305{
306 verboseLevel = value;
307}

References G4VDecayChannel::verboseLevel.

Referenced by G4Decay::DecayIt(), and G4MuonicAtomDecay::DecayIt().

Field Documentation

◆ applyARM

G4bool G4ECDecay::applyARM
private

Definition at line 64 of file G4ECDecay.hh.

Referenced by DecayIt(), and SetARM().

◆ daughterEx

const G4double G4NuclearDecay::daughterEx
privateinherited

Definition at line 71 of file G4NuclearDecay.hh.

Referenced by G4NuclearDecay::GetDaughterExcitation().

◆ daughters_name

G4String** G4VDecayChannel::daughters_name = nullptr
protectedinherited

◆ daughtersMutex

G4Mutex G4VDecayChannel::daughtersMutex
protectedinherited

◆ floatingLevel

const G4Ions::G4FloatLevelBase G4NuclearDecay::floatingLevel
privateinherited

Definition at line 72 of file G4NuclearDecay.hh.

Referenced by G4NuclearDecay::GetFloatingLevel().

◆ G4MT_daughters

G4ParticleDefinition** G4VDecayChannel::G4MT_daughters = nullptr
protectedinherited

◆ G4MT_daughters_mass

G4double* G4VDecayChannel::G4MT_daughters_mass = nullptr
protectedinherited

◆ G4MT_daughters_width

G4double* G4VDecayChannel::G4MT_daughters_width = nullptr
protectedinherited

◆ G4MT_parent

G4ParticleDefinition* G4VDecayChannel::G4MT_parent = nullptr
protectedinherited

◆ G4MT_parent_mass

G4double G4VDecayChannel::G4MT_parent_mass = 0.0
protectedinherited

◆ halflifeThreshold

G4double G4NuclearDecay::halflifeThreshold
privateinherited

◆ kinematics_name

G4String G4VDecayChannel::kinematics_name = ""
protectedinherited

◆ noName

const G4String G4VDecayChannel::noName = " "
staticprotectedinherited

Definition at line 170 of file G4VDecayChannel.hh.

Referenced by G4VDecayChannel::GetNoName().

◆ numberOfDaughters

G4int G4VDecayChannel::numberOfDaughters = 0
protectedinherited

◆ parent_name

G4String* G4VDecayChannel::parent_name = nullptr
protectedinherited

◆ parent_polarization

G4ThreeVector G4VDecayChannel::parent_polarization
protectedinherited

◆ parentMutex

G4Mutex G4VDecayChannel::parentMutex
protectedinherited

◆ particletable

G4ParticleTable* G4VDecayChannel::particletable = nullptr
protectedinherited

◆ PL1

G4double G4ECDecay::PL1
private

Definition at line 65 of file G4ECDecay.hh.

Referenced by DecayIt(), and DefineSubshellProbabilities().

◆ PL2

G4double G4ECDecay::PL2
private

Definition at line 65 of file G4ECDecay.hh.

Referenced by DecayIt(), and DefineSubshellProbabilities().

◆ PL2overPL1

const G4double G4ECDecay::PL2overPL1
staticprivate
Initial value:
= {
0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 1.8722e-04,
2.6438e-04, 3.5456e-04, 4.5790e-04, 6.1588e-04, 7.8944e-04, 9.8530e-04, 1.2030e-03,
1.4361e-03, 1.6886e-03, 1.9609e-03, 2.2641e-03, 2.5674e-03, 2.9019e-03, 3.2577e-03,
3.6338e-03, 4.0310e-03, 4.4541e-03, 4.8943e-03, 5.3620e-03, 5.8523e-03, 6.3650e-03,
6.9061e-03, 7.4607e-03, 8.0398e-03, 8.6417e-03, 9.2665e-03, 9.9150e-03, 1.0588e-02,
1.1284e-02, 1.2004e-02, 1.2744e-02, 1.3518e-02, 1.4312e-02, 1.5136e-02, 1.5981e-02,
1.6857e-02, 1.7764e-02, 1.8696e-02, 1.9682e-02, 2.0642e-02, 2.1661e-02, 2.2708e-02,
2.3788e-02, 2.4896e-02, 2.6036e-02, 2.7217e-02, 2.8409e-02, 2.9646e-02, 3.0917e-02,
3.2220e-02, 3.3561e-02, 3.4937e-02, 3.6353e-02, 3.7805e-02, 3.9297e-02, 4.0826e-02,
4.2399e-02, 4.4010e-02, 4.5668e-02, 4.7368e-02, 4.9115e-02, 5.0896e-02, 5.2744e-02,
5.4625e-02, 5.6565e-02, 5.8547e-02, 6.0593e-02, 6.2690e-02, 6.4844e-02, 6.7068e-02,
6.9336e-02, 7.1667e-02, 7.4075e-02, 7.6544e-02, 7.9085e-02, 8.1688e-02, 8.4371e-02,
8.7135e-02, 8.9995e-02, 9.2919e-02, 9.5949e-02, 9.9036e-02, 1.0226e-01, 1.0555e-01,
1.0899e-01, 1.1249e-01, 1.1613e-01, 1.1989e-01, 1.2379e-01, 1.2780e-01, 1.3196e-01,
1.3627e-01, 1.4071e-01}

Definition at line 68 of file G4ECDecay.hh.

Referenced by DefineSubshellProbabilities().

◆ PM1

G4double G4ECDecay::PM1
private

Definition at line 65 of file G4ECDecay.hh.

Referenced by DecayIt(), and DefineSubshellProbabilities().

◆ PM2

G4double G4ECDecay::PM2
private

Definition at line 65 of file G4ECDecay.hh.

Referenced by DecayIt(), and DefineSubshellProbabilities().

◆ PM2overPM1

const G4double G4ECDecay::PM2overPM1
staticprivate
Initial value:
= {
0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00,
0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00,
1.0210e-03, 1.2641e-03, 1.5231e-03, 1.7990e-03, 2.1938e-03, 2.5863e-03, 2.9621e-03,
3.3637e-03, 3.7909e-03, 4.2049e-03, 4.7021e-03, 5.1791e-03, 5.6766e-03, 6.1952e-03,
6.7045e-03, 7.2997e-03, 7.9438e-03, 8.6271e-03, 9.3294e-03, 1.0058e-02, 1.0813e-02,
1.1594e-02, 1.2408e-02, 1.3244e-02, 1.4118e-02, 1.5023e-02, 1.5962e-02, 1.6919e-02,
1.7910e-02, 1.8934e-02, 1.9986e-02, 2.1072e-02, 2.2186e-02, 2.3336e-02, 2.4524e-02,
2.5750e-02, 2.7006e-02, 2.8302e-02, 2.9629e-02, 3.0994e-02, 3.2399e-02, 3.3845e-02,
3.5328e-02, 3.6852e-02, 3.8414e-02, 4.0025e-02, 4.1673e-02, 4.3368e-02, 4.5123e-02,
4.6909e-02, 4.8767e-02, 5.0662e-02, 5.2612e-02, 5.4612e-02, 5.6662e-02, 5.8773e-02,
6.0930e-02, 6.3141e-02, 6.5413e-02, 6.7752e-02, 7.0139e-02, 7.2603e-02, 7.5127e-02,
7.7721e-02, 8.0408e-02, 8.3128e-02, 8.5949e-02, 8.8843e-02, 9.1824e-02, 9.4888e-02,
9.8025e-02, 1.0130e-01, 1.0463e-01, 1.0806e-01, 1.1159e-01, 1.1526e-01, 1.1900e-01,
1.2290e-01, 1.2688e-01, 1.3101e-01, 1.3528e-01, 1.3969e-01, 1.4425e-01, 1.4896e-01,
1.5384e-01, 1.5887e-01}

Definition at line 69 of file G4ECDecay.hh.

Referenced by DefineSubshellProbabilities().

◆ PN1

G4double G4ECDecay::PN1
private

Definition at line 65 of file G4ECDecay.hh.

Referenced by DecayIt(), and DefineSubshellProbabilities().

◆ PN2

G4double G4ECDecay::PN2
private

Definition at line 65 of file G4ECDecay.hh.

Referenced by DecayIt(), and DefineSubshellProbabilities().

◆ PN2overPN1

const G4double G4ECDecay::PN2overPN1
staticprivate
Initial value:
= {
0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00,
0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00,
0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00,
0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00,
0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00,
0.0000e+00, 9.6988e-03, 1.0797e-02, 1.1706e-02, 1.2603e-02, 1.3408e-02, 1.4352e-02,
1.5511e-02, 1.6579e-02, 1.7646e-02, 1.8731e-02, 1.9886e-02, 2.1069e-02, 2.2359e-02,
2.3710e-02, 2.5058e-02, 2.6438e-02, 2.7843e-02, 2.9283e-02, 3.0762e-02, 3.2275e-02,
3.3843e-02, 3.5377e-02, 3.6886e-02, 3.8502e-02, 4.0159e-02, 4.1867e-02, 4.3617e-02,
4.5470e-02, 4.7247e-02, 4.9138e-02, 5.1065e-02, 5.3049e-02, 5.5085e-02, 5.7173e-02,
5.9366e-02, 6.1800e-02, 6.3945e-02, 6.6333e-02, 6.8785e-02, 7.1303e-02, 7.3801e-02,
7.6538e-02, 7.9276e-02, 8.2070e-02, 8.4959e-02, 8.7940e-02, 9.0990e-02, 9.4124e-02,
9.7337e-02, 1.0069e-01, 1.0410e-01, 1.0761e-01, 1.1122e-01, 1.1499e-01, 1.1882e-01,
1.2282e-01, 1.2709e-01, 1.3114e-01, 1.3549e-01, 1.4001e-01, 1.4465e-01, 1.4946e-01,
1.5443e-01, 1.5954e-01}

Definition at line 70 of file G4ECDecay.hh.

Referenced by DefineSubshellProbabilities().

◆ rangeMass

G4double G4VDecayChannel::rangeMass = 2.5
protectedinherited

◆ rbranch

G4double G4VDecayChannel::rbranch = 0.0
protectedinherited

◆ theMode

const G4RadioactiveDecayMode G4NuclearDecay::theMode
protectedinherited

Definition at line 67 of file G4NuclearDecay.hh.

Referenced by DecayIt(), DumpNuclearInfo(), and G4NuclearDecay::GetDecayMode().

◆ transitionQ

const G4double G4ECDecay::transitionQ
private

Definition at line 63 of file G4ECDecay.hh.

Referenced by DecayIt(), and DumpNuclearInfo().

◆ verboseLevel

G4int G4VDecayChannel::verboseLevel = 1
protectedinherited

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