Belle II Software development
G4LongLivedNeutralTransportation Class Reference

Concrete class that does the geometrical transport. More...

#include <G4LongLivedNeutralTransportation.h>

Inheritance diagram for G4LongLivedNeutralTransportation:

Public Member Functions

 G4LongLivedNeutralTransportation (G4int verbosityLevel=1)
 Constructor.
 
 ~G4LongLivedNeutralTransportation ()
 Destructor.
 
G4double AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)
 G4VProcess::AlongStepGetPhysicalInteractionLength() implementation,.
 
G4VParticleChange * AlongStepDoIt (const G4Track &track, const G4Step &stepData)
 G4VProcess::AlongStepDoIt() implementation,.
 
G4VParticleChange * PostStepDoIt (const G4Track &track, const G4Step &stepData)
 G4VProcess::PostStepDoIt() implementation,.
 
G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *pForceCond)
 G4VProcess::PostStepGetPhysicalInteractionLength() implementation.
 
G4bool FieldExertedForce ()
 References fFieldExertedForce.
 
G4PropagatorInField * GetPropagatorInField ()
 Access fFieldPropagator, the assistant class that Propagate in a Field.
 
void SetPropagatorInField (G4PropagatorInField *pFieldPropagator)
 Access/set the assistant class that Propagate in a Field.
 
G4double GetThresholdWarningEnergy () const
 Access fThreshold_Warning_Energy.
 
G4double GetThresholdImportantEnergy () const
 Access fThreshold_Important_Energy.
 
G4int GetThresholdTrials () const
 Access fThresholdTrials.
 
void SetThresholdWarningEnergy (G4double newEnWarn)
 Set fThreshold_Warning_Energy.
 
void SetThresholdImportantEnergy (G4double newEnImp)
 Set fThreshold_Important_Energy.
 
void SetThresholdTrials (G4int newMaxTrials)
 Set fThresholdTrials.
 
void SetHighLooperThresholds ()
 Get/Set parameters for killing loopers: Above 'important' energy a 'looping' particle in field will NOT be abandoned, except after fThresholdTrials attempts.
 
void SetLowLooperThresholds ()
 Shortcut method - old values (meant for HEP)
 
void ReportLooperThresholds ()
 Set low thresholds - for low-E applications.
 
G4double GetMaxEnergyKilled () const
 Print values of looper thresholds.
 
G4double GetSumEnergyKilled () const
 Access fSumEnergyKilled.
 
void ResetKilledStatistics (G4int report=1)
 Statistics for tracks killed (currently due to looping in field)
 
void EnableShortStepOptimisation (G4bool optimise=true)
 Whether short steps < safety will avoid to call Navigator (if field=0)
 
G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 No operation in AtRestGPIL.
 
G4VParticleChange * AtRestDoIt (const G4Track &, const G4Step &)
 No operation in AtRestDoIt.
 
void StartTracking (G4Track *aTrack)
 Reset state for new (potentially resumed) track.
 
virtual void ProcessDescription (std::ostream &outStream) const
 G4LongLivedNeutralTransportation::ProcessDescription()
 
void PrintStatistics (std::ostream &outStr) const
 returns current logging info of the algorithm
 

Static Public Member Functions

static void SetSilenceLooperWarnings (G4bool val)
 Do not warn about 'looping' particles.
 
static G4bool GetSilenceLooperWarnings ()
 Do not throw exception about 'looping' particles.
 

Protected Member Functions

G4bool DoesGlobalFieldExist ()
 Checks whether a field exists for the "global" field manager.
 

Private Attributes

G4Navigator * fLinearNavigator
 Propagator used to transport the particle.
 
G4PropagatorInField * fFieldPropagator
 Propagator used to transport the particle.
 
G4ThreeVector fTransportEndPosition = G4ThreeVector(0.0, 0.0, 0.0)
 The particle's state after this Step, Store for DoIt.
 
G4ThreeVector fTransportEndMomentumDir = G4ThreeVector(0.0, 0.0, 0.0)
 The particle's state after this Step, Store for DoIt.
 
G4double fTransportEndKineticEnergy = 0.0
 The particle's state after this Step, Store for DoIt.
 
G4ThreeVector fTransportEndSpin = G4ThreeVector(0.0, 0.0, 0.0)
 The particle's state after this Step, Store for DoIt.
 
G4bool fMomentumChanged = true
 The particle's state after this Step, Store for DoIt.
 
G4bool fEndGlobalTimeComputed = false
 The particle's state after this Step, Store for DoIt.
 
G4double fCandidateEndGlobalTime = 0.0
 The particle's state after this Step, Store for DoIt.
 
G4bool fAnyFieldExists = false
 Flag for existing fields.
 
G4bool fNewTrack = true
 Flag from StartTracking.
 
G4bool fFirstStepInVolume = true
 Flag first step in a geom.
 
G4bool fLastStepInVolume = false
 Flag last step in a geom.
 
G4bool fGeometryLimitedStep = true
 Flag to determine whether a boundary was reached.
 
G4bool fFieldExertedForce = false
 During current step.
 
G4TouchableHandle fCurrentTouchableHandle
 Current touchable handle.
 
G4ThreeVector fPreviousSftOrigin
 Remember last safety origin.
 
G4double fPreviousSafety
 Remember last safety value.
 
G4ParticleChangeForTransport fParticleChange
 New ParticleChange.
 
G4double fEndPointDistance
 Endpoint distance.
 
G4double fThreshold_Warning_Energy = 1.0 * CLHEP::keV
 Warn above this energy about looping particle.
 
G4double fThreshold_Important_Energy = 1.0 * CLHEP::MeV
 Give a few trial above this E for looping particle.
 
G4int fThresholdTrials = 10
 Number of trials an important looper survives.
 
G4int fAbandonUnstableTrials = 0
 Number of trials after which to abandon.
 
G4int fNoLooperTrials = 0
 Counter for steps in which particle reports 'looping', if it is above 'Important' Energy.
 
G4double fSumEnergyKilled = 0.0
 Sum of abandoned looping tracks energies.
 
G4bool fShortStepOptimisation
 Whether to avoid calling G4Navigator for short step ( < safety) If using it, the safety estimate for endpoint will likely be smaller.
 
G4SafetyHelper * fpSafetyHelper
 To pass it the safety value obtained.
 

Static Private Attributes

static G4bool fUseMagneticMoment
 Flag take into account magnetic moment.
 
static G4bool fUseGravity
 Flag take into account gravity.
 
static G4bool fSilenceLooperWarnings = false
 Flag to Supress all 'looper' warnings.
 

Friends

class G4CoupledTransportation
 

Detailed Description

Concrete class that does the geometrical transport.

Definition at line 34 of file G4LongLivedNeutralTransportation.h.

Constructor & Destructor Documentation

◆ G4LongLivedNeutralTransportation()

G4LongLivedNeutralTransportation ( G4int verbosityLevel = 1)

Constructor.

Parameters
verbosityLevel

Definition at line 37 of file G4LongLivedNeutralTransportation.cc.

38 : G4VProcess(G4String("LongLivedNeutralTransportation"), fTransportation),
39 fFieldExertedForce(false),
40 fPreviousSftOrigin(0., 0., 0.),
41 fPreviousSafety(0.0),
43 fShortStepOptimisation(false) // Old default: true (=fast short steps)
44{
45 SetProcessSubType(static_cast<G4int>(TRANSPORTATION));
46 pParticleChange = &fParticleChange; // Required to conform to G4VProcess
47 SetVerboseLevel(verbosity);
48
49 G4TransportationManager* transportMgr ;
50
51 transportMgr = G4TransportationManager::GetTransportationManager() ;
52
53 fLinearNavigator = transportMgr->GetNavigatorForTracking() ;
54
55 fFieldPropagator = transportMgr->GetPropagatorInField() ;
56
57
58
59
61 // Use the old defaults: Warning = 100 MeV, Important = 250 MeV, No Trials = 10;
62
63 // Cannot determine whether a field exists here, as it would
64 // depend on the relative order of creating the detector's
65 // field and this process. That order is not guaranteed.
67 // This value must be updated using DoesAnyFieldExist() at least at the
68 // start of each Run -- for now this is at the Start of every Track. TODO
69
70 static G4ThreadLocal G4TouchableHandle* pNullTouchableHandle = 0;
71 pNullTouchableHandle = new G4TouchableHandle;
72 fCurrentTouchableHandle = *pNullTouchableHandle;
73 // Points to (G4VTouchable*) 0
74
75
76#ifdef G4VERBOSE
77 if (verboseLevel > 0) {
78 G4cout << " G4LongLivedNeutralTransportation constructor> set fShortStepOptimisation to ";
79 if (fShortStepOptimisation) { G4cout << "true" << G4endl; }
80 else { G4cout << "false" << G4endl; }
81 }
82#endif
83}
G4TouchableHandle fCurrentTouchableHandle
Current touchable handle.
void SetHighLooperThresholds()
Get/Set parameters for killing loopers: Above 'important' energy a 'looping' particle in field will N...
G4ThreeVector fPreviousSftOrigin
Remember last safety origin.
G4Navigator * fLinearNavigator
Propagator used to transport the particle.
G4bool fShortStepOptimisation
Whether to avoid calling G4Navigator for short step ( < safety) If using it, the safety estimate for ...
G4ParticleChangeForTransport fParticleChange
New ParticleChange.
G4bool DoesGlobalFieldExist()
Checks whether a field exists for the "global" field manager.
G4PropagatorInField * fFieldPropagator
Propagator used to transport the particle.

◆ ~G4LongLivedNeutralTransportation()

Destructor.

Definition at line 87 of file G4LongLivedNeutralTransportation.cc.

88{
89 if (fSumEnergyKilled > 0.0) {
90 PrintStatistics(G4cout);
91 }
92
93}
void PrintStatistics(std::ostream &outStr) const
returns current logging info of the algorithm
G4double fSumEnergyKilled
Sum of abandoned looping tracks energies.

Member Function Documentation

◆ AlongStepDoIt()

G4VParticleChange * AlongStepDoIt ( const G4Track & track,
const G4Step & stepData )

G4VProcess::AlongStepDoIt() implementation,.

Parameters
trackPropagating particle track reference
stepDataCurrent step reference

Definition at line 263 of file G4LongLivedNeutralTransportation.cc.

265{
266 // static G4ThreadLocal G4long noCallsASDI = 0;
267 // noCallsASDI++;
268
269 fParticleChange.Initialize(track) ;
270
271 // Code for specific process
272 //
273 fParticleChange.ProposePosition(fTransportEndPosition) ;
274 fParticleChange.ProposeMomentumDirection(fTransportEndMomentumDir) ;
276 fParticleChange.SetMomentumChanged(fMomentumChanged) ;
277
278 fParticleChange.ProposePolarization(fTransportEndSpin);
279
280 G4double deltaTime = 0.0 ;
281
282 // Calculate Lab Time of Flight (ONLY if field Equations used it!)
283 // G4double endTime = fCandidateEndGlobalTime;
284 // G4double delta_time = endTime - startTime;
285
286 G4double startTime = track.GetGlobalTime() ;
287
288 const G4DynamicParticle* pParticle = track.GetDynamicParticle() ;
289 const G4PrimaryParticle* primaryParticle = pParticle->GetPrimaryParticle();
290
291
293 // The time was not integrated .. make the best estimate possible
294 //
295
296 // use c*beta from primary particle, llp will not be slowed down anyways
297 // default was: stepData.GetPreStepPoint()->GetVelocity();, but assumes wrong values from particle definition
298 B2DEBUG(0, "Velocity calculated from Step data: " << stepData.GetPreStepPoint()->GetVelocity());
299
300 G4double initialVelocity = c_light * primaryParticle->GetTotalMomentum() / primaryParticle->GetTotalEnergy();
301 G4double stepLength = track.GetStepLength();
302
303 deltaTime = 0.0; // in case initialVelocity = 0
304 if (initialVelocity > 0.0) { deltaTime = stepLength / initialVelocity; }
305
306 fCandidateEndGlobalTime = startTime + deltaTime ;
307 fParticleChange.ProposeLocalTime(track.GetLocalTime() + deltaTime) ;
308 } else {
309 deltaTime = fCandidateEndGlobalTime - startTime ;
310 fParticleChange.ProposeGlobalTime(fCandidateEndGlobalTime) ;
311 }
312
313 // Now Correct by Lorentz factor to get delta "proper" Time
314 G4double restMass = primaryParticle->GetMass() ;
315 G4double deltaProperTime = deltaTime * (restMass / primaryParticle->GetTotalEnergy()) ;
316
317 fParticleChange.ProposeProperTime(track.GetProperTime() + deltaProperTime) ;
318 //fParticleChange.ProposeTrueStepLength( track.GetStepLength() ) ;
319
320 // Another (sometimes better way) is to use a user-limit maximum Step size
321 // to alleviate this problem ..
322
323 // Introduce smooth curved trajectories to particle-change
324 //
325 fParticleChange.SetPointerToVectorOfAuxiliaryPoints
326 (fFieldPropagator->GimmeTrajectoryVectorAndForgetIt());
327
328 return &fParticleChange ;
329}
G4double fCandidateEndGlobalTime
The particle's state after this Step, Store for DoIt.
G4ThreeVector fTransportEndMomentumDir
The particle's state after this Step, Store for DoIt.
G4ThreeVector fTransportEndSpin
The particle's state after this Step, Store for DoIt.
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 fMomentumChanged
The particle's state after this Step, Store for DoIt.
G4double fTransportEndKineticEnergy
The particle's state after this Step, Store for DoIt.

◆ AlongStepGetPhysicalInteractionLength()

G4double AlongStepGetPhysicalInteractionLength ( const G4Track & track,
G4double previousStepSize,
G4double currentMinimumStep,
G4double & currentSafety,
G4GPILSelection * selection )

G4VProcess::AlongStepGetPhysicalInteractionLength() implementation,.

Parameters
trackPropagating particle track reference
previousStepSizeThis argument of base function is ignored
currentMinimumStepCurrent minimum step size
currentSafetyReference to current step safety
selectionPointer for return value of GPILSelection, which is set to default value of CandidateForSelection
Returns
Next geometry step length

Definition at line 112 of file G4LongLivedNeutralTransportation.cc.

118{
119 G4double geometryStepLength = -1.0, newSafety = -1.0;
120
121 // Initial actions moved to StartTrack()
122 // --------------------------------------
123 // Note: in case another process changes touchable handle
124 // it will be necessary to add here (for all steps)
125 // fCurrentTouchableHandle = aTrack->GetTouchableHandle();
126
127 // GPILSelection is set to default value of CandidateForSelection
128 // It is a return value
129 //
130 *selection = CandidateForSelection ;
131
133 fLastStepInVolume = false;
134 fNewTrack = false;
135
136 fParticleChange.ProposeFirstStepInVolume(fFirstStepInVolume);
137
138 // Get initial Energy/Momentum of the track
139 //
140 const G4DynamicParticle* pParticle = track.GetDynamicParticle() ;
141 const G4PrimaryParticle* primaryParticle = pParticle->GetPrimaryParticle();
142 G4ThreeVector startMomentumDir = primaryParticle->GetMomentumDirection() ;
143 G4ThreeVector startPosition = track.GetPosition() ;
144
145 // The Step Point safety can be limited by other geometries and/or the
146 // assumptions of any process - it's not always the geometrical safety.
147 // We calculate the starting point's isotropic safety here.
148 //
149 G4ThreeVector OriginShift = startPosition - fPreviousSftOrigin ;
150 G4double MagSqShift = OriginShift.mag2() ;
151 if (MagSqShift >= sqr(fPreviousSafety)) {
152 currentSafety = 0.0 ;
153 } else {
154 currentSafety = fPreviousSafety - std::sqrt(MagSqShift) ;
155 }
156
157 // Is the particle charged or has it a magnetic moment?
158 //
159 G4double particleCharge = 0 ;
160
161 fGeometryLimitedStep = false ;
162
163 // There is no need to locate the current volume. It is Done elsewhere:
164 // On track construction
165 // By the tracking, after all AlongStepDoIts, in "Relocation"
166
167 G4double linearStepLength ;
168 if (fShortStepOptimisation && (currentMinimumStep <= currentSafety)) {
169 // The Step is guaranteed to be taken
170 //
171 geometryStepLength = currentMinimumStep ;
172 fGeometryLimitedStep = false ;
173 } else {
174 // Find whether the straight path intersects a volume
175 //
176 linearStepLength = fLinearNavigator->ComputeStep(startPosition,
177 startMomentumDir,
178 currentMinimumStep,
179 newSafety) ;
180 // Remember last safety origin & value.
181 //
182 fPreviousSftOrigin = startPosition ;
183 fPreviousSafety = newSafety ;
184
185
186 currentSafety = newSafety ;
187
188 fGeometryLimitedStep = (linearStepLength <= currentMinimumStep);
190 // The geometry limits the Step size (an intersection was found.)
191 geometryStepLength = linearStepLength ;
192 } else {
193 // The full Step is taken.
194 geometryStepLength = currentMinimumStep ;
195 }
196
197 fEndPointDistance = geometryStepLength ;
198
199 // Calculate final position
200 //
201 fTransportEndPosition = startPosition + geometryStepLength * startMomentumDir;
202
203 // Momentum direction, energy and polarisation are unchanged by transport
204 //
205 fTransportEndMomentumDir = startMomentumDir ;
206 fTransportEndKineticEnergy = primaryParticle->GetKineticEnergy() ;
207 fTransportEndSpin = track.GetPolarization();
208 fMomentumChanged = false ;
209 fEndGlobalTimeComputed = false ;
210 }
211
212
213 // If we are asked to go a step length of 0, and we are on a boundary
214 // then a boundary will also limit the step -> we must flag this.
215 //
216 if (currentMinimumStep == 0.0) {
217 if (currentSafety == 0.0) { fGeometryLimitedStep = true; }
218 }
219
220 // Update the safety starting from the end-point,
221 // if it will become negative at the end-point.
222 //
223 if (currentSafety < fEndPointDistance) {
224 if (particleCharge != 0.0) {
225 G4double endSafety =
227 currentSafety = endSafety ;
229 fPreviousSafety = currentSafety ;
230
231
232 // Because the Stepping Manager assumes it is from the start point,
233 // add the StepLength
234 //
235 currentSafety += fEndPointDistance ;
236
237#ifdef G4DEBUG_TRANSPORT
238 G4cout.precision(12) ;
239 G4cout << "***G4LongLivedNeutralTransportation::AlongStepGPIL ** " << G4endl ;
240 G4cout << " Called Navigator->ComputeSafety at " << fTransportEndPosition
241 << " and it returned safety= " << endSafety << G4endl ;
242 G4cout << " Adding endpoint distance " << fEndPointDistance
243 << " to obtain pseudo-safety= " << currentSafety << G4endl ;
244 } else {
245 G4cout << "***G4LongLivedNeutralTransportation::AlongStepGPIL ** " << G4endl ;
246 G4cout << " Avoiding call to ComputeSafety : " << G4endl;
247 G4cout << " charge = " << particleCharge << G4endl;
248 G4cout << " mag moment = " << magneticMoment << G4endl;
249#endif
250 }
251 }
252
253 fParticleChange.ProposeTrueStepLength(geometryStepLength) ;
254
255 return geometryStepLength ;
256}
G4bool fGeometryLimitedStep
Flag to determine whether a boundary was reached.

◆ AtRestDoIt()

G4VParticleChange * AtRestDoIt ( const G4Track & ,
const G4Step &  )
inline

No operation in AtRestDoIt.

Definition at line 175 of file G4LongLivedNeutralTransportation.h.

176 { return 0; }

◆ AtRestGetPhysicalInteractionLength()

G4double AtRestGetPhysicalInteractionLength ( const G4Track & ,
G4ForceCondition *  )
inline

No operation in AtRestGPIL.

Definition at line 168 of file G4LongLivedNeutralTransportation.h.

170 { return -1.0; }

◆ DoesGlobalFieldExist()

G4bool DoesGlobalFieldExist ( )
inlineprotected

Checks whether a field exists for the "global" field manager.

Definition at line 196 of file G4LongLivedNeutralTransportation.h.

197 {
198 G4TransportationManager* transportMgr;
199 transportMgr = G4TransportationManager::GetTransportationManager();
200
201 // fFieldExists= transportMgr->GetFieldManager()->DoesFieldExist();
202 // return fFieldExists;
203 return transportMgr->GetFieldManager()->DoesFieldExist();
204 }

◆ FieldExertedForce()

G4bool FieldExertedForce ( )
inline

References fFieldExertedForce.

Definition at line 107 of file G4LongLivedNeutralTransportation.h.

◆ GetMaxEnergyKilled()

G4double GetMaxEnergyKilled ( ) const
inline

Print values of looper thresholds.

Access fMaxEnergyKilled

◆ GetSilenceLooperWarnings()

G4bool GetSilenceLooperWarnings ( )
static

Do not throw exception about 'looping' particles.

Definition at line 519 of file G4LongLivedNeutralTransportation.cc.

520{
522}
static G4bool fSilenceLooperWarnings
Flag to Supress all 'looper' warnings.

◆ PostStepDoIt()

G4VParticleChange * PostStepDoIt ( const G4Track & track,
const G4Step & stepData )

G4VProcess::PostStepDoIt() implementation,.

Parameters
track
stepData

Definition at line 350 of file G4LongLivedNeutralTransportation.cc.

352{
353 G4TouchableHandle retCurrentTouchable ; // The one to return
354 G4bool isLastStep = false;
355
356 // Initialize ParticleChange (by setting all its members equal
357 // to corresponding members in G4Track)
358 // fParticleChange.Initialize(track) ; // To initialise TouchableChange
359
360 fParticleChange.ProposeTrackStatus(track.GetTrackStatus()) ;
361
362 // If the Step was determined by the volume boundary,
363 // logically relocate the particle
364
366 // fCurrentTouchable will now become the previous touchable,
367 // and what was the previous will be freed.
368 // (Needed because the preStepPoint can point to the previous touchable)
369
370 fLinearNavigator->SetGeometricallyLimitedStep() ;
372 LocateGlobalPointAndUpdateTouchableHandle(track.GetPosition(),
373 track.GetMomentumDirection(),
375 true) ;
376 // Check whether the particle is out of the world volume
377 // If so it has exited and must be killed.
378 //
379 if (fCurrentTouchableHandle->GetVolume() == 0) {
380 fParticleChange.ProposeTrackStatus(fStopAndKill) ;
381 }
382 retCurrentTouchable = fCurrentTouchableHandle ;
383 fParticleChange.SetTouchableHandle(fCurrentTouchableHandle) ;
384
385 // Update the Step flag which identifies the Last Step in a volume
387 isLastStep = fLinearNavigator->ExitedMotherVolume()
388 or fLinearNavigator->EnteredDaughterVolume() ;
389 else
390 isLastStep = fFieldPropagator->IsLastStepInVolume();
391 } else { // fGeometryLimitedStep is false
392 // This serves only to move the Navigator's location
393 //
394 fLinearNavigator->LocateGlobalPointWithinVolume(track.GetPosition()) ;
395
396 // The value of the track's current Touchable is retained.
397 // (and it must be correct because we must use it below to
398 // overwrite the (unset) one in particle change)
399 // It must be fCurrentTouchable too ??
400 //
401 fParticleChange.SetTouchableHandle(track.GetTouchableHandle()) ;
402 retCurrentTouchable = track.GetTouchableHandle() ;
403
404 isLastStep = false;
405 } // endif ( fGeometryLimitedStep )
406 fLastStepInVolume = isLastStep;
407
408 fParticleChange.ProposeFirstStepInVolume(fFirstStepInVolume);
409 fParticleChange.ProposeLastStepInVolume(isLastStep);
410
411 const G4VPhysicalVolume* pNewVol = retCurrentTouchable->GetVolume() ;
412 G4Material* pNewMaterial = nullptr ;
413 G4VSensitiveDetector* pNewSensitiveDetector = nullptr ;
414
415 if (pNewVol != 0) {
416 pNewMaterial = pNewVol->GetLogicalVolume()->GetMaterial();
417 pNewSensitiveDetector = pNewVol->GetLogicalVolume()->GetSensitiveDetector();
418 }
419
420 fParticleChange.SetMaterialInTouchable(pNewMaterial) ;
421 fParticleChange.SetSensitiveDetectorInTouchable(pNewSensitiveDetector) ;
422
423 const G4MaterialCutsCouple* pNewMaterialCutsCouple = 0;
424 if (pNewVol != 0) {
425 pNewMaterialCutsCouple = pNewVol->GetLogicalVolume()->GetMaterialCutsCouple();
426 }
427
428 if (pNewVol != 0 && pNewMaterialCutsCouple != 0
429 && pNewMaterialCutsCouple->GetMaterial() != pNewMaterial) {
430 // for parametrized volume
431 //
432 pNewMaterialCutsCouple =
433 G4ProductionCutsTable::GetProductionCutsTable()
434 ->GetMaterialCutsCouple(pNewMaterial,
435 pNewMaterialCutsCouple->GetProductionCuts());
436 }
437 fParticleChange.SetMaterialCutsCoupleInTouchable(pNewMaterialCutsCouple);
438
439 // temporarily until Get/Set Material of ParticleChange,
440 // and StepPoint can be made const.
441 // Set the touchable in ParticleChange
442 // this must always be done because the particle change always
443 // uses this value to overwrite the current touchable pointer.
444 //
445 fParticleChange.SetTouchableHandle(retCurrentTouchable) ;
446
447 return &fParticleChange ;
448}

◆ PostStepGetPhysicalInteractionLength()

G4double PostStepGetPhysicalInteractionLength ( const G4Track & track,
G4double previousStepSize,
G4ForceCondition * pForceCond )

G4VProcess::PostStepGetPhysicalInteractionLength() implementation.

Parameters
trackThis argument of base function is ignored
previousStepSizeThis argument of base function is ignored
pForceCondForce condition by default

Forces the PostStepDoIt action to be called, but does not limit the step

Definition at line 337 of file G4LongLivedNeutralTransportation.cc.

341{
342 fFieldExertedForce = false; // Not known
343 *pForceCond = Forced ;
344 return DBL_MAX ; // was kInfinity ; but convention now is DBL_MAX
345}

◆ PrintStatistics()

void PrintStatistics ( std::ostream & outStr) const

returns current logging info of the algorithm

Definition at line 98 of file G4LongLivedNeutralTransportation.cc.

99{
100 outStr << " G4LongLivedNeutralTransportation: Statistics for looping particles " << G4endl;
101 outStr << " No looping tracks found or killed. " << G4endl;
102
103}

◆ ProcessDescription()

void ProcessDescription ( std::ostream & outStream) const
virtual

G4LongLivedNeutralTransportation::ProcessDescription()

Parameters
outStreamOutput file with a description of process

Definition at line 562 of file G4LongLivedNeutralTransportation.cc.

566{
567 G4String indent = " "; // : "");
568 G4int oldPrec = outStream.precision(6);
569 // outStream << std::setprecision(6);
570 outStream << G4endl << indent << GetProcessName() << ": ";
571
572 outStream << " Parameters for looping particles: " << G4endl
573 << " warning-E = " << fThreshold_Warning_Energy / CLHEP::MeV << " MeV " << G4endl
574 << " important E = " << fThreshold_Important_Energy / CLHEP::MeV << " MeV " << G4endl
575 << " thresholdTrials " << fThresholdTrials << G4endl;
576 outStream.precision(oldPrec);
577}
G4int fThresholdTrials
Number of trials an important looper survives.
G4double fThreshold_Important_Energy
Give a few trial above this E for looping particle.
G4double fThreshold_Warning_Energy
Warn above this energy about looping particle.

◆ SetHighLooperThresholds()

void SetHighLooperThresholds ( )

Get/Set parameters for killing loopers: Above 'important' energy a 'looping' particle in field will NOT be abandoned, except after fThresholdTrials attempts.

Below Warning energy, no verbosity for looping particles is issued

Definition at line 527 of file G4LongLivedNeutralTransportation.cc.

528{
529 // Restores the old high values -- potentially appropriate for energy-frontier
530 // HEP experiments.
531 // Caution: All tracks with E < 100 MeV that are found to loop are
532 SetThresholdWarningEnergy(100.0 * CLHEP::MeV); // Warn above this energy
533 SetThresholdImportantEnergy(250.0 * CLHEP::MeV); // Extra trial above this En
534
535 G4int maxTrials = 10;
536 SetThresholdTrials(maxTrials);
537
538
539
540}
void SetThresholdWarningEnergy(G4double newEnWarn)
Set fThreshold_Warning_Energy.
void SetThresholdImportantEnergy(G4double newEnImp)
Set fThreshold_Important_Energy.
void SetThresholdTrials(G4int newMaxTrials)
Set fThresholdTrials.

◆ SetLowLooperThresholds()

void SetLowLooperThresholds ( )

Shortcut method - old values (meant for HEP)

Definition at line 543 of file G4LongLivedNeutralTransportation.cc.

544{
545 // These values were the default in Geant4 10.5 - beta
546 SetThresholdWarningEnergy(1.0 * CLHEP::keV); // Warn above this En
547 SetThresholdImportantEnergy(1.0 * CLHEP::MeV); // Extra trials above it
548
549 G4int maxTrials = 30; // A new value - was 10
550 SetThresholdTrials(maxTrials);
551
552}

◆ SetSilenceLooperWarnings()

void SetSilenceLooperWarnings ( G4bool val)
static

Do not warn about 'looping' particles.

Definition at line 511 of file G4LongLivedNeutralTransportation.cc.

512{
513 fSilenceLooperWarnings = val; // Flag to *Suppress* all 'looper' warnings
514 // G4CoupledTransportation::fSilenceLooperWarnings= val;
515}

◆ SetThresholdImportantEnergy()

void SetThresholdImportantEnergy ( G4double newEnImp)
inline

Set fThreshold_Important_Energy.

Definition at line 125 of file G4LongLivedNeutralTransportation.h.

◆ SetThresholdTrials()

void SetThresholdTrials ( G4int newMaxTrials)
inline

Set fThresholdTrials.

Definition at line 129 of file G4LongLivedNeutralTransportation.h.

◆ SetThresholdWarningEnergy()

void SetThresholdWarningEnergy ( G4double newEnWarn)
inline

Set fThreshold_Warning_Energy.

Definition at line 121 of file G4LongLivedNeutralTransportation.h.

◆ StartTracking()

void StartTracking ( G4Track * aTrack)

Reset state for new (potentially resumed) track.

Definition at line 457 of file G4LongLivedNeutralTransportation.cc.

458{
459 G4VProcess::StartTracking(aTrack);
460 fNewTrack = true;
461 fFirstStepInVolume = true;
462 fLastStepInVolume = false;
463
464 // The actions here are those that were taken in AlongStepGPIL
465 // when track.GetCurrentStepNumber()==1
466
467 // Whether field exists should be determined at run level -- TODO
469
470 // reset safety value and center
471 //
472 fPreviousSafety = 0.0 ;
473 fPreviousSftOrigin = G4ThreeVector(0., 0., 0.) ;
474
475 // reset looping counter -- for motion in field
476 fNoLooperTrials = 0;
477 // Must clear this state .. else it depends on last track's value
478 // --> a better solution would set this from state of suspended track TODO ?
479 // Was if( aTrack->GetCurrentStepNumber()==1 ) { .. }
480
481 // ChordFinder reset internal state
482 //
484 fFieldPropagator->ClearPropagatorState();
485 // Resets all state of field propagator class (ONLY) including safety
486 // values (in case of overlaps and to wipe for first track).
487 }
488
489 // Make sure to clear the chord finders of all fields (i.e. managers)
490 //
491 G4FieldManagerStore* fieldMgrStore = G4FieldManagerStore::GetInstance();
492 fieldMgrStore->ClearAllChordFindersState();
493
494 // Update the current touchable handle (from the track's)
495 //
496 fCurrentTouchableHandle = aTrack->GetTouchableHandle();
497
498 // Inform field propagator of new track
499 //
500 if (fFieldPropagator) fFieldPropagator->PrepareNewTrack();
501}
G4int fNoLooperTrials
Counter for steps in which particle reports 'looping', if it is above 'Important' Energy.

Friends And Related Symbol Documentation

◆ G4CoupledTransportation

friend class G4CoupledTransportation
friend

Definition at line 263 of file G4LongLivedNeutralTransportation.h.

Member Data Documentation

◆ fAbandonUnstableTrials

G4int fAbandonUnstableTrials = 0
private

Number of trials after which to abandon.

Definition at line 244 of file G4LongLivedNeutralTransportation.h.

◆ fAnyFieldExists

G4bool fAnyFieldExists = false
private

Flag for existing fields.

Definition at line 219 of file G4LongLivedNeutralTransportation.h.

◆ fCandidateEndGlobalTime

G4double fCandidateEndGlobalTime = 0.0
private

The particle's state after this Step, Store for DoIt.

Definition at line 217 of file G4LongLivedNeutralTransportation.h.

◆ fCurrentTouchableHandle

G4TouchableHandle fCurrentTouchableHandle
private

Current touchable handle.

Definition at line 228 of file G4LongLivedNeutralTransportation.h.

◆ fEndGlobalTimeComputed

G4bool fEndGlobalTimeComputed = false
private

The particle's state after this Step, Store for DoIt.

Definition at line 216 of file G4LongLivedNeutralTransportation.h.

◆ fEndPointDistance

G4double fEndPointDistance
private

Endpoint distance.

Definition at line 235 of file G4LongLivedNeutralTransportation.h.

◆ fFieldExertedForce

G4bool fFieldExertedForce = false
private

During current step.

Definition at line 226 of file G4LongLivedNeutralTransportation.h.

◆ fFieldPropagator

G4PropagatorInField* fFieldPropagator
private

Propagator used to transport the particle.

Definition at line 209 of file G4LongLivedNeutralTransportation.h.

◆ fFirstStepInVolume

G4bool fFirstStepInVolume = true
private

Flag first step in a geom.

volume

Definition at line 222 of file G4LongLivedNeutralTransportation.h.

◆ fGeometryLimitedStep

G4bool fGeometryLimitedStep = true
private

Flag to determine whether a boundary was reached.

Definition at line 224 of file G4LongLivedNeutralTransportation.h.

◆ fLastStepInVolume

G4bool fLastStepInVolume = false
private

Flag last step in a geom.

volume

Definition at line 223 of file G4LongLivedNeutralTransportation.h.

◆ fLinearNavigator

G4Navigator* fLinearNavigator
private

Propagator used to transport the particle.

Definition at line 208 of file G4LongLivedNeutralTransportation.h.

◆ fMomentumChanged

G4bool fMomentumChanged = true
private

The particle's state after this Step, Store for DoIt.

Definition at line 215 of file G4LongLivedNeutralTransportation.h.

◆ fNewTrack

G4bool fNewTrack = true
private

Flag from StartTracking.

Definition at line 221 of file G4LongLivedNeutralTransportation.h.

◆ fNoLooperTrials

G4int fNoLooperTrials = 0
private

Counter for steps in which particle reports 'looping', if it is above 'Important' Energy.

Definition at line 248 of file G4LongLivedNeutralTransportation.h.

◆ fParticleChange

G4ParticleChangeForTransport fParticleChange
private

New ParticleChange.

Definition at line 233 of file G4LongLivedNeutralTransportation.h.

◆ fPreviousSafety

G4double fPreviousSafety
private

Remember last safety value.

Definition at line 231 of file G4LongLivedNeutralTransportation.h.

◆ fPreviousSftOrigin

G4ThreeVector fPreviousSftOrigin
private

Remember last safety origin.

Definition at line 230 of file G4LongLivedNeutralTransportation.h.

◆ fpSafetyHelper

G4SafetyHelper* fpSafetyHelper
private

To pass it the safety value obtained.

Definition at line 259 of file G4LongLivedNeutralTransportation.h.

◆ fShortStepOptimisation

G4bool fShortStepOptimisation
private

Whether to avoid calling G4Navigator for short step ( < safety) If using it, the safety estimate for endpoint will likely be smaller.

Definition at line 257 of file G4LongLivedNeutralTransportation.h.

◆ fSilenceLooperWarnings

G4bool fSilenceLooperWarnings = false
staticprivate

Flag to Supress all 'looper' warnings.

Definition at line 266 of file G4LongLivedNeutralTransportation.h.

◆ fSumEnergyKilled

G4double fSumEnergyKilled = 0.0
private

Sum of abandoned looping tracks energies.

Definition at line 252 of file G4LongLivedNeutralTransportation.h.

◆ fThreshold_Important_Energy

G4double fThreshold_Important_Energy = 1.0 * CLHEP::MeV
private

Give a few trial above this E for looping particle.

Definition at line 240 of file G4LongLivedNeutralTransportation.h.

◆ fThreshold_Warning_Energy

G4double fThreshold_Warning_Energy = 1.0 * CLHEP::keV
private

Warn above this energy about looping particle.

Definition at line 239 of file G4LongLivedNeutralTransportation.h.

◆ fThresholdTrials

G4int fThresholdTrials = 10
private

Number of trials an important looper survives.

Definition at line 241 of file G4LongLivedNeutralTransportation.h.

◆ fTransportEndKineticEnergy

G4double fTransportEndKineticEnergy = 0.0
private

The particle's state after this Step, Store for DoIt.

Definition at line 213 of file G4LongLivedNeutralTransportation.h.

◆ fTransportEndMomentumDir

G4ThreeVector fTransportEndMomentumDir = G4ThreeVector(0.0, 0.0, 0.0)
private

The particle's state after this Step, Store for DoIt.

Definition at line 212 of file G4LongLivedNeutralTransportation.h.

◆ fTransportEndPosition

G4ThreeVector fTransportEndPosition = G4ThreeVector(0.0, 0.0, 0.0)
private

The particle's state after this Step, Store for DoIt.

Definition at line 211 of file G4LongLivedNeutralTransportation.h.

◆ fTransportEndSpin

G4ThreeVector fTransportEndSpin = G4ThreeVector(0.0, 0.0, 0.0)
private

The particle's state after this Step, Store for DoIt.

Definition at line 214 of file G4LongLivedNeutralTransportation.h.

◆ fUseGravity

G4bool fUseGravity
staticprivate

Flag take into account gravity.

Definition at line 265 of file G4LongLivedNeutralTransportation.h.

◆ fUseMagneticMoment

G4bool fUseMagneticMoment
staticprivate

Flag take into account magnetic moment.

Definition at line 264 of file G4LongLivedNeutralTransportation.h.


The documentation for this class was generated from the following files: