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

#include <G4ParticleChangeForLoss.hh>

+ Inheritance diagram for G4ParticleChangeForLoss:

Public Member Functions

 G4ParticleChangeForLoss ()
 
virtual ~G4ParticleChangeForLoss ()
 
G4StepUpdateStepForAlongStep (G4Step *Step)
 
G4StepUpdateStepForPostStep (G4Step *Step)
 
void InitializeForAlongStep (const G4Track &)
 
void InitializeForPostStep (const G4Track &)
 
void AddSecondary (G4DynamicParticle *aParticle)
 
G4double GetProposedCharge () const
 
void SetProposedCharge (G4double theCharge)
 
G4double GetCharge () const
 
void ProposeCharge (G4double finalCharge)
 
G4double GetProposedKineticEnergy () const
 
void SetProposedKineticEnergy (G4double proposedKinEnergy)
 
const G4ThreeVectorGetProposedMomentumDirection () const
 
void SetProposedMomentumDirection (const G4ThreeVector &dir)
 
const G4ThreeVectorGetMomentumDirection () const
 
void ProposeMomentumDirection (G4double Px, G4double Py, G4double Pz)
 
void ProposeMomentumDirection (const G4ThreeVector &Pfinal)
 
const G4ThreeVectorGetProposedPolarization () const
 
void ProposePolarization (const G4ThreeVector &dir)
 
void ProposePolarization (G4double Px, G4double Py, G4double Pz)
 
const G4TrackGetCurrentTrack () const
 
void SetLowEnergyLimit (G4double elimit)
 
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

 G4ParticleChangeForLoss (const G4ParticleChangeForLoss &right)
 
G4ParticleChangeForLossoperator= (const G4ParticleChangeForLoss &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
 

Additional Inherited Members

- 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
 
- Static Protected Attributes inherited from G4VParticleChange
static const G4double accuracyForWarning = 1.0e-9
 
static const G4double accuracyForException = 0.001
 

Detailed Description

Definition at line 60 of file G4ParticleChangeForLoss.hh.

Constructor & Destructor Documentation

◆ G4ParticleChangeForLoss() [1/2]

G4ParticleChangeForLoss::G4ParticleChangeForLoss ( )

Definition at line 52 of file G4ParticleChangeForLoss.cc.

53 : G4VParticleChange(), currentTrack(0), proposedKinEnergy(0.),
54 lowEnergyLimit(1.0*eV), currentCharge(0.)
55{
57 debugFlag = false;
58#ifdef G4VERBOSE
59 if (verboseLevel>2) {
60 G4cout << "G4ParticleChangeForLoss::G4ParticleChangeForLoss() " << G4endl;
61 }
62#endif
63}
@ NormalCondition
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
G4SteppingControl theSteppingControlFlag

◆ ~G4ParticleChangeForLoss()

G4ParticleChangeForLoss::~G4ParticleChangeForLoss ( )
virtual

Definition at line 65 of file G4ParticleChangeForLoss.cc.

66{
67#ifdef G4VERBOSE
68 if (verboseLevel>2) {
69 G4cout << "G4ParticleChangeForLoss::~G4ParticleChangeForLoss() " << G4endl;
70 }
71#endif
72}

◆ G4ParticleChangeForLoss() [2/2]

G4ParticleChangeForLoss::G4ParticleChangeForLoss ( const G4ParticleChangeForLoss right)
protected

Definition at line 74 of file G4ParticleChangeForLoss.cc.

76 : G4VParticleChange(right)
77{
78 if (verboseLevel>1) {
79 G4cout << "G4ParticleChangeForLoss:: copy constructor is called " << G4endl;
80 }
81 currentTrack = right.currentTrack;
82 proposedKinEnergy = right.proposedKinEnergy;
83 lowEnergyLimit = right.lowEnergyLimit;
84 currentCharge = right.currentCharge;
85 proposedMomentumDirection = right.proposedMomentumDirection;
86}

Member Function Documentation

◆ AddSecondary()

void G4ParticleChangeForLoss::AddSecondary ( G4DynamicParticle aParticle)
inline

Definition at line 260 of file G4ParticleChangeForLoss.hh.

261{
262 // create track
263 G4Track* aTrack = new G4Track(aParticle, currentTrack->GetGlobalTime(),
264 currentTrack->GetPosition());
265
266 // Touchable handle is copied to keep the pointer
267 aTrack->SetTouchableHandle(currentTrack->GetTouchableHandle());
268
269 // add a secondary
271}
const G4ThreeVector & GetPosition() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4double GetGlobalTime() const
const G4TouchableHandle & GetTouchableHandle() const
void AddSecondary(G4Track *aSecondary)

◆ CheckIt()

G4bool G4ParticleChangeForLoss::CheckIt ( const G4Track aTrack)
virtual

Reimplemented from G4VParticleChange.

Definition at line 160 of file G4ParticleChangeForLoss.cc.

161{
162 G4bool itsOK = true;
163 G4bool exitWithError = false;
164
165 G4double accuracy;
166
167 // Energy should not be lager than initial value
168 accuracy = ( proposedKinEnergy - aTrack.GetKineticEnergy())/MeV;
169 if (accuracy > accuracyForWarning) {
170 itsOK = false;
171 exitWithError = (accuracy > accuracyForException);
172#ifdef G4VERBOSE
173 G4cout << "G4ParticleChangeForLoss::CheckIt: ";
174 G4cout << "KinEnergy become larger than the initial value!"
175 << " Difference: " << accuracy << "[MeV] " <<G4endl;
177 << " E=" << aTrack.GetKineticEnergy()/MeV
178 << " pos=" << aTrack.GetPosition().x()/m
179 << ", " << aTrack.GetPosition().y()/m
180 << ", " << aTrack.GetPosition().z()/m
181 <<G4endl;
182#endif
183 }
184
185 // dump out information of this particle change
186#ifdef G4VERBOSE
187 if (!itsOK) DumpInfo();
188#endif
189
190 // Exit with error
191 if (exitWithError) {
192 G4Exception("G4ParticleChangeForLoss::CheckIt",
193 "TRACK004", EventMustBeAborted,
194 "energy was illegal");
195 }
196
197 //correction
198 if (!itsOK) {
199 proposedKinEnergy = aTrack.GetKineticEnergy();
200 }
201
202 itsOK = (itsOK) && G4VParticleChange::CheckIt(aTrack);
203 return itsOK;
204}
@ EventMustBeAborted
double G4double
Definition: G4Types.hh:64
bool G4bool
Definition: G4Types.hh:67
double z() const
double x() const
double y() const
const G4String & GetParticleName() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
virtual G4bool CheckIt(const G4Track &)
static const G4double accuracyForException
static const G4double accuracyForWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

◆ DumpInfo()

void G4ParticleChangeForLoss::DumpInfo ( ) const
virtual

Reimplemented from G4VParticleChange.

Definition at line 136 of file G4ParticleChangeForLoss.cc.

137{
138// use base-class DumpInfo
140
141 G4int oldprc = G4cout.precision(3);
142 G4cout << " Charge (eplus) : "
143 << std::setw(20) << currentCharge/eplus
144 << G4endl;
145 G4cout << " Kinetic Energy (MeV): "
146 << std::setw(20) << proposedKinEnergy/MeV
147 << G4endl;
148 G4cout << " Momentum Direct - x : "
149 << std::setw(20) << proposedMomentumDirection.x()
150 << G4endl;
151 G4cout << " Momentum Direct - y : "
152 << std::setw(20) << proposedMomentumDirection.y()
153 << G4endl;
154 G4cout << " Momentum Direct - z : "
155 << std::setw(20) << proposedMomentumDirection.z()
156 << G4endl;
157 G4cout.precision(oldprc);
158}
int G4int
Definition: G4Types.hh:66
virtual void DumpInfo() const

Referenced by CheckIt().

◆ GetCharge()

G4double G4ParticleChangeForLoss::GetCharge ( ) const
inline

Definition at line 160 of file G4ParticleChangeForLoss.hh.

161{
162 return currentCharge;
163}

◆ GetCurrentTrack()

const G4Track * G4ParticleChangeForLoss::GetCurrentTrack ( ) const
inline

Definition at line 207 of file G4ParticleChangeForLoss.hh.

208{
209 return currentTrack;
210}

Referenced by G4PolarizedMollerBhabhaModel::SampleSecondaries().

◆ GetMomentumDirection()

const G4ThreeVector & G4ParticleChangeForLoss::GetMomentumDirection ( ) const
inline

Definition at line 182 of file G4ParticleChangeForLoss.hh.

183{
184 return proposedMomentumDirection;
185}

◆ GetProposedCharge()

G4double G4ParticleChangeForLoss::GetProposedCharge ( ) const
inline

Definition at line 155 of file G4ParticleChangeForLoss.hh.

156{
157 return currentCharge;
158}

◆ GetProposedKineticEnergy()

G4double G4ParticleChangeForLoss::GetProposedKineticEnergy ( ) const
inline

Definition at line 145 of file G4ParticleChangeForLoss.hh.

146{
147 return proposedKinEnergy;
148}

Referenced by G4EmBiasingManager::ApplySecondaryBiasing(), and G4VEnergyLossProcess::PostStepDoIt().

◆ GetProposedMomentumDirection()

const G4ThreeVector & G4ParticleChangeForLoss::GetProposedMomentumDirection ( ) const
inline

Definition at line 176 of file G4ParticleChangeForLoss.hh.

177{
178 return proposedMomentumDirection;
179}

Referenced by G4EmBiasingManager::ApplySecondaryBiasing(), and G4ePolarizedBremsstrahlungModel::SampleSecondaries().

◆ GetProposedPolarization()

const G4ThreeVector & G4ParticleChangeForLoss::GetProposedPolarization ( ) const
inline

Definition at line 213 of file G4ParticleChangeForLoss.hh.

214{
215 return proposedPolarization;
216}

◆ InitializeForAlongStep()

void G4ParticleChangeForLoss::InitializeForAlongStep ( const G4Track track)
inline

Definition at line 232 of file G4ParticleChangeForLoss.hh.

233{
238 theParentWeight = track.GetWeight();
240 proposedKinEnergy = track.GetKineticEnergy();
241 currentCharge = track.GetDynamicParticle()->GetCharge();
242}
G4double GetCharge() const
G4TrackStatus GetTrackStatus() const
G4double GetWeight() const
const G4DynamicParticle * GetDynamicParticle() const
G4TrackStatus theStatusChange
G4double theNonIonizingEnergyDeposit
void InitializeSecondaries(const G4Track &)

Referenced by G4VEnergyLossProcess::AlongStepDoIt(), and G4NuclearStopping::AlongStepDoIt().

◆ InitializeForPostStep()

void G4ParticleChangeForLoss::InitializeForPostStep ( const G4Track track)
inline

Definition at line 244 of file G4ParticleChangeForLoss.hh.

245{
250 theParentWeight = track.GetWeight();
252 proposedKinEnergy = track.GetKineticEnergy();
253 currentCharge = track.GetDynamicParticle()->GetCharge();
254 proposedMomentumDirection = track.GetMomentumDirection();
255 proposedPolarization = track.GetPolarization();
256 currentTrack = &track;
257}
const G4ThreeVector & GetMomentumDirection() const
const G4ThreeVector & GetPolarization() const

Referenced by G4VEnergyLossProcess::PostStepDoIt().

◆ operator=()

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

Definition at line 89 of file G4ParticleChangeForLoss.cc.

91{
92#ifdef G4VERBOSE
93 if (verboseLevel>1) {
94 G4cout << "G4ParticleChangeForLoss:: assignment operator is called " << G4endl;
95 }
96#endif
97
98 if (this != &right) {
100#ifdef G4VERBOSE
101 if (verboseLevel>0) {
102 G4cout << "G4ParticleChangeForLoss: assignment operator Warning ";
103 G4cout << "theListOfSecondaries is not empty ";
104 }
105#endif
106 for (G4int index= 0; index<theNumberOfSecondaries; index++){
107 if ( (*theListOfSecondaries)[index] ) delete (*theListOfSecondaries)[index] ;
108 }
109 }
110 delete theListOfSecondaries;
113 for (G4int index = 0; index<theNumberOfSecondaries; index++){
114 G4Track* newTrack = new G4Track(*((*right.theListOfSecondaries)[index] ));
115 theListOfSecondaries->SetElement(index, newTrack); }
116
123
124 currentTrack = right.currentTrack;
125 proposedKinEnergy = right.proposedKinEnergy;
126 currentCharge = right.currentCharge;
127 proposedMomentumDirection = right.proposedMomentumDirection;
128 }
129 return *this;
130}
G4FastVector< G4Track, G4TrackFastVectorSize > G4TrackFastVector
void SetElement(G4int anIndex, Type *anElement)
Definition: G4FastVector.hh:76
G4TrackFastVector * theListOfSecondaries

◆ ProposeCharge()

void G4ParticleChangeForLoss::ProposeCharge ( G4double  finalCharge)
inline

Definition at line 170 of file G4ParticleChangeForLoss.hh.

171{
172 currentCharge = theCharge;
173}

◆ ProposeMomentumDirection() [1/2]

void G4ParticleChangeForLoss::ProposeMomentumDirection ( const G4ThreeVector Pfinal)
inline

Definition at line 188 of file G4ParticleChangeForLoss.hh.

189{
190 proposedMomentumDirection = dir;
191}

◆ ProposeMomentumDirection() [2/2]

void G4ParticleChangeForLoss::ProposeMomentumDirection ( G4double  Px,
G4double  Py,
G4double  Pz 
)
inline

◆ ProposePolarization() [1/2]

void G4ParticleChangeForLoss::ProposePolarization ( const G4ThreeVector dir)
inline

Definition at line 219 of file G4ParticleChangeForLoss.hh.

220{
221 proposedPolarization = dir;
222}

Referenced by G4ePolarizedBremsstrahlungModel::SampleSecondaries(), and G4PolarizedMollerBhabhaModel::SampleSecondaries().

◆ ProposePolarization() [2/2]

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

Definition at line 225 of file G4ParticleChangeForLoss.hh.

226{
227 proposedPolarization.setX(Px);
228 proposedPolarization.setY(Py);
229 proposedPolarization.setZ(Pz);
230}

◆ SetLowEnergyLimit()

void G4ParticleChangeForLoss::SetLowEnergyLimit ( G4double  elimit)
inline

Definition at line 273 of file G4ParticleChangeForLoss.hh.

274{
275 lowEnergyLimit = elimit;
276}

Referenced by G4VEnergyLossProcess::BuildPhysicsTable().

◆ SetProposedCharge()

void G4ParticleChangeForLoss::SetProposedCharge ( G4double  theCharge)
inline

Definition at line 165 of file G4ParticleChangeForLoss.hh.

166{
167 currentCharge = theCharge;
168}

Referenced by G4VEnergyLossProcess::AlongStepDoIt().

◆ SetProposedKineticEnergy()

◆ SetProposedMomentumDirection()

◆ UpdateStepForAlongStep()

G4Step * G4ParticleChangeForLoss::UpdateStepForAlongStep ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Definition at line 210 of file G4ParticleChangeForLoss.cc.

211{
212 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
213
214 // accumulate change of the kinetic energy
215 G4double kinEnergy = pPostStepPoint->GetKineticEnergy() +
216 (proposedKinEnergy - pStep->GetPreStepPoint()->GetKineticEnergy());
217
218 // update kinetic energy and charge
219 if (kinEnergy < lowEnergyLimit) {
220 theLocalEnergyDeposit += kinEnergy;
221 kinEnergy = 0.0;
222 } else {
223 pPostStepPoint->SetCharge( currentCharge );
224 }
225 pPostStepPoint->SetKineticEnergy( kinEnergy );
226 // calculate velocity
227 pStep->GetTrack()->SetKineticEnergy(kinEnergy);
228 pPostStepPoint->SetVelocity(pStep->GetTrack()->CalculateVelocity());
229 pStep->GetTrack()->SetKineticEnergy(pStep->GetPreStepPoint()->GetKineticEnergy());
230
232 pPostStepPoint->SetWeight( theParentWeight );
233 }
234
235 pStep->AddTotalEnergyDeposit( theLocalEnergyDeposit );
236 pStep->AddNonIonizingEnergyDeposit( theNonIonizingEnergyDeposit );
237 return pStep;
238}
void SetKineticEnergy(const G4double aValue)
void SetWeight(G4double aValue)
void SetCharge(G4double value)
void SetVelocity(G4double v)
G4double GetKineticEnergy() const

◆ UpdateStepForPostStep()

G4Step * G4ParticleChangeForLoss::UpdateStepForPostStep ( G4Step Step)
virtual

Reimplemented from G4VParticleChange.

Definition at line 240 of file G4ParticleChangeForLoss.cc.

241{
242 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
243 pPostStepPoint->SetCharge( currentCharge );
244 pPostStepPoint->SetMomentumDirection( proposedMomentumDirection );
245 pPostStepPoint->SetKineticEnergy( proposedKinEnergy );
246 pStep->GetTrack()->SetKineticEnergy( proposedKinEnergy );
247 pPostStepPoint->SetVelocity(pStep->GetTrack()->CalculateVelocity());
248 pPostStepPoint->SetPolarization( proposedPolarization );
249
251 pPostStepPoint->SetWeight( theParentWeight );
252 }
253
254 pStep->AddTotalEnergyDeposit( theLocalEnergyDeposit );
255 pStep->AddNonIonizingEnergyDeposit( theNonIonizingEnergyDeposit );
256 return pStep;
257}
void SetMomentumDirection(const G4ThreeVector &aValue)
void SetPolarization(const G4ThreeVector &aValue)

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