Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4ParticleChange Class Reference

#include <G4ParticleChange.hh>

Inheritance diagram for G4ParticleChange:
Collaboration diagram for G4ParticleChange:

Public Member Functions

 G4ParticleChange ()
 
virtual ~G4ParticleChange ()
 
G4bool operator== (const G4ParticleChange &right) const
 
G4bool operator!= (const G4ParticleChange &right) const
 
virtual G4StepUpdateStepForAlongStep (G4Step *Step)
 
virtual G4StepUpdateStepForAtRest (G4Step *Step)
 
virtual G4StepUpdateStepForPostStep (G4Step *Step)
 
virtual void Initialize (const G4Track &)
 
const G4ThreeVectorGetMomentumDirection () const
 
void ProposeMomentumDirection (G4double Px, G4double Py, G4double Pz)
 
void ProposeMomentumDirection (const G4ThreeVector &Pfinal)
 
const G4ThreeVectorGetPolarization () const
 
void ProposePolarization (G4double Px, G4double Py, G4double Pz)
 
void ProposePolarization (const G4ThreeVector &finalPoralization)
 
G4double GetEnergy () const
 
void ProposeEnergy (G4double finalEnergy)
 
G4double GetVelocity () const
 
void ProposeVelocity (G4double finalVelocity)
 
G4double GetProperTime () const
 
void ProposeProperTime (G4double finalProperTime)
 
const G4ThreeVectorGetPosition () const
 
void ProposePosition (G4double x, G4double y, G4double z)
 
void ProposePosition (const G4ThreeVector &finalPosition)
 
void ProposeGlobalTime (G4double t)
 
void ProposeLocalTime (G4double t)
 
G4double GetGlobalTime (G4double timeDelay=0.0) const
 
G4double GetLocalTime (G4double timeDelay=0.0) const
 
G4double GetMass () const
 
void ProposeMass (G4double finalMass)
 
G4double GetCharge () const
 
void ProposeCharge (G4double finalCharge)
 
G4double GetMagneticMoment () const
 
void ProposeMagneticMoment (G4double finalMagneticMoment)
 
G4ThreeVector GetGlobalPosition (const G4ThreeVector &displacement) const
 
G4ThreeVector CalcMomentum (G4double energy, G4ThreeVector direction, G4double mass) const
 
void AddSecondary (G4Track *aSecondary)
 
void AddSecondary (G4DynamicParticle *aSecondary, G4bool IsGoodForTracking=false)
 
void AddSecondary (G4DynamicParticle *aSecondary, G4ThreeVector position, G4bool IsGoodForTracking=false)
 
void AddSecondary (G4DynamicParticle *aSecondary, G4double time, G4bool IsGoodForTracking=false)
 
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
 
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

 G4ParticleChange (const G4ParticleChange &right)
 
G4ParticleChangeoperator= (const G4ParticleChange &right)
 
G4StepUpdateStepInfo (G4Step *Step)
 
- 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

G4ThreeVector theMomentumDirectionChange
 
G4ThreeVector thePolarizationChange
 
G4double theEnergyChange
 
G4double theVelocityChange
 
G4bool isVelocityChanged
 
G4ThreeVector thePositionChange
 
G4double theGlobalTime0
 
G4double theLocalTime0
 
G4double theTimeChange
 
G4double theProperTimeChange
 
G4double theMassChange
 
G4double theChargeChange
 
G4double theMagneticMomentChange
 
const G4TracktheCurrentTrack
 
- 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 77 of file G4ParticleChange.hh.

Constructor & Destructor Documentation

G4ParticleChange::G4ParticleChange ( )

Definition at line 54 of file G4ParticleChange.cc.

55  : G4VParticleChange(),
58  theEnergyChange(0.),
65 {
66 }
G4double theProperTimeChange
G4double theMagneticMomentChange
G4ThreeVector thePositionChange
G4ThreeVector thePolarizationChange
G4ThreeVector theMomentumDirectionChange
const G4Track * theCurrentTrack
G4ParticleChange::~G4ParticleChange ( )
virtual

Definition at line 68 of file G4ParticleChange.cc.

69 {
70 #ifdef G4VERBOSE
71  if (verboseLevel>2) {
72  G4cout << "G4ParticleChange::~G4ParticleChange() " << G4endl;
73  }
74 #endif
75 }
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
G4ParticleChange::G4ParticleChange ( const G4ParticleChange right)
protected

Definition at line 78 of file G4ParticleChange.cc.

79  : G4VParticleChange(right)
80 {
81  if (verboseLevel>1) {
82  G4cout << "G4ParticleChange:: copy constructor is called " << G4endl;
83  }
85 
95  isVelocityChanged = true;
99 }
G4double theProperTimeChange
G4double theMagneticMomentChange
G4ThreeVector thePositionChange
G4ThreeVector thePolarizationChange
G4GLOB_DLL std::ostream G4cout
G4ThreeVector theMomentumDirectionChange
const G4Track * theCurrentTrack
#define G4endl
Definition: G4ios.hh:61

Member Function Documentation

void G4ParticleChange::AddSecondary ( G4Track aSecondary)

Definition at line 218 of file G4ParticleChange.cc.

219 {
220  // add a secondary
222 }
void AddSecondary(G4Track *aSecondary)

Here is the call graph for this function:

Here is the caller graph for this function:

void G4ParticleChange::AddSecondary ( G4DynamicParticle aSecondary,
G4bool  IsGoodForTracking = false 
)

Definition at line 168 of file G4ParticleChange.cc.

170 {
171  // create track
172  G4Track* aTrack = new G4Track(aParticle, GetGlobalTime(), thePositionChange);
173 
174  // set IsGoodGorTrackingFlag
175  if (IsGoodForTracking) aTrack->SetGoodForTrackingFlag();
176 
177  // Touchable handle is copied to keep the pointer
179 
180  // add a secondary
182 }
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4double GetGlobalTime(G4double timeDelay=0.0) const
G4ThreeVector thePositionChange
void AddSecondary(G4Track *aSecondary)
const G4TouchableHandle & GetTouchableHandle() const
const G4Track * theCurrentTrack
void SetGoodForTrackingFlag(G4bool value=true)

Here is the call graph for this function:

void G4ParticleChange::AddSecondary ( G4DynamicParticle aSecondary,
G4ThreeVector  position,
G4bool  IsGoodForTracking = false 
)

Definition at line 184 of file G4ParticleChange.cc.

187 {
188  // create track
189  G4Track* aTrack = new G4Track(aParticle, GetGlobalTime(), newPosition);
190 
191  // set IsGoodGorTrackingFlag
192  if (IsGoodForTracking) aTrack->SetGoodForTrackingFlag();
193 
194  // Touchable is a temporary object, so you cannot keep the pointer
195  aTrack->SetTouchableHandle((G4VTouchable*)0);
196 
197  // add a secondary
199 }
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4double GetGlobalTime(G4double timeDelay=0.0) const
void AddSecondary(G4Track *aSecondary)
void SetGoodForTrackingFlag(G4bool value=true)

Here is the call graph for this function:

void G4ParticleChange::AddSecondary ( G4DynamicParticle aSecondary,
G4double  time,
G4bool  IsGoodForTracking = false 
)

Definition at line 201 of file G4ParticleChange.cc.

204 {
205  // create track
206  G4Track* aTrack = new G4Track(aParticle, newTime, thePositionChange);
207 
208  // set IsGoodGorTrackingFlag
209  if (IsGoodForTracking) aTrack->SetGoodForTrackingFlag();
210 
211  // Touchable handle is copied to keep the pointer
213 
214  // add a secondary
216 }
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4ThreeVector thePositionChange
void AddSecondary(G4Track *aSecondary)
const G4TouchableHandle & GetTouchableHandle() const
const G4Track * theCurrentTrack
void SetGoodForTrackingFlag(G4bool value=true)

Here is the call graph for this function:

G4ThreeVector G4ParticleChange::CalcMomentum ( G4double  energy,
G4ThreeVector  direction,
G4double  mass 
) const

Here is the caller graph for this function:

G4bool G4ParticleChange::CheckIt ( const G4Track aTrack)
virtual

Reimplemented from G4VParticleChange.

Definition at line 506 of file G4ParticleChange.cc.

507 {
508  G4bool exitWithError = false;
509  G4double accuracy;
510  static G4ThreadLocal G4int nError = 0;
511 #ifdef G4VERBOSE
512  const G4int maxError = 30;
513 #endif
514 
515  // No check in case of "fStopAndKill"
516  if (GetTrackStatus() == fStopAndKill ) return G4VParticleChange::CheckIt(aTrack);
517 
518  // MomentumDirection should be unit vector
519  G4bool itsOKforMomentum = true;
520  if ( theEnergyChange >0.) {
521  accuracy = std::fabs(theMomentumDirectionChange.mag2()-1.0);
522  if (accuracy > accuracyForWarning) {
523  itsOKforMomentum = false;
524  nError += 1;
525  exitWithError = exitWithError || (accuracy > accuracyForException);
526 #ifdef G4VERBOSE
527  if (nError < maxError) {
528  G4cout << " G4ParticleChange::CheckIt : ";
529  G4cout << "the Momentum Change is not unit vector !!"
530  << " Difference: " << accuracy << G4endl;
531  G4cout << aTrack.GetDefinition()->GetParticleName()
532  << " E=" << aTrack.GetKineticEnergy()/MeV
533  << " pos=" << aTrack.GetPosition().x()/m
534  << ", " << aTrack.GetPosition().y()/m
535  << ", " << aTrack.GetPosition().z()/m
536  <<G4endl;
537  }
538 #endif
539  }
540  }
541 
542  // Both global and proper time should not go back
543  G4bool itsOKforGlobalTime = true;
544  accuracy = (aTrack.GetLocalTime()- theTimeChange)/ns;
545  if (accuracy > accuracyForWarning) {
546  itsOKforGlobalTime = false;
547  nError += 1;
548  exitWithError = exitWithError || (accuracy > accuracyForException);
549 #ifdef G4VERBOSE
550  if (nError < maxError) {
551  G4cout << " G4ParticleChange::CheckIt : ";
552  G4cout << "the local time goes back !!"
553  << " Difference: " << accuracy << "[ns] " <<G4endl;
554  G4cout << aTrack.GetDefinition()->GetParticleName()
555  << " E=" << aTrack.GetKineticEnergy()/MeV
556  << " pos=" << aTrack.GetPosition().x()/m
557  << ", " << aTrack.GetPosition().y()/m
558  << ", " << aTrack.GetPosition().z()/m
559  << " global time=" << aTrack.GetGlobalTime()/ns
560  << " local time=" << aTrack.GetLocalTime()/ns
561  << " proper time=" << aTrack.GetProperTime()/ns
562  << G4endl;
563  }
564 #endif
565  }
566 
567  G4bool itsOKforProperTime = true;
568  accuracy = (aTrack.GetProperTime() - theProperTimeChange )/ns;
569  if (accuracy > accuracyForWarning) {
570  itsOKforProperTime = false;
571  nError += 1;
572  exitWithError = exitWithError || (accuracy > accuracyForException);
573 #ifdef G4VERBOSE
574  if (nError < maxError) {
575  G4cout << " G4ParticleChange::CheckIt : ";
576  G4cout << "the proper time goes back !!"
577  << " Difference: " << accuracy << "[ns] " <<G4endl;
578  G4cout << aTrack.GetDefinition()->GetParticleName()
579  << " E=" << aTrack.GetKineticEnergy()/MeV
580  << " pos=" << aTrack.GetPosition().x()/m
581  << ", " << aTrack.GetPosition().y()/m
582  << ", " << aTrack.GetPosition().z()/m
583  << " global time=" << aTrack.GetGlobalTime()/ns
584  << " local time=" << aTrack.GetLocalTime()/ns
585  << " proper time=" << aTrack.GetProperTime()/ns
586  <<G4endl;
587  }
588 #endif
589  }
590 
591  // Kinetic Energy should not be negative
592  G4bool itsOKforEnergy = true;
593  accuracy = -1.0*theEnergyChange/MeV;
594  if (accuracy > accuracyForWarning) {
595  itsOKforEnergy = false;
596  nError += 1;
597  exitWithError = exitWithError || (accuracy > accuracyForException);
598 #ifdef G4VERBOSE
599  if (nError < maxError) {
600  G4cout << " G4ParticleChange::CheckIt : ";
601  G4cout << "the kinetic energy is negative !!"
602  << " Difference: " << accuracy << "[MeV] " <<G4endl;
603  G4cout << aTrack.GetDefinition()->GetParticleName()
604  << " E=" << aTrack.GetKineticEnergy()/MeV
605  << " pos=" << aTrack.GetPosition().x()/m
606  << ", " << aTrack.GetPosition().y()/m
607  << ", " << aTrack.GetPosition().z()/m
608  <<G4endl;
609  }
610 #endif
611  }
612 
613  // Velocity should not be less than c_light
614  G4bool itsOKforVelocity = true;
615  if (theVelocityChange < 0.) {
616  itsOKforVelocity = false;
617  nError += 1;
618  exitWithError = true;
619 #ifdef G4VERBOSE
620  if (nError < maxError) {
621  G4cout << " G4ParticleChange::CheckIt : ";
622  G4cout << "the velocity is negative !!"
623  << " Velocity: " << theVelocityChange/c_light <<G4endl;
624  G4cout << aTrack.GetDefinition()->GetParticleName()
625  << " E=" << aTrack.GetKineticEnergy()/MeV
626  << " pos=" << aTrack.GetPosition().x()/m
627  << ", " << aTrack.GetPosition().y()/m
628  << ", " << aTrack.GetPosition().z()/m
629  <<G4endl;
630  }
631 #endif
632  }
633 
634  accuracy = theVelocityChange/c_light - 1.0;
635  if (accuracy > accuracyForWarning) {
636  itsOKforVelocity = false;
637  nError += 1;
638  exitWithError = exitWithError || (accuracy > accuracyForException);
639 #ifdef G4VERBOSE
640  if (nError < maxError) {
641  G4cout << " G4ParticleChange::CheckIt : ";
642  G4cout << "the velocity is greater than c_light !!" << G4endl;
643  G4cout << " Velocity: " << theVelocityChange/c_light <<G4endl;
644  G4cout << aTrack.GetDefinition()->GetParticleName()
645  << " E=" << aTrack.GetKineticEnergy()/MeV
646  << " pos=" << aTrack.GetPosition().x()/m
647  << ", " << aTrack.GetPosition().y()/m
648  << ", " << aTrack.GetPosition().z()/m
649  <<G4endl;
650  }
651 #endif
652  }
653 
654  G4bool itsOK = itsOKforMomentum && itsOKforEnergy && itsOKforVelocity && itsOKforProperTime && itsOKforGlobalTime;
655  // dump out information of this particle change
656 #ifdef G4VERBOSE
657  if (!itsOK) {
658  DumpInfo();
659  }
660 #endif
661 
662  // Exit with error
663  if (exitWithError) {
664  G4Exception("G4ParticleChange::CheckIt",
665  "TRACK003", EventMustBeAborted,
666  "momentum, energy, and/or time was illegal");
667  }
668  //correction
669  if (!itsOKforMomentum) {
672  }
673  if (!itsOKforGlobalTime) {
674  theTimeChange = aTrack.GetLocalTime();
675  }
676  if (!itsOKforProperTime) {
678  }
679  if (!itsOKforEnergy) {
680  theEnergyChange = 0.0;
681  }
682  if (!itsOKforVelocity) {
684  }
685 
686  itsOK = (itsOK) && G4VParticleChange::CheckIt(aTrack);
687  return itsOK;
688 }
G4ParticleDefinition * GetDefinition() const
G4double theProperTimeChange
G4double GetLocalTime() const
G4double GetProperTime() const
double x() const
const G4ThreeVector & GetPosition() const
#define G4ThreadLocal
Definition: tls.hh:89
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
double z() const
static const G4double accuracyForException
virtual G4bool CheckIt(const G4Track &)
G4double GetKineticEnergy() const
G4GLOB_DLL std::ostream G4cout
static constexpr double m
Definition: G4SIunits.hh:129
bool G4bool
Definition: G4Types.hh:79
G4double GetGlobalTime() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4ThreeVector theMomentumDirectionChange
virtual void DumpInfo() const
static constexpr double c_light
double y() const
double mag2() const
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
G4TrackStatus GetTrackStatus() const
double G4double
Definition: G4Types.hh:76
static const G4double accuracyForWarning
double mag() const
#define ns
Definition: xmlparse.cc:614

Here is the call graph for this function:

Here is the caller graph for this function:

void G4ParticleChange::DumpInfo ( ) const
virtual

Reimplemented from G4VParticleChange.

Reimplemented in G4ParticleChangeForTransport.

Definition at line 445 of file G4ParticleChange.cc.

446 {
447 // use base-class DumpInfo
449 
450  G4int oldprc = G4cout.precision(3);
451 
452  G4cout << " Mass (GeV) : "
453  << std::setw(20) << theMassChange/GeV
454  << G4endl;
455  G4cout << " Charge (eplus) : "
456  << std::setw(20) << theChargeChange/eplus
457  << G4endl;
458  G4cout << " MagneticMoment : "
459  << std::setw(20) << theMagneticMomentChange << G4endl;
460  G4cout << " : = " << std::setw(20)
462  << "*[e hbar]/[2 m]"
463  << G4endl;
464  G4cout << " Position - x (mm) : "
465  << std::setw(20) << thePositionChange.x()/mm
466  << G4endl;
467  G4cout << " Position - y (mm) : "
468  << std::setw(20) << thePositionChange.y()/mm
469  << G4endl;
470  G4cout << " Position - z (mm) : "
471  << std::setw(20) << thePositionChange.z()/mm
472  << G4endl;
473  G4cout << " Time (ns) : "
474  << std::setw(20) << theTimeChange/ns
475  << G4endl;
476  G4cout << " Proper Time (ns) : "
477  << std::setw(20) << theProperTimeChange/ns
478  << G4endl;
479  G4cout << " Momentum Direct - x : "
480  << std::setw(20) << theMomentumDirectionChange.x()
481  << G4endl;
482  G4cout << " Momentum Direct - y : "
483  << std::setw(20) << theMomentumDirectionChange.y()
484  << G4endl;
485  G4cout << " Momentum Direct - z : "
486  << std::setw(20) << theMomentumDirectionChange.z()
487  << G4endl;
488  G4cout << " Kinetic Energy (MeV): "
489  << std::setw(20) << theEnergyChange/MeV
490  << G4endl;
491  G4cout << " Velocity (/c): "
492  << std::setw(20) << theVelocityChange/c_light
493  << G4endl;
494  G4cout << " Polarization - x : "
495  << std::setw(20) << thePolarizationChange.x()
496  << G4endl;
497  G4cout << " Polarization - y : "
498  << std::setw(20) << thePolarizationChange.y()
499  << G4endl;
500  G4cout << " Polarization - z : "
501  << std::setw(20) << thePolarizationChange.z()
502  << G4endl;
503  G4cout.precision(oldprc);
504 }
G4double theProperTimeChange
static constexpr double mm
Definition: G4SIunits.hh:115
double x() const
virtual void DumpInfo() const
int G4int
Definition: G4Types.hh:78
double z() const
G4double theMagneticMomentChange
G4ThreeVector thePositionChange
G4ThreeVector thePolarizationChange
G4GLOB_DLL std::ostream G4cout
static constexpr double eplus
Definition: G4SIunits.hh:199
static constexpr double c_squared
G4ThreeVector theMomentumDirectionChange
static constexpr double c_light
double y() const
static constexpr double GeV
Definition: G4SIunits.hh:217
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
static constexpr double hbar_Planck
#define ns
Definition: xmlparse.cc:614

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4ParticleChange::GetCharge ( ) const
G4double G4ParticleChange::GetEnergy ( ) const

Here is the caller graph for this function:

G4ThreeVector G4ParticleChange::GetGlobalPosition ( const G4ThreeVector displacement) const
G4double G4ParticleChange::GetGlobalTime ( G4double  timeDelay = 0.0) const

Here is the caller graph for this function:

G4double G4ParticleChange::GetLocalTime ( G4double  timeDelay = 0.0) const
G4double G4ParticleChange::GetMagneticMoment ( ) const
G4double G4ParticleChange::GetMass ( ) const
const G4ThreeVector* G4ParticleChange::GetMomentumDirection ( ) const

Here is the caller graph for this function:

const G4ThreeVector* G4ParticleChange::GetPolarization ( ) const

Here is the caller graph for this function:

const G4ThreeVector* G4ParticleChange::GetPosition ( ) const
G4double G4ParticleChange::GetProperTime ( ) const
G4double G4ParticleChange::GetVelocity ( ) const
void G4ParticleChange::Initialize ( const G4Track track)
virtual

Reimplemented from G4VParticleChange.

Reimplemented in G4ParticleChangeForTransport.

Definition at line 228 of file G4ParticleChange.cc.

229 {
230  // use base class's method at first
232  theCurrentTrack= &track;
233 
234  // set Energy/Momentum etc. equal to those of the parent particle
235  const G4DynamicParticle* pParticle = track.GetDynamicParticle();
236  theEnergyChange = pParticle->GetKineticEnergy();
237  theVelocityChange = track.GetVelocity();
238  isVelocityChanged = false;
240  thePolarizationChange = pParticle->GetPolarization();
241  theProperTimeChange = pParticle->GetProperTime();
242 
243  // Set mass/charge/MagneticMoment of DynamicParticle
244  theMassChange = pParticle->GetMass();
245  theChargeChange = pParticle->GetCharge();
247 
248  // set Position equal to those of the parent track
249  thePositionChange = track.GetPosition();
250 
251  // set TimeChange equal to local time of the parent track
252  theTimeChange = track.GetLocalTime();
253 
254  // set initial Local/Global time of the parent track
255  theLocalTime0 = track.GetLocalTime();
256  theGlobalTime0 = track.GetGlobalTime();
257 
258 }
virtual void Initialize(const G4Track &)
G4double theProperTimeChange
G4double GetLocalTime() const
G4double GetKineticEnergy() const
G4double GetVelocity() const
const G4DynamicParticle * GetDynamicParticle() const
G4double GetProperTime() const
const G4ThreeVector & GetPosition() const
G4double theMagneticMomentChange
G4ThreeVector thePositionChange
G4ThreeVector thePolarizationChange
G4double GetMass() const
const G4ThreeVector & GetMomentumDirection() const
G4double GetCharge() const
G4double GetGlobalTime() const
G4ThreeVector theMomentumDirectionChange
const G4ThreeVector & GetPolarization() const
const G4Track * theCurrentTrack
G4double GetMagneticMoment() const

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 158 of file G4ParticleChange.cc.

159 {
160  return ((G4VParticleChange *)this != (G4VParticleChange *) &right);
161 }
G4ParticleChange & G4ParticleChange::operator= ( const G4ParticleChange right)
protected

Definition at line 102 of file G4ParticleChange.cc.

103 {
104 #ifdef G4VERBOSE
105  if (verboseLevel>1) {
106  G4cout << "G4ParticleChange:: assignment operator is called " << G4endl;
107  }
108 #endif
109  if (this != &right){
110  if (theNumberOfSecondaries>0) {
111 #ifdef G4VERBOSE
112  if (verboseLevel>0) {
113  G4cout << "G4ParticleChange: assignment operator Warning ";
114  G4cout << "theListOfSecondaries is not empty ";
115  }
116 #endif
117  for (G4int index= 0; index<theNumberOfSecondaries; index++){
118  if ( (*theListOfSecondaries)[index] ) delete (*theListOfSecondaries)[index] ;
119  }
120  }
121  delete theListOfSecondaries;
122 
124  theNumberOfSecondaries = right.theNumberOfSecondaries;
125  for (G4int index = 0; index<theNumberOfSecondaries; index++){
126  G4Track* newTrack = new G4Track(*((*right.theListOfSecondaries)[index] ));
127  theListOfSecondaries->SetElement(index, newTrack); }
128 
131 
141  isVelocityChanged = true;
145 
149  }
150  return *this;
151 }
G4double theProperTimeChange
void SetElement(G4int anIndex, Type *anElement)
Definition: G4FastVector.hh:76
G4TrackFastVector * theListOfSecondaries
int G4int
Definition: G4Types.hh:78
G4double theMagneticMomentChange
G4ThreeVector thePositionChange
G4ThreeVector thePolarizationChange
G4GLOB_DLL std::ostream G4cout
G4SteppingControl theSteppingControlFlag
G4ThreeVector theMomentumDirectionChange
G4FastVector< G4Track, G4TrackFastVectorSize > G4TrackFastVector
const G4Track * theCurrentTrack
#define G4endl
Definition: G4ios.hh:61
G4TrackStatus theStatusChange

Here is the call graph for this function:

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

Definition at line 153 of file G4ParticleChange.cc.

154 {
155  return ((G4VParticleChange *)this == (G4VParticleChange *) &right);
156 }
void G4ParticleChange::ProposeCharge ( G4double  finalCharge)
void G4ParticleChange::ProposeEnergy ( G4double  finalEnergy)

Here is the caller graph for this function:

void G4ParticleChange::ProposeGlobalTime ( G4double  t)

Here is the caller graph for this function:

void G4ParticleChange::ProposeLocalTime ( G4double  t)

Here is the caller graph for this function:

void G4ParticleChange::ProposeMagneticMoment ( G4double  finalMagneticMoment)
void G4ParticleChange::ProposeMass ( G4double  finalMass)
void G4ParticleChange::ProposeMomentumDirection ( G4double  Px,
G4double  Py,
G4double  Pz 
)

Here is the caller graph for this function:

void G4ParticleChange::ProposeMomentumDirection ( const G4ThreeVector Pfinal)
void G4ParticleChange::ProposePolarization ( G4double  Px,
G4double  Py,
G4double  Pz 
)

Here is the caller graph for this function:

void G4ParticleChange::ProposePolarization ( const G4ThreeVector finalPoralization)
void G4ParticleChange::ProposePosition ( G4double  x,
G4double  y,
G4double  z 
)

Here is the caller graph for this function:

void G4ParticleChange::ProposePosition ( const G4ThreeVector finalPosition)
void G4ParticleChange::ProposeProperTime ( G4double  finalProperTime)

Here is the caller graph for this function:

void G4ParticleChange::ProposeVelocity ( G4double  finalVelocity)

Here is the caller graph for this function:

G4Step * G4ParticleChange::UpdateStepForAlongStep ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Reimplemented in G4ParticleChangeForTransport.

Definition at line 264 of file G4ParticleChange.cc.

265 {
266  // A physics process always calculates the final state of the
267  // particle relative to the initial state at the beginning
268  // of the Step, i.e., based on information of G4Track (or
269  // equivalently the PreStepPoint).
270  // So, the differences (delta) between these two states have to be
271  // calculated and be accumulated in PostStepPoint.
272 
273  // Take note that the return type of GetMomentumDirectionChange is a
274  // pointer to G4ParticleMometum. Also it is a normalized
275  // momentum vector.
276 
277  G4StepPoint* pPreStepPoint = pStep->GetPreStepPoint();
278  G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
279  G4Track* pTrack = pStep->GetTrack();
280  G4double mass = theMassChange;
281 
282  // Set Mass/Charge/MagneticMoment
283  pPostStepPoint->SetMass(theMassChange);
284  pPostStepPoint->SetCharge(theChargeChange);
285  pPostStepPoint->SetMagneticMoment(theMagneticMomentChange);
286 
287  // calculate new kinetic energy
288  G4double preEnergy = pPreStepPoint->GetKineticEnergy();
289  G4double energy = pPostStepPoint->GetKineticEnergy()
290  + (theEnergyChange - preEnergy);
291 
292  // update kinetic energy and momentum direction
293  if (energy > 0.0) {
294  // calculate new momentum
295  G4ThreeVector pMomentum = pPostStepPoint->GetMomentum()
297  - pPreStepPoint->GetMomentum());
298  G4double tMomentum = pMomentum.mag();
299  G4ThreeVector direction(1.0,0.0,0.0);
300  if( tMomentum > 0. ){
301  G4double inv_Momentum= 1.0 / tMomentum;
302  direction= pMomentum * inv_Momentum;
303  }
304  pPostStepPoint->SetMomentumDirection(direction);
305  pPostStepPoint->SetKineticEnergy( energy );
306  } else {
307  // stop case
308  //pPostStepPoint->SetMomentumDirection(G4ThreeVector(1., 0., 0.));
309  pPostStepPoint->SetKineticEnergy(0.0);
310  }
311  // calculate velocity
312  if (!isVelocityChanged) {
313  if(energy > 0.0) {
314  pTrack->SetKineticEnergy(energy);
316  pTrack->SetKineticEnergy(preEnergy);
317  } else if(theMassChange > 0.0) {
318  theVelocityChange = 0.0;
319  }
320  }
321  pPostStepPoint->SetVelocity(theVelocityChange);
322 
323  // update polarization
324  pPostStepPoint->AddPolarization( thePolarizationChange
325  - pPreStepPoint->GetPolarization());
326 
327  // update position and time
328  pPostStepPoint->AddPosition( thePositionChange
329  - pPreStepPoint->GetPosition() );
330  pPostStepPoint->AddGlobalTime(theTimeChange - theLocalTime0);
331  pPostStepPoint->AddLocalTime( theTimeChange - theLocalTime0 );
332  pPostStepPoint->AddProperTime( theProperTimeChange
333  - pPreStepPoint->GetProperTime());
334 
336  pPostStepPoint->SetWeight( theParentWeight );
337  }
338 
339 #ifdef G4VERBOSE
340  G4Track* aTrack = pStep->GetTrack();
341  if (debugFlag) CheckIt(*aTrack);
342 #endif
343 
344  // Update the G4Step specific attributes
345  return UpdateStepInfo(pStep);
346 }
void AddGlobalTime(const G4double aValue)
G4double theProperTimeChange
void AddPosition(const G4ThreeVector &aValue)
G4Step * UpdateStepInfo(G4Step *Step)
void SetMagneticMoment(G4double value)
void SetWeight(G4double aValue)
G4ThreeVector GetMomentum() const
G4ThreeVector CalcMomentum(G4double energy, G4ThreeVector direction, G4double mass) const
void AddLocalTime(const G4double aValue)
void AddPolarization(const G4ThreeVector &aValue)
void SetMomentumDirection(const G4ThreeVector &aValue)
G4double theMagneticMomentChange
G4ThreeVector thePositionChange
G4ThreeVector thePolarizationChange
const G4ThreeVector & GetPosition() const
G4double CalculateVelocity() const
Definition: G4Track.cc:222
virtual G4bool CheckIt(const G4Track &)
G4ThreeVector theMomentumDirectionChange
void SetVelocity(G4double v)
G4double energy(const ThreeVector &p, const G4double m)
void SetCharge(G4double value)
G4double GetProperTime() const
void SetMass(G4double value)
void SetKineticEnergy(const G4double aValue)
G4double GetKineticEnergy() const
double G4double
Definition: G4Types.hh:76
void SetKineticEnergy(const G4double aValue)
double mag() const
const G4ThreeVector & GetPolarization() const
void AddProperTime(const G4double aValue)

Here is the call graph for this function:

G4Step * G4ParticleChange::UpdateStepForAtRest ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Reimplemented in G4ParticleChangeForTransport.

Definition at line 402 of file G4ParticleChange.cc.

403 {
404  // A physics process always calculates the final state of the particle
405 
406  G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
407 
408  // Set Mass/Charge
409  pPostStepPoint->SetMass(theMassChange);
410  pPostStepPoint->SetCharge(theChargeChange);
411  pPostStepPoint->SetMagneticMoment(theMagneticMomentChange);
412 
413  // update kinetic energy and momentum direction
415  pPostStepPoint->SetKineticEnergy( theEnergyChange );
416  if (!isVelocityChanged) theVelocityChange = pStep->GetTrack()->CalculateVelocity();
417  pPostStepPoint->SetVelocity(theVelocityChange);
418 
419  // update polarization
420  pPostStepPoint->SetPolarization( thePolarizationChange );
421 
422  // update position and time
423  pPostStepPoint->SetPosition( thePositionChange );
424  pPostStepPoint->AddGlobalTime(theTimeChange - theLocalTime0);
425  pPostStepPoint->SetLocalTime( theTimeChange );
426  pPostStepPoint->SetProperTime( theProperTimeChange );
427 
429  pPostStepPoint->SetWeight( theParentWeight );
430  }
431 
432 #ifdef G4VERBOSE
433  G4Track* aTrack = pStep->GetTrack();
434  if (debugFlag) CheckIt(*aTrack);
435 #endif
436 
437  // Update the G4Step specific attributes
438  return UpdateStepInfo(pStep);
439 }
void AddGlobalTime(const G4double aValue)
G4double theProperTimeChange
void SetPosition(const G4ThreeVector &aValue)
G4Step * UpdateStepInfo(G4Step *Step)
void SetMagneticMoment(G4double value)
void SetWeight(G4double aValue)
void SetMomentumDirection(const G4ThreeVector &aValue)
void SetLocalTime(const G4double aValue)
G4double theMagneticMomentChange
G4ThreeVector thePositionChange
G4ThreeVector thePolarizationChange
void SetPolarization(const G4ThreeVector &aValue)
void SetProperTime(const G4double aValue)
virtual G4bool CheckIt(const G4Track &)
G4ThreeVector theMomentumDirectionChange
void SetVelocity(G4double v)
void SetCharge(G4double value)
void SetMass(G4double value)
void SetKineticEnergy(const G4double aValue)

Here is the call graph for this function:

G4Step * G4ParticleChange::UpdateStepForPostStep ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Reimplemented in G4ParticleChangeForTransport.

Definition at line 348 of file G4ParticleChange.cc.

349 {
350  // A physics process always calculates the final state of the particle
351 
352  // Take note that the return type of GetMomentumChange is a
353  // pointer to G4ParticleMometum. Also it is a normalized
354  // momentum vector.
355 
356  G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
357  G4Track* pTrack = pStep->GetTrack();
358 
359  // Set Mass/Charge
360  pPostStepPoint->SetMass(theMassChange);
361  pPostStepPoint->SetCharge(theChargeChange);
362  pPostStepPoint->SetMagneticMoment(theMagneticMomentChange);
363 
364  // update kinetic energy and momentum direction
366  pPostStepPoint->SetKineticEnergy( theEnergyChange );
367 
368  // calculate velocity
370  if (!isVelocityChanged) {
371  if(theEnergyChange > 0.0) {
373  } else if(theMassChange > 0.0) {
374  theVelocityChange = 0.0;
375  }
376  }
377  pPostStepPoint->SetVelocity(theVelocityChange);
378 
379  // update polarization
380  pPostStepPoint->SetPolarization( thePolarizationChange );
381 
382  // update position and time
383  pPostStepPoint->SetPosition( thePositionChange );
384  pPostStepPoint->AddGlobalTime(theTimeChange - theLocalTime0);
385  pPostStepPoint->SetLocalTime( theTimeChange );
386  pPostStepPoint->SetProperTime( theProperTimeChange );
387 
389  pPostStepPoint->SetWeight( theParentWeight );
390  }
391 
392 #ifdef G4VERBOSE
393  G4Track* aTrack = pStep->GetTrack();
394  if (debugFlag) CheckIt(*aTrack);
395 #endif
396 
397  // Update the G4Step specific attributes
398  return UpdateStepInfo(pStep);
399 }
void AddGlobalTime(const G4double aValue)
G4double theProperTimeChange
void SetPosition(const G4ThreeVector &aValue)
G4Step * UpdateStepInfo(G4Step *Step)
void SetMagneticMoment(G4double value)
void SetWeight(G4double aValue)
void SetMomentumDirection(const G4ThreeVector &aValue)
void SetLocalTime(const G4double aValue)
G4double theMagneticMomentChange
G4ThreeVector thePositionChange
G4ThreeVector thePolarizationChange
void SetPolarization(const G4ThreeVector &aValue)
G4double CalculateVelocity() const
Definition: G4Track.cc:222
void SetProperTime(const G4double aValue)
virtual G4bool CheckIt(const G4Track &)
G4ThreeVector theMomentumDirectionChange
void SetVelocity(G4double v)
void SetCharge(G4double value)
void SetMass(G4double value)
void SetKineticEnergy(const G4double aValue)
void SetKineticEnergy(const G4double aValue)

Here is the call graph for this function:

G4Step* G4ParticleChange::UpdateStepInfo ( G4Step Step)
protected

Here is the caller graph for this function:

Member Data Documentation

G4bool G4ParticleChange::isVelocityChanged
protected

Definition at line 235 of file G4ParticleChange.hh.

G4double G4ParticleChange::theChargeChange
protected

Definition at line 255 of file G4ParticleChange.hh.

const G4Track* G4ParticleChange::theCurrentTrack
protected

Definition at line 261 of file G4ParticleChange.hh.

G4double G4ParticleChange::theEnergyChange
protected

Definition at line 231 of file G4ParticleChange.hh.

G4double G4ParticleChange::theGlobalTime0
protected

Definition at line 241 of file G4ParticleChange.hh.

G4double G4ParticleChange::theLocalTime0
protected

Definition at line 243 of file G4ParticleChange.hh.

G4double G4ParticleChange::theMagneticMomentChange
protected

Definition at line 258 of file G4ParticleChange.hh.

G4double G4ParticleChange::theMassChange
protected

Definition at line 252 of file G4ParticleChange.hh.

G4ThreeVector G4ParticleChange::theMomentumDirectionChange
protected

Definition at line 221 of file G4ParticleChange.hh.

G4ThreeVector G4ParticleChange::thePolarizationChange
protected

Definition at line 228 of file G4ParticleChange.hh.

G4ThreeVector G4ParticleChange::thePositionChange
protected

Definition at line 238 of file G4ParticleChange.hh.

G4double G4ParticleChange::theProperTimeChange
protected

Definition at line 249 of file G4ParticleChange.hh.

G4double G4ParticleChange::theTimeChange
protected

Definition at line 246 of file G4ParticleChange.hh.

G4double G4ParticleChange::theVelocityChange
protected

Definition at line 234 of file G4ParticleChange.hh.


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