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