Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Protected Member Functions | Protected Attributes
G4ParticleChangeForDecay Class Reference

#include <G4ParticleChangeForDecay.hh>

Inheritance diagram for G4ParticleChangeForDecay:
G4VParticleChange G4ParticleChangeForRadDecay

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 G4StepUpdateStepForAlongStep (G4Step *Step)
 
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
 
void SetVerboseLevel (G4int vLevel)
 
G4int GetVerboseLevel () const
 
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
 
G4double theLocalTime0
 
G4double theTimeChange
 
G4ThreeVector thePolarizationChange
 
- Protected Attributes inherited from G4VParticleChange
G4TrackFastVectortheListOfSecondaries
 
G4int theNumberOfSecondaries
 
G4int theSizeOftheListOfSecondaries
 
G4TrackStatus theStatusChange
 
G4SteppingControl theSteppingControlFlag
 
G4double theLocalEnergyDeposit
 
G4double theNonIonizingEnergyDeposit
 
G4double theTrueStepLength
 
G4bool theFirstStepInVolume
 
G4bool theLastStepInVolume
 
G4double theParentWeight
 
G4bool isParentWeightProposed
 
G4bool fSetSecondaryWeightByProcess
 
G4double theParentGlobalTime
 
G4int verboseLevel
 
G4bool debugFlag
 

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 54 of file G4ParticleChangeForDecay.hh.

Constructor & Destructor Documentation

G4ParticleChangeForDecay::G4ParticleChangeForDecay ( )

Definition at line 48 of file G4ParticleChangeForDecay.cc.

References G4cout, G4endl, and G4VParticleChange::verboseLevel.

49  : G4VParticleChange(),
51 {
52 #ifdef G4VERBOSE
53  if (verboseLevel>2) {
54  G4cout << "G4ParticleChangeForDecay::G4ParticleChangeForDecay() " << G4endl;
55  }
56 #endif
57 }
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
G4ParticleChangeForDecay::~G4ParticleChangeForDecay ( )
virtual

Definition at line 59 of file G4ParticleChangeForDecay.cc.

References G4cout, G4endl, and G4VParticleChange::verboseLevel.

60 {
61 #ifdef G4VERBOSE
62  if (verboseLevel>2) {
63  G4cout << "G4ParticleChangeForDecay::~G4ParticleChangeForDecay() " << G4endl;
64  }
65 #endif
66 }
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
G4ParticleChangeForDecay::G4ParticleChangeForDecay ( const G4ParticleChangeForDecay right)
protected

Member Function Documentation

G4bool G4ParticleChangeForDecay::CheckIt ( const G4Track aTrack)
virtual

Reimplemented from G4VParticleChange.

Definition at line 197 of file G4ParticleChangeForDecay.cc.

References G4VParticleChange::accuracyForException, G4VParticleChange::accuracyForWarning, G4VParticleChange::CheckIt(), DumpInfo(), EventMustBeAborted, G4cout, G4endl, G4Exception(), G4Track::GetDefinition(), G4Track::GetKineticEnergy(), G4Track::GetLocalTime(), G4ParticleDefinition::GetParticleName(), G4Track::GetPosition(), python.hepunit::m, python.hepunit::MeV, ns, theGlobalTime0, theLocalTime0, theTimeChange, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

Referenced by UpdateStepForAtRest().

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

Reimplemented from G4VParticleChange.

Definition at line 182 of file G4ParticleChangeForDecay.cc.

References G4VParticleChange::DumpInfo(), G4cout, G4endl, ns, theGlobalTime0, theLocalTime0, and theTimeChange.

Referenced by CheckIt().

183 {
184 // Show header
186 
187  G4int oldprc = G4cout.precision(3);
188  G4cout << " proposed local Time (ns) : "
189  << std::setw(20) << theTimeChange/ns << G4endl;
190  G4cout << " initial local Time (ns) : "
191  << std::setw(20) << theLocalTime0/ns << G4endl;
192  G4cout << " initial global Time (ns) : "
193  << std::setw(20) << theGlobalTime0/ns << G4endl;
194  G4cout.precision(oldprc);
195 }
virtual void DumpInfo() const
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
#define ns
Definition: xmlparse.cc:597
G4double G4ParticleChangeForDecay::GetGlobalTime ( G4double  timeDelay = 0.0) const
inline

Definition at line 132 of file G4ParticleChangeForDecay.hh.

References theGlobalTime0, theLocalTime0, and theTimeChange.

Referenced by UpdateStepForAtRest().

133 {
134  // Convert the time delay to the global time.
135  return theGlobalTime0 + (theTimeChange-theLocalTime0) + timeDelay;
136 }
G4double G4ParticleChangeForDecay::GetLocalTime ( G4double  timeDelay = 0.0) const
inline

Definition at line 145 of file G4ParticleChangeForDecay.hh.

References theTimeChange.

146 {
147  // Convert the time delay to the local time.
148  return theTimeChange + timeDelay;
149 }
const G4ThreeVector * G4ParticleChangeForDecay::GetPolarization ( void  ) const
inline

Definition at line 152 of file G4ParticleChangeForDecay.hh.

References thePolarizationChange.

153 {
154  return &thePolarizationChange;
155 }
void G4ParticleChangeForDecay::Initialize ( const G4Track track)
virtual

Reimplemented from G4VParticleChange.

Definition at line 124 of file G4ParticleChangeForDecay.cc.

References G4Track::GetDynamicParticle(), G4Track::GetGlobalTime(), G4Track::GetLocalTime(), G4DynamicParticle::GetPolarization(), G4VParticleChange::Initialize(), theGlobalTime0, theLocalTime0, thePolarizationChange, and theTimeChange.

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

125 {
126  // use base class's method at first
128 
129  const G4DynamicParticle* pParticle = track.GetDynamicParticle();
130 
131  // set TimeChange equal to local time of the parent track
132  theTimeChange = track.GetLocalTime();
133 
134  // set initial Local/Global time of the parent track
135  theLocalTime0 = track.GetLocalTime();
136  theGlobalTime0 = track.GetGlobalTime();
137 
138  // set the Polarization equal to those of the parent track
139  thePolarizationChange = pParticle->GetPolarization();
140 }
virtual void Initialize(const G4Track &)
G4double GetLocalTime() const
const G4DynamicParticle * GetDynamicParticle() const
G4double GetGlobalTime() const
const G4ThreeVector & GetPolarization() const
G4bool G4ParticleChangeForDecay::operator!= ( const G4ParticleChangeForDecay right) const

Definition at line 116 of file G4ParticleChangeForDecay.cc.

117 {
118  return ((G4VParticleChange *)this != (G4VParticleChange *) &right);
119 }
G4ParticleChangeForDecay & G4ParticleChangeForDecay::operator= ( const G4ParticleChangeForDecay right)
protected

Definition at line 77 of file G4ParticleChangeForDecay.cc.

References G4cout, G4FastVector< Type, N >::SetElement(), theGlobalTime0, G4VParticleChange::theListOfSecondaries, G4VParticleChange::theLocalEnergyDeposit, theLocalTime0, G4VParticleChange::theNumberOfSecondaries, thePolarizationChange, G4VParticleChange::theStatusChange, G4VParticleChange::theSteppingControlFlag, theTimeChange, G4VParticleChange::theTrueStepLength, and G4VParticleChange::verboseLevel.

78 {
79  if (this != &right){
80  if (theNumberOfSecondaries>0) {
81 #ifdef G4VERBOSE
82  if (verboseLevel>0) {
83  G4cout << "G4ParticleChangeForDecay: assignment operator Warning ";
84  G4cout << "theListOfSecondaries is not empty ";
85  }
86 #endif
87  for (G4int index= 0; index<theNumberOfSecondaries; index++){
88  if ( (*theListOfSecondaries)[index] ) delete (*theListOfSecondaries)[index] ;
89  }
90  }
91  delete theListOfSecondaries;
93  theNumberOfSecondaries = right.theNumberOfSecondaries;
94  for (G4int index = 0; index<theNumberOfSecondaries; index++){
95  G4Track* newTrack = new G4Track(*((*right.theListOfSecondaries)[index] ));
96  theListOfSecondaries->SetElement(index, newTrack); }
97 
102 
107  }
108  return *this;
109 }
void SetElement(G4int anIndex, Type *anElement)
Definition: G4FastVector.hh:76
G4TrackFastVector * theListOfSecondaries
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
G4SteppingControl theSteppingControlFlag
G4FastVector< G4Track, G4TrackFastVectorSize > G4TrackFastVector
G4TrackStatus theStatusChange
G4bool G4ParticleChangeForDecay::operator== ( const G4ParticleChangeForDecay right) const

Definition at line 111 of file G4ParticleChangeForDecay.cc.

112 {
113  return ((G4VParticleChange *)this == (G4VParticleChange *) &right);
114 }
void G4ParticleChangeForDecay::ProposeGlobalTime ( G4double  t)
inline
void G4ParticleChangeForDecay::ProposeLocalTime ( G4double  t)
inline

Definition at line 139 of file G4ParticleChangeForDecay.hh.

References theTimeChange.

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

140 {
141  theTimeChange = t;
142 }
void G4ParticleChangeForDecay::ProposePolarization ( G4double  Px,
G4double  Py,
G4double  Pz 
)
inline
void G4ParticleChangeForDecay::ProposePolarization ( const G4ThreeVector finalPoralization)
inline

Definition at line 158 of file G4ParticleChangeForDecay.hh.

References thePolarizationChange.

159 {
160  thePolarizationChange = finalPoralization;
161 }
G4Step * G4ParticleChangeForDecay::UpdateStepForAtRest ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Definition at line 157 of file G4ParticleChangeForDecay.cc.

References G4StepPoint::AddProperTime(), CheckIt(), G4VParticleChange::debugFlag, GetGlobalTime(), G4Step::GetPostStepPoint(), G4Step::GetTrack(), G4VParticleChange::isParentWeightProposed, G4StepPoint::SetGlobalTime(), G4StepPoint::SetLocalTime(), G4StepPoint::SetPolarization(), G4StepPoint::SetWeight(), theLocalTime0, G4VParticleChange::theParentWeight, thePolarizationChange, theTimeChange, and G4VParticleChange::UpdateStepInfo().

158 {
159  // A physics process always calculates the final state of the particle
160 
161  G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
162 
163  // update polarization
164  pPostStepPoint->SetPolarization( thePolarizationChange );
165 
166  // update time
167  pPostStepPoint->SetGlobalTime( GetGlobalTime() );
168  pPostStepPoint->SetLocalTime( theTimeChange );
169  pPostStepPoint->AddProperTime (theTimeChange-theLocalTime0);
170 
171 #ifdef G4VERBOSE
172  G4Track* aTrack = pStep->GetTrack();
173  if (debugFlag) CheckIt(*aTrack);
174 #endif
175 
176  if (isParentWeightProposed )pPostStepPoint->SetWeight( theParentWeight );
177 
178  // Update the G4Step specific attributes
179  return UpdateStepInfo(pStep);
180 }
G4double GetGlobalTime(G4double timeDelay=0.0) const
virtual G4bool CheckIt(const G4Track &)
void SetWeight(G4double aValue)
void SetGlobalTime(const G4double aValue)
void SetLocalTime(const G4double aValue)
void SetPolarization(const G4ThreeVector &aValue)
G4Step * UpdateStepInfo(G4Step *Step)
void AddProperTime(const G4double aValue)
G4Step * G4ParticleChangeForDecay::UpdateStepForPostStep ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Definition at line 146 of file G4ParticleChangeForDecay.cc.

References G4Step::GetPostStepPoint(), G4VParticleChange::isParentWeightProposed, G4StepPoint::SetWeight(), G4VParticleChange::theParentWeight, and G4VParticleChange::UpdateStepInfo().

147 {
149  pStep->GetPostStepPoint()->SetWeight( theParentWeight );
150  }
151 
152  // Update the G4Step specific attributes
153  return UpdateStepInfo(pStep);
154 }
G4Step * UpdateStepInfo(G4Step *Step)

Field Documentation

G4double G4ParticleChangeForDecay::theGlobalTime0
protected
G4double G4ParticleChangeForDecay::theLocalTime0
protected
G4ThreeVector G4ParticleChangeForDecay::thePolarizationChange
protected
G4double G4ParticleChangeForDecay::theTimeChange
protected

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