Belle II Software development
G4MonopoleTransportation.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9// modified from GEANT4 exoticphysics/monopole/*
10
11#include <simulation/monopoles/G4MonopoleTransportation.h>
12#include <simulation/monopoles/G4Monopole.h>
13
14#include <G4Navigator.hh>
15#include <G4PropagatorInField.hh>
16#include <G4Track.hh>
17#include <G4Step.hh>
18#include <G4ProductionCutsTable.hh>
19#include <G4ParticleTable.hh>
20#include <G4ChordFinder.hh>
21#include <G4SafetyHelper.hh>
22#include <G4FieldManagerStore.hh>
23#include <G4TransportationProcessType.hh>
24#include <CLHEP/Units/SystemOfUnits.h>
25#include <framework/logging/Logger.h>
26
27class G4VSensitiveDetector;
28
29using namespace std;
30using namespace Belle2;
31using namespace Belle2::Monopoles;
32using namespace CLHEP;
33
35 G4int verb)
36 : G4VProcess(G4String("MonopoleTransportation"), fTransportation),
37 fParticleDef(mpl),
38 fMagSetup(0),
41 fParticleIsLooping(false),
42 fPreviousSftOrigin(0., 0., 0.),
43 fPreviousSafety(0.0),
44 fThreshold_Trap_Energy(10 * MeV),
50 fShortStepOptimisation(false), // Old default: true (=fast short steps)
52{
53 verboseLevel = verb;
54
55 // set Process Sub Type
56 SetProcessSubType(TRANSPORTATION);
57
59
60 G4TransportationManager* transportMgr ;
61
62 transportMgr = G4TransportationManager::GetTransportationManager() ;
63
64 fLinearNavigator = transportMgr->GetNavigatorForTracking() ;
65
66 // fGlobalFieldMgr = transportMgr->GetFieldManager() ;
67
68 fFieldPropagator = transportMgr->GetPropagatorInField() ;
69
70 fpSafetyHelper = transportMgr->GetSafetyHelper(); // New
71
72 // Cannot determine whether a field exists here,
73 // because it would only work if the field manager has informed
74 // about the detector's field before this transportation process
75 // is constructed.
76 // Instead later the method DoesGlobalFieldExist() is called
77
78 static G4TouchableHandle nullTouchableHandle; // Points to (G4VTouchable*) 0
79 fCurrentTouchableHandle = nullTouchableHandle;
80
83}
84
86{
87 if ((verboseLevel > 0) && (fSumEnergyKilled > 0.0)) {
88 B2DEBUG(25, " G4MonopoleTransportation: Statistics for looping particles ");
89 B2DEBUG(25, " Sum of energy of loopers killed: " << fSumEnergyKilled);
90 B2DEBUG(25, " Max energy of loopers killed: " << fMaxEnergyKilled);
91 }
92}
93
94// Responsibilities:
95// Find whether the geometry limits the Step, and to what length
96// Calculate the new value of the safety and return it.
97// Store the final time, position and momentum.
98
100AlongStepGetPhysicalInteractionLength(const G4Track& track,
101 G4double, // previousStepSize
102 G4double currentMinimumStep,
103 G4double& currentSafety,
104 G4GPILSelection* selection)
105{
106 fMagSetup->SwitchChordFinder(1);
107 // change to monopole equation
108
109 G4double geometryStepLength, newSafety ;
110 fParticleIsLooping = false ;
111
112 // Initial actions moved to StartTrack()
113 // --------------------------------------
114 // Note: in case another process changes touchable handle
115 // it will be necessary to add here (for all steps)
116 // fCurrentTouchableHandle = aTrack->GetTouchableHandle();
117
118 // GPILSelection is set to defaule value of CandidateForSelection
119 // It is a return value
120 //
121 *selection = CandidateForSelection ;
122
123 // Get initial Energy/Momentum of the track
124 //
125 const G4DynamicParticle* pParticle = track.GetDynamicParticle() ;
126 G4ThreeVector startMomentumDir = pParticle->GetMomentumDirection() ;
127 G4ThreeVector startPosition = track.GetPosition() ;
128
129 // G4double theTime = track.GetGlobalTime() ;
130
131 // The Step Point safety can be limited by other geometries and/or the
132 // assumptions of any process - it's not always the geometrical safety.
133 // We calculate the starting point's isotropic safety here.
134 //
135 G4ThreeVector OriginShift = startPosition - fPreviousSftOrigin ;
136 G4double MagSqShift = OriginShift.mag2() ;
137 if (MagSqShift >= sqr(fPreviousSafety)) {
138 currentSafety = 0.0 ;
139 } else {
140 currentSafety = fPreviousSafety - std::sqrt(MagSqShift) ;
141 }
142
143 // Is the monopole charged ?
144 //
145 G4double particleMagneticCharge = fParticleDef->MagneticCharge() ;
146 G4double particleElectricCharge = pParticle->GetCharge();
147
148 fGeometryLimitedStep = false ;
149 // fEndGlobalTimeComputed = false ;
150
151 // There is no need to locate the current volume. It is Done elsewhere:
152 // On track construction
153 // By the tracking, after all AlongStepDoIts, in "Relocation"
154
155 // Check whether the particle have an (EM) field force exerting upon it
156 //
157 G4FieldManager* fieldMgr = nullptr;
158 G4bool fieldExertsForce = false ;
159
160 if ((particleMagneticCharge != 0.0)) {
161 fieldMgr = fFieldPropagator->FindAndSetFieldManager(track.GetVolume());
162 if (fieldMgr != 0) {
163 // Message the field Manager, to configure it for this track
164 fieldMgr->ConfigureForTrack(&track);
165 // Moved here, in order to allow a transition
166 // from a zero-field status (with fieldMgr->(field)0
167 // to a finite field status
168
169 // If the field manager has no field, there is no field !
170 fieldExertsForce = (fieldMgr->GetDetectorField() != 0);
171 }
172 }
173
174 // Choose the calculation of the transportation: Field or not
175 //
176 if (!fieldExertsForce) {
177 G4double linearStepLength ;
178 if (fShortStepOptimisation && (currentMinimumStep <= currentSafety)) {
179 // The Step is guaranteed to be taken
180 //
181 geometryStepLength = currentMinimumStep ;
182 fGeometryLimitedStep = false ;
183 } else {
184 // Find whether the straight path intersects a volume
185 //
186 linearStepLength = fLinearNavigator->ComputeStep(startPosition,
187 startMomentumDir,
188 currentMinimumStep,
189 newSafety) ;
190 // Remember last safety origin & value.
191 //
192 fPreviousSftOrigin = startPosition ;
193 fPreviousSafety = newSafety ;
194 // fpSafetyHelper->SetCurrentSafety( newSafety, startPosition);
195
196 // The safety at the initial point has been re-calculated:
197 //
198 currentSafety = newSafety ;
199
200 fGeometryLimitedStep = (linearStepLength <= currentMinimumStep);
202 // The geometry limits the Step size (an intersection was found.)
203 geometryStepLength = linearStepLength ;
204 } else {
205 // The full Step is taken.
206 geometryStepLength = currentMinimumStep ;
207 }
208 }
209 endpointDistance = geometryStepLength ;
210
211 // Calculate final position
212 //
213 fTransportEndPosition = startPosition + geometryStepLength * startMomentumDir ;
214
215 // Momentum direction, energy and polarisation are unchanged by transport
216 //
217 fTransportEndMomentumDir = startMomentumDir ;
218 fTransportEndKineticEnergy = track.GetKineticEnergy() ;
219 fTransportEndSpin = track.GetPolarization();
220 fParticleIsLooping = false ;
221 fMomentumChanged = false ;
222 fEndGlobalTimeComputed = false ;
223 } else { // A field exerts force
224 G4double momentumMagnitude = pParticle->GetTotalMomentum() ;
225 G4ThreeVector EndUnitMomentum ;
226 G4double restMass = fParticleDef->GetPDGMass() ;
227
228 G4ChargeState chargeState(particleElectricCharge, // The charge can change (dynamic)
229 fParticleDef->GetPDGSpin(),
230 0, // Magnetic moment: pParticleDef->GetMagneticMoment(),
231 0, // Electric Dipole moment - not in Particle Definition
232 particleMagneticCharge); // in units of the positron charge
233
234 G4EquationOfMotion* equationOfMotion =
235 (fFieldPropagator->GetChordFinder()->GetIntegrationDriver()->GetStepper())
236 ->GetEquationOfMotion();
237
238 equationOfMotion
239 ->SetChargeMomentumMass(chargeState,
240 momentumMagnitude,
241 restMass) ;
242 // SetChargeMomentumMass now passes both the electric and magnetic charge - in chargeState
243
244 G4ThreeVector spin = track.GetPolarization() ;
245 G4FieldTrack aFieldTrack = G4FieldTrack(startPosition,
246 track.GetMomentumDirection(),
247 0.0,
248 track.GetKineticEnergy(),
249 restMass,
250 track.GetVelocity(),
251 track.GetGlobalTime(), // Lab.
252 track.GetProperTime(), // Part.
253 &spin) ;
254 if (currentMinimumStep > 0) {
255 // Do the Transport in the field (non recti-linear)
256 //
257 G4double lengthAlongCurve = fFieldPropagator->ComputeStep(
258 aFieldTrack,
259 currentMinimumStep,
260 currentSafety,
261 track.GetVolume()) ;
262 fGeometryLimitedStep = lengthAlongCurve < currentMinimumStep;
264 geometryStepLength = lengthAlongCurve ;
265 } else {
266 geometryStepLength = currentMinimumStep ;
267 }
268 } else {
269 geometryStepLength = 0.0 ;
270 fGeometryLimitedStep = false ;
271 }
272
273 // Remember last safety origin & value.
274 //
275 fPreviousSftOrigin = startPosition ;
276 fPreviousSafety = currentSafety ;
277 // fpSafetyHelper->SetCurrentSafety( newSafety, startPosition);
278
279 // Get the End-Position and End-Momentum (Dir-ection)
280 //
281 fTransportEndPosition = aFieldTrack.GetPosition() ;
282
283 // Momentum: Magnitude and direction can be changed too now ...
284 //
285 fMomentumChanged = true ;
286 fTransportEndMomentumDir = aFieldTrack.GetMomentumDir() ;
287
288 fTransportEndKineticEnergy = aFieldTrack.GetKineticEnergy() ;
289
290 fCandidateEndGlobalTime = aFieldTrack.GetLabTimeOfFlight();
292
293 fTransportEndSpin = aFieldTrack.GetSpin();
294 fParticleIsLooping = fFieldPropagator->IsParticleLooping() ;
295 endpointDistance = (fTransportEndPosition - startPosition).mag() ;
296 }
297
298 // If we are asked to go a step length of 0, and we are on a boundary
299 // then a boundary will also limit the step -> we must flag this.
300 //
301 if (currentMinimumStep == 0.0) {
302 if (currentSafety == 0.0) fGeometryLimitedStep = true ;
303 }
304
305 // Update the safety starting from the end-point,
306 // if it will become negative at the end-point.
307 //
308 if (currentSafety < endpointDistance) {
309 // if( particleMagneticCharge == 0.0 )
310 // G4cout << " Avoiding call to ComputeSafety : charge = 0.0 " << G4endl;
311
312 if (particleMagneticCharge != 0.0) {
313
314 G4double endSafety =
316 currentSafety = endSafety ;
318 fPreviousSafety = currentSafety ;
319 fpSafetyHelper->SetCurrentSafety(currentSafety, fTransportEndPosition);
320
321 // Because the Stepping Manager assumes it is from the start point,
322 // add the StepLength
323 //
324 currentSafety += endpointDistance ;
325
326#ifdef G4DEBUG_TRANSPORT
327 G4cout.precision(12) ;
328 G4cout << "***G4MonopoleTransportation::AlongStepGPIL ** " << G4endl ;
329 G4cout << " Called Navigator->ComputeSafety at " << fTransportEndPosition
330 << " and it returned safety= " << endSafety << G4endl ;
331 G4cout << " Adding endpoint distance " << endpointDistance
332 << " to obtain pseudo-safety= " << currentSafety << G4endl ;
333#endif
334 }
335 }
336
337 fParticleChange.ProposeTrueStepLength(geometryStepLength) ;
338
339 fMagSetup->SwitchChordFinder(0);
340 // change back to usual equation
341
342 return geometryStepLength ;
343}
344
345//
346// Initialize ParticleChange (by setting all its members equal
347// to corresponding members in G4Track)
348
349G4VParticleChange* G4MonopoleTransportation::AlongStepDoIt(const G4Track& track,
350 const G4Step& stepData)
351{
352#ifdef G4VERBOSE
353 static G4int noCalls = 0;
354 noCalls++;
355#endif
356
357 static const G4ParticleDefinition* fOpticalPhoton =
358 G4ParticleTable::GetParticleTable()->FindParticle("opticalphoton");
359
360 fParticleChange.Initialize(track) ;
361
362 // Code for specific process
363 //
364 fParticleChange.ProposePosition(fTransportEndPosition) ;
365 fParticleChange.ProposeMomentumDirection(fTransportEndMomentumDir) ;
367 fParticleChange.SetMomentumChanged(fMomentumChanged) ;
368
369 fParticleChange.ProposePolarization(fTransportEndSpin);
370
371 G4double deltaTime = 0.0 ;
372
373 // Calculate Lab Time of Flight (ONLY if field Equations used it!)
374 // G4double endTime = fCandidateEndGlobalTime;
375 // G4double delta_time = endTime - startTime;
376
377 G4double startTime = track.GetGlobalTime() ;
378
380 // The time was not integrated .. make the best estimate possible
381 //
382 G4double finalVelocity = track.GetVelocity() ;
383 G4double initialVelocity = stepData.GetPreStepPoint()->GetVelocity() ;
384 G4double stepLength = track.GetStepLength() ;
385
386 deltaTime = 0.0; // in case initialVelocity = 0
387 const G4DynamicParticle* fpDynamicParticle = track.GetDynamicParticle();
388 if (fpDynamicParticle->GetDefinition() == fOpticalPhoton) {
389 // A photon is in the medium of the final point
390 // during the step, so it has the final velocity.
391 deltaTime = stepLength / finalVelocity ;
392 } else if (finalVelocity > 0.0) {
393 G4double meanInverseVelocity ;
394 // deltaTime = stepLength/finalVelocity ;
395 meanInverseVelocity = 0.5
396 * (1.0 / initialVelocity + 1.0 / finalVelocity) ;
397 deltaTime = stepLength * meanInverseVelocity ;
398 } else if (initialVelocity > 0.0) {
399 deltaTime = stepLength / initialVelocity ;
400 }
401 fCandidateEndGlobalTime = startTime + deltaTime ;
402 } else {
403 deltaTime = fCandidateEndGlobalTime - startTime ;
404 }
405
406 fParticleChange.ProposeGlobalTime(fCandidateEndGlobalTime) ;
407
408 // Now Correct by Lorentz factor to get "proper" deltaTime
409
410 G4double restMass = track.GetDynamicParticle()->GetMass() ;
411 G4double deltaProperTime = deltaTime * (restMass / track.GetTotalEnergy()) ;
412
413 fParticleChange.ProposeProperTime(track.GetProperTime() + deltaProperTime) ;
414 //fParticleChange. ProposeTrueStepLength( track.GetStepLength() ) ;
415
416 // If the particle is caught looping or is stuck (in very difficult
417 // boundaries) in a magnetic field (doing many steps)
418 // THEN this kills it ...
419 //
420 if (fParticleIsLooping) {
421 G4double endEnergy = fTransportEndKineticEnergy;
422
423 if ((endEnergy < fThreshold_Important_Energy)
425 // Kill the looping particle
426 //
427 fParticleChange.ProposeTrackStatus(fStopAndKill) ;
428
429 // 'Bare' statistics
430 fSumEnergyKilled += endEnergy;
431 if (endEnergy > fMaxEnergyKilled) { fMaxEnergyKilled = endEnergy; }
432
433#ifdef G4VERBOSE
434 if ((verboseLevel > 1) ||
435 (endEnergy > fThreshold_Warning_Energy)) {
436 B2DEBUG(25, " G4MonopoleTransportation is killing track that is looping or stuck ");
437 B2DEBUG(25, " This track has " << track.GetKineticEnergy() / MeV << " MeV energy.");
438 B2DEBUG(25, " Number of trials = " << fNoLooperTrials
439 << " No of calls to AlongStepDoIt = " << noCalls);
440 }
441#endif
442 fNoLooperTrials = 0;
443 } else {
445#ifdef G4VERBOSE
446 if ((verboseLevel > 2)) {
447 B2DEBUG(25, " G4MonopoleTransportation::AlongStepDoIt(): Particle looping - "
448 << " Number of trials = " << fNoLooperTrials
449 << " No of calls to = " << noCalls);
450 }
451#endif
452 }
453 } else {
454 fNoLooperTrials = 0;
455 }
456
457 // Another (sometimes better way) is to use a user-limit maximum Step size
458 // to alleviate this problem ..
459
460 // Introduce smooth curved trajectories to particle-change
461 //
462 fParticleChange.SetPointerToVectorOfAuxiliaryPoints
463 (fFieldPropagator->GimmeTrajectoryVectorAndForgetIt());
464
465 //Monopole bounding
467 fParticleChange.ProposeTrackStatus(fStopAndKill);
468 B2DEBUG(150, "Monopole bound in " << fCurrentTouchableHandle->GetVolume()->GetName());
469 B2DEBUG(150, "with energy " << fTransportEndKineticEnergy << " MeV");
470 }
471
472 return &fParticleChange ;
473}
474
475// This ensures that the PostStep action is always called,
476// so that it can do the relocation if it is needed.
477//
478
481 G4double, // previousStepSize
482 G4ForceCondition* pForceCond)
483{
484 *pForceCond = Forced ;
485 return DBL_MAX ; // was kInfinity ; but convention now is DBL_MAX
486}
487
488G4VParticleChange* G4MonopoleTransportation::PostStepDoIt(const G4Track& track,
489 const G4Step&)
490{
491 G4TouchableHandle retCurrentTouchable ; // The one to return
492
493 // Initialize ParticleChange (by setting all its members equal
494 // to corresponding members in G4Track)
495 // fParticleChange.Initialize(track) ; // To initialise TouchableChange
496
497 fParticleChange.ProposeTrackStatus(track.GetTrackStatus()) ;
498
499 // If the Step was determined by the volume boundary,
500 // logically relocate the particle
501
503 // fCurrentTouchable will now become the previous touchable,
504 // and what was the previous will be freed.
505 // (Needed because the preStepPoint can point to the previous touchable)
506
507 fLinearNavigator->SetGeometricallyLimitedStep() ;
509 LocateGlobalPointAndUpdateTouchableHandle(track.GetPosition(),
510 track.GetMomentumDirection(),
512 true) ;
513 // Check whether the particle is out of the world volume
514 // If so it has exited and must be killed.
515 //
516 if (fCurrentTouchableHandle->GetVolume() == 0) {
517 fParticleChange.ProposeTrackStatus(fStopAndKill) ;
518 }
519 retCurrentTouchable = fCurrentTouchableHandle ;
520 fParticleChange.SetTouchableHandle(fCurrentTouchableHandle) ;
521 } else { // fGeometryLimitedStep is false
522 // This serves only to move the Navigator's location
523 //
524 fLinearNavigator->LocateGlobalPointWithinVolume(track.GetPosition()) ;
525
526 // The value of the track's current Touchable is retained.
527 // (and it must be correct because we must use it below to
528 // overwrite the (unset) one in particle change)
529 // It must be fCurrentTouchable too ??
530 //
531 fParticleChange.SetTouchableHandle(track.GetTouchableHandle()) ;
532 retCurrentTouchable = track.GetTouchableHandle() ;
533 } // endif ( fGeometryLimitedStep )
534
535 const G4VPhysicalVolume* pNewVol = retCurrentTouchable->GetVolume() ;
536 G4Material* pNewMaterial = nullptr ;
537 G4VSensitiveDetector* pNewSensitiveDetector = nullptr ;
538
539 if (pNewVol != 0) {
540 pNewMaterial = pNewVol->GetLogicalVolume()->GetMaterial();
541 pNewSensitiveDetector = pNewVol->GetLogicalVolume()->GetSensitiveDetector();
542 }
543
544 fParticleChange.SetMaterialInTouchable(
545 pNewMaterial) ;
546 fParticleChange.SetSensitiveDetectorInTouchable(
547 pNewSensitiveDetector) ;
548
549 const G4MaterialCutsCouple* pNewMaterialCutsCouple = 0;
550 if (pNewVol != 0) {
551 pNewMaterialCutsCouple = pNewVol->GetLogicalVolume()->GetMaterialCutsCouple();
552 }
553
554 if (pNewVol != 0 && pNewMaterialCutsCouple != 0 &&
555 pNewMaterialCutsCouple->GetMaterial() != pNewMaterial) {
556 // for parametrized volume
557 //
558 pNewMaterialCutsCouple =
559 G4ProductionCutsTable::GetProductionCutsTable()
560 ->GetMaterialCutsCouple(pNewMaterial,
561 pNewMaterialCutsCouple->GetProductionCuts());
562 }
563 fParticleChange.SetMaterialCutsCoupleInTouchable(pNewMaterialCutsCouple);
564
565 // temporarily until Get/Set Material of ParticleChange,
566 // and StepPoint can be made const.
567 // Set the touchable in ParticleChange
568 // this must always be done because the particle change always
569 // uses this value to overwrite the current touchable pointer.
570 //
571 fParticleChange.SetTouchableHandle(retCurrentTouchable) ;
572
573 return &fParticleChange ;
574}
575
576// New method takes over the responsibility to reset the state
577// of G4MonopoleTransportation object at the start of a new track
578// or the resumption of a suspended track.
579
580void
582{
583 G4VProcess::StartTracking(aTrack);
584
585// The actions here are those that were taken in AlongStepGPIL
586// when track.GetCurrentStepNumber()==1
587
588 // reset safety value and center
589 //
590 fPreviousSafety = 0.0 ;
591 fPreviousSftOrigin = G4ThreeVector(0., 0., 0.) ;
592
593 // reset looping counter -- for motion in field
594 fNoLooperTrials = 0;
595 // Must clear this state .. else it depends on last track's value
596 // --> a better solution would set this from state of suspended track TODO ?
597 // Was if( aTrack->GetCurrentStepNumber()==1 ) { .. }
598
599 // ChordFinder reset internal state
600 //
601 if (DoesGlobalFieldExist()) {
602 fFieldPropagator->ClearPropagatorState();
603 // Resets all state of field propagator class (ONLY)
604 // including safety values (in case of overlaps and to wipe for first track).
605
606 // G4ChordFinder* chordF= fFieldPropagator->GetChordFinder();
607 // if( chordF ) chordF->ResetStepEstimate();
608 }
609
610 // Make sure to clear the chord finders of all fields (ie managers)
611 static G4FieldManagerStore* fieldMgrStore = G4FieldManagerStore::GetInstance();
612 fieldMgrStore->ClearAllChordFindersState();
613
614 // Update the current touchable handle (from the track's)
615 //
616 fCurrentTouchableHandle = aTrack->GetTouchableHandle();
617}
static G4MonopoleFieldSetup * GetMonopoleFieldSetup()
Returns G4MonopoleFieldSetup instance.
G4TouchableHandle fCurrentTouchableHandle
Current touchable handle.
virtual G4VParticleChange * PostStepDoIt(const G4Track &track, const G4Step &stepData)
G4VProcess::PostStepDoIt() implementation.
G4double fCandidateEndGlobalTime
The particle's state after this Step, Store for DoIt.
virtual G4double AlongStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)
G4VProcess::AlongStepGetPhysicalInteractionLength() implementation.
G4int fThresholdTrials
Nubmer of trials for looping particles.
G4double fMaxEnergyKilled
Max of abandoned looping tracks energies.
G4bool fGeometryLimitedStep
Flag to determine whether a boundary was reached.
G4ThreeVector fTransportEndMomentumDir
The particle's state after this Step, Store for DoIt.
G4double fPreviousSafety
Remember last safety value.
virtual G4VParticleChange * AlongStepDoIt(const G4Track &track, const G4Step &stepData)
G4VProcess::AlongStepDoIt() implementation, Proposes changes during step to fParticleChange of this c...
G4SafetyHelper * fpSafetyHelper
To pass it the safety value obtained.
G4double fSumEnergyKilled
Sum of abandoned looping tracks energies.
G4ThreeVector fPreviousSftOrigin
Remember last safety origin.
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *pForceCond)
G4VProcess::PostStepGetPhysicalInteractionLength() implementation.
G4ThreeVector fTransportEndSpin
The particle's state after this Step, Store for DoIt.
G4Navigator * fLinearNavigator
Propagator used to transport the particle.
G4double fThreshold_Important_Energy
Hesitate above this about looping particle for a certain no of trials.
G4ThreeVector fTransportEndPosition
The particle's state after this Step, Store for DoIt.
G4bool fEndGlobalTimeComputed
The particle's state after this Step, Store for DoIt.
G4bool fShortStepOptimisation
Whether to avoid calling G4Navigator for short step ( < safety) If using it, the safety estimate for ...
G4int fNoLooperTrials
Counter for steps in which particle reports 'looping', if it is above 'Important' Energy.
G4bool fMomentumChanged
The particle's state after this Step, Store for DoIt.
virtual void StartTracking(G4Track *aTrack)
Reset state for new (potentially resumed) track.
G4ParticleChangeForTransport fParticleChange
New ParticleChange.
G4MonopoleTransportation(const G4Monopole *mpl, G4int verb=1)
Constructor.
G4double fThreshold_Warning_Energy
Warn above this energy about looping particle.
G4bool DoesGlobalFieldExist()
Checks whether a field exists for the "global" field manager.
const G4Monopole * fParticleDef
Monopole definition for charge and mass reference.
G4bool fParticleIsLooping
Is the monopole stuck in looping.
G4double fThreshold_Trap_Energy
Assume monopoles below this can bound to material.
G4PropagatorInField * fFieldPropagator
Propagator used to transport the particle.
G4double fTransportEndKineticEnergy
The particle's state after this Step, Store for DoIt.
G4MonopoleFieldSetup * fMagSetup
Monpole field setup.
A class to hold monopole description as a particle.
Definition G4Monopole.h:31
Abstract base class for different kinds of events.
STL namespace.