Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/geometry/navigation/src/G4PathFinder.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 ]

Diff markup

Differences between /geometry/navigation/src/G4PathFinder.cc (Version 11.3.0) and /geometry/navigation/src/G4PathFinder.cc (Version 9.5)


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