64 fMax_loop_count(1000),
65 fUseSafetyForOptimisation(true),
66 fZeroStepThreshold( 0.0 ),
67 fDetectorFieldMgr(detectorFieldMgr),
68 fpTrajectoryFilter( 0 ),
69 fNavigator(theNavigator),
70 fCurrentFieldMgr(detectorFieldMgr),
74 fParticleIsLooping(false),
78 fFirstStepInVolume(true),
79 fLastStepInVolume(true),
97 G4cout <<
" PiF: Zero Step Threshold set to "
100 G4cout <<
" PiF: Value of kCarTolerance = "
172 G4cout <<
"G4PropagatorInField::ComputeStep() called" <<
G4endl;
173 G4cout <<
" Starting FT: " << pFieldTrack;
174 G4cout <<
" Requested length = " << CurrentProposedStepLength <<
G4endl;
186 G4double TruePathLength = CurrentProposedStepLength;
190 G4bool first_substep =
true;
219 G4double trialProposedStep = 1.e2 * ( 10.0 *
cm +
221 GetSolid()->DistanceToOut(StartPointA, VelocityUnit) );
222 CurrentProposedStepLength=
std::min( trialProposedStep,
229 if( epsilon < epsilonMin ) epsilon = epsilonMin;
230 if( epsilon > epsilonMax ) epsilon = epsilonMax;
252 decreaseFactor= 0.25;
261 decreaseFactor = 0.35;
265 decreaseFactor= 0.75;
269 stepTrial *= decreaseFactor;
272 G4cerr <<
" G4PropagatorInField::ComputeStep(): " <<
G4endl
273 <<
" Decreasing step - in volume " << pPhysVol;
275 G4cerr <<
" with name " << pPhysVol->GetName();
278 stepTrial, pFieldTrack);
280 if( stepTrial == 0.0 )
282 std::ostringstream message;
283 message <<
"Particle abandoned due to lack of progress in field."
285 <<
" Properties : " << pFieldTrack << G4endl
286 <<
" Attempting a zero step = " << stepTrial << G4endl
287 <<
" while attempting to progress after " <<
fNoZeroStep
288 <<
" trial steps. Will abandon step.";
289 G4Exception(
"G4PropagatorInField::ComputeStep()",
"GeomNav1002",
294 if( stepTrial < CurrentProposedStepLength )
295 CurrentProposedStepLength = stepTrial;
299 G4int do_loop_count = 0;
309 G4cout <<
" PiF: Calling Nav/Locate Global Point within-Volume "
317 h_TrialStepSize = CurrentProposedStepLength - StepTaken;
338 NewSafety, LinearStepLength,
339 InterSectionPointE );
343 if( first_substep ) {
344 currentSafety = NewSafety;
353 G4bool recalculatedEndPt=
false;
356 EstimateIntersectionPoint( SubStepStartState, CurrentState,
357 InterSectionPointE, IntersectPointVelct_G,
359 intersects = intersects && found_intersection;
360 if( found_intersection ) {
362 StepTaken = TruePathLength = IntersectPointVelct_G.
GetCurveLength()
366 if( recalculatedEndPt ){
367 CurrentState= IntersectPointVelct_G;
373 StepTaken += s_length_taken;
379 first_substep =
false;
384 CurrentState, CurrentProposedStepLength,
385 NewSafety, do_loop_count, pPhysVol );
389 G4cout <<
" G4PropagatorInField::ComputeStep(): " <<
G4endl
390 <<
" Difficult track - taking many sub steps." <<
G4endl;
392 printStatus( SubStepStartState, CurrentState, CurrentProposedStepLength,
393 NewSafety, do_loop_count, pPhysVol );
399 }
while( (!intersects )
417 TruePathLength = StepTaken;
434 std::ostringstream message;
435 message <<
"Curve length mis-match between original state "
436 <<
"and proposed endpoint of propagation." <<
G4endl
437 <<
" The curve length of the endpoint should be: "
439 <<
" and it is instead: "
441 <<
" A difference of: "
444 <<
" Original state = " << OriginalState <<
G4endl
451 if( TruePathLength+
kCarTolerance >= CurrentProposedStepLength )
476 return TruePathLength;
501 if( ((stepNo == 0) && (verboseLevel <3)) || (verboseLevel >= 3) )
503 oldprec =
G4cout.precision(4);
504 G4cout << std::setw( 6) <<
" "
505 << std::setw( 25) <<
" Current Position and Direction" <<
" "
507 G4cout << std::setw( 5) <<
"Step#"
508 << std::setw(10) <<
" s " <<
" "
509 << std::setw(10) <<
"X(mm)" <<
" "
510 << std::setw(10) <<
"Y(mm)" <<
" "
511 << std::setw(10) <<
"Z(mm)" <<
" "
512 << std::setw( 7) <<
" N_x " <<
" "
513 << std::setw( 7) <<
" N_y " <<
" "
514 << std::setw( 7) <<
" N_z " <<
" " ;
515 G4cout << std::setw( 7) <<
" Delta|N|" <<
" "
516 << std::setw( 9) <<
"StepLen" <<
" "
517 << std::setw(12) <<
"StartSafety" <<
" "
518 << std::setw( 9) <<
"PhsStep" <<
" ";
520 {
G4cout << std::setw(18) <<
"NextVolume" <<
" "; }
521 G4cout.precision(oldprec);
524 if((stepNo == 0) && (verboseLevel <=3))
528 printStatus( StartFT, StartFT, -1.0, safety, -1, startVolume);
530 if( verboseLevel <= 3 )
533 {
G4cout << std::setw( 4) << stepNo <<
" "; }
535 {
G4cout << std::setw( 5) <<
"Start" ; }
536 oldprec =
G4cout.precision(8);
539 G4cout << std::setw(10) << CurrentPosition.x() <<
" "
540 << std::setw(10) << CurrentPosition.y() <<
" "
541 << std::setw(10) << CurrentPosition.z() <<
" ";
543 G4cout << std::setw( 7) << CurrentUnitVelocity.x() <<
" "
544 << std::setw( 7) << CurrentUnitVelocity.y() <<
" "
545 << std::setw( 7) << CurrentUnitVelocity.z() <<
" ";
549 G4cout << std::setw( 9) << step_len <<
" ";
550 G4cout << std::setw(12) << safety <<
" ";
551 if( requestStep != -1.0 )
552 {
G4cout << std::setw( 9) << requestStep <<
" "; }
554 {
G4cout << std::setw( 9) <<
"Init/NotKnown" <<
" "; }
555 if( startVolume != 0)
556 {
G4cout << std::setw(12) << startVolume->
GetName() <<
" "; }
557 G4cout.precision(oldprec);
564 G4cout <<
"Step taken was " << step_len
565 <<
" out of PhysicalStep = " << requestStep <<
G4endl;
567 G4cout <<
"Chord length = " << (CurrentPosition-StartPosition).mag()
585 G4cout <<
" " << std::setw(12) <<
" PiF: NoZeroStep "
586 <<
" " << std::setw(20) <<
" CurrentProposed len "
587 <<
" " << std::setw(18) <<
" Full_curvelen_last"
588 <<
" " << std::setw(18) <<
" last proposed len "
589 <<
" " << std::setw(18) <<
" decrease factor "
590 <<
" " << std::setw(15) <<
" step trial "
594 <<
" " << std::setw(20) << CurrentProposedStepLength
597 <<
" " << std::setw(18) << decreaseFactor
598 <<
" " << std::setw(15) << stepTrial
600 G4cout.precision( iprec );
613 std::vector<G4ThreeVector>*
646 0.0,0.0,0.0,0.0,0.0);
660 if( pCurrentPhysicalVolume)
670 if( pRegionFieldMgr )
671 currentFieldMgr= pRegionFieldMgr;
677 currentFieldMgr = localFieldMgr;
685 return currentFieldMgr;
707 std::ostringstream message;
708 message <<
" Killing looping particle "
710 <<
" after " << count <<
" field substeps "
711 <<
" totaling " << StepTaken /
mm <<
" mm " ;
715 G4cout <<
" in unknown or null volume. " ;
717 G4Exception(
"G4PropagatorInField::ComputeStep()",
"GeomNav1002",
726 std::ostringstream message;
727 message <<
"Particle is stuck; it will be killed." <<
G4endl
728 <<
" Zero progress for " << noZeroSteps <<
" attempted steps."
730 <<
" Proposed Step is " << proposedStep
731 <<
" but Step Taken is "<< lastTriedStep <<
G4endl;
733 message <<
" in volume " << physVol->
GetName() ;
735 message <<
" in unknown or null volume. " ;
G4FieldManager * fCurrentFieldMgr
void SetEpsilonStep(G4double newEps)
G4VCurvedTrajectoryFilter * fpTrajectoryFilter
static const G4double kInfinity
G4double GetCurveLength() const
std::vector< G4ThreeVector > * GimmeTrajectoryVectorAndForgetIt() const
CLHEP::Hep3Vector G4ThreeVector
void SetVerboseLevel(G4int newLevel)
virtual void TakeIntermediatePoint(G4ThreeVector newPoint)=0
void PrintStepLengthDiagnostic(G4double currentProposedStepLength, G4double decreaseFactor, G4double stepTrial, const G4FieldTrack &aFieldTrack)
G4double fLast_ProposedStepLength
void RefreshIntersectionLocator()
G4double GetSurfaceTolerance() const
const G4ThreeVector & GetMomentumDir() const
G4double GetDeltaOneStep() const
std::vector< G4ThreeVector > * GimmeThePointsAndForgetThem()
G4PropagatorInField(G4Navigator *theNavigator, G4FieldManager *detectorFieldMgr, G4VIntersectionLocator *vLocator=0)
void SetChordFinderFor(G4ChordFinder *fCFinder)
G4Region * GetRegion() const
G4ThreeVector fPreviousSftOrigin
G4double fFull_CurveLen_of_LastAttempt
G4bool fUseSafetyForOptimisation
G4int SetVerboseLevel(G4int verbose)
G4FieldManager * GetFieldManager() const
G4ThreeVector GetPosition() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetName() const
static const double meter
void ReportLoopingParticle(G4int count, double StepTaken, G4VPhysicalVolume *pPhysVol)
G4double GetMaximumEpsilonStep() const
G4double fZeroStepThreshold
G4double fLargestAcceptableStep
void CreateNewTrajectorySegment()
G4FieldManager * FindAndSetFieldManager(G4VPhysicalVolume *pCurrentPhysVol)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
G4FieldManager * GetFieldManager() const
G4double AdvanceChordLimited(G4FieldTrack &yCurrent, G4double stepInitial, G4double epsStep_Relative, const G4ThreeVector latestSafetyOrigin, G4double lasestSafetyRadius)
static const double micrometer
G4double ComputeStep(G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4double &pNewSafety, G4VPhysicalVolume *pPhysVol=0)
G4int fSevereActionThreshold_NoZeroSteps
G4double GetDeltaIntersection() const
G4LogicalVolume * GetLogicalVolume() const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
void ReportStuckParticle(G4int noZeroSteps, G4double proposedStep, G4double lastTriedStep, G4VPhysicalVolume *physVol)
G4ChordFinder * GetChordFinder()
G4bool fFirstStepInVolume
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4FieldManager * fDetectorFieldMgr
G4bool IntersectChord(const G4ThreeVector &StartPointA, const G4ThreeVector &EndPointB, G4double &NewSafety, G4double &LinearStepLength, G4ThreeVector &IntersectionPoint)
void SetTrajectoryFilter(G4VCurvedTrajectoryFilter *filter)
static const double millimeter
void ClearPropagatorState()
G4bool fParticleIsLooping
void SetSafetyParametersFor(G4bool UseSafety)
G4int fAbandonThreshold_NoZeroSteps
G4FieldTrack End_PointAndTangent
G4ThreeVector GetMomentum() const
G4MagInt_Driver * GetIntegrationDriver()
G4int fActionThreshold_NoZeroSteps
G4VPhysicalVolume * GetWorldVolume() const
void SetEpsilonStepFor(G4double EpsilonStep)
static G4GeometryTolerance * GetInstance()
virtual void LocateGlobalPointWithinVolume(const G4ThreeVector &position)
G4VIntersectionLocator * fIntersectionLocator
G4GLOB_DLL std::ostream G4cerr
void SetDeltaIntersectionFor(G4double deltaIntersection)
void printStatus(const G4FieldTrack &startFT, const G4FieldTrack ¤tFT, G4double requestStep, G4double safety, G4int step, G4VPhysicalVolume *startVolume)
G4double GetMinimumEpsilonStep() const