Geant4  10.00.p03
G4PathFinder.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4PathFinder.cc 81596 2014-06-03 14:08:49Z gcosmo $
28 // GEANT4 tag $ Name: $
29 //
30 // class G4PathFinder Implementation
31 //
32 // Original author: John Apostolakis, April 2006
33 //
34 // --------------------------------------------------------------------
35 
36 #include <iomanip>
37 
38 #include "G4PathFinder.hh"
39 
40 #include "G4SystemOfUnits.hh"
41 #include "G4GeometryTolerance.hh"
42 #include "G4Navigator.hh"
43 #include "G4PropagatorInField.hh"
45 #include "G4MultiNavigator.hh"
46 #include "G4SafetyHelper.hh"
47 
48 // Initialise the static instance of the singleton
49 //
51 
52 // ----------------------------------------------------------------------------
53 // GetInstance()
54 //
55 // Retrieve the static instance of the singleton
56 //
58 {
59  if( ! fpPathFinder )
60  {
62  }
63  return fpPathFinder;
64 }
65 
66 // ----------------------------------------------------------------------------
67 // Constructor
68 //
70  : fEndState( G4ThreeVector(), G4ThreeVector(), 0., 0., 0., 0., 0.),
71  fFieldExertedForce(false),
72  fRelocatedPoint(true),
73  fLastStepNo(-1), fCurrentStepNo(-1),
74  fVerboseLevel(0)
75 {
77 
80 
82 
85  fLastLocatedPosition= Big3Vector;
86  fSafetyLocation= Big3Vector;
87  fPreSafetyLocation= Big3Vector;
88  fPreStepLocation= Big3Vector;
89 
90  fPreSafetyMinValue= -1.0;
91  fMinSafety_PreStepPt= -1.0;
93  fMinStep= -1.0;
94  fTrueMinStep= -1.0;
95  fPreStepCenterRenewed= false;
96  fNewTrack= false;
98 
99  for( G4int num=0; num< fMaxNav; ++num )
100  {
101  fpNavigator[num] = 0;
102  fLimitTruth[num] = false;
104  fCurrentStepSize[num] = -1.0;
105  fLocatedVolume[num] = 0;
106  fPreSafetyValues[num]= -1.0;
107  fCurrentPreStepSafety[num] = -1.0;
108  fNewSafetyComputed[num]= -1.0;
109  }
110 }
111 
112 // ----------------------------------------------------------------------------
113 // Destructor
114 //
116 {
117  delete fpMultiNavigator;
118  if (fpPathFinder) { delete fpPathFinder; fpPathFinder=0; }
119 }
120 
121 // ----------------------------------------------------------------------------
122 //
123 void
125 {
126  G4Navigator *navigatorForPropagation=0, *massNavigator=0;
127 
128  massNavigator= fpTransportManager->GetNavigatorForTracking();
129  if( enableChoice )
130  {
131  navigatorForPropagation= fpMultiNavigator;
132 
133  // Enable SafetyHelper to use PF
134  //
136  }
137  else
138  {
139  navigatorForPropagation= massNavigator;
140 
141  // Disable SafetyHelper to use PF
142  //
144  }
145  fpFieldPropagator->SetNavigatorForPropagating(navigatorForPropagation);
146 }
147 
148 // ----------------------------------------------------------------------------
149 //
150 G4double
151 G4PathFinder::ComputeStep( const G4FieldTrack &InitialFieldTrack,
152  G4double proposedStepLength,
153  G4int navigatorNo,
154  G4int stepNo, // find next step
155  G4double &pNewSafety, // for this geom
156  ELimited &limitedStep,
157  G4FieldTrack &EndState,
158  G4VPhysicalVolume* currentVolume)
159 {
160  G4double possibleStep= -1.0;
161 
162 #ifdef G4DEBUG_PATHFINDER
163  if( fVerboseLevel > 2 )
164  {
165  G4cout << " -------------------------" << G4endl;
166  G4cout << " G4PathFinder::ComputeStep - entered " << G4endl;
167  G4cout << " - stepNo = " << std::setw(4) << stepNo << " "
168  << " navigatorId = " << std::setw(2) << navigatorNo << " "
169  << " proposed step len = " << proposedStepLength << " " << G4endl;
170  G4cout << " PF::ComputeStep requested step "
171  << " from " << InitialFieldTrack.GetPosition()
172  << " dir " << InitialFieldTrack.GetMomentumDirection() << G4endl;
173  }
174 #endif
175 #ifdef G4VERBOSE
176  if( navigatorNo >= fNoActiveNavigators )
177  {
178  std::ostringstream message;
179  message << "Bad Navigator ID !" << G4endl
180  << " Requested Navigator ID = " << navigatorNo << G4endl
181  << " Number of active navigators = " << fNoActiveNavigators;
182  G4Exception("G4PathFinder::ComputeStep()", "GeomNav0002",
183  FatalException, message);
184  }
185 #endif
186 
187  if( fNewTrack || (stepNo != fLastStepNo) )
188  {
189  // This is a new track or a new step, so we must make the step
190  // ( else we can simply retrieve its results for this Navigator Id )
191 
192  G4FieldTrack currentState= InitialFieldTrack;
193 
194  fCurrentStepNo = stepNo;
195 
196  // Check whether a process shifted the position
197  // since the last step -- by physics processes
198  //
199  G4ThreeVector newPosition = InitialFieldTrack.GetPosition();
200  G4ThreeVector moveVector= newPosition - fLastLocatedPosition;
201  G4double moveLenSq= moveVector.mag2();
202  if( moveLenSq > kCarTolerance * kCarTolerance )
203  {
204  G4ThreeVector newDirection = InitialFieldTrack.GetMomentumDirection();
205 #ifdef G4DEBUG_PATHFINDER
206  if( fVerboseLevel > 2 )
207  {
208  G4double moveLen= std::sqrt( moveLenSq );
209  G4cout << " G4PathFinder::ComputeStep : Point moved since last step "
210  << " -- at step # = " << stepNo << G4endl
211  << " by " << moveLen << " to " << newPosition << G4endl;
212  }
213 #endif
214  MovePoint(); // Unintentional changed -- ????
215 
216  // Relocate to cope with this move -- else could abort !?
217  //
218  Locate( newPosition, newDirection );
219  }
220 
221  // Check whether the particle have an (EM) field force exerting upon it
222  //
223  G4double particleCharge= currentState.GetCharge();
224 
225  G4FieldManager* fieldMgr=0;
226  G4bool fieldExertsForce = false ;
227  if( (particleCharge != 0.0) )
228  {
229  fieldMgr= fpFieldPropagator->FindAndSetFieldManager( currentVolume );
230 
231  // Protect for case where field manager has no field (= field is zero)
232  //
233  fieldExertsForce = (fieldMgr != 0)
234  && (fieldMgr->GetDetectorField() != 0);
235  }
236  fFieldExertedForce = fieldExertsForce; // Store for use in later calls
237  // referring to this 'step'.
238 
239  fNoGeometriesLimiting= -1; // At start of track, no process limited step
240  if( fieldExertsForce )
241  {
242  DoNextCurvedStep( currentState, proposedStepLength, currentVolume );
243  //--------------
244  }else{
245  DoNextLinearStep( currentState, proposedStepLength );
246  //--------------
247  }
248  fLastStepNo= stepNo;
249 
250 #ifdef G4DEBUG_PATHFINDER
251  if ( (fNoGeometriesLimiting < 0)
253  {
254  std::ostringstream message;
255  message << "Number of geometries limiting the step not set." << G4endl
256  << " Number of geometries limiting step = "
258  G4Exception("G4PathFinder::ComputeStep()",
259  "GeomNav0002", FatalException, message);
260  }
261 #endif
262  }
263 #ifdef G4DEBUG_PATHFINDER
264  else
265  {
266  if( proposedStepLength < fTrueMinStep ) // For 2nd+ geometry
267  {
268  std::ostringstream message;
269  message << "Problem in step size request." << G4endl
270  << " Error can be caused by incorrect process ordering."
271  << " Being requested to make a step which is shorter"
272  << " than the minimum Step " << G4endl
273  << " already computed for any Navigator/geometry during"
274  << " this tracking-step: " << G4endl
275  << " This can happen due to an error in process ordering."
276  << G4endl
277  << " Check that all physics processes are registered"
278  << G4endl
279  << " before all processes with a navigator/geometry."
280  << G4endl
281  << " If using pre-packaged physics list and/or"
282  << G4endl
283  << " functionality, please report this error."
284  << G4endl << G4endl
285  << " Additional information for problem: " << G4endl
286  << " Steps request/proposed = " << proposedStepLength
287  << G4endl
288  << " MinimumStep (true) = " << fTrueMinStep
289  << G4endl
290  << " MinimumStep (navraw) = " << fMinStep
291  << G4endl
292  << " Navigator raw return value" << G4endl
293  << " Requested step now = " << proposedStepLength
294  << G4endl
295  << " Difference min-req = "
296  << fTrueMinStep-proposedStepLength << G4endl
297  << " -- Step info> stepNo= " << stepNo
298  << " last= " << fLastStepNo
299  << " newTr= " << fNewTrack << G4endl;
300  G4Exception("G4PathFinder::ComputeStep()",
301  "GeomNav0003", FatalException, message);
302  }
303  else
304  {
305  // This is neither a new track nor a new step -- just another
306  // client accessing information for the current track, step
307  // We will simply retrieve the results of the synchronous
308  // stepping for this Navigator Id below.
309  //
310  if( fVerboseLevel > 1 )
311  {
312  G4cout << " G4P::CS -> Not calling DoNextLinearStep: "
313  << " stepNo= " << stepNo << " last= " << fLastStepNo
314  << " new= " << fNewTrack << " Step already done" << G4endl;
315  }
316  }
317  }
318 #endif
319 
320  fNewTrack= false;
321 
322  // Prepare the information to return
323 
324  pNewSafety = fCurrentPreStepSafety[ navigatorNo ];
325  limitedStep = fLimitedStep[ navigatorNo ];
326  fRelocatedPoint= false;
327 
328  possibleStep= std::min(proposedStepLength, fCurrentStepSize[ navigatorNo ]);
329  EndState = fEndState; // now corrected for smaller step, if needed
330 
331 #ifdef G4DEBUG_PATHFINDER
332  if( fVerboseLevel > 0 )
333  {
334  G4cout << " G4PathFinder::ComputeStep returns "
335  << fCurrentStepSize[ navigatorNo ]
336  << " for Navigator " << navigatorNo
337  << " Limited step = " << limitedStep
338  << " Safety(mm) = " << pNewSafety / mm
339  << G4endl;
340  }
341 #endif
342 
343  return possibleStep;
344 }
345 
346 // ----------------------------------------------------------------------
347 
348 void
350  const G4ThreeVector& direction,
351  G4VPhysicalVolume* massStartVol)
352 {
353  // Key purposes:
354  // - Check and cache set of active navigators
355  // - Reset state for new track
356 
357  G4int num=0;
358 
360  // Switch PropagatorInField to use MultiNavigator
361 
363  // Reinitialise state of safety helper -- avoid problems with overlaps
364 
365  fNewTrack= true;
366  this->MovePoint(); // Signal further that the last status is wiped
367 
368  // Message the G4NavigatorPanel / Dispatcher to find active navigators
369  //
370  std::vector<G4Navigator*>::iterator pNavigatorIter;
371 
372  fNoActiveNavigators= fpTransportManager-> GetNoActiveNavigators();
374  {
375  std::ostringstream message;
376  message << "Too many active Navigators / worlds." << G4endl
377  << " Transportation Manager has "
378  << fNoActiveNavigators << " active navigators." << G4endl
379  << " This is more than the number allowed = "
380  << fMaxNav << " !";
381  G4Exception("G4PathFinder::PrepareNewTrack()", "GeomNav0002",
382  FatalException, message);
383  }
384 
386  //------------------------------------
387 
389  for( num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num )
390  {
391  // Keep information in C-array ... for creating touchables - at least
392 
393  fpNavigator[num] = *pNavigatorIter;
394  fLimitTruth[num] = false;
395  fLimitedStep[num] = kDoNot;
396  fCurrentStepSize[num] = 0.0;
397  fLocatedVolume[num] = 0;
398  }
399  fNoGeometriesLimiting= 0; // At start of track, no process limited step
400 
401  // In case of one geometry, the tracking will have done the locating!!
402 
403  if( fNoActiveNavigators > 1 )
404  {
405  Locate( position, direction, false );
406  }
407  else
408  {
409  // Update state -- depending on the tracking's call to Mass Navigator
410 
412  fLocatedVolume[0]= massStartVol; // This information must be given
413  // by transportation
414  fLimitedStep[0] = kDoNot;
415  fCurrentStepSize[0] = 0.0;
416  }
417 
418  // Reset Safety Information -- as in case of overlaps this can cause
419  // inconsistencies ...
420  //
422 
423  for( num=0; num< fNoActiveNavigators; ++num )
424  {
425  fPreSafetyValues[num]= 0.0;
426  fNewSafetyComputed[num]= 0.0;
427  fCurrentPreStepSafety[num] = 0.0;
428  }
429 
430  // The first location for each Navigator must be non-relative
431  // or else call ResetStackAndState() for each Navigator
432 
433  fRelocatedPoint= false;
434 }
435 
436 
438  // Signal end of tracking of current track.
439  // Reset TransportationManager to use 'ordinary' Navigator
440  // Reset internal state, if needed
441 {
442  EnableParallelNavigation(false); // Else it will be continue to be used
443 }
444 
445 void G4PathFinder::ReportMove( const G4ThreeVector& OldVector,
446  const G4ThreeVector& NewVector,
447  const G4String& Quantity ) const
448 {
449  G4ThreeVector moveVec = ( NewVector - OldVector );
450 
451  G4int prc= G4cerr.precision(12);
452  std::ostringstream message;
453  message << "Endpoint moved between value returned by ComputeStep()"
454  << " and call to Locate(). " << G4endl
455  << " Change of " << Quantity << " is "
456  << moveVec.mag() / mm << " mm long" << G4endl
457  << " and its vector is "
458  << (1.0/mm) * moveVec << " mm " << G4endl
459  << " Endpoint of ComputeStep() was " << OldVector << G4endl
460  << " and current position to locate is " << NewVector;
461  G4Exception("G4PathFinder::ReportMove()", "GeomNav1002",
462  JustWarning, message);
463  G4cerr.precision(prc);
464 }
465 
466 void
468  const G4ThreeVector& direction,
469  G4bool relative)
470 {
471  // Locate the point in each geometry
472 
473  std::vector<G4Navigator*>::iterator pNavIter=
475 
476  G4ThreeVector lastEndPosition= fEndState.GetPosition();
477  G4ThreeVector moveVec = (position - lastEndPosition );
478  G4double moveLenSq= moveVec.mag2();
479  if( (!fNewTrack) && (!fRelocatedPoint)
480  && ( moveLenSq> 10*kCarTolerance*kCarTolerance ) )
481  {
482  ReportMove( position, lastEndPosition, "Position" );
483  }
485 
486 #ifdef G4DEBUG_PATHFINDER
487  if( fVerboseLevel > 2 )
488  {
489  G4cout << G4endl;
490  G4cout << " G4PathFinder::Locate : entered " << G4endl;
491  G4cout << " -------------------- -------" << G4endl;
492  G4cout << " Locating at position " << position
493  << " with direction " << direction
494  << " relative= " << relative << G4endl;
495  if ( (fVerboseLevel > 1) || ( moveLenSq > 0.0) )
496  {
497  G4cout << " lastEndPosition = " << lastEndPosition
498  << " moveVec = " << moveVec
499  << " newTr = " << fNewTrack
500  << " relocated = " << fRelocatedPoint << G4endl;
501  }
502 
503  G4cout << " Located at " << position ;
504  if( fNoActiveNavigators > 1 ) { G4cout << G4endl; }
505  }
506 #endif
507 
508  for ( G4int num=0; num< fNoActiveNavigators ; ++pNavIter,++num )
509  {
510  // ... who limited the step ....
511 
512  if( fLimitTruth[num] ) { (*pNavIter)->SetGeometricallyLimitedStep(); }
513 
514  G4VPhysicalVolume *pLocated=
515  (*pNavIter)->LocateGlobalPointAndSetup( position, &direction,
516  relative,
517  false);
518  // Set the state related to the location
519  //
520  fLocatedVolume[num] = pLocated;
521 
522  // Clear state related to the step
523  //
524  fLimitedStep[num] = kDoNot;
525  fCurrentStepSize[num] = 0.0;
526 
527 #ifdef G4DEBUG_PATHFINDER
528  if( fVerboseLevel > 2 )
529  {
530  G4cout << " - In world " << num << " geomLimStep= " << fLimitTruth[num]
531  << " gives volume= " << pLocated ;
532  if( pLocated )
533  {
534  G4cout << " name = '" << pLocated->GetName() << "'";
535  G4cout << " - CopyNo= " << pLocated->GetCopyNo();
536  }
537  G4cout << G4endl;
538  }
539 #endif
540  }
541 
542  fRelocatedPoint= false;
543 }
544 
546 {
547  // Locate the point in each geometry
548 
549  std::vector<G4Navigator*>::iterator pNavIter=
551 
552 #ifdef G4DEBUG_PATHFINDER
553 
554  // Check that this relocation does not violate safety
555  // - at endpoint (computed from start point) AND
556  // - at last safety location (likely just called)
557 
558  G4ThreeVector lastEndPosition= fEndState.GetPosition();
559 
560  // Calculate end-point safety ...
561  //
562  G4double DistanceStartEnd= (lastEndPosition - fPreStepLocation).mag();
563  G4double endPointSafety_raw = fMinSafety_PreStepPt - DistanceStartEnd;
564  G4double endPointSafety_Est1 = std::max( 0.0, endPointSafety_raw );
565 
566  // ... and check move from endpoint against this endpoint safety
567  //
568  G4ThreeVector moveVecEndPos = position - lastEndPosition;
569  G4double moveLenEndPosSq = moveVecEndPos.mag2();
570 
571  // Check that move from endpoint of last step is within safety
572  // -- or check against last location or relocation ??
573  //
574  G4ThreeVector moveVecSafety= position - fSafetyLocation;
575  G4double moveLenSafSq= moveVecSafety.mag2();
576 
577  G4double distCheckEnd_sq= ( moveLenEndPosSq - endPointSafety_Est1
578  *endPointSafety_Est1 );
579  G4double distCheckSaf_sq= ( moveLenSafSq - fMinSafety_atSafLocation
581 
582  G4bool longMoveEnd = distCheckEnd_sq > 0.0;
583  G4bool longMoveSaf = distCheckSaf_sq > 0.0;
584 
585  G4double revisedSafety= 0.0;
586 
587  if( (!fNewTrack) && ( longMoveEnd && longMoveSaf ) )
588  {
589  // Recompute ComputeSafety for end position
590  //
591  revisedSafety= ComputeSafety(lastEndPosition);
592 
593  const G4double kRadTolerance =
595  const G4double cErrorTolerance=1e-12;
596  // Maximum relative error from roundoff of arithmetic
597 
598  G4double distCheckRevisedEnd= moveLenEndPosSq-revisedSafety*revisedSafety;
599 
600  G4bool longMoveRevisedEnd= ( distCheckRevisedEnd > 0. ) ;
601 
602  G4double moveMinusSafety= 0.0;
603  G4double moveLenEndPosition= std::sqrt( moveLenEndPosSq );
604  moveMinusSafety = moveLenEndPosition - revisedSafety;
605 
606  if ( longMoveRevisedEnd && (moveMinusSafety > 0.0 )
607  && ( revisedSafety > 0.0 ) )
608  {
609  // Take into account possibility of roundoff error causing
610  // this apparent move further than safety
611 
612  if( fVerboseLevel > 0 )
613  {
614  G4cout << " G4PF:Relocate> Ratio to revised safety is "
615  << std::fabs(moveMinusSafety)/revisedSafety << G4endl;
616  }
617 
618  G4double absMoveMinusSafety= std::fabs(moveMinusSafety);
619  G4bool smallRatio= absMoveMinusSafety < kRadTolerance * revisedSafety ;
620  G4double maxCoordPos = std::max(
621  std::max( std::fabs(position.x()),
622  std::fabs(position.y())),
623  std::fabs(position.z()) );
624  G4bool smallValue= absMoveMinusSafety < cErrorTolerance * maxCoordPos;
625  if( ! (smallRatio || smallValue) )
626  {
627  G4cout << " G4PF:Relocate> Ratio to revised safety is "
628  << std::fabs(moveMinusSafety)/revisedSafety << G4endl;
629  G4cout << " Difference of move and safety is not very small."
630  << G4endl;
631  }
632  else
633  {
634  moveMinusSafety = 0.0;
635  longMoveRevisedEnd = false; // Numerical issue -- not too long!
636 
637  G4cout << " Difference of move & safety is very small in magnitude, "
638  << absMoveMinusSafety << G4endl;
639  if( smallRatio )
640  {
641  G4cout << " ratio to safety " << revisedSafety
642  << " is " << absMoveMinusSafety / revisedSafety
643  << "smaller than " << kRadTolerance << " of safety ";
644  }
645  else
646  {
647  G4cout << " as fraction " << absMoveMinusSafety / maxCoordPos
648  << " of position vector max-coord " << maxCoordPos
649  << " smaller than " << cErrorTolerance ;
650  }
651  G4cout << " -- reset moveMinusSafety to "
652  << moveMinusSafety << G4endl;
653  }
654  }
655 
656  if ( longMoveEnd && longMoveSaf
657  && longMoveRevisedEnd && (moveMinusSafety>0.0) )
658  {
659  G4int oldPrec= G4cout.precision(9);
660  std::ostringstream message;
661  message << "ReLocation is further than end-safety value." << G4endl
662  << " Moved from last endpoint by " << moveLenEndPosition
663  << " compared to end safety (from preStep point) = "
664  << endPointSafety_Est1 << G4endl
665  << " --> last PreSafety Location was " << fPreSafetyLocation
666  << G4endl
667  << " safety value = " << fPreSafetyMinValue << G4endl
668  << " --> last PreStep Location was " << fPreStepLocation
669  << G4endl
670  << " safety value = " << fMinSafety_PreStepPt << G4endl
671  << " --> last EndStep Location was " << lastEndPosition
672  << G4endl
673  << " safety value = " << endPointSafety_Est1
674  << " raw-value = " << endPointSafety_raw << G4endl
675  << " --> Calling again at this endpoint, we get "
676  << revisedSafety << " as safety value." << G4endl
677  << " --> last position for safety " << fSafetyLocation
678  << G4endl
679  << " its safety value = " << fMinSafety_atSafLocation
680  << G4endl
681  << " move from safety location = "
682  << std::sqrt(moveLenSafSq) << G4endl
683  << " again= " << moveVecSafety.mag() << G4endl
684  << " safety - Move-from-end= "
685  << revisedSafety - moveLenEndPosition
686  << " (negative is Bad.)" << G4endl
687  << " Debug: distCheckRevisedEnd = "
688  << distCheckRevisedEnd;
689  ReportMove( lastEndPosition, position, "Position" );
690  G4Exception("G4PathFinder::ReLocate", "GeomNav0003",
691  FatalException, message);
692  G4cout.precision(oldPrec);
693  }
694  }
695 
696  if( fVerboseLevel > 2 )
697  {
698  G4cout << G4endl;
699  G4cout << " G4PathFinder::ReLocate : entered " << G4endl;
700  G4cout << " ---------------------- -------" << G4endl;
701  G4cout << " *Re*Locating at position " << position << G4endl;
702  // << " with direction " << direction
703  // << " relative= " << relative << G4endl;
704  if ( (fVerboseLevel > -1) || ( moveLenEndPosSq > 0.0) )
705  {
706  G4cout << " lastEndPosition = " << lastEndPosition
707  << " moveVec from step-end = " << moveVecEndPos
708  << " is new Track = " << fNewTrack
709  << " relocated = " << fRelocatedPoint << G4endl;
710  }
711  }
712 #endif // G4DEBUG_PATHFINDER
713 
714  for ( G4int num=0; num< fNoActiveNavigators ; ++pNavIter,++num )
715  {
716  // ... none limited the step
717 
718  (*pNavIter)->LocateGlobalPointWithinVolume( position );
719 
720  // Clear state related to the step
721  //
722  fLimitedStep[num] = kDoNot;
723  fCurrentStepSize[num] = 0.0;
724  fLimitTruth[num] = false;
725  }
726 
728  fRelocatedPoint= false;
729 
730 #ifdef G4DEBUG_PATHFINDER
731  if( fVerboseLevel > 2 )
732  {
733  G4cout << " G4PathFinder::ReLocate : exiting "
734  << " at position " << fLastLocatedPosition << G4endl << G4endl;
735  }
736 #endif
737 }
738 
739 // -----------------------------------------------------------------------------
740 
742 {
743  // Recompute safety for the relevant point
744 
745  G4double minSafety= kInfinity;
746 
747  std::vector<G4Navigator*>::iterator pNavigatorIter;
749 
750  for( G4int num=0; num<fNoActiveNavigators; ++pNavigatorIter,++num )
751  {
752  G4double safety = (*pNavigatorIter)->ComputeSafety( position,true );
753  if( safety < minSafety ) { minSafety = safety; }
754  fNewSafetyComputed[num]= safety;
755  }
756 
758  fMinSafety_atSafLocation = minSafety;
759 
760 #ifdef G4DEBUG_PATHFINDER
761  if( fVerboseLevel > 1 )
762  {
763  G4cout << " G4PathFinder::ComputeSafety - returns "
764  << minSafety << " at location " << position << G4endl;
765  }
766 #endif
767  return minSafety;
768 }
769 
770 
771 // -----------------------------------------------------------------------------
772 
775 {
776 #ifdef G4DEBUG_PATHFINDER
777  if( fVerboseLevel > 2 )
778  {
779  G4cout << "G4PathFinder::CreateTouchableHandle : navId = "
780  << navId << " -- " << GetNavigator(navId) << G4endl;
781  }
782 #endif
783 
784  G4TouchableHistory* touchHist;
785  touchHist= GetNavigator(navId) -> CreateTouchableHistory();
786 
787  G4VPhysicalVolume* locatedVolume= fLocatedVolume[navId];
788  if( locatedVolume == 0 )
789  {
790  // Workaround to ensure that the touchable is fixed !! // TODO: fix
791 
792  touchHist->UpdateYourself( locatedVolume, touchHist->GetHistory() );
793  }
794 
795 #ifdef G4DEBUG_PATHFINDER
796  if( fVerboseLevel > 2 )
797  {
798  G4String VolumeName("None");
799  if( locatedVolume ) { VolumeName= locatedVolume->GetName(); }
800  G4cout << " Touchable History created at address " << touchHist
801  << " volume = " << locatedVolume << " name= " << VolumeName
802  << G4endl;
803  }
804 #endif
805 
806  return G4TouchableHandle(touchHist);
807 }
808 
809 G4double
811  G4double proposedStepLength )
812 {
813  std::vector<G4Navigator*>::iterator pNavigatorIter;
814  G4double safety= 0.0, step=0.0;
815  G4double minSafety= kInfinity, minStep;
816 
817  const G4int IdTransport= 0; // Id of Mass Navigator !!
818  G4int num=0;
819 
820 #ifdef G4DEBUG_PATHFINDER
821  if( fVerboseLevel > 2 )
822  {
823  G4cout << " G4PathFinder::DoNextLinearStep : entered " << G4endl;
824  G4cout << " Input field track= " << initialState << G4endl;
825  G4cout << " Requested step= " << proposedStepLength << G4endl;
826  }
827 #endif
828 
829  G4ThreeVector initialPosition= initialState.GetPosition();
830  G4ThreeVector initialDirection= initialState.GetMomentumDirection();
831 
832  G4ThreeVector OriginShift = initialPosition - fPreSafetyLocation;
833  G4double MagSqShift = OriginShift.mag2() ;
834  G4double MagShift; // Only given value if it larger than minimum safety
835 
836  // Potential optimisation using Maximum Value of safety!
837  // if( MagSqShift >= sqr(fPreSafetyMaxValue ) ){
838  // MagShift= kInfinity; // Not a useful value -- all will not use/ignore
839  // else
840  // MagShift= std::sqrt(MagSqShift) ;
841 
842  MagShift= std::sqrt(MagSqShift) ;
843 
844 #ifdef G4PATHFINDER_OPTIMISATION
845 
846  G4double fullSafety; // For all geometries, for prestep point
847 
848  if( MagSqShift >= sqr(fPreSafetyMinValue ) )
849  {
850  fullSafety = 0.0 ;
851  }
852  else
853  {
854  fullSafety = fPreSafetyMinValue - MagShift;
855  }
856  if( proposedStepLength < fullSafety )
857  {
858  // Move is smaller than all safeties
859  // -> so we do not have to move the safety center
860 
861  fPreStepCenterRenewed= false;
862 
863  for( num=0; num< fNoActiveNavigators; ++num )
864  {
866  safety = std::max( 0.0, fPreSafetyValues[num] - MagShift);
867  minSafety= std::min( safety, minSafety );
868  fCurrentPreStepSafety[num]= safety;
869  }
870  minStep= kInfinity;
871 
872 #ifdef G4DEBUG_PATHFINDER
873  if( fVerboseLevel > 2 )
874  {
875  G4cout << "G4PathFinder::DoNextLinearStep : Quick Stepping. " << G4endl
876  << " proposedStepLength " << proposedStepLength
877  << " < (full) safety = " << fullSafety
878  << " at " << initialPosition
879  << G4endl;
880  }
881 #endif
882  }
883  else
884 #endif // End of G4PATHFINDER_OPTIMISATION 1
885  {
886  // Move is larger than at least one of the safeties
887  // -> so we must move the safety center!
888 
889  fPreStepCenterRenewed= true;
890  pNavigatorIter= fpTransportManager-> GetActiveNavigatorsIterator();
891 
892  minStep= kInfinity; // Not proposedStepLength;
893 
894  for( num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num )
895  {
896  safety = std::max( 0.0, fPreSafetyValues[num] - MagShift);
897 
898 #ifdef G4PATHFINDER_OPTIMISATION
899  if( proposedStepLength <= safety ) // Should be just < safety ?
900  {
901  // The Step is guaranteed to be taken
902 
903  step= kInfinity; // ComputeStep Would return this
904 
905 #ifdef G4DEBUG_PATHFINDER
906  G4cout.precision(8);
907  G4cout << "PathFinder::ComputeStep> small proposed step = "
908  << proposedStepLength
909  << " <= safety = " << safety << " for nav " << num
910  << " Step fully taken. " << G4endl;
911 #endif
912  }
913  else
914 #endif // End of G4PATHFINDER_OPTIMISATION 2
915  {
916 #ifdef G4DEBUG_PATHFINDER
917  G4double previousSafety= safety;
918 #endif
919  step= (*pNavigatorIter)->ComputeStep( initialPosition,
920  initialDirection,
921  proposedStepLength,
922  safety );
923  minStep = std::min( step, minStep);
924 
925  // TODO: consider whether/how to reduce the proposed step
926  // to the latest minStep value - to reduce calculations
927 
928 #ifdef G4DEBUG_PATHFINDER
929  if( fVerboseLevel > 0)
930  {
931  G4cout.precision(8);
932  G4cout << "PathFinder::ComputeStep> long proposed step = "
933  << proposedStepLength
934  << " > safety = " << previousSafety
935  << " for nav " << num
936  << " . New safety = " << safety << " step= " << step
937  << G4endl;
938  }
939 #endif
940  }
941  fCurrentStepSize[num] = step;
942 
943  // Save safety value, must be done for all geometries "together"
944  // (even if not recomputed using call to ComputeStep)
945  // since they share the fPreSafetyLocation
946 
947  fPreSafetyValues[num]= safety;
948  fCurrentPreStepSafety[num]= safety;
949 
950  minSafety= std::min( safety, minSafety );
951 
952 #ifdef G4DEBUG_PATHFINDER
953  if( fVerboseLevel > 2 )
954  {
955  G4cout << "G4PathFinder::DoNextLinearStep : Navigator ["
956  << num << "] -- step size " << step << G4endl;
957  }
958 #endif
959  }
960 
961  // Only change these when safety is recalculated
962  // it is good/relevant only for safety calculations
963 
964  fPreSafetyLocation= initialPosition;
965  fPreSafetyMinValue= minSafety;
966  } // end of else for if( proposedStepLength <= fullSafety)
967 
968  // For use in Relocation, need PreStep point location, min-safety
969  //
970  fPreStepLocation= initialPosition;
971  fMinSafety_PreStepPt= minSafety;
972 
973  fMinStep= minStep;
974 
975  if( fMinStep == kInfinity )
976  {
977  minStep = proposedStepLength; // Use this below for endpoint !!
978  }
979  fTrueMinStep = minStep;
980 
981  // Set the EndState
982 
984 
985  fEndState= initialState;
986  endPosition= initialPosition + minStep * initialDirection ;
987 
988 #ifdef G4DEBUG_PATHFINDER
989  if( fVerboseLevel > 1 )
990  {
991  G4cout << "G4PathFinder::DoNextLinearStep : "
992  << " initialPosition = " << initialPosition
993  << " and endPosition = " << endPosition<< G4endl;
994  }
995 #endif
996 
997  fEndState.SetPosition( endPosition );
998  fEndState.SetProperTimeOfFlight( -1.000 ); // Not defined YET
999 
1000  if( fNoActiveNavigators == 1 )
1001  {
1002  G4bool transportLimited = (fMinStep!= kInfinity);
1003  fLimitTruth[IdTransport] = transportLimited;
1004  fLimitedStep[IdTransport] = transportLimited ? kUnique : kDoNot;
1005 
1006  // Set fNoGeometriesLimiting - as WhichLimited does
1007  fNoGeometriesLimiting = transportLimited ? 1 : 0;
1008  }
1009  else
1010  {
1011  WhichLimited();
1012  }
1013 
1014 #ifdef G4DEBUG_PATHFINDER
1015  if( fVerboseLevel > 2 )
1016  {
1017  G4cout << " G4PathFinder::DoNextLinearStep : exits returning "
1018  << minStep << G4endl;
1019  G4cout << " Endpoint values = " << fEndState << G4endl;
1020  G4cout << G4endl;
1021  }
1022 #endif
1023 
1024  return minStep;
1025 }
1026 
1028 {
1029  // Flag which processes limited the step
1030 
1031  G4int num=-1, last=-1;
1032  G4int noLimited=0;
1033  ELimited shared= kSharedOther;
1034 
1035  const G4int IdTransport= 0; // Id of Mass Navigator !!
1036 
1037  // Assume that [IdTransport] is Mass / Transport
1038  //
1039  G4bool transportLimited = (fCurrentStepSize[IdTransport] == fMinStep)
1040  && ( fMinStep!= kInfinity) ;
1041 
1042  if( transportLimited ) {
1043  shared= kSharedTransport;
1044  }
1045 
1046  for ( num= 0; num < fNoActiveNavigators; num++ )
1047  {
1048  G4bool limitedStep;
1049 
1050  G4double step= fCurrentStepSize[num];
1051 
1052  limitedStep = ( std::fabs(step - fMinStep) < kCarTolerance )
1053  && ( step != kInfinity);
1054 
1055  fLimitTruth[ num ] = limitedStep;
1056  if( limitedStep )
1057  {
1058  noLimited++;
1059  fLimitedStep[num] = shared;
1060  last= num;
1061  }
1062  else
1063  {
1064  fLimitedStep[num] = kDoNot;
1065  }
1066  }
1067  fNoGeometriesLimiting= noLimited; // Save # processes limiting step
1068 
1069  if( (last > -1) && (noLimited == 1 ) )
1070  {
1071  fLimitedStep[ last ] = kUnique;
1072  }
1073 
1074 #ifdef G4DEBUG_PATHFINDER
1075  if( fVerboseLevel > 1 )
1076  {
1077  PrintLimited(); // --> for tracing
1078  if( fVerboseLevel > 4 ) {
1079  G4cout << " G4PathFinder::WhichLimited - exiting. " << G4endl;
1080  }
1081  }
1082 #endif
1083 }
1084 
1086 {
1087  // Report results -- for checking
1088 
1089  G4cout << "G4PathFinder::PrintLimited reports: " ;
1090  G4cout << " Minimum step (true)= " << fTrueMinStep
1091  << " reported min = " << fMinStep
1092  << G4endl;
1093  if( (fCurrentStepNo <= 2) || (fVerboseLevel>=2) )
1094  {
1095  G4cout << std::setw(5) << " Step#" << " "
1096  << std::setw(5) << " NavId" << " "
1097  << std::setw(12) << " step-size " << " "
1098  << std::setw(12) << " raw-size " << " "
1099  << std::setw(12) << " pre-safety " << " "
1100  << std::setw(15) << " Limited / flag" << " "
1101  << std::setw(15) << " World " << " "
1102  << G4endl;
1103  }
1104  G4int num;
1105  for ( num= 0; num < fNoActiveNavigators; num++ )
1106  {
1107  G4double rawStep = fCurrentStepSize[num];
1108  G4double stepLen = fCurrentStepSize[num];
1109  if( stepLen > fTrueMinStep )
1110  {
1111  stepLen = fTrueMinStep; // did not limit (went as far as asked)
1112  }
1113  G4int oldPrec= G4cout.precision(9);
1114 
1115  G4cout << std::setw(5) << fCurrentStepNo << " "
1116  << std::setw(5) << num << " "
1117  << std::setw(12) << stepLen << " "
1118  << std::setw(12) << rawStep << " "
1119  << std::setw(12) << fCurrentPreStepSafety[num] << " "
1120  << std::setw(5) << (fLimitTruth[num] ? "YES" : " NO") << " ";
1121  G4String limitedStr= LimitedString(fLimitedStep[num]);
1122  G4cout << " " << std::setw(15) << limitedStr << " ";
1123  G4cout.precision(oldPrec);
1124 
1125  G4Navigator *pNav= GetNavigator( num );
1126  G4String WorldName( "Not-Set" );
1127  if (pNav)
1128  {
1129  G4VPhysicalVolume *pWorld= pNav->GetWorldVolume();
1130  if( pWorld )
1131  {
1132  WorldName = pWorld->GetName();
1133  }
1134  }
1135  G4cout << " " << WorldName ;
1136  G4cout << G4endl;
1137  }
1138 
1139  if( fVerboseLevel > 4 )
1140  {
1141  G4cout << " G4PathFinder::PrintLimited - exiting. " << G4endl;
1142  }
1143 }
1144 
1145 G4double
1147  G4double proposedStepLength,
1148  G4VPhysicalVolume* pCurrentPhysicalVolume )
1149 {
1150  const G4double toleratedRelativeError= 1.0e-10;
1151  G4double minStep= kInfinity, newSafety=0.0;
1152  G4int numNav;
1153  G4FieldTrack fieldTrack= initialState;
1154  G4ThreeVector startPoint= initialState.GetPosition();
1155 
1156 #ifdef G4DEBUG_PATHFINDER
1157  G4int prc= G4cout.precision(9);
1158  if( fVerboseLevel > 2 )
1159  {
1160  G4cout << " G4PathFinder::DoNextCurvedStep ****** " << G4endl;
1161  G4cout << " Initial value of field track is " << fieldTrack
1162  << " and proposed step= " << proposedStepLength << G4endl;
1163  }
1164 #endif
1165 
1166  fPreStepCenterRenewed= true; // Always update PreSafety with PreStep point
1167 
1168  if( fNoActiveNavigators > 1 )
1169  {
1170  // Calculate the safety values before making the step
1171 
1172  G4double minSafety= kInfinity, safety;
1173  for( numNav=0; numNav < fNoActiveNavigators; ++numNav )
1174  {
1175  safety= fpNavigator[numNav]->ComputeSafety( startPoint, false );
1176  fPreSafetyValues[numNav]= safety;
1177  fCurrentPreStepSafety[numNav]= safety;
1178  minSafety = std::min( safety, minSafety );
1179  }
1180 
1181  // Save safety value, related position
1182 
1183  fPreSafetyLocation= startPoint;
1184  fPreSafetyMinValue= minSafety;
1185  fPreStepLocation= startPoint;
1186  fMinSafety_PreStepPt= minSafety;
1187  }
1188 
1189  // Allow Propagator In Field to do the hard work, calling G4MultiNavigator
1190  //
1191  minStep= fpFieldPropagator->ComputeStep( fieldTrack,
1192  proposedStepLength,
1193  newSafety,
1194  pCurrentPhysicalVolume );
1195 
1196  // fieldTrack now contains the endpoint information
1197  //
1198  fEndState= fieldTrack;
1199  fMinStep= minStep;
1200  fTrueMinStep = std::min( minStep, proposedStepLength );
1201 
1202  if( fNoActiveNavigators== 1 )
1203  {
1204  // Update the 'PreSafety' sphere - as any ComputeStep was called
1205  // (must be done anyway in field)
1206 
1207  fPreSafetyValues[0]= newSafety;
1208  fPreSafetyLocation= startPoint;
1209  fPreSafetyMinValue= newSafety;
1210 
1211  // Update the current 'PreStep' point's values - mandatory
1212  //
1213  fCurrentPreStepSafety[0]= newSafety;
1214  fPreStepLocation= startPoint;
1215  fMinSafety_PreStepPt= newSafety;
1216  }
1217 
1218 #ifdef G4DEBUG_PATHFINDER
1219  if( fVerboseLevel > 2 )
1220  {
1221  G4cout << "G4PathFinder::DoNextCurvedStep : " << G4endl
1222  << " initialState = " << initialState << G4endl
1223  << " and endState = " << fEndState << G4endl;
1224  G4cout << "G4PathFinder::DoNextCurvedStep : "
1225  << " minStep = " << minStep
1226  << " proposedStepLength " << proposedStepLength
1227  << " safety = " << newSafety << G4endl;
1228  }
1229 #endif
1230  G4double currentStepSize; // = 0.0;
1231  if( minStep < proposedStepLength ) // if == , then a boundary found at end ??
1232  {
1233  // Recover the remaining information from MultiNavigator
1234  // especially regarding which Navigator limited the step
1235 
1236  G4int noLimited= 0; // No geometries limiting step
1237  for( numNav=0; numNav < fNoActiveNavigators; ++numNav )
1238  {
1239  G4double finalStep, lastPreSafety=0.0, minStepLast;
1240  ELimited didLimit;
1241  G4bool limited;
1242 
1243  finalStep= fpMultiNavigator->ObtainFinalStep( numNav, lastPreSafety,
1244  minStepLast, didLimit );
1245 
1246  // Calculate the step for this geometry, using the
1247  // final step (the only one which can differ.)
1248 
1249  currentStepSize = fTrueMinStep;
1250  G4double diffStep= 0.0;
1251  if( (minStepLast != kInfinity) )
1252  {
1253  diffStep = (finalStep-minStepLast);
1254  if ( std::abs(diffStep) <= toleratedRelativeError * finalStep )
1255  {
1256  diffStep = 0.0;
1257  }
1258  currentStepSize += diffStep;
1259  }
1260  fCurrentStepSize[numNav] = currentStepSize;
1261 
1262  // TODO: could refine the way to obtain safeties for > 1 geometries
1263  // - for pre step safety
1264  // notify MultiNavigator about new set of sub-steps
1265  // allow it to return this value in ObtainFinalStep
1266  // instead of lastPreSafety (or as well?)
1267  // - for final step start (available)
1268  // get final Step start from MultiNavigator
1269  // and corresponding safety values
1270  // and/or ALSO calculate ComputeSafety at endpoint
1271  // endSafety= fpNavigator[numNav]->ComputeSafety( endPoint );
1272 
1273  fLimitedStep[numNav] = didLimit;
1274  fLimitTruth[numNav] = limited = (didLimit != kDoNot );
1275  if( limited ) { noLimited++; }
1276 
1277 #ifdef G4DEBUG_PATHFINDER
1278  G4bool StepError= (currentStepSize < 0)
1279  || ( (minStepLast != kInfinity) && (diffStep < 0) ) ;
1280  if( StepError || (fVerboseLevel > 2) )
1281  {
1282  G4String limitedString= LimitedString( fLimitedStep[numNav] );
1283 
1284  G4cout << " G4PathFinder::ComputeStep. Geometry " << numNav
1285  << " step= " << fCurrentStepSize[numNav]
1286  << " from final-step= " << finalStep
1287  << " fTrueMinStep= " << fTrueMinStep
1288  << " minStepLast= " << minStepLast
1289  << " limited = " << (fLimitTruth[numNav] ? "YES" : " NO")
1290  << " ";
1291  G4cout << " status = " << limitedString << " #= " << didLimit
1292  << G4endl;
1293 
1294  if( StepError )
1295  {
1296  std::ostringstream message;
1297  message << "Incorrect calculation of step size for one navigator"
1298  << G4endl
1299  << " currentStepSize = " << currentStepSize
1300  << ", diffStep= " << diffStep << G4endl
1301  << "ERROR in computing step size for this navigator.";
1302  G4Exception("G4PathFinder::DoNextCurvedStep",
1303  "GeomNav0003", FatalException, message);
1304  }
1305  }
1306 #endif
1307  } // for num Navigators
1308 
1309  fNoGeometriesLimiting= noLimited; // Save # processes limiting step
1310  }
1311  else if ( (minStep == proposedStepLength)
1312  || (minStep == kInfinity)
1313  || ( std::abs(minStep-proposedStepLength)
1314  < toleratedRelativeError * proposedStepLength ) )
1315  {
1316  // In case the step was not limited, use default responses
1317  // --> all Navigators
1318  // Also avoid problems in case of PathFinder using safety to optimise
1319  // - it is possible that the Navigators were not called
1320  // if the safety was already satisfactory.
1321  // (In that case calling ObtainFinalStep gives invalid results.)
1322 
1323  currentStepSize= minStep;
1324  for( numNav=0; numNav < fNoActiveNavigators; ++numNav )
1325  {
1326  fCurrentStepSize[numNav] = minStep;
1327  // Safety for endpoint ?? // Can eventuall improve it -- see TODO above
1328  fLimitedStep[numNav] = kDoNot;
1329  fLimitTruth[numNav] = false;
1330  }
1331  fNoGeometriesLimiting= 0; // Save # processes limiting step
1332  }
1333  else // (minStep > proposedStepLength) and not (minStep == kInfinity)
1334  {
1335  std::ostringstream message;
1336  message << "Incorrect calculation of step size for one navigator." << G4endl
1337  << " currentStepSize = " << minStep << " is larger than "
1338  << " proposed StepSize = " << proposedStepLength << ".";
1339  G4Exception("G4PathFinder::DoNextCurvedStep()",
1340  "GeomNav0003", FatalException, message);
1341  }
1342 
1343 #ifdef G4DEBUG_PATHFINDER
1344  if( fVerboseLevel > 2 )
1345  {
1346  G4cout << " Exiting G4PathFinder::DoNextCurvedStep " << G4endl;
1347  PrintLimited();
1348  }
1349  G4cout.precision(prc);
1350 #endif
1351 
1352  return minStep;
1353 }
1354 
1356 {
1357  static G4String StrDoNot("DoNot"),
1358  StrUnique("Unique"),
1359  StrUndefined("Undefined"),
1360  StrSharedTransport("SharedTransport"),
1361  StrSharedOther("SharedOther");
1362 
1363  G4String* limitedStr;
1364  switch ( lim )
1365  {
1366  case kDoNot: limitedStr= &StrDoNot; break;
1367  case kUnique: limitedStr = &StrUnique; break;
1368  case kSharedTransport: limitedStr= &StrSharedTransport; break;
1369  case kSharedOther: limitedStr = &StrSharedOther; break;
1370  default: limitedStr = &StrUndefined; break;
1371  }
1372  return *limitedStr;
1373 }
1374 
1376 {
1379  for( G4int nav=0; nav < fNoActiveNavigators; ++nav )
1380  {
1382  }
1383 }
void PrepareNewTrack(const G4ThreeVector &position, const G4ThreeVector &direction, G4VPhysicalVolume *massStartVol=0)
G4bool fLimitTruth[fMaxNav]
G4String & LimitedString(ELimited lim)
static G4PathFinder * GetInstance()
Definition: G4PathFinder.cc:57
ELimited fLimitedStep[fMaxNav]
void SetPosition(G4ThreeVector nPos)
G4MultiNavigator * fpMultiNavigator
G4double fNewSafetyComputed[fMaxNav]
G4TransportationManager * fpTransportManager
void Locate(const G4ThreeVector &position, const G4ThreeVector &direction, G4bool relativeSearch=true)
void SetNavigatorForPropagating(G4Navigator *SimpleOrMultiNavigator)
G4SafetyHelper * GetSafetyHelper() const
void InitialiseHelper()
void EnableParallelNavigation(G4bool enableChoice=true)
static const G4int fMaxNav
G4double fMinSafety_PreStepPt
static const G4double kInfinity
Definition: geomdefs.hh:42
CLHEP::Hep3Vector G4ThreeVector
void ReLocate(const G4ThreeVector &position)
ELimited
G4double GetSurfaceTolerance() const
G4Navigator * GetNavigatorForTracking() const
G4int fVerboseLevel
G4VPhysicalVolume * fLocatedVolume[fMaxNav]
G4TouchableHandle CreateTouchableHandle(G4int navId) const
void WhichLimited()
G4int fCurrentStepNo
G4ThreeVector fPreStepLocation
#define G4ThreadLocal
Definition: tls.hh:52
G4double fMinSafety_atSafLocation
G4FieldTrack fEndState
int G4int
Definition: G4Types.hh:78
void MovePoint()
G4double DoNextCurvedStep(const G4FieldTrack &FieldTrack, G4double proposedStepLength, G4VPhysicalVolume *pCurrentPhysVolume)
G4ReferenceCountedHandle< G4VTouchable > G4TouchableHandle
G4int fNoActiveNavigators
G4double fCurrentPreStepSafety[fMaxNav]
static G4ThreadLocal G4PathFinder * fpPathFinder
G4ThreeVector GetPosition() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetName() const
G4double GetRadialTolerance() const
bool G4bool
Definition: G4Types.hh:79
G4Navigator * fpNavigator[fMaxNav]
void EnableParallelNavigation(G4bool parallel)
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
G4double DoNextLinearStep(const G4FieldTrack &FieldTrack, G4double proposedStepLength)
G4bool fRelocatedPoint
G4double ComputeSafety(const G4ThreeVector &globalPoint)
G4FieldManager * FindAndSetFieldManager(G4VPhysicalVolume *pCurrentPhysVol)
G4double kCarTolerance
void PushPostSafetyToPreSafety()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static G4TransportationManager * GetTransportationManager()
G4double GetCharge() const
G4double ComputeStep(G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4double &pNewSafety, G4VPhysicalVolume *pPhysVol=0)
G4double fPreSafetyValues[fMaxNav]
G4PropagatorInField * fpFieldPropagator
G4double fCurrentStepSize[fMaxNav]
G4ThreeVector fLastLocatedPosition
T max(const T t1, const T t2)
brief Return the largest of the two arguments
std::vector< G4Navigator * >::iterator GetActiveNavigatorsIterator()
G4double fTrueMinStep
int position
Definition: filter.cc:7
G4double ObtainFinalStep(G4int navigatorId, G4double &pNewSafety, G4double &minStepLast, ELimited &limitedStep)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4ThreeVector fPreSafetyLocation
#define G4endl
Definition: G4ios.hh:61
void UpdateYourself(G4VPhysicalVolume *pPhysVol, const G4NavigationHistory *history=0)
static char * endPosition
const G4NavigationHistory * GetHistory() const
G4bool fNewTrack
T sqr(const T &x)
Definition: templates.hh:145
double G4double
Definition: G4Types.hh:76
void PrintLimited()
const G4Field * GetDetectorField() const
virtual G4double ComputeSafety(const G4ThreeVector &globalpoint, const G4double pProposedMaxLength=DBL_MAX, const G4bool keepState=true)
G4bool fPreStepCenterRenewed
G4PropagatorInField * GetPropagatorInField() const
G4bool fFieldExertedForce
void ReportMove(const G4ThreeVector &OldV, const G4ThreeVector &NewV, const G4String &Quantity) const
static const double mm
Definition: G4SIunits.hh:102
G4int fNoGeometriesLimiting
G4VPhysicalVolume * GetWorldVolume() const
G4Navigator * GetNavigator(G4int n) const
G4ThreeVector fSafetyLocation
G4double fPreSafetyMinValue
G4double fMinStep
static G4GeometryTolerance * GetInstance()
G4GLOB_DLL std::ostream G4cerr
G4ThreeVector GetMomentumDirection() const
void SetProperTimeOfFlight(G4double tofProper)