64 fTransportEndPosition(0.0, 0.0, 0.0),
65 fTransportEndMomentumDir(0.0, 0.0, 0.0),
66 fTransportEndKineticEnergy(0.0),
67 fTransportEndSpin(0.0, 0.0, 0.0),
68 fMomentumChanged(false),
69 fEndGlobalTimeComputed(false),
70 fCandidateEndGlobalTime(0.0),
71 fParticleIsLooping( false ),
74 fPreviousMassSafety( 0.0 ),
75 fPreviousFullSafety( 0.0 ),
76 fMassGeometryLimitedStep( false ),
77 fAnyGeometryLimitedStep( false ),
78 fEndpointDistance( -1.0 ),
79 fSumEnergyKilled( 0.0 ), fMaxEnergyKilled( 0.0 ),
80 fFirstStepInMassVolume( true ),
81 fFirstStepInAnyVolume( true )
96 G4cout <<
" G4CoupledTransportation constructor: ----- " <<
G4endl;
98 G4cout <<
" Navigator Id obtained in G4CoupledTransportation constructor "
100 G4cout <<
" Reports First/Last in "
102 <<
" geometry " <<
G4endl;
141 outStr <<
" G4CoupledTransportation: Statistics for looping particles "
143 outStr <<
" Sum of energy of loopers killed: "
145 outStr <<
" Max energy of loopers killed: "
150 outStr <<
" Sum of energy of loopers 'saved': "
152 outStr <<
" Sum of energy of unstable loopers 'saved': "
195#ifdef G4DEBUG_TRANSPORT
196 G4cout <<
" CoupledTransport::AlongStep GPIL: "
215#ifdef G4DEBUG_TRANSPORT
218 G4cout <<
"G4CoupledTransportation::AlongStepGPIL> called in volume "
230 startMassSafety = 0.0;
231 startFullSafety= 0.0;
237 G4double mag_shift= std::sqrt(MagSqShift);
264 G4bool fieldExertsForce = false ;
266 const G4Field* ptrField=
nullptr;
269 G4bool eligibleEM = (particleCharge != 0.0)
273 if( (fieldMgr!=
nullptr) && (eligibleEM||eligibleGrav) )
285 if( ptrField !=
nullptr)
287 fieldExertsForce = eligibleEM
293 if( fieldExertsForce )
300 if( equationOfMotion )
302 equationOfMotion->SetChargeMomentumMass( chargeState,
306#ifdef G4DEBUG_TRANSPORT
309 G4cerr <<
" ERROR in G4CoupledTransportation> "
310 <<
"Cannot find valid Equation of motion: " <<
G4endl;
311 <<
" Unable to pass Charge, Momentum and Mass " <<
G4endl;
334 if( currentMinimumStep > 0 )
359#ifdef G4DEBUG_TRANSPORT
362 std::ostringstream message;
363 message <<
" ERROR in determining geometries limiting the step" <<
G4endl;
366 message <<
"Incompatible conditions - by which geometries was it limited ?";
367 G4Exception(
"G4CoupledTransportation::AlongStepGetPhysicalInteractionLength()",
372 geometryStepLength =
std::min( lengthAlongCurve, currentMinimumStep);
385#ifdef G4DEBUG_TRANSPORT
388 G4cout <<
"G4Transport:CompStep> "
389 <<
" called the pathfinder for a new step at " << startPosition
390 <<
" and obtained step = " << lengthAlongCurve <<
G4endl;
391 G4cout <<
" New safety (preStep) = " << newMassSafety
392 <<
" versus precalculated = " << startMassSafety <<
G4endl;
397 startMassSafety = newMassSafety ;
398 startFullSafety = newFullSafety ;
406 geometryStepLength = lengthAlongCurve= 0.0 ;
415 endTrackState= aFieldTrack;
419 if( startMassSafety == 0.0 )
428 if( !fieldExertsForce )
438#ifdef G4DEBUG_TRANSPORT
441 G4cout <<
" G4CT::CS End Position = "
443 G4cout <<
" G4CT::CS End Direction = "
472 G4double absEdiff = std::fabs(startEnergy- endEnergy);
496 safetyProposal= startFullSafety;
502 && ( particleCharge != 0.0 ) )
530#ifdef G4DEBUG_TRANSPORT
532 G4cout <<
"***CoupledTransportation::AlongStepGPIL ** " <<
G4endl ;
534 <<
" give safety values: Mass= " << endMassSafety
535 <<
" All= " << endFullSafety <<
G4endl ;
537 <<
" to obtain pseudo-safety= " << safetyProposal <<
G4endl ;
543 G4cout <<
"***CoupledTransportation::AlongStepGPIL ** " <<
G4endl ;
544 G4cout <<
" Quick Safety estimate at endpoint "
546 <<
" gives safety endpoint value = "
548 <<
" using start-point value " << startFullSafety
554 proposedSafetyForStart= safetyProposal;
557 return geometryStepLength ;
567 const char *methodName=
"AlongStepDoIt";
598 if( finalVelocity > 0.0 ) { finalInverseVel= 1.0 / finalVelocity; }
600 if( initialVelocity > 0.0 ) { initialInverseVel= 1.0 / initialVelocity; }
603 if (finalVelocity > 0.0)
607 * (initialInverseVel + finalInverseVel);
608 deltaTime = stepLength * meanInverseVelocity ;
612 deltaTime = stepLength * initialInverseVel ;
646 if( candidateForEnd && stable )
648 const G4int electronPDG= 11;
681 noCallsCT_ASDI, methodName );
699 G4cout <<
" ** G4CoupledTransportation::AlongStepDoIt():"
700 <<
" Particle is looping but is saved ..." <<
G4endl
705 <<
" Total no of calls to this method (all tracks) = "
706 << noCallsCT_ASDI <<
G4endl;
752 <<
"**************************************************************"
754 G4cerr <<
"Endpoint has moved between value expected from TransportEndPosition "
755 <<
" and value from Track in PostStepDoIt. " <<
G4endl
756 <<
"Change of " << Quantity <<
" is " << moveVec.
mag() /
mm
758 <<
" and its vector is " << (1.0/
mm) * moveVec <<
" mm " <<
G4endl
759 <<
"Endpoint of ComputeStep was " << OldVector
760 <<
" and current position to locate is " << NewVector <<
G4endl;
788#ifdef G4DEBUG_TRANSPORT
793 "End of Step Position" );
794 G4cerr <<
" Problem in G4CoupledTransportation::PostStepDoIt " <<
G4endl;
802 G4cout <<
" Calling PathFinder::Locate() from "
803 <<
" G4CoupledTransportation::PostStepDoIt " <<
G4endl;
823#ifdef G4DEBUG_TRANSPORT
826 G4cout <<
"G4CoupledTransportation::PostStepDoIt --- fNavigatorId = "
832 G4cout <<
"CHECK !!!!!!!!!!! fCurrentTouchableHandle->GetVolume() = "
851#ifdef G4DEBUG_TRANSPORT
854 G4cout <<
"G4CoupledTransportation::PostStepDoIt -- "
856 <<
" must be false " <<
G4endl;
875#ifdef G4DEBUG_NAVIGATION
876 G4cout <<
" CoupledTransport::AlongStep GPIL: "
910 if( pNewMaterialCutsCouple!=0
911 && pNewMaterialCutsCouple->
GetMaterial()!=pNewMaterial )
991#ifdef G4DEBUG_TRANSPORT
1021 moduloFactor= 10, no_large_ediff= 0;
1023 if( std::fabs(startEnergy- endEnergy) >
perThousand * endEnergy )
1026 if( (no_large_ediff% warnModulo) == 0 )
1029 std::ostringstream message;
1030 message <<
"Energy change in Step is above 1^-3 relative value. "
1032 <<
" Relative change in 'tracking' step = "
1033 << std::setw(15) << (endEnergy-startEnergy)/startEnergy
1035 <<
" Starting E= " << std::setw(12) << startEnergy /
MeV
1037 <<
" Ending E= " << std::setw(12) << endEnergy /
MeV
1039 <<
"Energy has been corrected -- however, review"
1040 <<
" field propagation parameters for accuracy." <<
G4endl;
1042 || (no_large_ediff == warnModulo * moduloFactor) )
1044 message <<
"These include EpsilonStepMax(/Min) in G4FieldManager,"
1046 <<
"which determine fractional error per step for integrated quantities."
1048 <<
"Note also the influence of the permitted number of integration steps."
1051 message <<
"Bad 'endpoint'. Energy change detected and corrected."
1053 <<
"Has occurred already " << no_large_ediff <<
" times.";
1054 G4Exception(
"G4CoupledTransportation::AlongStepGetPIL()",
1056 if( no_large_ediff == warnModulo * moduloFactor )
1058 warnModulo *= moduloFactor;
1111 G4int maxTrials = 10;
1124 G4int maxTrials = 30;
1135 const char* message=
"Logger object missing from G4CoupledTransportation";
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
#define fPreviousSftOrigin
static constexpr double perMillion
static constexpr double mm
static constexpr double MeV
static constexpr double perThousand
CLHEP::Hep3Vector G4ThreeVector
G4ReferenceCountedHandle< G4VTouchable > G4TouchableHandle
G4GLOB_DLL std::ostream G4cerr
G4GLOB_DLL std::ostream G4cout
unsigned long fNumLoopersKilled_NonElectron
G4double fMaxEnergyKilled_NonElectron
G4double fSumEnergyKilled
G4ThreeVector fTransportEndSpin
G4double fSumEnergyKilled_NonElectron
void SetLowLooperThresholds()
G4double fPreviousFullSafety
void PushThresholdsToLogger()
G4Navigator * fMassNavigator
static G4bool GetSilenceLooperWarnings()
static G4bool fSilenceLooperWarnings
void ReportLooperThresholds()
void SetHighLooperThresholds()
G4VParticleChange * AlongStepDoIt(const G4Track &track, const G4Step &stepData)
G4VParticleChange * PostStepDoIt(const G4Track &track, const G4Step &stepData)
G4bool fFirstStepInAnyVolume
G4double fSumEnerSqKilled
G4double fSumEnerSqKilled_NonElectron
void SetThresholdTrials(G4int newMaxTrials)
G4TransportationLogger * fpLogger
static G4bool EnableGravity(G4bool useGravity)
static G4bool fSignifyStepInAnyVolume
G4bool fAnyGeometryLimitedStep
G4double fTransportEndKineticEnergy
G4double fMaxEnergyKilled
G4bool fFirstStepInMassVolume
void ReportMove(G4ThreeVector OldVector, G4ThreeVector NewVector, const G4String &Quantity)
G4double fCandidateEndGlobalTime
static void SetSilenceLooperWarnings(G4bool val)
void ReportInexactEnergy(G4double startEnergy, G4double endEnergy)
static G4bool EnableMagneticMoment(G4bool useMoment=true)
G4double fSumEnergyUnstableSaved
G4bool fParticleIsLooping
G4SafetyHelper * fpSafetyHelper
void PrintStatistics(std::ostream &outStr) const
G4bool fMassGeometryLimitedStep
G4int fMaxEnergyKilledPDG
G4TouchableHandle fCurrentTouchableHandle
G4double fEndpointDistance
G4int fMaxEnergyKilled_NonElecPDG
G4ThreeVector fTransportEndMomentumDir
unsigned long fNumLoopersKilled
G4ThreeVector fTransportEndPosition
G4CoupledTransportation(G4int verbosityLevel=0)
G4bool DoesAnyFieldExist()
G4ThreeVector fPreviousSftOrigin
G4double PostStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4ForceCondition *pForceCond)
G4double fPreviousMassSafety
static G4bool fUseGravity
static G4bool fUseMagneticMoment
G4bool fEndGlobalTimeComputed
G4PropagatorInField * fFieldPropagator
void SetThresholdImportantEnergy(G4double newEnImp)
G4double fThreshold_Important_Energy
void SetThresholdWarningEnergy(G4double newEnWarn)
void StartTracking(G4Track *aTrack)
G4double AlongStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double ¤tSafety, G4GPILSelection *selection)
~G4CoupledTransportation()
G4double fThreshold_Warning_Energy
void ReportMissingLogger(const char *methodName)
G4PathFinder * fPathFinder
G4ParticleChangeForTransport fParticleChange
const G4ThreeVector & GetMomentumDirection() const
G4double GetCharge() const
G4ParticleDefinition * GetDefinition() const
G4double GetMagneticMoment() const
G4double GetTotalMomentum() const
static G4FieldManagerStore * GetInstance()
void ClearAllChordFindersState()
G4bool DoesFieldChangeEnergy() const
virtual void ConfigureForTrack(const G4Track *)
const G4Field * GetDetectorField() const
const G4ThreeVector & GetMomentumDir() const
G4double GetKineticEnergy() const
G4ThreeVector GetPosition() const
G4ThreeVector GetSpin() const
G4double GetLabTimeOfFlight() const
G4bool IsGravityActive() const
G4VSensitiveDetector * GetSensitiveDetector() const
G4Material * GetMaterial() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
const G4Material * GetMaterial() const
G4ProductionCuts * GetProductionCuts() const
void SetPointerToVectorOfAuxiliaryPoints(std::vector< G4ThreeVector > *vec)
void SetMaterialInTouchable(G4Material *fMaterial)
void SetTouchableHandle(const G4TouchableHandle &fTouchable)
void SetMaterialCutsCoupleInTouchable(const G4MaterialCutsCouple *fMaterialCutsCouple)
virtual void Initialize(const G4Track &)
void SetSensitiveDetectorInTouchable(G4VSensitiveDetector *fSensitiveDetector)
void SetMomentumChanged(G4bool b)
void ProposePolarization(G4double Px, G4double Py, G4double Pz)
void ProposePosition(G4double x, G4double y, G4double z)
void ProposeLocalTime(G4double t)
void ProposeProperTime(G4double finalProperTime)
void ProposeEnergy(G4double finalEnergy)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
void ProposeGlobalTime(G4double t)
G4double GetPDGMagneticMoment() const
G4bool GetPDGStable() const
G4int GetPDGEncoding() const
G4double GetPDGSpin() const
G4double ComputeSafety(const G4ThreeVector &globalPoint)
void ReLocate(const G4ThreeVector &position)
unsigned int GetNumberGeometriesLimitingStep() const
G4double ObtainSafety(G4int navId, G4ThreeVector &globalCenterPoint)
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
void Locate(const G4ThreeVector &position, const G4ThreeVector &direction, G4bool relativeSearch=true)
G4double GetCurrentSafety() const
void PrepareNewTrack(const G4ThreeVector &position, const G4ThreeVector &direction, G4VPhysicalVolume *massStartVol=nullptr)
static G4PathFinder * GetInstance()
G4TouchableHandle CreateTouchableHandle(G4int navId) const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
static G4ProductionCutsTable * GetProductionCutsTable()
G4FieldManager * FindAndSetFieldManager(G4VPhysicalVolume *pCurrentPhysVol)
G4FieldManager * GetCurrentFieldManager()
G4bool IsParticleLooping() const
G4ChordFinder * GetChordFinder()
G4EquationOfMotion * GetCurrentEquationOfMotion()
void ClearPropagatorState()
std::vector< G4ThreeVector > * GimmeTrajectoryVectorAndForgetIt() const
void SetCurrentSafety(G4double val, const G4ThreeVector &pos)
G4double GetVelocity() const
G4StepPoint * GetPreStepPoint() const
G4TrackStatus GetTrackStatus() const
G4double GetVelocity() const
G4VPhysicalVolume * GetVolume() const
const G4ThreeVector & GetPosition() const
G4double GetGlobalTime() const
G4double GetProperTime() const
G4int GetCurrentStepNumber() const
G4double GetLocalTime() const
const G4DynamicParticle * GetDynamicParticle() const
const G4TouchableHandle & GetTouchableHandle() const
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
const G4ThreeVector & GetPolarization() const
G4double GetStepLength() const
G4double GetTotalEnergy() const
void ReportLoopingTrack(const G4Track &track, const G4Step &stepInfo, G4int numTrials, G4long noCalls, const char *methodName) const
void ReportLooperThresholds(const char *className)
static G4TransportationManager * GetTransportationManager()
G4PropagatorInField * GetPropagatorInField() const
G4SafetyHelper * GetSafetyHelper() const
G4Navigator * GetNavigatorForTracking() const
G4int ActivateNavigator(G4Navigator *aNavigator)
static G4bool fUseGravity
static G4bool fUseMagneticMoment
void ProposeTrackStatus(G4TrackStatus status)
void ProposeLastStepInVolume(G4bool flag)
void ProposeFirstStepInVolume(G4bool flag)
void ProposeTrueStepLength(G4double truePathLength)
G4LogicalVolume * GetLogicalVolume() const
const G4String & GetName() const
void SetVerboseLevel(G4int value)
void SetProcessSubType(G4int)
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
static constexpr double keV
static constexpr double MeV
T max(const T t1, const T t2)
brief Return the largest of the two arguments
T min(const T t1, const T t2)
brief Return the smallest of the two arguments