11#include <simulation/monopoles/G4MonopoleTransportation.h>
12#include <simulation/monopoles/G4Monopole.h>
14#include <G4Navigator.hh>
15#include <G4PropagatorInField.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>
27class G4VSensitiveDetector;
31using namespace Belle2::Monopoles;
36 : G4VProcess(G4String(
"MonopoleTransportation"), fTransportation),
56 SetProcessSubType(TRANSPORTATION);
60 G4TransportationManager* transportMgr ;
62 transportMgr = G4TransportationManager::GetTransportationManager() ;
78 static G4TouchableHandle nullTouchableHandle;
88 B2DEBUG(25,
" G4MonopoleTransportation: Statistics for looping particles ");
102 G4double currentMinimumStep,
103 G4double& currentSafety,
104 G4GPILSelection* selection)
109 G4double geometryStepLength, newSafety ;
121 *selection = CandidateForSelection ;
125 const G4DynamicParticle* pParticle = track.GetDynamicParticle() ;
126 G4ThreeVector startMomentumDir = pParticle->GetMomentumDirection() ;
127 G4ThreeVector startPosition = track.GetPosition() ;
136 G4double MagSqShift = OriginShift.mag2() ;
138 currentSafety = 0.0 ;
145 G4double particleMagneticCharge =
fParticleDef->MagneticCharge() ;
146 G4double particleElectricCharge = pParticle->GetCharge();
157 G4FieldManager* fieldMgr =
nullptr;
158 G4bool fieldExertsForce = false ;
160 if ((particleMagneticCharge != 0.0)) {
164 fieldMgr->ConfigureForTrack(&track);
170 fieldExertsForce = (fieldMgr->GetDetectorField() != 0);
176 if (!fieldExertsForce) {
177 G4double linearStepLength ;
181 geometryStepLength = currentMinimumStep ;
198 currentSafety = newSafety ;
203 geometryStepLength = linearStepLength ;
206 geometryStepLength = currentMinimumStep ;
224 G4double momentumMagnitude = pParticle->GetTotalMomentum() ;
225 G4ThreeVector EndUnitMomentum ;
228 G4ChargeState chargeState(particleElectricCharge,
232 particleMagneticCharge);
234 G4EquationOfMotion* equationOfMotion =
236 ->GetEquationOfMotion();
239 ->SetChargeMomentumMass(chargeState,
244 G4ThreeVector spin = track.GetPolarization() ;
245 G4FieldTrack aFieldTrack = G4FieldTrack(startPosition,
246 track.GetMomentumDirection(),
248 track.GetKineticEnergy(),
251 track.GetGlobalTime(),
252 track.GetProperTime(),
254 if (currentMinimumStep > 0) {
264 geometryStepLength = lengthAlongCurve ;
266 geometryStepLength = currentMinimumStep ;
269 geometryStepLength = 0.0 ;
301 if (currentMinimumStep == 0.0) {
312 if (particleMagneticCharge != 0.0) {
316 currentSafety = endSafety ;
326#ifdef G4DEBUG_TRANSPORT
327 G4cout.precision(12) ;
328 G4cout <<
"***G4MonopoleTransportation::AlongStepGPIL ** " << G4endl ;
330 <<
" and it returned safety= " << endSafety << G4endl ;
332 <<
" to obtain pseudo-safety= " << currentSafety << G4endl ;
342 return geometryStepLength ;
350 const G4Step& stepData)
353 static G4int noCalls = 0;
357 static const G4ParticleDefinition* fOpticalPhoton =
358 G4ParticleTable::GetParticleTable()->FindParticle(
"opticalphoton");
371 G4double deltaTime = 0.0 ;
377 G4double startTime = track.GetGlobalTime() ;
382 G4double finalVelocity = track.GetVelocity() ;
383 G4double initialVelocity = stepData.GetPreStepPoint()->GetVelocity() ;
384 G4double stepLength = track.GetStepLength() ;
387 const G4DynamicParticle* fpDynamicParticle = track.GetDynamicParticle();
388 if (fpDynamicParticle->GetDefinition() == fOpticalPhoton) {
391 deltaTime = stepLength / finalVelocity ;
392 }
else if (finalVelocity > 0.0) {
393 G4double meanInverseVelocity ;
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 ;
410 G4double restMass = track.GetDynamicParticle()->GetMass() ;
411 G4double deltaProperTime = deltaTime * (restMass / track.GetTotalEnergy()) ;
413 fParticleChange.ProposeProperTime(track.GetProperTime() + deltaProperTime) ;
434 if ((verboseLevel > 1) ||
436 B2DEBUG(25,
" G4MonopoleTransportation is killing track that is looping or stuck ");
437 B2DEBUG(25,
" This track has " << track.GetKineticEnergy() / MeV <<
" MeV energy.");
439 <<
" No of calls to AlongStepDoIt = " << noCalls);
446 if ((verboseLevel > 2)) {
447 B2DEBUG(25,
" G4MonopoleTransportation::AlongStepDoIt(): Particle looping - "
449 <<
" No of calls to = " << noCalls);
482 G4ForceCondition* pForceCond)
484 *pForceCond = Forced ;
491 G4TouchableHandle retCurrentTouchable ;
509 LocateGlobalPointAndUpdateTouchableHandle(track.GetPosition(),
510 track.GetMomentumDirection(),
532 retCurrentTouchable = track.GetTouchableHandle() ;
535 const G4VPhysicalVolume* pNewVol = retCurrentTouchable->GetVolume() ;
536 G4Material* pNewMaterial = nullptr ;
537 G4VSensitiveDetector* pNewSensitiveDetector = nullptr ;
540 pNewMaterial = pNewVol->GetLogicalVolume()->GetMaterial();
541 pNewSensitiveDetector = pNewVol->GetLogicalVolume()->GetSensitiveDetector();
547 pNewSensitiveDetector) ;
549 const G4MaterialCutsCouple* pNewMaterialCutsCouple = 0;
551 pNewMaterialCutsCouple = pNewVol->GetLogicalVolume()->GetMaterialCutsCouple();
554 if (pNewVol != 0 && pNewMaterialCutsCouple != 0 &&
555 pNewMaterialCutsCouple->GetMaterial() != pNewMaterial) {
558 pNewMaterialCutsCouple =
559 G4ProductionCutsTable::GetProductionCutsTable()
560 ->GetMaterialCutsCouple(pNewMaterial,
561 pNewMaterialCutsCouple->GetProductionCuts());
563 fParticleChange.SetMaterialCutsCoupleInTouchable(pNewMaterialCutsCouple);
583 G4VProcess::StartTracking(aTrack);
611 static G4FieldManagerStore* fieldMgrStore = G4FieldManagerStore::GetInstance();
612 fieldMgrStore->ClearAllChordFindersState();
static G4MonopoleFieldSetup * GetMonopoleFieldSetup()
Returns G4MonopoleFieldSetup instance.
G4TouchableHandle fCurrentTouchableHandle
Current touchable handle.
G4double endpointDistance
Endpint distance.
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 ¤tSafety, 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.
~G4MonopoleTransportation()
Destructor.
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.
Abstract base class for different kinds of events.