Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/dna/management/src/G4ITPathFinder.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

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