Geant4  10.02.p03
G4Decay Class Reference

#include <G4Decay.hh>

Inheritance diagram for G4Decay:
Collaboration diagram for G4Decay:

Public Member Functions

 G4Decay (const G4String &processName="Decay")
 
virtual ~G4Decay ()
 
virtual G4VParticleChange * PostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual G4VParticleChange * AtRestDoIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual void StartTracking (G4Track *)
 
virtual void EndTracking ()
 
void SetExtDecayer (G4VExtDecayer *)
 
const G4VExtDecayerGetExtDecayer () const
 
G4double GetRemainderLifeTime () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
- Public Member Functions inherited from G4VRestDiscreteProcess
 G4VRestDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VRestDiscreteProcess (G4VRestDiscreteProcess &)
 
virtual ~G4VRestDiscreteProcess ()
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4VParticleChange * AlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

virtual G4VParticleChange * DecayIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual void DaughterPolarization (const G4Track &aTrack, G4DecayProducts *products)
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4double GetMeanLifeTime (const G4Track &aTrack, G4ForceCondition *condition)
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Protected Attributes

G4int verboseLevel
 
const G4double HighestValue
 
G4double fRemainderLifeTime
 
G4ParticleChangeForDecay fParticleChangeForDecay
 
G4VExtDecayerpExtDecayer
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChange * pParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Private Member Functions

 G4Decay (const G4Decay &right)
 
G4Decayoperator= (const G4Decay &right)
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 

Detailed Description

Definition at line 63 of file G4Decay.hh.

Constructor & Destructor Documentation

◆ G4Decay() [1/2]

G4Decay::G4Decay ( const G4String processName = "Decay")

Definition at line 63 of file G4Decay.cc.

64  :G4VRestDiscreteProcess(processName, fDecay),
65  verboseLevel(1),
66  HighestValue(20.0),
67  fRemainderLifeTime(-1.0),
68  pExtDecayer(0)
69 {
70  // set Process Sub Type
71  SetProcessSubType(static_cast<int>(DECAY));
72 
73 #ifdef G4VERBOSE
74  if (GetVerboseLevel()>1) {
75  G4cout << "G4Decay constructor " << " Name:" << processName << G4endl;
76  }
77 #endif
78 
80 }
const G4double HighestValue
Definition: G4Decay.hh:172
G4VExtDecayer * pExtDecayer
Definition: G4Decay.hh:181
G4double fRemainderLifeTime
Definition: G4Decay.hh:175
G4int verboseLevel
Definition: G4Decay.hh:164
G4GLOB_DLL std::ostream G4cout
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
G4int GetVerboseLevel() const
Definition: G4Decay.hh:188
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
#define G4endl
Definition: G4ios.hh:61
G4ParticleChangeForDecay fParticleChangeForDecay
Definition: G4Decay.hh:178
Here is the call graph for this function:

◆ ~G4Decay()

G4Decay::~G4Decay ( )
virtual

Definition at line 82 of file G4Decay.cc.

83 {
84  if (pExtDecayer) {
85  delete pExtDecayer;
86  }
87 }
G4VExtDecayer * pExtDecayer
Definition: G4Decay.hh:181

◆ G4Decay() [2/2]

G4Decay::G4Decay ( const G4Decay right)
private

Member Function Documentation

◆ AtRestDoIt()

G4VParticleChange * G4Decay::AtRestDoIt ( const G4Track &  aTrack,
const G4Step &  aStep 
)
inlinevirtual

Reimplemented from G4VRestDiscreteProcess.

Reimplemented in G4DecayWithSpin.

Definition at line 191 of file G4Decay.hh.

195 {
196  return DecayIt(aTrack, aStep);
197 }
virtual G4VParticleChange * DecayIt(const G4Track &aTrack, const G4Step &aStep)
Definition: G4Decay.cc:181
Here is the call graph for this function:

◆ AtRestGetPhysicalInteractionLength()

G4double G4Decay::AtRestGetPhysicalInteractionLength ( const G4Track &  track,
G4ForceCondition *  condition 
)
virtual

Reimplemented from G4VRestDiscreteProcess.

Definition at line 452 of file G4Decay.cc.

456 {
457  // condition is set to "Not Forced"
458  *condition = NotForced;
459 
460  G4double pTime = track.GetDynamicParticle()->GetPreAssignedDecayProperTime();
461  if (pTime >= 0.) {
462  fRemainderLifeTime = pTime - track.GetProperTime();
464  } else {
467  }
468  return fRemainderLifeTime;
469 }
G4double condition(const G4ErrorSymMatrix &m)
virtual G4double GetMeanLifeTime(const G4Track &aTrack, G4ForceCondition *condition)
Definition: G4Decay.cc:101
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
G4double fRemainderLifeTime
Definition: G4Decay.hh:175
#define DBL_MIN
Definition: templates.hh:75
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:

◆ BuildPhysicsTable()

void G4Decay::BuildPhysicsTable ( const G4ParticleDefinition )
virtual

Reimplemented from G4VProcess.

Definition at line 176 of file G4Decay.cc.

177 {
178  return;
179 }

◆ DaughterPolarization()

void G4Decay::DaughterPolarization ( const G4Track &  aTrack,
G4DecayProducts products 
)
protectedvirtual

Reimplemented in G4PionDecayMakeSpin.

Definition at line 362 of file G4Decay.cc.

363 {
364 }
Here is the caller graph for this function:

◆ DecayIt()

G4VParticleChange * G4Decay::DecayIt ( const G4Track &  aTrack,
const G4Step &  aStep 
)
protectedvirtual

Definition at line 181 of file G4Decay.cc.

182 {
183  // The DecayIt() method returns by pointer a particle-change object.
184  // Units are expressed in GEANT4 internal units.
185 
186  // Initialize ParticleChange
187  // all members of G4VParticleChange are set to equal to
188  // corresponding member in G4Track
189  fParticleChangeForDecay.Initialize(aTrack);
190 
191  // get particle
192  const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
193  const G4ParticleDefinition* aParticleDef = aParticle->GetDefinition();
194 
195  // check if the particle is stable
196  if (aParticleDef->GetPDGStable()) return &fParticleChangeForDecay ;
197 
198 
199  //check if thePreAssignedDecayProducts exists
200  const G4DecayProducts* o_products = (aParticle->GetPreAssignedDecayProducts());
201  G4bool isPreAssigned = (o_products != 0);
202  G4DecayProducts* products = 0;
203 
204  // decay table
205  G4DecayTable *decaytable = aParticleDef->GetDecayTable();
206 
207  // check if external decayer exists
208  G4bool isExtDecayer = (decaytable == 0) && (pExtDecayer !=0);
209 
210  // Error due to NO Decay Table
211  if ( (decaytable == 0) && !isExtDecayer &&!isPreAssigned ){
212  if (GetVerboseLevel()>0) {
213  G4cout << "G4Decay::DoIt : decay table not defined for ";
214  G4cout << aParticle->GetDefinition()->GetParticleName()<< G4endl;
215  }
216  G4Exception( "G4Decay::DecayIt ",
217  "DECAY101",JustWarning,
218  "Decay table is not defined");
219 
220  fParticleChangeForDecay.SetNumberOfSecondaries(0);
221  // Kill the parent particle
222  fParticleChangeForDecay.ProposeTrackStatus( fStopAndKill ) ;
223  fParticleChangeForDecay.ProposeLocalEnergyDeposit(0.0);
224 
226  return &fParticleChangeForDecay ;
227  }
228 
229  if (isPreAssigned) {
230  // copy decay products
231  products = new G4DecayProducts(*o_products);
232  } else if ( isExtDecayer ) {
233  // decay according to external decayer
234  products = pExtDecayer->ImportDecayProducts(aTrack);
235  } else {
236  // Decay according to decay table.
237  // Keep trying to choose a candidate decay channel if the dynamic mass
238  // of the decaying particle is below the sum of the PDG masses of the
239  // candidate daughter particles.
240  // This is needed because the decay table used in Geant4 is based on
241  // the assumption of nominal PDG masses, but a wide resonance can have
242  // a dynamic masses well below its nominal PDG masses, and therefore
243  // some of its decay channels can be below the kinematical threshold.
244  // Note that, for simplicity, we ignore here the possibility that
245  // one or more of the candidate daughter particles can be, in turn,
246  // wide resonance. However, if this is the case, and the channel is
247  // accepted, then the masses of the resonance daughter particles will
248  // be sampled by taking into account their widths.
249  G4VDecayChannel* decaychannel = 0;
250  G4double massParent = aParticle->GetMass();
251  decaychannel = decaytable->SelectADecayChannel(massParent);
252  if ( decaychannel ==0) {
253  // decay channel not found
254  G4Exception("G4Decay::DoIt", "DECAY003", FatalException,
255  " can not determine decay channel ");
256  } else {
257  // execute DecayIt()
258 #ifdef G4VERBOSE
259  G4int temp = decaychannel->GetVerboseLevel();
260  if (GetVerboseLevel()>1) {
261  G4cout << "G4Decay::DoIt : selected decay channel addr:"
262  << decaychannel <<G4endl;
263  decaychannel->SetVerboseLevel(GetVerboseLevel());
264  }
265 #endif
266  products = decaychannel->DecayIt(aParticle->GetMass());
267 #ifdef G4VERBOSE
268  if (GetVerboseLevel()>1) {
269  decaychannel->SetVerboseLevel(temp);
270  }
271 #endif
272 #ifdef G4VERBOSE
273  if (GetVerboseLevel()>2) {
274  if (! products->IsChecked() ) products->DumpInfo();
275  }
276 #endif
277  }
278  }
279 
280  // get parent particle information ...................................
281  G4double ParentEnergy = aParticle->GetTotalEnergy();
282  G4double ParentMass = aParticle->GetMass();
283  if (ParentEnergy < ParentMass) {
284  if (GetVerboseLevel()>0) {
285  G4cout << "G4Decay::DoIt : Total Energy is less than its mass" << G4endl;
286  G4cout << " Particle: " << aParticle->GetDefinition()->GetParticleName();
287  G4cout << " Energy:" << ParentEnergy/MeV << "[MeV]";
288  G4cout << " Mass:" << ParentMass/MeV << "[MeV]";
289  G4cout << G4endl;
290  }
291  G4Exception( "G4Decay::DecayIt ",
292  "DECAY102",JustWarning,
293  "Total Energy is less than its mass");
294  ParentEnergy = ParentMass;
295  }
296 
297  G4ThreeVector ParentDirection(aParticle->GetMomentumDirection());
298 
299  //boost all decay products to laboratory frame
300  G4double energyDeposit = 0.0;
301  G4double finalGlobalTime = aTrack.GetGlobalTime();
302  G4double finalLocalTime = aTrack.GetLocalTime();
303  if (aTrack.GetTrackStatus() == fStopButAlive ){
304  // AtRest case
305  finalGlobalTime += fRemainderLifeTime;
306  finalLocalTime += fRemainderLifeTime;
307  energyDeposit += aParticle->GetKineticEnergy();
308  if (isPreAssigned) products->Boost( ParentEnergy, ParentDirection);
309  } else {
310  // PostStep case
311  if (!isExtDecayer) products->Boost( ParentEnergy, ParentDirection);
312  }
313  // set polarization for daughter particles
314  DaughterPolarization(aTrack, products);
315 
316 
317  //add products in fParticleChangeForDecay
318  G4int numberOfSecondaries = products->entries();
319  fParticleChangeForDecay.SetNumberOfSecondaries(numberOfSecondaries);
320 #ifdef G4VERBOSE
321  if (GetVerboseLevel()>1) {
322  G4cout << "G4Decay::DoIt : Decay vertex :";
323  G4cout << " Time: " << finalGlobalTime/ns << "[ns]";
324  G4cout << " X:" << (aTrack.GetPosition()).x() /cm << "[cm]";
325  G4cout << " Y:" << (aTrack.GetPosition()).y() /cm << "[cm]";
326  G4cout << " Z:" << (aTrack.GetPosition()).z() /cm << "[cm]";
327  G4cout << G4endl;
328  G4cout << "G4Decay::DoIt : decay products in Lab. Frame" << G4endl;
329  products->DumpInfo();
330  }
331 #endif
332  G4int index;
333  G4ThreeVector currentPosition;
334  const G4TouchableHandle thand = aTrack.GetTouchableHandle();
335  for (index=0; index < numberOfSecondaries; index++)
336  {
337  // get current position of the track
338  currentPosition = aTrack.GetPosition();
339  // create a new track object
340  G4Track* secondary = new G4Track( products->PopProducts(),
341  finalGlobalTime ,
342  currentPosition );
343  // switch on good for tracking flag
344  secondary->SetGoodForTrackingFlag();
345  secondary->SetTouchableHandle(thand);
346  // add the secondary track in the List
347  fParticleChangeForDecay.AddSecondary(secondary);
348  }
349  delete products;
350 
351  // Kill the parent particle
352  fParticleChangeForDecay.ProposeTrackStatus( fStopAndKill ) ;
353  fParticleChangeForDecay.ProposeLocalEnergyDeposit(energyDeposit);
354  fParticleChangeForDecay.ProposeLocalTime( finalLocalTime );
355 
356  // Clear NumberOfInteractionLengthLeft
358 
359  return &fParticleChangeForDecay ;
360 }
static const double cm
Definition: G4SIunits.hh:118
G4double GetMass() const
static const double MeV
Definition: G4SIunits.hh:211
Int_t index
virtual G4DecayProducts * ImportDecayProducts(const G4Track &aTrack)=0
G4VExtDecayer * pExtDecayer
Definition: G4Decay.hh:181
void Boost(G4double totalEnergy, const G4ThreeVector &momentumDirection)
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:447
G4bool GetPDGStable() const
int G4int
Definition: G4Types.hh:78
G4double GetTotalEnergy() const
const G4DecayProducts * GetPreAssignedDecayProducts() const
G4bool IsChecked() const
G4double fRemainderLifeTime
Definition: G4Decay.hh:175
Double_t y
G4double GetKineticEnergy() const
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
void SetVerboseLevel(G4int value)
bool G4bool
Definition: G4Types.hh:79
G4VDecayChannel * SelectADecayChannel(G4double parentMass=-1.)
Definition: G4DecayTable.cc:81
G4DecayTable * GetDecayTable() const
G4int entries() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4int GetVerboseLevel() const
Definition: G4Decay.hh:188
const G4ThreeVector & GetMomentumDirection() const
G4DynamicParticle * PopProducts()
void DumpInfo() const
virtual void DaughterPolarization(const G4Track &aTrack, G4DecayProducts *products)
Definition: G4Decay.cc:362
G4ParticleDefinition * GetDefinition() const
#define G4endl
Definition: G4ios.hh:61
G4int GetVerboseLevel() const
virtual G4DecayProducts * DecayIt(G4double parentMass=-1.0)=0
double G4double
Definition: G4Types.hh:76
#define ns
Definition: xmlparse.cc:614
G4ParticleChangeForDecay fParticleChangeForDecay
Definition: G4Decay.hh:178
Here is the call graph for this function:
Here is the caller graph for this function:

◆ EndTracking()

void G4Decay::EndTracking ( )
virtual

Reimplemented from G4VProcess.

Definition at line 376 of file G4Decay.cc.

377 {
378  // Clear NumberOfInteractionLengthLeft
380 
382 }
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:447
G4double currentInteractionLength
Definition: G4VProcess.hh:297
Here is the call graph for this function:

◆ GetExtDecayer()

const G4VExtDecayer * G4Decay::GetExtDecayer ( ) const
inline

Definition at line 201 of file G4Decay.hh.

202 {
203  return pExtDecayer;
204 }
G4VExtDecayer * pExtDecayer
Definition: G4Decay.hh:181

◆ GetMeanFreePath()

G4double G4Decay::GetMeanFreePath ( const G4Track &  aTrack,
G4double  previousStepSize,
G4ForceCondition *  condition 
)
protectedvirtual

Implements G4VRestDiscreteProcess.

Definition at line 130 of file G4Decay.cc.

131 {
132  // get particle
133  const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
134  const G4ParticleDefinition* aParticleDef = aParticle->GetDefinition();
135  G4double aMass = aParticle->GetMass();
136  G4double aLife = aParticleDef->GetPDGLifeTime();
137 
138 
139  // returns the mean free path in GEANT4 internal units
140  G4double pathlength;
141  G4double aCtau = c_light * aLife;
142 
143  // check if the particle is stable?
144  if (aParticleDef->GetPDGStable()) {
145  pathlength = DBL_MAX;
146 
147  //check if the particle has very short life time ?
148  } else if (aCtau < DBL_MIN) {
149  pathlength = DBL_MIN;
150 
151  } else {
152  //calculate the mean free path
153  // by using normalized kinetic energy (= Ekin/mass)
154  G4double rKineticEnergy = aParticle->GetKineticEnergy()/aMass;
155  if ( rKineticEnergy > HighestValue) {
156  // gamma >> 1
157  pathlength = ( rKineticEnergy + 1.0)* aCtau;
158  } else if ( rKineticEnergy < DBL_MIN ) {
159  // too slow particle
160 #ifdef G4VERBOSE
161  if (GetVerboseLevel()>1) {
162  G4cout << "G4Decay::GetMeanFreePath() !!particle stops!!";
163  G4cout << aParticleDef->GetParticleName() << G4endl;
164  G4cout << "KineticEnergy:" << aParticle->GetKineticEnergy()/GeV <<"[GeV]";
165  }
166 #endif
167  pathlength = DBL_MIN;
168  } else {
169  // beta <1
170  pathlength = (aParticle->GetTotalMomentum())/aMass*aCtau ;
171  }
172  }
173  return pathlength;
174 }
G4double GetMass() const
const G4double HighestValue
Definition: G4Decay.hh:172
G4double GetTotalMomentum() const
G4double GetPDGLifeTime() const
G4bool GetPDGStable() const
G4double GetKineticEnergy() const
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
static const double GeV
Definition: G4SIunits.hh:214
G4int GetVerboseLevel() const
Definition: G4Decay.hh:188
#define DBL_MIN
Definition: templates.hh:75
G4ParticleDefinition * GetDefinition() const
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
float c_light
Definition: hepunit.py:257
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetMeanLifeTime()

G4double G4Decay::GetMeanLifeTime ( const G4Track &  aTrack,
G4ForceCondition *  condition 
)
protectedvirtual

Implements G4VRestDiscreteProcess.

Definition at line 101 of file G4Decay.cc.

103 {
104  // returns the mean free path in GEANT4 internal units
105  G4double meanlife;
106 
107  // get particle
108  const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
109  const G4ParticleDefinition* aParticleDef = aParticle->GetDefinition();
110  G4double aLife = aParticleDef->GetPDGLifeTime();
111 
112  // check if the particle is stable?
113  if (aParticleDef->GetPDGStable()) {
114  //1000000 times the life time of the universe
115  meanlife = 1e24 * s;
116 
117  } else {
118  meanlife = aLife;
119  }
120 
121 #ifdef G4VERBOSE
122  if (GetVerboseLevel()>1) {
123  G4cout << "mean life time: "<< meanlife/ns << "[ns]" << G4endl;
124  }
125 #endif
126 
127  return meanlife;
128 }
G4double GetPDGLifeTime() const
G4bool GetPDGStable() const
static const double s
Definition: G4SIunits.hh:168
G4GLOB_DLL std::ostream G4cout
G4int GetVerboseLevel() const
Definition: G4Decay.hh:188
G4ParticleDefinition * GetDefinition() const
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
#define ns
Definition: xmlparse.cc:614
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetRemainderLifeTime()

G4double G4Decay::GetRemainderLifeTime ( ) const
inline

Definition at line 207 of file G4Decay.hh.

208 {
209  return fRemainderLifeTime;
210 }
G4double fRemainderLifeTime
Definition: G4Decay.hh:175

◆ GetVerboseLevel()

G4int G4Decay::GetVerboseLevel ( ) const
inline

Definition at line 188 of file G4Decay.hh.

188 { return verboseLevel; }
G4int verboseLevel
Definition: G4Decay.hh:164
Here is the caller graph for this function:

◆ IsApplicable()

G4bool G4Decay::IsApplicable ( const G4ParticleDefinition aParticleType)
virtual

Reimplemented from G4VProcess.

Definition at line 89 of file G4Decay.cc.

90 {
91  // check if the particle is stable?
92  if (aParticleType.GetPDGLifeTime() <0.0) {
93  return false;
94  } else if (aParticleType.GetPDGMass() <= 0.0*MeV) {
95  return false;
96  } else {
97  return true;
98  }
99 }
static const double MeV
Definition: G4SIunits.hh:211
G4double GetPDGLifeTime() const
Here is the call graph for this function:
Here is the caller graph for this function:

◆ operator=()

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

◆ PostStepDoIt()

G4VParticleChange * G4Decay::PostStepDoIt ( const G4Track &  aTrack,
const G4Step &  aStep 
)
virtual

Reimplemented from G4VRestDiscreteProcess.

Reimplemented in G4DecayWithSpin.

Definition at line 482 of file G4Decay.cc.

486 {
487  if ( (aTrack.GetTrackStatus() == fStopButAlive ) ||
488  (aTrack.GetTrackStatus() == fStopAndKill ) ){
489  fParticleChangeForDecay.Initialize(aTrack);
490  return &fParticleChangeForDecay;
491  } else {
492  return DecayIt(aTrack, aStep);
493  }
494 }
virtual G4VParticleChange * DecayIt(const G4Track &aTrack, const G4Step &aStep)
Definition: G4Decay.cc:181
G4ParticleChangeForDecay fParticleChangeForDecay
Definition: G4Decay.hh:178
Here is the call graph for this function:

◆ PostStepGetPhysicalInteractionLength()

G4double G4Decay::PostStepGetPhysicalInteractionLength ( const G4Track &  track,
G4double  previousStepSize,
G4ForceCondition *  condition 
)
virtual

Reimplemented from G4VRestDiscreteProcess.

Definition at line 385 of file G4Decay.cc.

390 {
391  // condition is set to "Not Forced"
392  *condition = NotForced;
393 
394  // pre-assigned Decay time
395  G4double pTime = track.GetDynamicParticle()->GetPreAssignedDecayProperTime();
396  G4double aLife = track.GetDynamicParticle()->GetDefinition()->GetPDGLifeTime();
397 
398  if (pTime < 0.) {
399  // normal case
400  if ( previousStepSize > 0.0){
401  // subtract NumberOfInteractionLengthLeft
402  SubtractNumberOfInteractionLengthLeft(previousStepSize);
405  }
407  }
408  // get mean free path
410 
411 #ifdef G4VERBOSE
412  if ((currentInteractionLength <=0.0) || (verboseLevel>2)){
413  G4cout << "G4Decay::PostStepGetPhysicalInteractionLength " << G4endl;
414  track.GetDynamicParticle()->DumpInfo();
415  G4cout << " in Material " << track.GetMaterial()->GetName() <<G4endl;
416  G4cout << "MeanFreePath = " << currentInteractionLength/cm << "[cm]" <<G4endl;
417  }
418 #endif
419 
420  G4double value;
423  //fRemainderLifeTime = theNumberOfInteractionLengthLeft*aLife;
424  } else {
425  value = DBL_MAX;
426  }
427 
428  return value;
429 
430  } else {
431  //pre-assigned Decay time case
432  // reminder proper time
433  fRemainderLifeTime = pTime - track.GetProperTime();
435 
436  G4double rvalue=0.0;
437  // use pre-assigned Decay time to determine PIL
438  if (aLife>0.0) {
439  // ordinary particle
440  rvalue = (fRemainderLifeTime/aLife)*GetMeanFreePath(track, previousStepSize, condition);
441  } else {
442  // shortlived particle
443  rvalue = c_light * fRemainderLifeTime;
444  // by using normalized kinetic energy (= Ekin/mass)
445  G4double aMass = track.GetDynamicParticle()->GetMass();
446  rvalue *= track.GetDynamicParticle()->GetTotalMomentum()/aMass;
447  }
448  return rvalue;
449  }
450 }
static const double cm
Definition: G4SIunits.hh:118
G4double condition(const G4ErrorSymMatrix &m)
virtual G4double GetMeanFreePath(const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
Definition: G4Decay.cc:130
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
G4double fRemainderLifeTime
Definition: G4Decay.hh:175
G4int verboseLevel
Definition: G4Decay.hh:164
G4GLOB_DLL std::ostream G4cout
G4double currentInteractionLength
Definition: G4VProcess.hh:297
static const double perMillion
Definition: G4SIunits.hh:331
#define DBL_MIN
Definition: templates.hh:75
#define G4endl
Definition: G4ios.hh:61
void SubtractNumberOfInteractionLengthLeft(G4double previousStepSize)
Definition: G4VProcess.hh:544
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
float c_light
Definition: hepunit.py:257
Here is the call graph for this function:

◆ SetExtDecayer()

void G4Decay::SetExtDecayer ( G4VExtDecayer val)

Definition at line 472 of file G4Decay.cc.

473 {
474  pExtDecayer = val;
475 
476  // set Process Sub Type
477  if ( pExtDecayer !=0 ) {
478  SetProcessSubType(static_cast<int>(DECAY_External));
479  }
480 }
G4VExtDecayer * pExtDecayer
Definition: G4Decay.hh:181
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
Here is the call graph for this function:
Here is the caller graph for this function:

◆ SetVerboseLevel()

void G4Decay::SetVerboseLevel ( G4int  value)
inline

Definition at line 185 of file G4Decay.hh.

185 { verboseLevel = value; }
G4int verboseLevel
Definition: G4Decay.hh:164

◆ StartTracking()

void G4Decay::StartTracking ( G4Track *  )
virtual

Reimplemented from G4VProcess.

Definition at line 368 of file G4Decay.cc.

369 {
372 
373  fRemainderLifeTime = -1.0;
374 }
virtual void ResetNumberOfInteractionLengthLeft()
Definition: G4VProcess.cc:95
G4double fRemainderLifeTime
Definition: G4Decay.hh:175
G4double currentInteractionLength
Definition: G4VProcess.hh:297
Here is the call graph for this function:

Member Data Documentation

◆ fParticleChangeForDecay

G4ParticleChangeForDecay G4Decay::fParticleChangeForDecay
protected

Definition at line 178 of file G4Decay.hh.

◆ fRemainderLifeTime

G4double G4Decay::fRemainderLifeTime
protected

Definition at line 175 of file G4Decay.hh.

◆ HighestValue

const G4double G4Decay::HighestValue
protected

Definition at line 172 of file G4Decay.hh.

◆ pExtDecayer

G4VExtDecayer* G4Decay::pExtDecayer
protected

Definition at line 181 of file G4Decay.hh.

◆ verboseLevel

G4int G4Decay::verboseLevel
protected

Definition at line 164 of file G4Decay.hh.


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