#include <G4ParticleChangeForMSC.hh>
Inheritance diagram for G4ParticleChangeForMSC:
Public Member Functions | |
G4ParticleChangeForMSC () | |
virtual | ~G4ParticleChangeForMSC () |
virtual G4Step * | UpdateStepForAlongStep (G4Step *Step) |
virtual G4Step * | UpdateStepForPostStep (G4Step *Step) |
virtual void | Initialize (const G4Track &) |
void | ProposeMomentumDirection (const G4ThreeVector &Pfinal) |
void | ProposeMomentumDirection (G4double Px, G4double Py, G4double Pz) |
const G4ThreeVector * | GetMomentumDirection () const |
const G4ThreeVector * | GetProposedMomentumDirection () const |
void | SetProposedMomentumDirection (const G4ThreeVector &Pfinal) |
const G4ThreeVector * | GetPosition () const |
void | ProposePosition (const G4ThreeVector &finalPosition) |
const G4ThreeVector * | GetProposedPosition () const |
void | SetProposedPosition (const G4ThreeVector &finalPosition) |
virtual void | DumpInfo () const |
virtual G4bool | CheckIt (const G4Track &) |
Protected Member Functions | |
G4ParticleChangeForMSC (const G4ParticleChangeForMSC &right) | |
G4ParticleChangeForMSC & | operator= (const G4ParticleChangeForMSC &right) |
Definition at line 54 of file G4ParticleChangeForMSC.hh.
G4ParticleChangeForMSC::G4ParticleChangeForMSC | ( | ) |
Definition at line 47 of file G4ParticleChangeForMSC.cc.
References G4cout, G4endl, and G4VParticleChange::verboseLevel.
00048 : G4VParticleChange() 00049 { 00050 #ifdef G4VERBOSE 00051 if (verboseLevel>2) { 00052 G4cout << "G4ParticleChangeForMSC::G4ParticleChangeForMSC() " << G4endl; 00053 } 00054 #endif 00055 }
G4ParticleChangeForMSC::~G4ParticleChangeForMSC | ( | ) | [virtual] |
Definition at line 57 of file G4ParticleChangeForMSC.cc.
References G4cout, G4endl, and G4VParticleChange::verboseLevel.
00058 { 00059 #ifdef G4VERBOSE 00060 if (verboseLevel>2) { 00061 G4cout << "G4ParticleChangeForMSC::~G4ParticleChangeForMSC() " << G4endl; 00062 } 00063 #endif 00064 }
G4ParticleChangeForMSC::G4ParticleChangeForMSC | ( | const G4ParticleChangeForMSC & | right | ) | [protected] |
Definition at line 67 of file G4ParticleChangeForMSC.cc.
References G4cout, G4endl, theMomentumDirection, thePosition, and G4VParticleChange::verboseLevel.
00068 : G4VParticleChange(right) 00069 { 00070 if (verboseLevel>1) { 00071 G4cout << "G4ParticleChangeForMSC:: copy constructor is called " << G4endl; 00072 } 00073 theMomentumDirection = right.theMomentumDirection; 00074 thePosition = right.thePosition; 00075 }
Reimplemented from G4VParticleChange.
Definition at line 169 of file G4ParticleChangeForMSC.cc.
References G4VParticleChange::accuracyForException, G4VParticleChange::accuracyForWarning, G4VParticleChange::CheckIt(), DumpInfo(), EventMustBeAborted, G4cout, G4endl, G4Exception(), G4Track::GetDefinition(), G4Track::GetKineticEnergy(), G4ParticleDefinition::GetParticleName(), and G4Track::GetPosition().
00170 { 00171 G4bool itsOK = true; 00172 G4bool exitWithError = false; 00173 00174 G4double accuracy; 00175 00176 // check 00177 00178 // MomentumDirection should be unit vector 00179 accuracy = std::fabs(theMomentumDirection.mag2()-1.0); 00180 if (accuracy > accuracyForWarning) { 00181 itsOK = false; 00182 exitWithError = (accuracy > accuracyForException); 00183 #ifdef G4VERBOSE 00184 G4cout << " G4ParticleChangeForMSC::CheckIt : "; 00185 G4cout << "the Momentum Change is not unit vector !!" 00186 << " Difference: " << accuracy << G4endl; 00187 G4cout << aTrack.GetDefinition()->GetParticleName() 00188 << " E=" << aTrack.GetKineticEnergy()/MeV 00189 << " pos=" << aTrack.GetPosition().x()/m 00190 << ", " << aTrack.GetPosition().y()/m 00191 << ", " << aTrack.GetPosition().z()/m 00192 <<G4endl; 00193 #endif 00194 } 00195 00196 // dump out information of this particle change 00197 #ifdef G4VERBOSE 00198 if (!itsOK) DumpInfo(); 00199 #endif 00200 00201 // Exit with error 00202 if (exitWithError) { 00203 G4Exception("G4ParticleChangeForMSC::CheckIt", 00204 "300", 00205 EventMustBeAborted, 00206 "momentum direction was illegal"); 00207 } 00208 //correction 00209 if (!itsOK) { 00210 G4double vmag = theMomentumDirection.mag(); 00211 theMomentumDirection = (1./vmag)*theMomentumDirection; 00212 } 00213 00214 itsOK = (itsOK) && G4VParticleChange::CheckIt(aTrack); 00215 return itsOK; 00216 }
void G4ParticleChangeForMSC::DumpInfo | ( | ) | const [virtual] |
Reimplemented from G4VParticleChange.
Definition at line 141 of file G4ParticleChangeForMSC.cc.
References G4VParticleChange::DumpInfo(), G4cout, and G4endl.
Referenced by CheckIt().
00142 { 00143 // use base-class DumpInfo 00144 G4VParticleChange::DumpInfo(); 00145 00146 G4int oldprc = G4cout.precision(3); 00147 G4cout << " Position - x (mm) : " 00148 << std::setw(20) << thePosition.x()/mm 00149 << G4endl; 00150 G4cout << " Position - y (mm) : " 00151 << std::setw(20) << thePosition.y()/mm 00152 << G4endl; 00153 G4cout << " Position - z (mm) : " 00154 << std::setw(20) << thePosition.z()/mm 00155 << G4endl; 00156 G4cout << " Momentum Direct - x : " 00157 << std::setw(20) << theMomentumDirection.x() 00158 << G4endl; 00159 G4cout << " Momentum Direct - y : " 00160 << std::setw(20) << theMomentumDirection.y() 00161 << G4endl; 00162 G4cout << " Momentum Direct - z : " 00163 << std::setw(20) << theMomentumDirection.z() 00164 << G4endl; 00165 G4cout.precision(oldprc); 00166 }
const G4ThreeVector * G4ParticleChangeForMSC::GetMomentumDirection | ( | ) | const [inline] |
const G4ThreeVector * G4ParticleChangeForMSC::GetPosition | ( | ) | const [inline] |
const G4ThreeVector * G4ParticleChangeForMSC::GetProposedMomentumDirection | ( | ) | const [inline] |
const G4ThreeVector * G4ParticleChangeForMSC::GetProposedPosition | ( | ) | const [inline] |
void G4ParticleChangeForMSC::Initialize | ( | const G4Track & | ) | [inline, virtual] |
Reimplemented from G4VParticleChange.
Definition at line 89 of file G4ParticleChangeForMSC.icc.
References G4Track::GetMomentumDirection(), G4Track::GetPosition(), G4Track::GetTrackStatus(), and G4VParticleChange::theStatusChange.
Referenced by G4VMultipleScattering::PostStepDoIt().
00090 { 00091 theStatusChange = track.GetTrackStatus(); 00092 theMomentumDirection = track.GetMomentumDirection(); 00093 thePosition = track.GetPosition(); 00094 }
G4ParticleChangeForMSC & G4ParticleChangeForMSC::operator= | ( | const G4ParticleChangeForMSC & | right | ) | [protected] |
Definition at line 78 of file G4ParticleChangeForMSC.cc.
References G4cout, G4endl, G4VParticleChange::theListOfSecondaries, G4VParticleChange::theLocalEnergyDeposit, theMomentumDirection, G4VParticleChange::theNumberOfSecondaries, thePosition, G4VParticleChange::theSizeOftheListOfSecondaries, G4VParticleChange::theStatusChange, G4VParticleChange::theSteppingControlFlag, G4VParticleChange::theTrueStepLength, and G4VParticleChange::verboseLevel.
00080 { 00081 if (verboseLevel>1) { 00082 G4cout << "G4ParticleChangeForMSC:: assignment operator is called " << G4endl; 00083 } 00084 if (this != &right) 00085 { 00086 theListOfSecondaries = right.theListOfSecondaries; 00087 theSizeOftheListOfSecondaries = right.theSizeOftheListOfSecondaries; 00088 theNumberOfSecondaries = right.theNumberOfSecondaries; 00089 theStatusChange = right.theStatusChange; 00090 theLocalEnergyDeposit = right.theLocalEnergyDeposit; 00091 theSteppingControlFlag = right.theSteppingControlFlag; 00092 theTrueStepLength = right.theTrueStepLength; 00093 00094 theMomentumDirection = right.theMomentumDirection; 00095 thePosition = right.thePosition; 00096 } 00097 return *this; 00098 }
void G4ParticleChangeForMSC::ProposeMomentumDirection | ( | G4double | Px, | |
G4double | Py, | |||
G4double | Pz | |||
) | [inline] |
Definition at line 45 of file G4ParticleChangeForMSC.icc.
00046 { 00047 theMomentumDirection.setX(Px); 00048 theMomentumDirection.setY(Py); 00049 theMomentumDirection.setZ(Pz); 00050 }
void G4ParticleChangeForMSC::ProposeMomentumDirection | ( | const G4ThreeVector & | Pfinal | ) | [inline] |
Definition at line 33 of file G4ParticleChangeForMSC.icc.
Referenced by G4VMultipleScattering::AlongStepDoIt(), G4WentzelVIRelModel::SampleScattering(), G4WentzelVIModel::SampleScattering(), G4UrbanMscModel96::SampleScattering(), G4UrbanMscModel95::SampleScattering(), G4UrbanMscModel93::SampleScattering(), G4UrbanMscModel92::SampleScattering(), G4UrbanMscModel90::SampleScattering(), and G4GoudsmitSaundersonMscModel::SampleScattering().
void G4ParticleChangeForMSC::ProposePosition | ( | const G4ThreeVector & | finalPosition | ) | [inline] |
Definition at line 83 of file G4ParticleChangeForMSC.icc.
Referenced by G4VMultipleScattering::AlongStepDoIt(), and G4VMultipleScattering::PostStepDoIt().
void G4ParticleChangeForMSC::SetProposedMomentumDirection | ( | const G4ThreeVector & | Pfinal | ) | [inline] |
void G4ParticleChangeForMSC::SetProposedPosition | ( | const G4ThreeVector & | finalPosition | ) | [inline] |
Reimplemented from G4VParticleChange.
Definition at line 104 of file G4ParticleChangeForMSC.cc.
References G4Step::GetPostStepPoint(), G4Step::GetTrack(), G4Track::GetTrackStatus(), G4StepPoint::SetMomentumDirection(), G4StepPoint::SetPosition(), G4Step::SetStepLength(), G4VParticleChange::theStatusChange, and G4VParticleChange::theTrueStepLength.
00105 { 00106 // Update the G4Step specific attributes 00107 pStep->SetStepLength(theTrueStepLength); 00108 theStatusChange = pStep->GetTrack()->GetTrackStatus(); 00109 00110 // Multiple scattering calculates the final state of the particle 00111 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint(); 00112 00113 // update momentum direction 00114 pPostStepPoint->SetMomentumDirection(theMomentumDirection); 00115 00116 // update position 00117 pPostStepPoint->SetPosition( thePosition ); 00118 return pStep; 00119 }
Reimplemented from G4VParticleChange.
Definition at line 121 of file G4ParticleChangeForMSC.cc.
References G4Step::GetPostStepPoint(), G4StepPoint::SetMomentumDirection(), and G4StepPoint::SetPosition().
00122 { 00123 // Multiple scattering calculates the final state of the particle 00124 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint(); 00125 00126 // update momentum direction 00127 pPostStepPoint->SetMomentumDirection(theMomentumDirection); 00128 00129 // update position 00130 pPostStepPoint->SetPosition( thePosition ); 00131 00132 // Update the G4Step specific attributes 00133 return pStep; 00134 }