Public Member Functions | |
Bug1021 (G4double _minStep) | |
G4VProcess * | clone () |
void | StartTracking (G4Track *track) |
virtual G4double | PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) |
G4VParticleChange * | PostStepDoIt (const G4Track &track, const G4Step &stepData) |
virtual G4double | AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection) |
virtual G4VParticleChange * | AlongStepDoIt (const G4Track &track, const G4Step &stepData) |
virtual G4double | AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition) |
virtual G4VParticleChange * | AtRestDoIt (const G4Track &track, const G4Step &stepData) |
Private Attributes | |
G4double | minStep |
bool | reflect |
ParticleChange | change |
Classes | |
class | ParticleChange |
Bug1021::Bug1021 | ( | G4double | _minStep | ) | [inline] |
References minStep, and reflect.
Referenced by clone().
00065 : 00066 G4VProcess("Bug1021",fUserDefined), change() { 00067 minStep = _minStep; 00068 reflect = false; 00069 G4ParticleTable::G4PTblDicIterator *theParticleIterator = 00070 G4ParticleTable::GetParticleTable()->GetIterator(); 00071 theParticleIterator->reset(); 00072 while((*theParticleIterator)()) { 00073 G4ParticleDefinition *pd = theParticleIterator->value(); 00074 G4ProcessManager *pmgr = pd->GetProcessManager(); 00075 if(!pmgr) continue; 00076 pmgr->AddProcess(this,-1,-1,3); 00077 pmgr->SetProcessOrdering(this,idxPostStep); 00078 pmgr->SetProcessOrdering(this,idxAtRest); 00079 } 00080 }
void Bug1021::StartTracking | ( | G4Track * | track | ) | [inline] |
virtual G4double Bug1021::PostStepGetPhysicalInteractionLength | ( | const G4Track & | track, | |
G4double | previousStepSize, | |||
G4ForceCondition * | condition | |||
) | [inline, virtual] |
References E, BLGlobalField::GetFieldValue(), BLGlobalField::getObject(), minStep, and reflect.
00089 { *condition = NotForced; 00090 G4double q = track.GetDefinition()->GetPDGCharge(); 00091 if(q == 0.0) return DBL_MAX; 00092 const G4ThreeVector &pos = track.GetPosition(); 00093 G4double point[4], field[6]; 00094 point[0] = pos[0]; 00095 point[1] = pos[1]; 00096 point[2] = pos[2]; 00097 point[3] = track.GetGlobalTime(); 00098 BLGlobalField::getObject()->GetFieldValue(point,field); 00099 G4ThreeVector E(field[3],field[4],field[5]); 00100 G4double Etot = E.mag(); 00101 if(Etot < 100*volt/meter) return DBL_MAX; 00102 G4ThreeVector dir = track.GetMomentumDirection(); 00103 if((q>0.0?1.0:-1.0)*dir.dot(E.unit()) > -0.99) return DBL_MAX; 00104 double distanceToStop = track.GetKineticEnergy()/Etot; 00105 double step = distanceToStop / 2.0; 00106 if(distanceToStop < minStep) { 00107 step = 0.0; 00108 reflect = true; 00109 } 00110 return step; 00111 }
G4VParticleChange* Bug1021::PostStepDoIt | ( | const G4Track & | track, | |
const G4Step & | stepData | |||
) | [inline] |
References change, minStep, and reflect.
00114 { change.Initialize(track); 00115 if(reflect) { 00116 G4ThreeVector dir=track.GetMomentumDirection(); 00117 //@@ need to reflect dir in plane perpendicular to E 00118 // but dir || E, so this is good enough 00119 dir[0] = -dir[0]; 00120 dir[1] = -dir[1]; 00121 dir[2] = -dir[2]; 00122 //@@ approximation for elapsed time 00123 G4double time = track.GetGlobalTime(); 00124 time += minStep/track.GetVelocity(); 00125 change.ProposeGlobalTime(time); 00126 change.ProposeMomentumDirection(dir); 00127 reflect = false; 00128 } 00129 return &change; 00130 }
virtual G4double Bug1021::AlongStepGetPhysicalInteractionLength | ( | const G4Track & | track, | |
G4double | previousStepSize, | |||
G4double | currentMinimumStep, | |||
G4double & | proposedSafety, | |||
G4GPILSelection * | selection | |||
) | [inline, virtual] |
virtual G4VParticleChange* Bug1021::AlongStepDoIt | ( | const G4Track & | track, | |
const G4Step & | stepData | |||
) | [inline, virtual] |
virtual G4double Bug1021::AtRestGetPhysicalInteractionLength | ( | const G4Track & | track, | |
G4ForceCondition * | condition | |||
) | [inline, virtual] |
virtual G4VParticleChange* Bug1021::AtRestDoIt | ( | const G4Track & | track, | |
const G4Step & | stepData | |||
) | [inline, virtual] |
G4double Bug1021::minStep [private] |
Referenced by Bug1021(), PostStepDoIt(), and PostStepGetPhysicalInteractionLength().
bool Bug1021::reflect [private] |
Referenced by Bug1021(), PostStepDoIt(), PostStepGetPhysicalInteractionLength(), and StartTracking().
ParticleChange Bug1021::change [private] |
Referenced by AlongStepDoIt(), AtRestDoIt(), and PostStepDoIt().