Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4Decay Class Reference

#include <G4Decay.hh>

+ Inheritance diagram for G4Decay:

Public Member Functions

 G4Decay (const G4String &processName="Decay")
 
virtual ~G4Decay ()
 
virtual G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual G4VParticleChangeAtRestDoIt (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 PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4VParticleChangeAlongStepDoIt (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
 
virtual G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)=0
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)=0
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)=0
 
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 G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
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 StartTracking (G4Track *)
 
virtual void EndTracking ()
 
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
 

Protected Member Functions

virtual G4VParticleChangeDecayIt (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)
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
 
virtual G4double GetMeanLifeTime (const G4Track &aTrack, G4ForceCondition *condition)=0
 
- 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
 
G4VParticleChangepParticleChange
 
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
 

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

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),
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}
@ fDecay
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
G4VExtDecayer * pExtDecayer
Definition: G4Decay.hh:181
G4double fRemainderLifeTime
Definition: G4Decay.hh:175
G4ParticleChangeForDecay fParticleChangeForDecay
Definition: G4Decay.hh:178
G4int GetVerboseLevel() const
Definition: G4Decay.hh:188
const G4double HighestValue
Definition: G4Decay.hh:172
G4int verboseLevel
Definition: G4Decay.hh:164
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:403
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283

◆ ~G4Decay()

G4Decay::~G4Decay ( )
virtual

Definition at line 82 of file G4Decay.cc.

83{
84 if (pExtDecayer) {
85 delete pExtDecayer;
86 }
87}

Member Function Documentation

◆ AtRestDoIt()

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

Reimplemented from G4VRestDiscreteProcess.

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

◆ AtRestGetPhysicalInteractionLength()

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

Reimplemented from G4VRestDiscreteProcess.

Definition at line 438 of file G4Decay.cc.

442{
443 // condition is set to "Not Forced"
445
447 if (pTime >= 0.) {
448 fRemainderLifeTime = pTime - track.GetProperTime();
450 } else {
453 }
454 return fRemainderLifeTime;
455}
G4double condition(const G4ErrorSymMatrix &m)
@ NotForced
double G4double
Definition: G4Types.hh:64
virtual G4double GetMeanLifeTime(const G4Track &aTrack, G4ForceCondition *condition)
Definition: G4Decay.cc:101
G4double GetPreAssignedDecayProperTime() const
G4double GetProperTime() const
const G4DynamicParticle * GetDynamicParticle() const
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
#define DBL_MIN
Definition: templates.hh:75

◆ BuildPhysicsTable()

void G4Decay::BuildPhysicsTable ( const G4ParticleDefinition )
virtual

Reimplemented from G4VProcess.

Definition at line 175 of file G4Decay.cc.

176{
177 return;
178}

◆ DaughterPolarization()

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

Reimplemented in G4PionDecayMakeSpin.

Definition at line 348 of file G4Decay.cc.

349{
350}

Referenced by DecayIt().

◆ DecayIt()

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

Reimplemented in G4DecayWithSpin.

Definition at line 180 of file G4Decay.cc.

181{
182 // The DecayIt() method returns by pointer a particle-change object.
183 // Units are expressed in GEANT4 internal units.
184
185 // Initialize ParticleChange
186 // all members of G4VParticleChange are set to equal to
187 // corresponding member in G4Track
189
190 // get particle
191 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
192 const G4ParticleDefinition* aParticleDef = aParticle->GetDefinition();
193
194 // check if the particle is stable
195 if (aParticleDef->GetPDGStable()) return &fParticleChangeForDecay ;
196
197
198 //check if thePreAssignedDecayProducts exists
199 const G4DecayProducts* o_products = (aParticle->GetPreAssignedDecayProducts());
200 G4bool isPreAssigned = (o_products != 0);
201 G4DecayProducts* products = 0;
202
203 // decay table
204 G4DecayTable *decaytable = aParticleDef->GetDecayTable();
205
206 // check if external decayer exists
207 G4bool isExtDecayer = (decaytable == 0) && (pExtDecayer !=0);
208
209 // Error due to NO Decay Table
210 if ( (decaytable == 0) && !isExtDecayer &&!isPreAssigned ){
211 if (GetVerboseLevel()>0) {
212 G4cout << "G4Decay::DoIt : decay table not defined for ";
213 G4cout << aParticle->GetDefinition()->GetParticleName()<< G4endl;
214 }
215 G4Exception( "G4Decay::DecayIt ",
216 "DECAY101",JustWarning,
217 "Decay table is not defined");
218
220 // Kill the parent particle
223
226 }
227
228 if (isPreAssigned) {
229 // copy decay products
230 products = new G4DecayProducts(*o_products);
231 } else if ( isExtDecayer ) {
232 // decay according to external decayer
233 products = pExtDecayer->ImportDecayProducts(aTrack);
234 } else {
235 // decay acoording to decay table
236 // choose a decay channel
237 G4VDecayChannel *decaychannel = decaytable->SelectADecayChannel();
238 if (decaychannel == 0 ){
239 // decay channel not found
240 G4Exception("G4Decay::DoIt", "DECAY003", FatalException,
241 " can not determine decay channel ");
242 } else {
243 // execute DecayIt()
244#ifdef G4VERBOSE
245 G4int temp = decaychannel->GetVerboseLevel();
246 if (GetVerboseLevel()>1) {
247 G4cout << "G4Decay::DoIt : selected decay channel addr:" << decaychannel <<G4endl;
248 decaychannel->SetVerboseLevel(GetVerboseLevel());
249 }
250#endif
251 products = decaychannel->DecayIt(aParticle->GetMass());
252#ifdef G4VERBOSE
253 if (GetVerboseLevel()>1) {
254 decaychannel->SetVerboseLevel(temp);
255 }
256#endif
257#ifdef G4VERBOSE
258 if (GetVerboseLevel()>2) {
259 if (! products->IsChecked() ) products->DumpInfo();
260 }
261#endif
262 }
263 }
264
265 // get parent particle information ...................................
266 G4double ParentEnergy = aParticle->GetTotalEnergy();
267 G4double ParentMass = aParticle->GetMass();
268 if (ParentEnergy < ParentMass) {
269 if (GetVerboseLevel()>0) {
270 G4cout << "G4Decay::DoIt : Total Energy is less than its mass" << G4endl;
271 G4cout << " Particle: " << aParticle->GetDefinition()->GetParticleName();
272 G4cout << " Energy:" << ParentEnergy/MeV << "[MeV]";
273 G4cout << " Mass:" << ParentMass/MeV << "[MeV]";
274 G4cout << G4endl;
275 }
276 G4Exception( "G4Decay::DecayIt ",
277 "DECAY102",JustWarning,
278 "Total Energy is less than its mass");
279 ParentEnergy = ParentMass;
280 }
281
282 G4ThreeVector ParentDirection(aParticle->GetMomentumDirection());
283
284 //boost all decay products to laboratory frame
285 G4double energyDeposit = 0.0;
286 G4double finalGlobalTime = aTrack.GetGlobalTime();
287 G4double finalLocalTime = aTrack.GetLocalTime();
288 if (aTrack.GetTrackStatus() == fStopButAlive ){
289 // AtRest case
290 finalGlobalTime += fRemainderLifeTime;
291 finalLocalTime += fRemainderLifeTime;
292 energyDeposit += aParticle->GetKineticEnergy();
293 if (isPreAssigned) products->Boost( ParentEnergy, ParentDirection);
294 } else {
295 // PostStep case
296 if (!isExtDecayer) products->Boost( ParentEnergy, ParentDirection);
297 }
298
299 // set polarization for daughter particles
300 DaughterPolarization(aTrack, products);
301
302
303 //add products in fParticleChangeForDecay
304 G4int numberOfSecondaries = products->entries();
306#ifdef G4VERBOSE
307 if (GetVerboseLevel()>1) {
308 G4cout << "G4Decay::DoIt : Decay vertex :";
309 G4cout << " Time: " << finalGlobalTime/ns << "[ns]";
310 G4cout << " X:" << (aTrack.GetPosition()).x() /cm << "[cm]";
311 G4cout << " Y:" << (aTrack.GetPosition()).y() /cm << "[cm]";
312 G4cout << " Z:" << (aTrack.GetPosition()).z() /cm << "[cm]";
313 G4cout << G4endl;
314 G4cout << "G4Decay::DoIt : decay products in Lab. Frame" << G4endl;
315 products->DumpInfo();
316 }
317#endif
318 G4int index;
319 G4ThreeVector currentPosition;
320 const G4TouchableHandle thand = aTrack.GetTouchableHandle();
321 for (index=0; index < numberOfSecondaries; index++)
322 {
323 // get current position of the track
324 currentPosition = aTrack.GetPosition();
325 // create a new track object
326 G4Track* secondary = new G4Track( products->PopProducts(),
327 finalGlobalTime ,
328 currentPosition );
329 // switch on good for tracking flag
330 secondary->SetGoodForTrackingFlag();
331 secondary->SetTouchableHandle(thand);
332 // add the secondary track in the List
334 }
335 delete products;
336
337 // Kill the parent particle
341
342 // Clear NumberOfInteractionLengthLeft
344
346}
@ JustWarning
@ FatalException
@ fStopAndKill
@ fStopButAlive
int G4int
Definition: G4Types.hh:66
bool G4bool
Definition: G4Types.hh:67
void DumpInfo() const
G4int entries() const
G4DynamicParticle * PopProducts()
G4bool IsChecked() const
void Boost(G4double totalEnergy, const G4ThreeVector &momentumDirection)
G4VDecayChannel * SelectADecayChannel()
Definition: G4DecayTable.cc:81
virtual void DaughterPolarization(const G4Track &aTrack, G4DecayProducts *products)
Definition: G4Decay.cc:348
G4double GetMass() const
const G4ThreeVector & GetMomentumDirection() const
const G4DecayProducts * GetPreAssignedDecayProducts() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
G4double GetTotalEnergy() const
virtual void Initialize(const G4Track &)
G4DecayTable * GetDecayTable() const
const G4String & GetParticleName() const
G4TrackStatus GetTrackStatus() const
const G4ThreeVector & GetPosition() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4double GetGlobalTime() const
G4double GetLocalTime() const
const G4TouchableHandle & GetTouchableHandle() const
void SetGoodForTrackingFlag(G4bool value=true)
void SetVerboseLevel(G4int value)
G4int GetVerboseLevel() const
virtual G4DecayProducts * DecayIt(G4double parentMass=-1.0)=0
virtual G4DecayProducts * ImportDecayProducts(const G4Track &aTrack)=0
void ProposeTrackStatus(G4TrackStatus status)
void AddSecondary(G4Track *aSecondary)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetNumberOfSecondaries(G4int totSecondaries)
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:418
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define ns
Definition: xmlparse.cc:597

Referenced by AtRestDoIt(), G4DecayWithSpin::DecayIt(), and PostStepDoIt().

◆ EndTracking()

void G4Decay::EndTracking ( )
virtual

Reimplemented from G4VProcess.

Definition at line 362 of file G4Decay.cc.

363{
364 // Clear NumberOfInteractionLengthLeft
366
368}
G4double currentInteractionLength
Definition: G4VProcess.hh:297

◆ GetExtDecayer()

const G4VExtDecayer * G4Decay::GetExtDecayer ( ) const
inline

Definition at line 209 of file G4Decay.hh.

210{
211 return pExtDecayer;
212}

◆ GetMeanFreePath()

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

Implements G4VRestDiscreteProcess.

Definition at line 129 of file G4Decay.cc.

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

Referenced by PostStepGetPhysicalInteractionLength().

◆ 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 meanlife = DBL_MAX;
115
116 } else {
117 meanlife = aLife;
118 }
119
120#ifdef G4VERBOSE
121 if (GetVerboseLevel()>1) {
122 G4cout << "mean life time: "<< meanlife/ns << "[ns]" << G4endl;
123 }
124#endif
125
126 return meanlife;
127}

Referenced by AtRestGetPhysicalInteractionLength().

◆ GetRemainderLifeTime()

G4double G4Decay::GetRemainderLifeTime ( ) const
inline

Definition at line 215 of file G4Decay.hh.

216{
217 return fRemainderLifeTime;
218}

◆ GetVerboseLevel()

G4int G4Decay::GetVerboseLevel ( ) const
inline

Definition at line 188 of file G4Decay.hh.

188{ return verboseLevel; }

Referenced by DecayIt(), G4Decay(), GetMeanFreePath(), and GetMeanLifeTime().

◆ 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}

Referenced by G4DecayPhysics::ConstructProcess().

◆ PostStepDoIt()

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

Reimplemented from G4VRestDiscreteProcess.

Definition at line 200 of file G4Decay.hh.

204{
205 return DecayIt(aTrack, aStep);
206}

◆ PostStepGetPhysicalInteractionLength()

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

Reimplemented from G4VRestDiscreteProcess.

Definition at line 371 of file G4Decay.cc.

376{
377
378 // condition is set to "Not Forced"
380
381 // pre-assigned Decay time
384
385 if (pTime < 0.) {
386 // normal case
387 if ( previousStepSize > 0.0){
388 // subtract NumberOfInteractionLengthLeft
392 }
394 }
395 // get mean free path
396 currentInteractionLength = GetMeanFreePath(track, previousStepSize, condition);
397
398#ifdef G4VERBOSE
399 if ((currentInteractionLength <=0.0) || (verboseLevel>2)){
400 G4cout << "G4Decay::PostStepGetPhysicalInteractionLength " << G4endl;
401 track.GetDynamicParticle()->DumpInfo();
402 G4cout << " in Material " << track.GetMaterial()->GetName() <<G4endl;
403 G4cout << "MeanFreePath = " << currentInteractionLength/cm << "[cm]" <<G4endl;
404 }
405#endif
406
407 G4double value;
410 } else {
411 value = DBL_MAX;
412 }
413
414 return value;
415
416 } else {
417 //pre-assigned Decay time case
418 // reminder proper time
419 fRemainderLifeTime = pTime - track.GetProperTime();
421
422 G4double rvalue=0.0;
423 // use pre-assigned Decay time to determine PIL
424 if (aLife>0.0) {
425 // ordinary particle
426 rvalue = (fRemainderLifeTime/aLife)*GetMeanFreePath(track, previousStepSize, condition);
427 } else {
428 // shortlived particle
429 rvalue = c_light * fRemainderLifeTime;
430 // by using normalized kinetic energy (= Ekin/mass)
431 G4double aMass = track.GetDynamicParticle()->GetMass();
432 rvalue *= track.GetDynamicParticle()->GetTotalMomentum()/aMass;
433 }
434 return rvalue;
435 }
436}
virtual G4double GetMeanFreePath(const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
Definition: G4Decay.cc:129
void DumpInfo(G4int mode=0) const
const G4String & GetName() const
Definition: G4Material.hh:177
G4Material * GetMaterial() const
void SubtractNumberOfInteractionLengthLeft(G4double previousStepSize)
Definition: G4VProcess.cc:98

◆ SetExtDecayer()

void G4Decay::SetExtDecayer ( G4VExtDecayer val)

Definition at line 458 of file G4Decay.cc.

459{
460 pExtDecayer = val;
461
462 // set Process Sub Type
463 if ( pExtDecayer !=0 ) {
464 SetProcessSubType(static_cast<int>(DECAY_External));
465 }
466}
@ DECAY_External

◆ SetVerboseLevel()

void G4Decay::SetVerboseLevel ( G4int  value)
inline

Definition at line 185 of file G4Decay.hh.

185{ verboseLevel = value; }

◆ StartTracking()

void G4Decay::StartTracking ( G4Track )
virtual

Reimplemented from G4VProcess.

Definition at line 354 of file G4Decay.cc.

355{
358
359 fRemainderLifeTime = -1.0;
360}
virtual void ResetNumberOfInteractionLengthLeft()
Definition: G4VProcess.cc:92

Member Data Documentation

◆ fParticleChangeForDecay

G4ParticleChangeForDecay G4Decay::fParticleChangeForDecay
protected

Definition at line 178 of file G4Decay.hh.

Referenced by DecayIt(), and G4Decay().

◆ fRemainderLifeTime

◆ HighestValue

const G4double G4Decay::HighestValue
protected

Definition at line 172 of file G4Decay.hh.

Referenced by GetMeanFreePath().

◆ pExtDecayer

G4VExtDecayer* G4Decay::pExtDecayer
protected

Definition at line 181 of file G4Decay.hh.

Referenced by DecayIt(), GetExtDecayer(), SetExtDecayer(), and ~G4Decay().

◆ verboseLevel

G4int G4Decay::verboseLevel
protected

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