Belle II Software development
G4LongLivedNeutralTransportation.h
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 G4Transportation class
10
11#pragma once
12
13#include "G4VProcess.hh"
14#include "G4TransportationManager.hh"
15#include "G4PropagatorInField.hh"
16#include "G4ParticleChangeForTransport.hh"
17
18class G4SafetyHelper;
19class G4CoupledTransportation;
20class G4FieldManager;
21class G4Navigator;
22class G4Track;
23class G4Step;
24
25namespace Belle2 {
30
34 class G4LongLivedNeutralTransportation : public G4VProcess {
35
36 public:
37
42 G4LongLivedNeutralTransportation(G4int verbosityLevel = 1);
43
48
60 const G4Track& track,
61 G4double previousStepSize,
62 G4double currentMinimumStep,
63 G4double& currentSafety,
64 G4GPILSelection* selection
65 );
66
73 G4VParticleChange* AlongStepDoIt(
74 const G4Track& track,
75 const G4Step& stepData
76 );
77
84 G4VParticleChange* PostStepDoIt(
85 const G4Track& track,
86 const G4Step& stepData
87 );
88
100 const G4Track& track,
101 G4double previousStepSize,
102 G4ForceCondition* pForceCond
103 );
104
105
106
107 inline G4bool FieldExertedForce() { return fFieldExertedForce; }
109
110 G4PropagatorInField* GetPropagatorInField();
112 void SetPropagatorInField(G4PropagatorInField* pFieldPropagator);
114
115 inline G4double GetThresholdWarningEnergy() const;
116 inline G4double GetThresholdImportantEnergy() const;
117 inline G4int GetThresholdTrials() const;
118
119
120
121 inline void SetThresholdWarningEnergy(G4double newEnWarn)
122 {
123 fThreshold_Warning_Energy = newEnWarn;
124 }
125 inline void SetThresholdImportantEnergy(G4double newEnImp)
126 {
128 }
129 inline void SetThresholdTrials(G4int newMaxTrials)
130 {
131 fThresholdTrials = newMaxTrials;
132 }
133
134
140
147
148 inline G4double GetMaxEnergyKilled() const;
150 inline G4double GetSumEnergyKilled() const;
152 inline void ResetKilledStatistics(G4int report = 1);
154
155 inline void EnableShortStepOptimisation(G4bool optimise = true);
157
158 static void SetSilenceLooperWarnings(G4bool val);
160 static G4bool GetSilenceLooperWarnings();
162
163 public: // without description
164
168 G4double AtRestGetPhysicalInteractionLength(const G4Track&,
169 G4ForceCondition*)
170 { return -1.0; }
171
175 G4VParticleChange* AtRestDoIt(const G4Track&, const G4Step&)
176 { return 0; }
177
178 void StartTracking(G4Track* aTrack);
180
186 virtual void ProcessDescription(std::ostream& outStream) const;
187
188 void PrintStatistics(std::ostream& outStr) const;
190
191 protected:
192
196 inline G4bool DoesGlobalFieldExist()
197 {
198 G4TransportationManager* transportMgr;
199 transportMgr = G4TransportationManager::GetTransportationManager();
200
201 // fFieldExists= transportMgr->GetFieldManager()->DoesFieldExist();
202 // return fFieldExists;
203 return transportMgr->GetFieldManager()->DoesFieldExist();
204 }
205
206 private:
207
208 G4Navigator* fLinearNavigator;
209 G4PropagatorInField* fFieldPropagator;
210
211 G4ThreeVector fTransportEndPosition = G4ThreeVector(0.0, 0.0, 0.0);
212 G4ThreeVector fTransportEndMomentumDir = G4ThreeVector(0.0, 0.0, 0.0);
214 G4ThreeVector fTransportEndSpin = G4ThreeVector(0.0, 0.0, 0.0);
215 G4bool fMomentumChanged = true;
216 G4bool fEndGlobalTimeComputed = false;
217 G4double fCandidateEndGlobalTime = 0.0;
218
219 G4bool fAnyFieldExists = false;
220
221 G4bool fNewTrack = true;
222 G4bool fFirstStepInVolume = true;
223 G4bool fLastStepInVolume = false;
224 G4bool fGeometryLimitedStep = true;
225
226 G4bool fFieldExertedForce = false;
227
228 G4TouchableHandle fCurrentTouchableHandle;
229
230 G4ThreeVector fPreviousSftOrigin;
232
233 G4ParticleChangeForTransport fParticleChange;
234
236
237 // Thresholds for looping particles:
238 //
239 G4double fThreshold_Warning_Energy = 1.0 * CLHEP::keV;
240 G4double fThreshold_Important_Energy = 1.0 * CLHEP::MeV;
241 G4int fThresholdTrials = 10;
242 // Above 'important' energy a 'looping' particle in field will
243 // *NOT* be abandoned, except after fThresholdTrials attempts.
245 // unstable loopers ( 0 = never )
246 // Counter for steps in which particle reports 'looping',
247 // ( Used if it is above 'Important' Energy. )
248 G4int fNoLooperTrials = 0;
249
250 // Statistics for tracks abandoned due to looping, only used as flag, no looping neutral particles
251 //
252 G4double fSumEnergyKilled = 0.0;
258
259 G4SafetyHelper* fpSafetyHelper;
260
261 private:
262
263 friend class G4CoupledTransportation;
264 static G4bool fUseMagneticMoment;
265 static G4bool fUseGravity;
267
268 };
270}
G4TouchableHandle fCurrentTouchableHandle
Current touchable handle.
G4int GetThresholdTrials() const
Access fThresholdTrials.
G4double GetThresholdImportantEnergy() const
Access fThreshold_Important_Energy.
void PrintStatistics(std::ostream &outStr) const
returns current logging info of the algorithm
G4VParticleChange * PostStepDoIt(const G4Track &track, const G4Step &stepData)
G4VProcess::PostStepDoIt() implementation,.
G4double AtRestGetPhysicalInteractionLength(const G4Track &, G4ForceCondition *)
No operation in AtRestGPIL.
G4double fCandidateEndGlobalTime
The particle's state after this Step, Store for DoIt.
G4double AlongStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)
G4VProcess::AlongStepGetPhysicalInteractionLength() implementation,.
G4int fThresholdTrials
Number of trials an important looper survives.
G4PropagatorInField * GetPropagatorInField()
Access fFieldPropagator, the assistant class that Propagate in a Field.
G4double GetMaxEnergyKilled() const
Print values of looper thresholds.
void SetThresholdWarningEnergy(G4double newEnWarn)
Set fThreshold_Warning_Energy.
void SetThresholdImportantEnergy(G4double newEnImp)
Set fThreshold_Important_Energy.
virtual void ProcessDescription(std::ostream &outStream) const
G4LongLivedNeutralTransportation::ProcessDescription()
G4bool fGeometryLimitedStep
Flag to determine whether a boundary was reached.
G4ThreeVector fTransportEndMomentumDir
The particle's state after this Step, Store for DoIt.
void SetThresholdTrials(G4int newMaxTrials)
Set fThresholdTrials.
G4VParticleChange * AlongStepDoIt(const G4Track &track, const G4Step &stepData)
G4VProcess::AlongStepDoIt() implementation,.
G4int fAbandonUnstableTrials
Number of trials after which to abandon.
G4SafetyHelper * fpSafetyHelper
To pass it the safety value obtained.
void SetPropagatorInField(G4PropagatorInField *pFieldPropagator)
Access/set the assistant class that Propagate in a Field.
void SetHighLooperThresholds()
Get/Set parameters for killing loopers: Above 'important' energy a 'looping' particle in field will N...
G4double fSumEnergyKilled
Sum of abandoned looping tracks energies.
G4ThreeVector fPreviousSftOrigin
Remember last safety origin.
static G4bool fUseGravity
Flag take into account gravity.
void ReportLooperThresholds()
Set low thresholds - for low-E applications.
static void SetSilenceLooperWarnings(G4bool val)
Do not warn about 'looping' particles.
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.
static G4bool fSilenceLooperWarnings
Flag to Supress all 'looper' warnings.
G4double fThreshold_Important_Energy
Give a few trial above this E for looping particle.
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.
G4LongLivedNeutralTransportation(G4int verbosityLevel=1)
Constructor.
void StartTracking(G4Track *aTrack)
Reset state for new (potentially resumed) track.
static G4bool fUseMagneticMoment
Flag take into account magnetic moment.
G4ParticleChangeForTransport fParticleChange
New ParticleChange.
void SetLowLooperThresholds()
Shortcut method - old values (meant for HEP)
G4double fThreshold_Warning_Energy
Warn above this energy about looping particle.
void EnableShortStepOptimisation(G4bool optimise=true)
Whether short steps < safety will avoid to call Navigator (if field=0)
G4bool FieldExertedForce()
References fFieldExertedForce.
G4bool DoesGlobalFieldExist()
Checks whether a field exists for the "global" field manager.
G4VParticleChange * AtRestDoIt(const G4Track &, const G4Step &)
No operation in AtRestDoIt.
static G4bool GetSilenceLooperWarnings()
Do not throw exception about 'looping' particles.
G4double GetSumEnergyKilled() const
Access fSumEnergyKilled.
void ResetKilledStatistics(G4int report=1)
Statistics for tracks killed (currently due to looping in field)
G4double GetThresholdWarningEnergy() const
Access fThreshold_Warning_Energy.
G4PropagatorInField * fFieldPropagator
Propagator used to transport the particle.
G4double fTransportEndKineticEnergy
The particle's state after this Step, Store for DoIt.
Abstract base class for different kinds of events.