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

#include <G4ParticleChangeForDecay.hh>

+ Inheritance diagram for G4ParticleChangeForDecay:

Public Member Functions

 G4ParticleChangeForDecay ()
 
virtual ~G4ParticleChangeForDecay ()
 
G4bool operator== (const G4ParticleChangeForDecay &right) const
 
G4bool operator!= (const G4ParticleChangeForDecay &right) const
 
virtual G4StepUpdateStepForAtRest (G4Step *Step)
 
virtual G4StepUpdateStepForPostStep (G4Step *Step)
 
virtual void Initialize (const G4Track &)
 
void ProposeGlobalTime (G4double t)
 
void ProposeLocalTime (G4double t)
 
G4double GetGlobalTime (G4double timeDelay=0.0) const
 
G4double GetLocalTime (G4double timeDelay=0.0) const
 
const G4ThreeVectorGetPolarization () const
 
void ProposePolarization (G4double Px, G4double Py, G4double Pz)
 
void ProposePolarization (const G4ThreeVector &finalPoralization)
 
virtual void DumpInfo () const
 
virtual G4bool CheckIt (const G4Track &)
 
- Public Member Functions inherited from G4VParticleChange
 G4VParticleChange ()
 
virtual ~G4VParticleChange ()
 
G4bool operator== (const G4VParticleChange &right) const
 
G4bool operator!= (const G4VParticleChange &right) const
 
virtual G4StepUpdateStepForAtRest (G4Step *Step)
 
virtual G4StepUpdateStepForAlongStep (G4Step *Step)
 
virtual G4StepUpdateStepForPostStep (G4Step *Step)
 
virtual void Initialize (const G4Track &)
 
G4double GetTrueStepLength () const
 
void ProposeTrueStepLength (G4double truePathLength)
 
G4double GetLocalEnergyDeposit () const
 
void ProposeLocalEnergyDeposit (G4double anEnergyPart)
 
G4double GetNonIonizingEnergyDeposit () const
 
void ProposeNonIonizingEnergyDeposit (G4double anEnergyPart)
 
G4TrackStatus GetTrackStatus () const
 
void ProposeTrackStatus (G4TrackStatus status)
 
G4SteppingControl GetSteppingControl () const
 
void ProposeSteppingControl (G4SteppingControl StepControlFlag)
 
G4bool GetFirstStepInVolume () const
 
G4bool GetLastStepInVolume () const
 
void ProposeFirstStepInVolume (G4bool flag)
 
void ProposeLastStepInVolume (G4bool flag)
 
void Clear ()
 
void SetNumberOfSecondaries (G4int totSecondaries)
 
G4int GetNumberOfSecondaries () const
 
G4TrackGetSecondary (G4int anIndex) const
 
void AddSecondary (G4Track *aSecondary)
 
G4double GetWeight () const
 
G4double GetParentWeight () const
 
void ProposeWeight (G4double finalWeight)
 
void ProposeParentWeight (G4double finalWeight)
 
void SetSecondaryWeightByProcess (G4bool)
 
G4bool IsSecondaryWeightSetByProcess () const
 
void SetParentWeightByProcess (G4bool)
 
G4bool IsParentWeightSetByProcess () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int vLevel)
 
G4int GetVerboseLevel () const
 
virtual G4bool CheckIt (const G4Track &)
 
void ClearDebugFlag ()
 
void SetDebugFlag ()
 
G4bool GetDebugFlag () const
 

Protected Member Functions

 G4ParticleChangeForDecay (const G4ParticleChangeForDecay &right)
 
G4ParticleChangeForDecayoperator= (const G4ParticleChangeForDecay &right)
 
- Protected Member Functions inherited from G4VParticleChange
 G4VParticleChange (const G4VParticleChange &right)
 
G4VParticleChangeoperator= (const G4VParticleChange &right)
 
G4StepUpdateStepInfo (G4Step *Step)
 
void InitializeTrueStepLength (const G4Track &)
 
void InitializeLocalEnergyDeposit (const G4Track &)
 
void InitializeSteppingControl (const G4Track &)
 
void InitializeParentWeight (const G4Track &)
 
void InitializeParentGlobalTime (const G4Track &)
 
void InitializeStatusChange (const G4Track &)
 
void InitializeSecondaries (const G4Track &)
 
void InitializeStepInVolumeFlags (const G4Track &)
 
G4bool CheckSecondary (G4Track &)
 
G4double GetAccuracyForWarning () const
 
G4double GetAccuracyForException () const
 

Protected Attributes

G4double theGlobalTime0 = 0.0
 
G4double theLocalTime0 = 0.0
 
G4double theTimeChange = 0.0
 
G4ThreeVector thePolarizationChange
 
- Protected Attributes inherited from G4VParticleChange
G4TrackFastVectortheListOfSecondaries = nullptr
 
G4TrackStatus theStatusChange = fAlive
 
G4SteppingControl theSteppingControlFlag = NormalCondition
 
G4double theLocalEnergyDeposit = 0.0
 
G4double theNonIonizingEnergyDeposit = 0.0
 
G4double theTrueStepLength = 0.0
 
G4double theParentWeight = 1.0
 
G4double theParentGlobalTime = 0.0
 
G4int theNumberOfSecondaries = 0
 
G4int theSizeOftheListOfSecondaries = 0
 
G4int verboseLevel = 1
 
G4bool theFirstStepInVolume = false
 
G4bool theLastStepInVolume = false
 
G4bool isParentWeightProposed = false
 
G4bool fSetSecondaryWeightByProcess = false
 
G4bool debugFlag = false
 

Additional Inherited Members

- Static Protected Attributes inherited from G4VParticleChange
static const G4double accuracyForWarning = 1.0e-9
 
static const G4double accuracyForException = 0.001
 

Detailed Description

Definition at line 47 of file G4ParticleChangeForDecay.hh.

Constructor & Destructor Documentation

◆ G4ParticleChangeForDecay() [1/2]

G4ParticleChangeForDecay::G4ParticleChangeForDecay ( )

Definition at line 40 of file G4ParticleChangeForDecay.cc.

◆ ~G4ParticleChangeForDecay()

G4ParticleChangeForDecay::~G4ParticleChangeForDecay ( )
virtual

Definition at line 46 of file G4ParticleChangeForDecay.cc.

47{
48}

◆ G4ParticleChangeForDecay() [2/2]

G4ParticleChangeForDecay::G4ParticleChangeForDecay ( const G4ParticleChangeForDecay right)
protected

Member Function Documentation

◆ CheckIt()

G4bool G4ParticleChangeForDecay::CheckIt ( const G4Track aTrack)
virtual

Reimplemented from G4VParticleChange.

Definition at line 190 of file G4ParticleChangeForDecay.cc.

191{
192 G4bool exitWithError = false;
193
194 G4double accuracy;
195
196 // local time should not go back
197 G4bool itsOK = true;
198 accuracy = -1.0 * (theTimeChange - theLocalTime0) / ns;
199 if(accuracy > accuracyForWarning)
200 {
201 itsOK = false;
202 exitWithError = (accuracy > accuracyForException);
203#ifdef G4VERBOSE
204 G4cout << " G4ParticleChangeForDecay::CheckIt : ";
205 G4cout << "the local time goes back !!"
206 << " Difference: " << accuracy << "[ns] " << G4endl;
207 G4cout << "initial local time " << theLocalTime0 / ns << "[ns] "
208 << "initial global time " << theGlobalTime0 / ns << "[ns] "
209 << G4endl;
211 << " E=" << aTrack.GetKineticEnergy() / MeV
212 << " pos=" << aTrack.GetPosition().x() / m << ", "
213 << aTrack.GetPosition().y() / m << ", "
214 << aTrack.GetPosition().z() / m << G4endl;
215#endif
216 }
217
218 // dump out information of this particle change
219#ifdef G4VERBOSE
220 if(!itsOK) { DumpInfo(); }
221#endif
222
223 // exit with error
224 if(exitWithError)
225 {
226 G4Exception("G4ParticleChangeForDecay::CheckIt()", "TRACK005",
227 EventMustBeAborted, "time was illegal");
228 }
229
230 // correction
231 if(!itsOK)
232 {
233 theTimeChange = aTrack.GetLocalTime();
234 }
235
236 itsOK = (itsOK) && G4VParticleChange::CheckIt(aTrack);
237 return itsOK;
238}
@ EventMustBeAborted
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
double z() const
double x() const
double y() const
const G4String & GetParticleName() const
const G4ThreeVector & GetPosition() const
G4double GetLocalTime() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
virtual G4bool CheckIt(const G4Track &)
static const G4double accuracyForException
static const G4double accuracyForWarning
#define ns
Definition: xmlparse.cc:614

Referenced by UpdateStepForAtRest().

◆ DumpInfo()

void G4ParticleChangeForDecay::DumpInfo ( ) const
virtual

Reimplemented from G4VParticleChange.

Definition at line 174 of file G4ParticleChangeForDecay.cc.

175{
176 // Show header
178
179 G4int oldprc = G4cout.precision(3);
180 G4cout << " proposed local Time (ns) : " << std::setw(20)
181 << theTimeChange / ns << G4endl;
182 G4cout << " initial local Time (ns) : " << std::setw(20)
183 << theLocalTime0 / ns << G4endl;
184 G4cout << " initial global Time (ns) : " << std::setw(20)
185 << theGlobalTime0 / ns << G4endl;
186 G4cout.precision(oldprc);
187}
int G4int
Definition: G4Types.hh:85
virtual void DumpInfo() const

Referenced by CheckIt().

◆ GetGlobalTime()

G4double G4ParticleChangeForDecay::GetGlobalTime ( G4double  timeDelay = 0.0) const
inline

Definition at line 123 of file G4ParticleChangeForDecay.hh.

124{
125 // Convert the time delay to the global time.
126 return theGlobalTime0 + (theTimeChange - theLocalTime0) + timeDelay;
127}

Referenced by UpdateStepForAtRest().

◆ GetLocalTime()

G4double G4ParticleChangeForDecay::GetLocalTime ( G4double  timeDelay = 0.0) const
inline

Definition at line 136 of file G4ParticleChangeForDecay.hh.

137{
138 // Convert the time delay to the local time.
139 return theTimeChange + timeDelay;
140}

◆ GetPolarization()

const G4ThreeVector * G4ParticleChangeForDecay::GetPolarization ( ) const
inline

Definition at line 143 of file G4ParticleChangeForDecay.hh.

144{
145 return &thePolarizationChange;
146}

◆ Initialize()

void G4ParticleChangeForDecay::Initialize ( const G4Track track)
virtual

Reimplemented from G4VParticleChange.

Definition at line 112 of file G4ParticleChangeForDecay.cc.

113{
114 // use base class's method at first
116
117 const G4DynamicParticle* pParticle = track.GetDynamicParticle();
118
119 // set TimeChange equal to local time of the parent track
120 theTimeChange = track.GetLocalTime();
121
122 // set initial Local/Global time of the parent track
123 theLocalTime0 = track.GetLocalTime();
125
126 // set the Polarization equal to those of the parent track
128}
const G4ThreeVector & GetPolarization() const
G4double GetGlobalTime() const
const G4DynamicParticle * GetDynamicParticle() const
virtual void Initialize(const G4Track &)

Referenced by G4Decay::DecayIt(), G4UnknownDecay::DecayIt(), G4Radioactivation::DecayIt(), G4RadioactiveDecay::DecayIt(), G4RadioactiveDecayBase::DecayIt(), G4Decay::PostStepDoIt(), and G4DecayWithSpin::PostStepDoIt().

◆ operator!=()

G4bool G4ParticleChangeForDecay::operator!= ( const G4ParticleChangeForDecay right) const

Definition at line 105 of file G4ParticleChangeForDecay.cc.

107{
108 return ((G4VParticleChange*) this != (G4VParticleChange*) &right);
109}

◆ operator=()

G4ParticleChangeForDecay & G4ParticleChangeForDecay::operator= ( const G4ParticleChangeForDecay right)
protected

Definition at line 63 of file G4ParticleChangeForDecay.cc.

64{
65 if(this != &right)
66 {
68 {
69 for(G4int index = 0; index < theNumberOfSecondaries; ++index)
70 {
71 if((*theListOfSecondaries)[index])
72 delete(*theListOfSecondaries)[index];
73 }
74 }
78 for(G4int index = 0; index < theNumberOfSecondaries; ++index)
79 {
80 G4Track* newTrack = new G4Track(*((*right.theListOfSecondaries)[index]));
81 theListOfSecondaries->SetElement(index, newTrack);
82 }
83
88
93 }
94 return *this;
95}
G4FastVector< G4Track, G4TrackFastVectorSize > G4TrackFastVector
void SetElement(G4int anIndex, Type *anElement)
Definition: G4FastVector.hh:72
G4TrackStatus theStatusChange
G4TrackFastVector * theListOfSecondaries
G4SteppingControl theSteppingControlFlag

◆ operator==()

G4bool G4ParticleChangeForDecay::operator== ( const G4ParticleChangeForDecay right) const

Definition at line 98 of file G4ParticleChangeForDecay.cc.

100{
101 return ((G4VParticleChange*) this == (G4VParticleChange*) &right);
102}

◆ ProposeGlobalTime()

void G4ParticleChangeForDecay::ProposeGlobalTime ( G4double  t)
inline

Definition at line 117 of file G4ParticleChangeForDecay.hh.

118{
120}

Referenced by G4UnknownDecay::DecayIt().

◆ ProposeLocalTime()

void G4ParticleChangeForDecay::ProposeLocalTime ( G4double  t)
inline

◆ ProposePolarization() [1/2]

void G4ParticleChangeForDecay::ProposePolarization ( const G4ThreeVector finalPoralization)
inline

Definition at line 149 of file G4ParticleChangeForDecay.hh.

151{
152 thePolarizationChange = finalPoralization;
153}

◆ ProposePolarization() [2/2]

void G4ParticleChangeForDecay::ProposePolarization ( G4double  Px,
G4double  Py,
G4double  Pz 
)
inline

◆ UpdateStepForAtRest()

G4Step * G4ParticleChangeForDecay::UpdateStepForAtRest ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Definition at line 147 of file G4ParticleChangeForDecay.cc.

148{
149 // A physics process always calculates the final state of the particle
150
151 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
152
153 // update polarization
154 pPostStepPoint->SetPolarization(thePolarizationChange);
155
156 // update time
157 pPostStepPoint->SetGlobalTime(GetGlobalTime());
158 pPostStepPoint->SetLocalTime(theTimeChange);
159 pPostStepPoint->AddProperTime(theTimeChange - theLocalTime0);
160
161#ifdef G4VERBOSE
162 G4Track* aTrack = pStep->GetTrack();
163 if(debugFlag) { CheckIt(*aTrack); }
164#endif
165
167 pPostStepPoint->SetWeight(theParentWeight);
168
169 // Update the G4Step specific attributes
170 return UpdateStepInfo(pStep);
171}
G4double GetGlobalTime(G4double timeDelay=0.0) const
virtual G4bool CheckIt(const G4Track &)
void SetLocalTime(const G4double aValue)
void SetWeight(G4double aValue)
void AddProperTime(const G4double aValue)
void SetGlobalTime(const G4double aValue)
void SetPolarization(const G4ThreeVector &aValue)
G4Step * UpdateStepInfo(G4Step *Step)

◆ UpdateStepForPostStep()

G4Step * G4ParticleChangeForDecay::UpdateStepForPostStep ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Definition at line 131 of file G4ParticleChangeForDecay.cc.

132{
133 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
134
136 {
137 pPostStepPoint->SetWeight(theParentWeight);
138 }
139 // update polarization
140 pPostStepPoint->SetPolarization(thePolarizationChange);
141
142 // update the G4Step specific attributes
143 return UpdateStepInfo(pStep);
144}

Member Data Documentation

◆ theGlobalTime0

G4double G4ParticleChangeForDecay::theGlobalTime0 = 0.0
protected

◆ theLocalTime0

G4double G4ParticleChangeForDecay::theLocalTime0 = 0.0
protected

◆ thePolarizationChange

G4ThreeVector G4ParticleChangeForDecay::thePolarizationChange
protected

◆ theTimeChange


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