Geant4  10.02.p02
G4ITNavigator2.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4ITNavigator2.cc 79045 2014-02-12 13:57:26Z matkara $
27 //
28 // class G4ITNavigator2 Implementation
29 //
30 // Original author: Paul Kent, July 95/96
31 //
32 // G4ITNavigator2 is a duplicate version of G4Navigator starting from Geant4.9.5
33 // initially written by Paul Kent and colleagues.
34 // The only difference resides in the way the information is saved and managed
35 //
36 // History:
37 // - Created. Paul Kent, Jul 95/96
38 // - Zero step protections J.A. / G.C., Nov 2004
39 // - Added check mode G. Cosmo, Mar 2004
40 // - Made Navigator Abstract G. Cosmo, Nov 2003
41 // - G4ITNavigator2 created M.K., Nov 2012
42 // --------------------------------------------------------------------
43 
44 #include <iomanip>
45 
46 //#include "G4ITNavigator2.hh"
47 #include "G4ITNavigator.hh"
48 #include "G4ios.hh"
49 #include "G4SystemOfUnits.hh"
50 #include "G4GeometryTolerance.hh"
51 #include "G4VPhysicalVolume.hh"
52 
53 //#define G4DEBUG_NAVIGATION 1
54 #include "G4VoxelSafety.hh"
55 #include "Randomize.hh"
56 #include "G4VoxelLimits.hh"
57 
58 #define fHistory fpNavigatorState->fHistory
59 #define fLastTriedStepComputation fpNavigatorState->fLastTriedStepComputation
60 #define fEnteredDaughter fpNavigatorState->fEnteredDaughter
61 #define fExitedMother fpNavigatorState->fExitedMother
62 #define fWasLimitedByGeometry fpNavigatorState->fWasLimitedByGeometry
63 #define fStepEndPoint fpNavigatorState->fStepEndPoint
64 #define fLastStepEndPointLocal fpNavigatorState->fLastStepEndPointLocal
65 #define fPushed fpNavigatorState->fPushed
66 #define fLastTriedStepComputation fpNavigatorState->fLastTriedStepComputation
67 #define fEntering fpNavigatorState->fEntering
68 #define fExiting fpNavigatorState->fExiting
69 #define fBlockedPhysicalVolume fpNavigatorState->fBlockedPhysicalVolume
70 #define fBlockedReplicaNo fpNavigatorState->fBlockedReplicaNo
71 #define fLastLocatedPointLocal fpNavigatorState->fLastLocatedPointLocal
72 #define fLocatedOutsideWorld fpNavigatorState->fLocatedOutsideWorld
73 #define fValidExitNormal fpNavigatorState->fValidExitNormal
74 #define fExitNormal fpNavigatorState->fExitNormal
75 #define fGrandMotherExitNormal fpNavigatorState->fGrandMotherExitNormal
76 #define fChangedGrandMotherRefFrame fpNavigatorState->fChangedGrandMotherRefFrame
77 #define fExitNormalGlobalFrame fpNavigatorState->fExitNormalGlobalFrame
78 #define fCalculatedExitNormal fpNavigatorState->fCalculatedExitNormal
79 #define fLastStepWasZero fpNavigatorState->fLastStepWasZero
80 #define fLocatedOnEdge fpNavigatorState->fLocatedOnEdge
81 #define fNumberZeroSteps fpNavigatorState->fNumberZeroSteps
82 #define fPreviousSftOrigin fpNavigatorState->fPreviousSftOrigin
83 #define fPreviousSafety fpNavigatorState->fPreviousSafety
84 
85 //------------------------------------------------------------------------------
86 
87 namespace BoundingBox
88 {
90 {
93 };
94 }
95 
96 // ********************************************************************
97 // Constructor
98 // ********************************************************************
99 //
101  : fVerbose(0),
102  fTopPhysical(0), fCheck(false),
103  fWarnPush(true)
104 {
105  fActive= false;
110 
111  fpNavigatorState = 0;
112 // fSaveState = 0;
114 
115  // this->SetVerboseLevel(3);
116  // this->CheckMode(true);
117 }
118 
119 // ********************************************************************
120 // Destructor
121 // ********************************************************************
122 //
124 { delete fpVoxelSafety; }
125 
126 // ********************************************************************
127 // ResetHierarchyAndLocate
128 // ********************************************************************
129 //
132  const G4ThreeVector &direction,
133  const G4TouchableHistory &h)
134 {
135 // ResetState();
136  fHistory = *h.GetHistory();
137  SetupHierarchy();
138  fLastTriedStepComputation= false; // Redundant, but best
139  return LocateGlobalPointAndSetup(p, &direction, true, false);
140 }
141 
142 // ********************************************************************
143 // LocateGlobalPointAndSetup
144 //
145 // Locate the point in the hierarchy return 0 if outside
146 // The direction is required
147 // - if on an edge shared by more than two surfaces
148 // (to resolve likely looping in tracking)
149 // - at initial location of a particle
150 // (to resolve potential ambiguity at boundary)
151 //
152 // Flags on exit: (comments to be completed)
153 // fEntering - True if entering `daughter' volume (or replica)
154 // whether daughter of last mother directly
155 // or daughter of that volume's ancestor.
156 // ********************************************************************
157 //
160  const G4ThreeVector* pGlobalDirection,
161  const G4bool relativeSearch,
162  const G4bool ignoreDirection )
163 {
165  G4bool notKnownContained=true, noResult;
166  G4VPhysicalVolume *targetPhysical;
167  G4LogicalVolume *targetLogical;
168  G4VSolid *targetSolid=0;
169  G4ThreeVector localPoint, globalDirection;
170  EInside insideCode;
171 
172  G4bool considerDirection = (!ignoreDirection) || fLocatedOnEdge;
173 
175  fChangedGrandMotherRefFrame= false; // For local exit normal
176 
177  if( considerDirection && pGlobalDirection != 0 )
178  {
179  globalDirection=*pGlobalDirection;
180  }
181 
182 
183 #ifdef G4VERBOSE
184  if( fVerbose > 2 )
185  {
186  static int nCalls = 0;
187  G4int oldcoutPrec = G4cout.precision(8);
188  G4cout << "*** G4ITNavigator2::LocateGlobalPointAndSetup: ***" << G4endl;
189  G4cout << " Called with arguments: " << G4endl
190  << " Globalpoint = " << globalPoint << G4endl
191  << " RelativeSearch = " << relativeSearch << G4endl;
192  if( fVerbose == 4 )
193  {
194  G4cout << " ----- Upon entering:" << G4endl;
195  PrintState();
196  }
197  G4cout.precision(oldcoutPrec);
198  nCalls++;
199 // if(nCalls == 2) abort();
200  }
201 #endif
202 
203  G4int noLevelsExited=0 ;
204  G4int noLevelsEntered= 0;
205 
206  if ( !relativeSearch )
207  {
209  }
210  else
211  {
212  if ( fWasLimitedByGeometry )
213  {
214  fWasLimitedByGeometry = false;
215  fEnteredDaughter = fEntering; // Remember
216  fExitedMother = fExiting; // Remember
217  if ( fExiting )
218  {
219  noLevelsExited++; // count this first level entered too
220 
221  if ( fHistory.GetDepth() )
222  {
223  fBlockedPhysicalVolume = fHistory.GetTopVolume();
224  fBlockedReplicaNo = fHistory.GetTopReplicaNo();
225  fHistory.BackLevel();
226  }
227  else
228  {
229  fLastLocatedPointLocal = localPoint;
230  fLocatedOutsideWorld = true;
231  fBlockedPhysicalVolume = 0; // to be sure
232  fBlockedReplicaNo = -1;
233  fEntering = false; // No longer
234  fEnteredDaughter = false;
235  fExitedMother = true; // ??
236 
237  return 0; // Have exited world volume
238  }
239  // A fix for the case where a volume is "entered" at an edge
240  // and a coincident surface exists outside it.
241  // - This stops it from exiting further volumes and cycling
242  // - However ReplicaNavigator treats this case itself
243  //
244  // assert( fBlockedPhysicalVolume!=0 );
245 
246  // Expect to be on edge => on surface
247  //
249  {
250  fExiting= false;
251  // Consider effect on Exit Normal !?
252  }
253  }
254  else
255  if ( fEntering )
256  {
257  // assert( fBlockedPhysicalVolume!=0 );
258 
259  noLevelsEntered++; // count the first level entered too
260 
262  {
263  case kNormal:
265  fBlockedPhysicalVolume->GetCopyNo());
266  break;
267  case kReplica:
273  break;
274  case kParameterised:
275  if( fBlockedPhysicalVolume->GetRegularStructureId() == 0 )
276  {
277  G4VSolid *pSolid;
278  G4VPVParameterisation *pParam;
279  G4TouchableHistory parentTouchable( fHistory );
280  pParam = fBlockedPhysicalVolume->GetParameterisation();
281  pSolid = pParam->ComputeSolid(fBlockedReplicaNo,
283  pSolid->ComputeDimensions(pParam, fBlockedReplicaNo,
290  //
291  // Set the correct solid and material in Logical Volume
292  //
293  G4LogicalVolume *pLogical;
294  pLogical = fBlockedPhysicalVolume->GetLogicalVolume();
295  pLogical->SetSolid( pSolid );
296  pLogical->UpdateMaterial(pParam ->
297  ComputeMaterial(fBlockedReplicaNo,
299  &parentTouchable));
300  }
301  break;
302  }
303  fEntering = false;
305  localPoint = fHistory.GetTopTransform().TransformPoint(globalPoint);
306  notKnownContained = false;
307  }
308  }
309  else
310  {
312  fEntering = false;
313  fEnteredDaughter = false; // Full Step was not taken, did not enter
314  fExiting = false;
315  fExitedMother = false; // Full Step was not taken, did not exit
316  }
317  }
318  //
319  // Search from top of history up through geometry until
320  // containing volume found:
321  // If on
322  // o OUTSIDE - Back up level, not/no longer exiting volumes
323  // o SURFACE and EXITING - Back up level, setting new blocking no.s
324  // else
325  // o containing volume found
326  //
327 
328  while (notKnownContained)
329  {
330  if ( fHistory.GetTopVolumeType()!=kReplica )
331  {
332  targetSolid = fHistory.GetTopVolume()->GetLogicalVolume()->GetSolid();
333  localPoint = fHistory.GetTopTransform().TransformPoint(globalPoint);
334  insideCode = targetSolid->Inside(localPoint);
335 #ifdef G4VERBOSE
336  if(( fVerbose == 1 ) && ( fCheck ))
337  {
338  G4String solidResponse = "-kInside-";
339  if (insideCode == kOutside)
340  solidResponse = "-kOutside-";
341  else if (insideCode == kSurface)
342  solidResponse = "-kSurface-";
343  G4cout << "*** G4ITNavigator2::LocateGlobalPointAndSetup(): ***" << G4endl
344  << " Invoked Inside() for solid: " << targetSolid->GetName()
345  << ". Solid replied: " << solidResponse << G4endl
346  << " For local point p: " << localPoint << G4endl;
347  }
348 #endif
349  }
350  else
351  {
352  insideCode = freplicaNav.BackLocate(fHistory, globalPoint, localPoint,
353  fExiting, notKnownContained);
354  // !CARE! if notKnownContained returns false then the point is within
355  // the containing placement volume of the replica(s). If insidecode
356  // will result in the history being backed up one level, then the
357  // local point returned is the point in the system of this new level
358  }
359 
360 
361  if ( insideCode==kOutside )
362  {
363  noLevelsExited++;
364  if ( fHistory.GetDepth() )
365  {
366  fBlockedPhysicalVolume = fHistory.GetTopVolume();
367  fBlockedReplicaNo = fHistory.GetTopReplicaNo();
368  fHistory.BackLevel();
369  fExiting = false;
370 
371  if( noLevelsExited > 1 )
372  {
373  // The first transformation was done by the sub-navigator
374  //
375  const G4RotationMatrix* mRot = fBlockedPhysicalVolume->GetRotation();
376  if( mRot )
377  {
378  fGrandMotherExitNormal *= (*mRot).inverse();
380  }
381  }
382  }
383  else
384  {
385  fLastLocatedPointLocal = localPoint;
386  fLocatedOutsideWorld = true;
387  // No extra transformation for ExitNormal - is in frame of Top Volume
388  return 0; // Have exited world volume
389  }
390  }
391  else
392  if ( insideCode==kSurface )
393  {
394  G4bool isExiting = fExiting;
395  if( (!fExiting)&&considerDirection )
396  {
397  // Figure out whether we are exiting this level's volume
398  // by using the direction
399  //
400  G4bool directionExiting = false;
401  G4ThreeVector localDirection =
402  fHistory.GetTopTransform().TransformAxis(globalDirection);
403 
404  // Make sure localPoint in correct reference frame
405  // ( Was it already correct ? How ? )
406  //
407  localPoint= fHistory.GetTopTransform().TransformPoint(globalPoint);
408  if ( fHistory.GetTopVolumeType()!=kReplica )
409  {
410  G4ThreeVector normal = targetSolid->SurfaceNormal(localPoint);
411  directionExiting = normal.dot(localDirection) > 0.0;
412  isExiting = isExiting || directionExiting;
413  }
414  }
415  if( isExiting )
416  {
417  noLevelsExited++;
418  if ( fHistory.GetDepth() )
419  {
420  fBlockedPhysicalVolume = fHistory.GetTopVolume();
421  fBlockedReplicaNo = fHistory.GetTopReplicaNo();
422  fHistory.BackLevel();
423  //
424  // Still on surface but exited volume not necessarily convex
425  //
426  fValidExitNormal = false;
427 
428  if( noLevelsExited > 1 )
429  {
430  // The first transformation was done by the sub-navigator
431  //
432  const G4RotationMatrix* mRot =
433  fBlockedPhysicalVolume->GetRotation();
434  if( mRot )
435  {
436  fGrandMotherExitNormal *= (*mRot).inverse();
438  }
439  }
440  }
441  else
442  {
443  fLastLocatedPointLocal = localPoint;
444  fLocatedOutsideWorld = true;
445  // No extra transformation for ExitNormal, is in frame of Top Vol
446  return 0; // Have exited world volume
447  }
448  }
449  else
450  {
451  notKnownContained=false;
452  }
453  }
454  else
455  {
456  notKnownContained=false;
457  }
458  } // END while (notKnownContained)
459  //
460  // Search downwards until deepest containing volume found,
461  // blocking fBlockedPhysicalVolume/BlockedReplicaNum
462  //
463  // 3 Cases:
464  //
465  // o Parameterised daughters
466  // =>Must be one G4PVParameterised daughter & voxels
467  // o Positioned daughters & voxels
468  // o Positioned daughters & no voxels
469 
470  noResult = true; // noResult should be renamed to
471  // something like enteredLevel, as that is its meaning.
472  do
473  {
474  // Determine `type' of current mother volume
475  //
476  targetPhysical = fHistory.GetTopVolume();
477  if (!targetPhysical) { break; }
478  targetLogical = targetPhysical->GetLogicalVolume();
479  switch( CharacteriseDaughters(targetLogical) )
480  {
481  case kNormal:
482  if ( targetLogical->GetVoxelHeader() ) // use optimised navigation
483  {
484  noResult = fvoxelNav.LevelLocate(fHistory,
487  globalPoint,
488  pGlobalDirection,
489  considerDirection,
490  localPoint);
491  }
492  else // do not use optimised navigation
493  {
494  noResult = fnormalNav.LevelLocate(fHistory,
497  globalPoint,
498  pGlobalDirection,
499  considerDirection,
500  localPoint);
501  }
502  break;
503  case kReplica:
504  noResult = freplicaNav.LevelLocate(fHistory,
507  globalPoint,
508  pGlobalDirection,
509  considerDirection,
510  localPoint);
511  break;
512  case kParameterised:
513  if( GetDaughtersRegularStructureId(targetLogical) != 1 )
514  {
515  noResult = fparamNav.LevelLocate(fHistory,
518  globalPoint,
519  pGlobalDirection,
520  considerDirection,
521  localPoint);
522  }
523  else // Regular structure
524  {
525  noResult = fregularNav.LevelLocate(fHistory,
528  globalPoint,
529  pGlobalDirection,
530  considerDirection,
531  localPoint);
532  }
533  break;
534  }
535 
536  // LevelLocate returns true if it finds a daughter volume
537  // in which globalPoint is inside (or on the surface).
538 
539  if ( noResult )
540  {
541  noLevelsEntered++;
542 
543  // Entering a daughter after ascending
544  //
545  // The blocked volume is no longer valid - it was for another level
546  //
548  fBlockedReplicaNo = -1;
549 
550  // fEntering should be false -- else blockedVolume is assumed good.
551  // fEnteredDaughter is used for ExitNormal
552  //
553  fEntering = false;
554  fEnteredDaughter = true;
555 
556  if( fExitedMother )
557  {
558  G4VPhysicalVolume* enteredPhysical = fHistory.GetTopVolume();
559  const G4RotationMatrix* mRot = enteredPhysical->GetRotation();
560  if( mRot )
561  {
562  // Go deeper, i.e. move 'down' in the hierarchy
563  // Apply direct rotation, not inverse
564  //
565  fGrandMotherExitNormal *= (*mRot);
567  }
568  }
569 
570 #ifdef G4DEBUG_NAVIGATION
571  if( fVerbose > 2 )
572  {
573  G4VPhysicalVolume* enteredPhysical = fHistory.GetTopVolume();
574  G4cout << "*** G4ITNavigator2::LocateGlobalPointAndSetup() ***" << G4endl;
575  G4cout << " Entering volume: " << enteredPhysical->GetName()
576  << G4endl;
577  }
578 #endif
579  }
580  } while (noResult);
581 
582  fLastLocatedPointLocal = localPoint;
583 
584 #ifdef G4VERBOSE
585  if( fVerbose >= 4 )
586  {
587  G4int oldcoutPrec = G4cout.precision(8);
588  G4String curPhysVol_Name("None");
589  if (targetPhysical) { curPhysVol_Name = targetPhysical->GetName(); }
590  G4cout << " Return value = new volume = " << curPhysVol_Name << G4endl;
591  G4cout << " ----- Upon exiting:" << G4endl;
592  PrintState();
593  if( fVerbose >= 5 )
594  {
595  G4cout << "Upon exiting LocateGlobalPointAndSetup():" << G4endl;
596  G4cout << " History = " << G4endl << fHistory << G4endl << G4endl;
597  }
598  G4cout.precision(oldcoutPrec);
599  }
600 #endif
601 
602  fLocatedOutsideWorld= false;
603 
604  return targetPhysical;
605 }
606 
607 // ********************************************************************
608 // LocateGlobalPointWithinVolume
609 //
610 // -> the state information of this Navigator and its subNavigators
611 // is updated in order to start the next step at pGlobalpoint
612 // -> no check is performed whether pGlobalpoint is inside the
613 // original volume (this must be the case).
614 //
615 // Note: a direction could be added to the arguments, to aid in future
616 // optional checking (via the old code below, flagged by OLD_LOCATE).
617 // [ This would be done only in verbose mode ]
618 // ********************************************************************
619 //
620 void
622 {
624 
625 #ifdef G4DEBUG_NAVIGATION
626  // Check: Either step was not limited by a boundary
627  // or else the full step is no longer being taken
628  assert( !fWasLimitedByGeometry );
629 #endif
630 
633  fChangedGrandMotherRefFrame= false; // Frame for Exit Normal
634 
635 #ifdef G4DEBUG_NAVIGATION
636  if( fVerbose > 2 )
637  {
638  G4cout << "Entering LocateGlobalWithinVolume(): History = " << G4endl;
639  G4cout << fHistory << G4endl;
640  }
641 #endif
642 
643  // For the case of Voxel (or Parameterised) volume the respective
644  // Navigator must be messaged to update its voxel information etc
645 
646  // Update the state of the Sub Navigators
647  // - in particular any voxel information they store/cache
648  //
649  G4VPhysicalVolume* motherPhysical = fHistory.GetTopVolume();
650  G4LogicalVolume* motherLogical = motherPhysical->GetLogicalVolume();
651  G4SmartVoxelHeader* pVoxelHeader = motherLogical->GetVoxelHeader();
652 
653  if ( fHistory.GetTopVolumeType()!=kReplica )
654  {
655  switch( CharacteriseDaughters(motherLogical) )
656  {
657  case kNormal:
658  if ( pVoxelHeader )
659  {
661  }
662  break;
663  case kParameterised:
664  if( GetDaughtersRegularStructureId(motherLogical) != 1 )
665  {
666  // Resets state & returns voxel node
667  //
669  }
670  break;
671  case kReplica:
672  G4Exception("G4ITNavigator2::LocateGlobalPointWithinVolume()",
673  "GeomNav0001", FatalException,
674  "Not applicable for replicated volumes.");
675  break;
676  }
677  }
678 
679  // Reset the state variables
680  // - which would have been affected
681  // by the 'equivalent' call to LocateGlobalPointAndSetup
682  // - who's values have been invalidated by the 'move'.
683  //
685  fBlockedReplicaNo = -1;
686  fEntering = false;
687  fEnteredDaughter = false; // Boundary not encountered, did not enter
688  fExiting = false;
689  fExitedMother = false; // Boundary not encountered, did not exit
690 }
691 
692 // !>
693 
695 {
696  if(fpNavigatorState == 0)
697  {
698  G4ExceptionDescription exceptionDescription;
699  exceptionDescription << "The navigator state is NULL. ";
700  exceptionDescription << "Either NewNavigatorStateAndLocate was not called ";
701  exceptionDescription << "or the provided navigator state was already NULL.";
702 
703  G4Exception("G4ITNavigator::CheckNavigatorStateIsValid",
704  "NavigatorStateNotValid",FatalException,exceptionDescription);
705  return;
706  }
707 }
708 
710 {
711  return fpNavigatorState;
712 }
713 
715 {
716  fpNavigatorState = (G4NavigatorState*) navState;
717  if(fpNavigatorState) SetupHierarchy(); // ?
718 // fpSaveState = (G4SaveNavigatorState*) navState;
719 // if(navState) RestoreSavedState();
720 }
721 
723 {
724  fpNavigatorState = 0;
725 }
726 
728 {
730  if(fTopPhysical == 0)
731  {
732  G4ExceptionDescription exceptionDescription;
733  exceptionDescription << "No World Volume";
734 
735  G4Exception("G4ITNavigator::NewNavigatorState",
736  "NoWorldVolume",FatalException,exceptionDescription);
737  return;
738  }
739 
740  fHistory.SetFirstEntry(fTopPhysical );
741  SetupHierarchy();
742 }
743 
745 {
747  if(fTopPhysical == 0)
748  {
749  G4ExceptionDescription exceptionDescription;
750  exceptionDescription << "No World Volume";
751 
752  G4Exception("G4ITNavigator::NewNavigatorState",
753  "NoWorldVolume",FatalException,exceptionDescription);
754  return;
755  }
756 
757  fHistory = *h.GetHistory();
758  fLastTriedStepComputation= false; // Redundant, but best
759  SetupHierarchy();
760 }
761 
763  const G4ThreeVector &direction)
764 {
766 
767  if(fTopPhysical == 0)
768  {
769  G4ExceptionDescription exceptionDescription;
770  exceptionDescription << "No World Volume";
771 
772  G4Exception("G4ITNavigator::NewNavigatorStateAndLocate",
773  "NoWorldVolume",FatalException,exceptionDescription);
774  return 0;
775  }
776 
777  fHistory.SetFirstEntry(fTopPhysical );
778  SetupHierarchy();
779  return LocateGlobalPointAndSetup(p, &direction, false, false);
780 }
781 
782 // ********************************************************************
783 // SetSavedState
784 //
785 // Save the state, in case this is a parasitic call
786 // Save fValidExitNormal, fExitNormal, fExiting, fEntering,
787 // fBlockedPhysicalVolume, fBlockedReplicaNo, fLastStepWasZero;
788 // ********************************************************************
789 //
790 /*
791 void G4ITNavigator2::SetSavedState()
792 {
793  fSaveState = *fpNavigatorState;
794 }
795 */
796 
797 /*
798 void G4ITNavigator2::SetSavedState()
799 {
800  // !>
801  // This check can be avoid if instead, at every first step of a track,
802  // the IT tracking uses NewNavigatorSate
803  // The normal tracking would just call once NewNavigatorState() before tracking
804 
805 // if(fpSaveState == 0)
806 // fpSaveState = new G4SaveNavigatorState;
807  // <!
808 
809  // fSaveExitNormal = fExitNormal;
810  fpSaveState->sExitNormal = fExitNormal;
811  fpSaveState->sValidExitNormal = fValidExitNormal;
812  fpSaveState->sExiting = fExiting;
813  fpSaveState->sEntering = fEntering;
814 
815  fpSaveState->spBlockedPhysicalVolume = fBlockedPhysicalVolume;
816  fpSaveState->sBlockedReplicaNo = fBlockedReplicaNo,
817 
818  fpSaveState->sLastStepWasZero = fLastStepWasZero;
819 
820  // !>
821  fpSaveState->sPreviousSftOrigin = fPreviousSftOrigin;
822  fpSaveState->sPreviousSafety = fPreviousSafety;
823  fpSaveState->sNumberZeroSteps = fNumberZeroSteps;
824  fpSaveState->sLocatedOnEdge = fLocatedOnEdge;
825  fpSaveState->sWasLimitedByGeometry= fWasLimitedByGeometry;
826  fpSaveState->sPushed=fPushed;
827  fpSaveState->sNumberZeroSteps=fNumberZeroSteps;
828  fpSaveState->sEnteredDaughter = fEnteredDaughter;
829  fpSaveState->sExitedMother = fExitedMother;
830 
831  fpSaveState->sLastLocatedPointLocal = fLastLocatedPointLocal;
832  fpSaveState->sLocatedOutsideWorld = fLocatedOutsideWorld;
833  // <!
834 }
835 */
836 // ********************************************************************
837 // RestoreSavedState
838 //
839 // Restore the state (in Compute Step), in case this is a parasitic call
840 // ********************************************************************
841 //
842 /*
843 void G4ITNavigator2::RestoreSavedState()
844 {
845  *fpNavigatorState = fSaveState;
846 }
847 */
848 /*
849 void G4ITNavigator2::RestoreSavedState()
850 {
851  fExitNormal = fpSaveState->sExitNormal;
852  fValidExitNormal = fpSaveState->sValidExitNormal;
853  fExiting = fpSaveState->sExiting;
854  fEntering = fpSaveState->sEntering;
855 
856  fBlockedPhysicalVolume = fpSaveState->spBlockedPhysicalVolume;
857  fBlockedReplicaNo = fpSaveState->sBlockedReplicaNo,
858 
859  fLastStepWasZero = fpSaveState->sLastStepWasZero;
860 
861  // !>
862  fPreviousSftOrigin = fpSaveState->sPreviousSftOrigin ;
863  fPreviousSafety = fpSaveState->sPreviousSafety ;
864  fNumberZeroSteps = fpSaveState->sNumberZeroSteps ;
865  fLocatedOnEdge = fpSaveState->sLocatedOnEdge ;
866  fWasLimitedByGeometry = fpSaveState->sWasLimitedByGeometry;
867  fPushed = fpSaveState->sPushed;
868  fNumberZeroSteps = fpSaveState->sNumberZeroSteps;
869  fEnteredDaughter= fpSaveState->sEnteredDaughter ;
870  fExitedMother = fpSaveState->sExitedMother ;
871 
872  fLastLocatedPointLocal = fpSaveState->sLastLocatedPointLocal ;
873  fLocatedOutsideWorld = fpSaveState->sLocatedOutsideWorld;
874  // <!
875 }
876 */
877 // <!
878 
879 // ********************************************************************
880 // ComputeStep
881 //
882 // Computes the next geometric Step: intersections with current
883 // mother and `daughter' volumes.
884 //
885 // NOTE:
886 //
887 // Flags on entry:
888 // --------------
889 // fValidExitNormal - Normal of exited volume is valid (convex, not a
890 // coincident boundary)
891 // fExitNormal - Surface normal of exited volume
892 // fExiting - True if have exited solid
893 //
894 // fBlockedPhysicalVolume - Ptr to exited volume (or 0)
895 // fBlockedReplicaNo - Replication no of exited volume
896 // fLastStepWasZero - True if last Step size was zero.
897 //
898 // Flags on exit:
899 // -------------
900 // fValidExitNormal - True if surface normal of exited volume is valid
901 // fExitNormal - Surface normal of exited volume rotated to mothers
902 // reference system
903 // fExiting - True if exiting mother
904 // fEntering - True if entering `daughter' volume (or replica)
905 // fBlockedPhysicalVolume - Ptr to candidate (entered) volume
906 // fBlockedReplicaNo - Replication no of candidate (entered) volume
907 // fLastStepWasZero - True if this Step size was zero.
908 // ********************************************************************
909 //
911  const G4ThreeVector &pDirection,
912  const G4double pCurrentProposedStepLength,
913  G4double &pNewSafety)
914 {
916  G4ThreeVector localDirection = ComputeLocalAxis(pDirection);
918  G4VPhysicalVolume *motherPhysical = fHistory.GetTopVolume();
919  G4LogicalVolume *motherLogical = motherPhysical->GetLogicalVolume();
920 
921  // All state relating to exiting normals must be reset
922  //
924  // Reset value - to erase its memory
926  // Reset - used for local exit normal
927  fGrandMotherExitNormal= G4ThreeVector( 0., 0., 0.);
928  fCalculatedExitNormal = false;
929  // Reset for new step
930 
931  static G4ThreadLocal G4int sNavCScalls=0;
932  sNavCScalls++;
933 
935 
936 #ifdef G4VERBOSE
937  if( fVerbose > 0 )
938  {
939  G4cout << "*** G4ITNavigator2::ComputeStep: ***" << G4endl;
940  G4cout << " Volume = " << motherPhysical->GetName()
941  << " - Proposed step length = " << pCurrentProposedStepLength
942  << G4endl;
943 #ifdef G4DEBUG_NAVIGATION
944  if( fVerbose >= 2 )
945  {
946  G4cout << " Called with the arguments: " << G4endl
947  << " Globalpoint = " << std::setw(25) << pGlobalpoint << G4endl
948  << " Direction = " << std::setw(25) << pDirection << G4endl;
949  if( fVerbose >= 4 )
950  {
951  G4cout << " ---- Upon entering : State" << G4endl;
952  PrintState();
953  }
954  }
955 #endif
956  }
957 #endif
958 
959  G4ThreeVector newLocalPoint = ComputeLocalPoint(pGlobalpoint);
960  if( newLocalPoint != fLastLocatedPointLocal )
961  {
962  // Check whether the relocation is within safety
963  //
964  G4ThreeVector oldLocalPoint = fLastLocatedPointLocal;
965  G4double moveLenSq = (newLocalPoint-oldLocalPoint).mag2();
966 
967  if ( moveLenSq >= kCarTolerance*kCarTolerance )
968  {
969 #ifdef G4VERBOSE
970  ComputeStepLog(pGlobalpoint, moveLenSq);
971 #endif
972  // Relocate the point within the same volume
973  //
974  LocateGlobalPointWithinVolume( pGlobalpoint );
975  fLastTriedStepComputation= true; // Ensure that this is set again !!
976  }
977  }
978  if ( fHistory.GetTopVolumeType()!=kReplica )
979  {
980  switch( CharacteriseDaughters(motherLogical) )
981  {
982  case kNormal:
983  if ( motherLogical->GetVoxelHeader() )
984  {
986  localDirection,
987  pCurrentProposedStepLength,
988  pNewSafety,
989  fHistory,
991  fExitNormal,
992  fExiting,
993  fEntering,
996 
997  }
998  else
999  {
1000  if( motherPhysical->GetRegularStructureId() == 0 )
1001  {
1003  localDirection,
1004  pCurrentProposedStepLength,
1005  pNewSafety,
1006  fHistory,
1008  fExitNormal,
1009  fExiting,
1010  fEntering,
1013  }
1014  else // Regular (non-voxelised) structure
1015  {
1016  LocateGlobalPointAndSetup( pGlobalpoint, &pDirection, true, true );
1017  fLastTriedStepComputation= true; // Ensure that this is set again !!
1018  //
1019  // if physical process limits the step, the voxel will not be the
1020  // one given by ComputeStepSkippingEqualMaterials() and the local
1021  // point will be wrongly calculated.
1022 
1023  // There is a problem: when msc limits the step and the point is
1024  // assigned wrongly to phantom in previous step (while it is out
1025  // of the container volume). Then LocateGlobalPointAndSetup() has
1026  // reset the history topvolume to world.
1027  //
1028  if(fHistory.GetTopVolume()->GetRegularStructureId() == 0 )
1029  {
1030  G4Exception("G4ITNavigator2::ComputeStep()",
1031  "GeomNav1001", JustWarning,
1032  "Point is relocated in voxels, while it should be outside!");
1034  localDirection,
1035  pCurrentProposedStepLength,
1036  pNewSafety,
1037  fHistory,
1039  fExitNormal,
1040  fExiting,
1041  fEntering,
1044  }
1045  else
1046  {
1047  Step = fregularNav.
1048  ComputeStepSkippingEqualMaterials(fLastLocatedPointLocal,
1049  localDirection,
1050  pCurrentProposedStepLength,
1051  pNewSafety,
1052  fHistory,
1054  fExitNormal,
1055  fExiting,
1056  fEntering,
1059  motherPhysical);
1060  }
1061  }
1062  }
1063  break;
1064  case kParameterised:
1065  if( GetDaughtersRegularStructureId(motherLogical) != 1 )
1066  {
1068  localDirection,
1069  pCurrentProposedStepLength,
1070  pNewSafety,
1071  fHistory,
1073  fExitNormal,
1074  fExiting,
1075  fEntering,
1078  }
1079  else // Regular structure
1080  {
1082  localDirection,
1083  pCurrentProposedStepLength,
1084  pNewSafety,
1085  fHistory,
1087  fExitNormal,
1088  fExiting,
1089  fEntering,
1092  }
1093  break;
1094  case kReplica:
1095  G4Exception("G4ITNavigator2::ComputeStep()", "GeomNav0001",
1096  FatalException, "Not applicable for replicated volumes.");
1097  break;
1098  }
1099  }
1100  else
1101  {
1102  // In the case of a replica, it must handle the exiting
1103  // edge/corner problem by itself
1104  //
1105  G4bool exitingReplica = fExitedMother;
1106  G4bool calculatedExitNormal;
1107  Step = freplicaNav.ComputeStep(pGlobalpoint,
1108  pDirection,
1110  localDirection,
1111  pCurrentProposedStepLength,
1112  pNewSafety,
1113  fHistory,
1115  calculatedExitNormal,
1116  fExitNormal,
1117  exitingReplica,
1118  fEntering,
1121  fExiting= exitingReplica;
1122  fCalculatedExitNormal= calculatedExitNormal;
1123  }
1124 
1125 
1126 // G4cout << " !!!! Step = " << Step << G4endl;
1127 
1128  // Remember last safety origin & value.
1129  //
1130  fPreviousSftOrigin = pGlobalpoint;
1131  fPreviousSafety = pNewSafety;
1132 
1133  // Count zero steps - one can occur due to changing momentum at a boundary
1134  // - one, two (or a few) can occur at common edges between
1135  // volumes
1136  // - more than two is likely a problem in the geometry
1137  // description or the Navigation
1138 
1139  // Rule of thumb: likely at an Edge if two consecutive steps are zero,
1140  // because at least two candidate volumes must have been
1141  // checked
1142  //
1143  fLocatedOnEdge = fLastStepWasZero && (Step==0.0);
1144  fLastStepWasZero = (Step==0.0);
1145  if (fPushed) { fPushed = fLastStepWasZero; }
1146 
1147  // Handle large number of consecutive zero steps
1148  //
1149  if ( fLastStepWasZero )
1150  {
1151  fNumberZeroSteps++;
1152 #ifdef G4DEBUG_NAVIGATION
1153  if( fNumberZeroSteps > 1 )
1154  {
1155  G4cout << "G4ITNavigator2::ComputeStep(): another zero step, # "
1156  << fNumberZeroSteps
1157  << " at " << pGlobalpoint
1158  << " in volume " << motherPhysical->GetName()
1159  << " nav-comp-step calls # " << sNavCScalls
1160  << G4endl;
1161  }
1162 #endif
1164  {
1165  // Act to recover this stuck track. Pushing it along direction
1166  //
1167  Step += 100*kCarTolerance;
1168 #ifdef G4VERBOSE
1169  if ((!fPushed) && (fWarnPush))
1170  {
1171  std::ostringstream message;
1172  message << "Track stuck or not moving." << G4endl
1173  << " Track stuck, not moving for "
1174  << fNumberZeroSteps << " steps" << G4endl
1175  << " in volume -" << motherPhysical->GetName()
1176  << "- at point " << pGlobalpoint << G4endl
1177  << " direction: " << pDirection << "." << G4endl
1178  << " Potential geometry or navigation problem !"
1179  << G4endl
1180  << " Trying pushing it of " << Step << " mm ...";
1181  G4Exception("G4ITNavigator2::ComputeStep()", "GeomNav1002",
1182  JustWarning, message, "Potential overlap in geometry!");
1183  }
1184 #endif
1185  fPushed = true;
1186  }
1188  {
1189  // Must kill this stuck track
1190  //
1191  std::ostringstream message;
1192  message << "Stuck Track: potential geometry or navigation problem."
1193  << G4endl
1194  << " Track stuck, not moving for "
1195  << fNumberZeroSteps << " steps" << G4endl
1196  << " in volume -" << motherPhysical->GetName()
1197  << "- at point " << pGlobalpoint << G4endl
1198  << " direction: " << pDirection << ".";
1199  motherPhysical->CheckOverlaps(5000, false);
1200  G4Exception("G4ITNavigator2::ComputeStep()", "GeomNav0003",
1201  EventMustBeAborted, message);
1202  }
1203  }
1204  else
1205  {
1206  if (!fPushed) fNumberZeroSteps = 0;
1207  }
1208 
1209  fEnteredDaughter = fEntering; // I expect to enter a volume in this Step
1211 
1212  fStepEndPoint = pGlobalpoint
1213  + std::min(Step,pCurrentProposedStepLength) * pDirection;
1214  fLastStepEndPointLocal = fLastLocatedPointLocal + Step * localDirection;
1215 
1216  if( fExiting )
1217  {
1218 #ifdef G4DEBUG_NAVIGATION
1219  if( fVerbose > 2 )
1220  {
1221  G4cout << " At G4Nav CompStep End - if(exiting) - fExiting= " << fExiting
1222  << " fValidExitNormal = " << fValidExitNormal << G4endl;
1223  G4cout << " fExitNormal= " << fExitNormal << G4endl;
1224  }
1225 #endif
1226 
1228  {
1229  if ( fHistory.GetTopVolumeType()!=kReplica )
1230  {
1231  // Convention: fExitNormal is in the 'grand-mother' coordinate system
1232  //
1234  fCalculatedExitNormal= true;
1235  }
1236  else
1237  {
1239  }
1240  }
1241  else
1242  {
1243  // We must calculate the normal anyway (in order to have it if requested)
1244  //
1245  G4ThreeVector finalLocalPoint =
1246  fLastLocatedPointLocal + localDirection*Step;
1247 
1248  if ( fHistory.GetTopVolumeType()!=kReplica )
1249  {
1250  // Find normal in the 'mother' coordinate system
1251  //
1252  G4ThreeVector exitNormalMotherFrame=
1253  motherLogical->GetSolid()->SurfaceNormal(finalLocalPoint);
1254 
1255  // Transform it to the 'grand-mother' coordinate system
1256  //
1257  const G4RotationMatrix* mRot = motherPhysical->GetRotation();
1258  if( mRot )
1259  {
1261  fGrandMotherExitNormal = (*mRot).inverse() * exitNormalMotherFrame;
1262  }
1263  else
1264  {
1265  fGrandMotherExitNormal = exitNormalMotherFrame;
1266  }
1267 
1268  // Do not set fValidExitNormal -- this signifies
1269  // that the solid is convex!
1270  //
1271  fCalculatedExitNormal= true;
1272  }
1273  else
1274  {
1275  fCalculatedExitNormal = false;
1276  //
1277  // Nothing can be done at this stage currently - to solve this
1278  // Replica Navigation must have calculated the normal for this case
1279  // already.
1280  // Cases: mother is not convex, and exit is at previous replica level
1281 
1282 #ifdef G4DEBUG_NAVIGATION
1284 
1285  desc << "Problem in ComputeStep: Replica Navigation did not provide"
1286  << " valid exit Normal. " << G4endl;
1287  desc << " Do not know how calculate it in this case." << G4endl;
1288  desc << " Location = " << finalLocalPoint << G4endl;
1289  desc << " Volume name = " << motherPhysical->GetName()
1290  << " copy/replica No = " << motherPhysical->GetCopyNo() << G4endl;
1291  G4Exception("G4ITNavigator2::ComputeStep()", "GeomNav0003",
1292  JustWarning, desc, "Normal not available for exiting.");
1293 #endif
1294  }
1295  }
1296 
1297  // Now transform it to the global reference frame !!
1298  //
1300  {
1301  G4int depth= fHistory.GetDepth();
1302  if( depth > 0 )
1303  {
1304  G4AffineTransform GrandMotherToGlobalTransf =
1305  fHistory.GetTransform(depth-1).Inverse();
1307  GrandMotherToGlobalTransf.TransformAxis( fGrandMotherExitNormal );
1308  }
1309  else
1310  {
1312  }
1313  }
1314  else
1315  {
1316  fExitNormalGlobalFrame= G4ThreeVector( 0., 0., 0.);
1317  }
1318  }
1319 
1320  if( (Step == pCurrentProposedStepLength) && (!fExiting) && (!fEntering) )
1321  {
1322  // This if Step is not really limited by the geometry.
1323  // The Navigator is obliged to return "infinity"
1324  //
1325  Step = kInfinity;
1326  }
1327 
1328 #ifdef G4VERBOSE
1329  if( fVerbose > 1 )
1330  {
1331  if( fVerbose >= 4 )
1332  {
1333  G4cout << " ----- Upon exiting :" << G4endl;
1334  PrintState();
1335  }
1336  G4cout << " Returned step= " << Step;
1337  if( fVerbose > 5 ) G4cout << G4endl;
1338  if( Step == kInfinity )
1339  {
1340  G4cout << " Requested step= " << pCurrentProposedStepLength ;
1341  if( fVerbose > 5) G4cout << G4endl;
1342  }
1343  G4cout << " Safety = " << pNewSafety << G4endl;
1344  }
1345 #endif
1346 
1347  return Step;
1348 }
1349 
1350 // ********************************************************************
1351 // CheckNextStep
1352 //
1353 // Compute the step without altering the navigator state
1354 // ********************************************************************
1355 //
1357  const G4ThreeVector& pDirection,
1358  const G4double pCurrentProposedStepLength,
1359  G4double& pNewSafety)
1360 {
1362  G4double step;
1363 
1364  // Save the state, for this parasitic call
1365  //
1366  //SetSavedState();
1367 // G4SaveNavigatorState savedState(fpNavigatorState);
1368  G4NavigatorState savedState(*fpNavigatorState);
1369 
1370  step = ComputeStep ( pGlobalpoint,
1371  pDirection,
1372  pCurrentProposedStepLength,
1373  pNewSafety );
1374 
1375  // If a parasitic call, then attempt to restore the key parts of the state
1376  //
1377  *fpNavigatorState = savedState;
1378  //RestoreSavedState();
1379  // NOTE: the state of the current subnavigator is NOT restored.
1380  // ***> TODO: restore subnavigator state
1381  // if( last_located) Need Position of last location
1382  // if( last_computed step) Need Endposition of last step
1383 
1384 
1385  return step;
1386 }
1387 
1388 // ********************************************************************
1389 // ResetState
1390 //
1391 // Resets stack and minimum of navigator state `machine'
1392 // ********************************************************************
1393 //
1395 {
1396  fWasLimitedByGeometry = false;
1397  fEntering = false;
1398  fExiting = false;
1399  fLocatedOnEdge = false;
1400  fLastStepWasZero = false;
1401  fEnteredDaughter = false;
1402  fExitedMother = false;
1403  fPushed = false;
1404 
1405  fValidExitNormal = false;
1407  fCalculatedExitNormal = false;
1408 
1409  fExitNormal = G4ThreeVector(0,0,0);
1412 
1414  fPreviousSafety = 0.0;
1415 
1416  fNumberZeroSteps = 0;
1417 
1419  fBlockedReplicaNo = -1;
1420 
1422  fLocatedOutsideWorld = false;
1423 }
1424 
1425 // ********************************************************************
1426 // SetupHierarchy
1427 //
1428 // Renavigates & resets hierarchy described by current history
1429 // o Reset volumes
1430 // o Recompute transforms and/or solids of replicated/parameterised volumes
1431 // ********************************************************************
1432 //
1434 {
1435  G4int i;
1436  const G4int cdepth = fHistory.GetDepth();
1437  G4VPhysicalVolume *current;
1438  G4VSolid *pSolid;
1439  G4VPVParameterisation *pParam;
1440 
1441  for ( i=1; i<=cdepth; i++ )
1442  {
1443  current = fHistory.GetVolume(i);
1444  switch ( fHistory.GetVolumeType(i) )
1445  {
1446  case kNormal:
1447  break;
1448  case kReplica:
1449  freplicaNav.ComputeTransformation(fHistory.GetReplicaNo(i), current);
1450  break;
1451  case kParameterised:
1452  G4int replicaNo;
1453  pParam = current->GetParameterisation();
1454  replicaNo = fHistory.GetReplicaNo(i);
1455  pSolid = pParam->ComputeSolid(replicaNo, current);
1456 
1457  // Set up dimensions & transform in solid/physical volume
1458  //
1459  pSolid->ComputeDimensions(pParam, replicaNo, current);
1460  pParam->ComputeTransformation(replicaNo, current);
1461 
1462  G4TouchableHistory *pTouchable= 0;
1463  if( pParam->IsNested() )
1464  {
1465  pTouchable= new G4TouchableHistory( fHistory );
1466  pTouchable->MoveUpHistory(); // Move up to the parent level
1467  // Adequate only if Nested at the Branch level (last)
1468  // To extend to other cases:
1469  // pTouchable->MoveUpHistory(cdepth-i-1);
1470  // Move to the parent level of *Current* level
1471  // Could replace this line and constructor with a revised
1472  // c-tor for History(levels to drop)
1473  }
1474  // Set up the correct solid and material in Logical Volume
1475  //
1476  G4LogicalVolume *pLogical = current->GetLogicalVolume();
1477  pLogical->SetSolid( pSolid );
1478  pLogical->UpdateMaterial( pParam ->
1479  ComputeMaterial(replicaNo, current, pTouchable) );
1480  delete pTouchable;
1481  break;
1482  }
1483  }
1484 }
1485 
1486 // ********************************************************************
1487 // GetLocalExitNormal
1488 //
1489 // Obtains the Normal vector to a surface (in local coordinates)
1490 // pointing out of previous volume and into current volume
1491 // ********************************************************************
1492 //
1494 {
1496  G4ThreeVector ExitNormal(0.,0.,0.);
1497  G4VSolid *currentSolid=0;
1498  G4LogicalVolume *candidateLogical;
1500  {
1501  // use fLastLocatedPointLocal and next candidate volume
1502  //
1503  G4ThreeVector nextSolidExitNormal(0.,0.,0.);
1504 
1505  if( fEntering && (fBlockedPhysicalVolume!=0) )
1506  {
1507  candidateLogical= fBlockedPhysicalVolume->GetLogicalVolume();
1508  if( candidateLogical )
1509  {
1510  // fLastStepEndPointLocal is in the coordinates of the mother
1511  // we need it in the daughter's coordinate system.
1512 
1513  // The following code should also work in case of Replica
1514  {
1515  // First transform fLastLocatedPointLocal to the new daughter
1516  // coordinates
1517  //
1518  G4AffineTransform MotherToDaughterTransform=
1522  G4ThreeVector daughterPointOwnLocal=
1523  MotherToDaughterTransform.TransformPoint( fLastStepEndPointLocal );
1524 
1525  // OK if it is a parameterised volume
1526  //
1527  EInside inSideIt;
1528  G4bool onSurface;
1529  G4double safety= -1.0;
1530  currentSolid= candidateLogical->GetSolid();
1531  inSideIt = currentSolid->Inside(daughterPointOwnLocal);
1532  onSurface = (inSideIt == kSurface);
1533  if( ! onSurface )
1534  {
1535  if( inSideIt == kOutside )
1536  {
1537  safety = (currentSolid->DistanceToIn(daughterPointOwnLocal));
1538  onSurface = safety < 100.0 * kCarTolerance;
1539  }
1540  else if (inSideIt == kInside )
1541  {
1542  safety = (currentSolid->DistanceToOut(daughterPointOwnLocal));
1543  onSurface = safety < 100.0 * kCarTolerance;
1544  }
1545  }
1546 
1547  if( onSurface )
1548  {
1549  nextSolidExitNormal =
1550  currentSolid->SurfaceNormal(daughterPointOwnLocal);
1551 
1552  // Entering the solid ==> opposite
1553  //
1554  ExitNormal = -nextSolidExitNormal;
1555  fCalculatedExitNormal= true;
1556  }
1557  else
1558  {
1559 #ifdef G4VERBOSE
1560  if(( fVerbose == 1 ) && ( fCheck ))
1561  {
1562  std::ostringstream message;
1563  message << "Point not on surface ! " << G4endl
1564  << " Point = "
1565  << daughterPointOwnLocal << G4endl
1566  << " Physical volume = "
1567  << fBlockedPhysicalVolume->GetName() << G4endl
1568  << " Logical volume = "
1569  << candidateLogical->GetName() << G4endl
1570  << " Solid = " << currentSolid->GetName()
1571  << " Type = "
1572  << currentSolid->GetEntityType() << G4endl
1573  << *currentSolid << G4endl;
1574  if( inSideIt == kOutside )
1575  {
1576  message << "Point is Outside. " << G4endl
1577  << " Safety (from outside) = " << safety << G4endl;
1578  }
1579  else // if( inSideIt == kInside )
1580  {
1581  message << "Point is Inside. " << G4endl
1582  << " Safety (from inside) = " << safety << G4endl;
1583  }
1584  G4Exception("G4ITNavigator2::GetLocalExitNormal()", "GeomNav1001",
1585  JustWarning, message);
1586  }
1587 #endif
1588  }
1589  *valid = onSurface; // was =true;
1590  }
1591  }
1592  }
1593  else if ( fExiting )
1594  {
1595  ExitNormal = fGrandMotherExitNormal;
1596  *valid = true;
1597  fCalculatedExitNormal= true; // Should be true already
1598  }
1599  else // i.e. ( fBlockedPhysicalVolume == 0 )
1600  {
1601  *valid = false;
1602  G4Exception("G4ITNavigator2::GetLocalExitNormal()",
1603  "GeomNav0003", JustWarning,
1604  "Incorrect call to GetLocalSurfaceNormal." );
1605  }
1606  }
1607  else // ( ! fLastTriedStepComputation ) ie. last call was to Locate
1608  {
1609  if ( EnteredDaughterVolume() )
1610  {
1611  G4VSolid* daughterSolid =fHistory.GetTopVolume()->GetLogicalVolume()
1612  ->GetSolid();
1613  ExitNormal= -(daughterSolid->SurfaceNormal(fLastLocatedPointLocal));
1614  if( std::fabs(ExitNormal.mag2()-1.0 ) > CLHEP::perMillion )
1615  {
1617  desc << " Parameters of solid: " << *daughterSolid
1618  << " Point for surface = " << fLastLocatedPointLocal << std::endl;
1619  G4Exception("G4ITNavigator2::GetLocalExitNormal()",
1620  "GeomNav0003", FatalException, desc,
1621  "Surface Normal returned by Solid is not a Unit Vector." );
1622  }
1623  fCalculatedExitNormal= true;
1624  *valid = true;
1625  }
1626  else
1627  {
1628  if( fExitedMother )
1629  {
1630  ExitNormal = fGrandMotherExitNormal;
1631  *valid = true;
1632  fCalculatedExitNormal= true;
1633  }
1634  else // We are not at a boundary. ExitNormal remains (0,0,0)
1635  {
1636  *valid = false;
1637  fCalculatedExitNormal= false;
1638  G4ExceptionDescription message;
1639  message << "Function called when *NOT* at a Boundary." << G4endl;
1640  G4Exception("G4ITNavigator2::GetLocalExitNormal()",
1641  "GeomNav0003", JustWarning, message);
1642  }
1643  }
1644  }
1645  return ExitNormal;
1646 }
1647 
1648 // ********************************************************************
1649 // GetMotherToDaughterTransform
1650 //
1651 // Obtains the mother to daughter affine transformation
1652 // ********************************************************************
1653 //
1656  G4int enteringReplicaNo,
1657  EVolume enteringVolumeType )
1658 {
1660  switch (enteringVolumeType)
1661  {
1662  case kNormal: // Nothing is needed to prepare the transformation
1663  break; // It is stored already in the physical volume (placement)
1664  case kReplica: // Sets the transform in the Replica - tbc
1665  G4Exception("G4ITNavigator2::GetMotherToDaughterTransform()",
1666  "GeomNav0001", FatalException,
1667  "Method NOT Implemented yet for replica volumes.");
1668  break;
1669  case kParameterised:
1670  if( pEnteringPhysVol->GetRegularStructureId() == 0 )
1671  {
1672  G4VPVParameterisation *pParam =
1673  pEnteringPhysVol->GetParameterisation();
1674  G4VSolid* pSolid =
1675  pParam->ComputeSolid(enteringReplicaNo, pEnteringPhysVol);
1676  pSolid->ComputeDimensions(pParam, enteringReplicaNo, pEnteringPhysVol);
1677 
1678  // Sets the transform in the Parameterisation
1679  //
1680  pParam->ComputeTransformation(enteringReplicaNo, pEnteringPhysVol);
1681 
1682  // Set the correct solid and material in Logical Volume
1683  //
1684  G4LogicalVolume* pLogical = pEnteringPhysVol->GetLogicalVolume();
1685  pLogical->SetSolid( pSolid );
1686  }
1687  break;
1688  }
1689  return G4AffineTransform(pEnteringPhysVol->GetRotation(),
1690  pEnteringPhysVol->GetTranslation()).Invert();
1691 }
1692 
1693 // ********************************************************************
1694 // GetLocalExitNormalAndCheck
1695 //
1696 // Obtains the Normal vector to a surface (in local coordinates)
1697 // pointing out of previous volume and into current volume, and
1698 // checks the current point against expected 'local' value.
1699 // ********************************************************************
1700 //
1703 #ifdef G4DEBUG_NAVIGATION
1704  const G4ThreeVector& ExpectedBoundaryPointGlobal,
1705 #else
1706  const G4ThreeVector&,
1707 #endif
1708  G4bool* pValid)
1709 {
1711 #ifdef G4DEBUG_NAVIGATION
1712  // Check Current point against expected 'local' value
1713  //
1715  {
1716  G4ThreeVector ExpectedBoundaryPointLocal;
1717 
1718  const G4AffineTransform& GlobalToLocal= GetGlobalToLocalTransform();
1719  ExpectedBoundaryPointLocal =
1720  GlobalToLocal.TransformPoint( ExpectedBoundaryPointGlobal );
1721 
1722  // Add here: Comparison against expected position,
1723  // i.e. the endpoint of ComputeStep
1724  }
1725 #endif
1726 
1727  return GetLocalExitNormal( pValid);
1728 }
1729 
1730 // ********************************************************************
1731 // GetGlobalExitNormal
1732 //
1733 // Obtains the Normal vector to a surface (in global coordinates)
1734 // pointing out of previous volume and into current volume
1735 // ********************************************************************
1736 //
1739  G4bool* pNormalCalculated)
1740 {
1742  G4bool validNormal;
1743  G4ThreeVector localNormal, globalNormal;
1744  G4bool usingStored;
1745 
1746  usingStored=
1748  ( ( fLastTriedStepComputation && fExiting) // Just calculated it
1749  || // No locate in between
1751  && (IntersectPointGlobal-fStepEndPoint).mag2()
1752  < (10.0*kCarTolerance*kCarTolerance)
1753  ) // Calculated it 'just' before & then called locate
1754  // but it did not move position
1755  );
1756 
1757  if( usingStored )
1758  {
1759  // This was computed in ComputeStep -- and only on arrival at boundary
1760  //
1761  globalNormal = fExitNormalGlobalFrame;
1762  G4double normMag2 = globalNormal.mag2();
1763  if( std::fabs ( normMag2 - 1.0 ) < perMillion ) // Value is good
1764  {
1765  *pNormalCalculated = true; // ComputeStep always computes it if Exiting
1766  // (fExiting==true)
1767  }
1768  else
1769  {
1770  G4ExceptionDescription message;
1771 
1772  message << " ERROR> Expected normal-global-frame to valid (unit vector) "
1773  << " - but |normal| = " << std::sqrt(normMag2)
1774  << " - and |normal|^ = " << normMag2
1775  << " which differs from 1.0 by " << normMag2 - 1.0 << G4endl
1776  << " n = " << fExitNormalGlobalFrame << G4endl;
1777  message << "============================================================"
1778  << G4endl;
1779  G4int oldVerbose = fVerbose;
1780  fVerbose=4;
1781  message << " State of Navigator: " << G4endl;
1782  message << *this << G4endl;
1783  fVerbose = oldVerbose;
1784  message << "============================================================"
1785  << G4endl;
1786 
1787  G4Exception("G4ITNavigator2::GetGlobalExitNormal()",
1788  "GeomNav0003",JustWarning, message,
1789  "Value obtained from stored global-normal is not a unit vector.");
1790 
1791  // (Re)Compute it now -- as either it was not computed, or it is wrong.
1792  //
1793  localNormal = GetLocalExitNormalAndCheck(IntersectPointGlobal,
1794  &validNormal);
1795  *pNormalCalculated = fCalculatedExitNormal;
1796 
1797  G4AffineTransform localToGlobal = GetLocalToGlobalTransform();
1798  globalNormal = localToGlobal.TransformAxis( localNormal );
1799  }
1800  }
1801  else
1802  {
1803  localNormal = GetLocalExitNormalAndCheck(IntersectPointGlobal,&validNormal);
1804  *pNormalCalculated = fCalculatedExitNormal;
1805 
1806 #ifdef G4DEBUG_NAVIGATION
1807  usingStored= false;
1808 
1809  if( (!validNormal) && !fCalculatedExitNormal)
1810  {
1812  edN << " Calculated = " << fCalculatedExitNormal << G4endl;
1813  edN << " Entering= " << fEntering << G4endl;
1814  G4int oldVerbose= this->GetVerboseLevel();
1815  this->SetVerboseLevel(4);
1816  edN << " State of Navigator: " << G4endl;
1817  edN << *this << G4endl;
1818  this->SetVerboseLevel( oldVerbose );
1819 
1820  G4Exception("G4ITNavigator2::GetGlobalExitNormal()",
1821  "GeomNav0003", JustWarning, edN,
1822  "LocalExitNormalAndCheck() did not calculate Normal.");
1823  }
1824 #endif
1825 
1826  G4double localMag2= localNormal.mag2();
1827  if( validNormal && (std::fabs(localMag2-1.0)) > CLHEP::perMillion )
1828  {
1830 
1831  edN << "G4ITNavigator2::GetGlobalExitNormal: "
1832  << " Using Local Normal - from call to GetLocalExitNormalAndCheck. "
1833  << G4endl
1834  << " Local Exit Normal : " << " || = " << std::sqrt(localMag2)
1835  << " vec = " << localNormal << G4endl
1836  << " Global Exit Normal : " << " || = " << globalNormal.mag()
1837  << " vec = " << globalNormal << G4endl;
1838  edN << " Calculated It = " << fCalculatedExitNormal << G4endl;
1839 
1840  G4Exception("G4ITNavigator2::GetGlobalExitNormal()",
1841  "GeomNav0003",JustWarning, edN,
1842  "Value obtained from new local *solid* is incorrect.");
1843  localNormal = localNormal.unit(); // Should we correct it ??
1844  }
1845  G4AffineTransform localToGlobal = GetLocalToGlobalTransform();
1846  globalNormal = localToGlobal.TransformAxis( localNormal );
1847  }
1848 
1849 #ifdef G4DEBUG_NAVIGATION
1850  if( usingStored )
1851  {
1852  G4ThreeVector globalNormAgn;
1853 
1854  localNormal= GetLocalExitNormalAndCheck(IntersectPointGlobal, &validNormal);
1855 
1856  G4AffineTransform localToGlobal = GetLocalToGlobalTransform();
1857  globalNormAgn = localToGlobal.TransformAxis( localNormal );
1858 
1859  // Check the value computed against fExitNormalGlobalFrame
1860  G4ThreeVector diffNorm = globalNormAgn - fExitNormalGlobalFrame;
1861  if( diffNorm.mag2() > perMillion*CLHEP::perMillion)
1862  {
1863  G4ExceptionDescription edDfn;
1864  edDfn << "Found difference in normals in case of exiting mother "
1865  << "- when Get is called after ComputingStep " << G4endl;
1866  edDfn << " Magnitude of diff = " << diffNorm.mag() << G4endl;
1867  edDfn << " Normal stored (Global) = " << fExitNormalGlobalFrame
1868  << G4endl;
1869  edDfn << " Global Computed from Local = " << globalNormAgn << G4endl;
1870  G4Exception("G4ITNavigator::GetGlobalExitNormal()", "GeomNav0003",
1871  JustWarning, edDfn);
1872  }
1873  }
1874 #endif
1875 
1876  return globalNormal;
1877 }
1878 
1879 // To make the new Voxel Safety the default, uncomment the next line
1880 #define G4NEW_SAFETY 1
1881 
1882 // ********************************************************************
1883 // ComputeSafety
1884 //
1885 // It assumes that it will be
1886 // i) called at the Point in the same volume as the EndPoint of the
1887 // ComputeStep.
1888 // ii) after (or at the end of) ComputeStep OR after the relocation.
1889 // ********************************************************************
1890 //
1892  const G4double pMaxLength,
1893  const G4bool keepState)
1894 {
1896  G4double newSafety = 0.0;
1897 
1898 #ifdef G4DEBUG_NAVIGATION
1899  G4int oldcoutPrec = G4cout.precision(8);
1900  if( fVerbose > 0 )
1901  {
1902  G4cout << "*** G4ITNavigator2::ComputeSafety: ***" << G4endl
1903  << " Called at point: " << pGlobalpoint << G4endl;
1904 
1905  G4VPhysicalVolume *motherPhysical = fHistory.GetTopVolume();
1906  G4cout << " Volume = " << motherPhysical->GetName()
1907  << " - Maximum length = " << pMaxLength << G4endl;
1908  if( fVerbose >= 4 )
1909  {
1910  G4cout << " ----- Upon entering Compute Safety:" << G4endl;
1911  PrintState();
1912  }
1913  }
1914 #endif
1915 
1916  G4SaveNavigatorState* savedState(0);
1917 
1918  G4double distEndpointSq = (pGlobalpoint-fStepEndPoint).mag2();
1919  G4bool stayedOnEndpoint = distEndpointSq < kCarTolerance*kCarTolerance;
1920  G4bool endpointOnSurface = fEnteredDaughter || fExitedMother;
1921 
1922  if( endpointOnSurface && stayedOnEndpoint )
1923  {
1924 #ifdef G4DEBUG_NAVIGATION
1925  if( fVerbose >= 2 )
1926  {
1927  G4cout << " G4Navigator::ComputeSafety() finds that point - "
1928  << pGlobalpoint << " - is on surface " << G4endl;
1929  if( fEnteredDaughter ) { G4cout << " entered new daughter volume"; }
1930  if( fExitedMother ) { G4cout << " and exited previous volume."; }
1931  G4cout << G4endl;
1932  G4cout << " EndPoint was = " << fStepEndPoint << G4endl;
1933  }
1934 #endif
1935  newSafety = 0.0;
1936  // return newSafety;
1937  }
1938  else // if( !(endpointOnSurface && stayedOnEndpoint) )
1939  {
1940 
1941  if (keepState)
1942  {
1943 // SetSavedState();
1944  savedState = new G4SaveNavigatorState(fpNavigatorState);
1945  }
1946 
1947  // Pseudo-relocate to this point (updates voxel information only)
1948  //
1949  LocateGlobalPointWithinVolume( pGlobalpoint );
1950  // --->> DANGER: Side effects on sub-navigator voxel information <<---
1951  // Could be replaced again by 'granular' calls to sub-navigator
1952  // locates (similar side-effects, but faster.
1953  // Solutions:
1954  // 1) Re-locate (to where?)
1955  // 2) Insure that the methods using (G4ComputeStep?)
1956  // does a relocation (if information is disturbed only ?)
1957 
1958 #ifdef G4DEBUG_NAVIGATION
1959  if( fVerbose >= 2 )
1960  {
1961  G4cout << " G4ITNavigator2::ComputeSafety() relocates-in-volume to point: "
1962  << pGlobalpoint << G4endl;
1963  }
1964 #endif
1965  G4VPhysicalVolume *motherPhysical = fHistory.GetTopVolume();
1966  G4LogicalVolume *motherLogical = motherPhysical->GetLogicalVolume();
1967  G4SmartVoxelHeader* pVoxelHeader = motherLogical->GetVoxelHeader();
1968  G4ThreeVector localPoint = ComputeLocalPoint(pGlobalpoint);
1969 
1970  if ( fHistory.GetTopVolumeType()!=kReplica )
1971  {
1972  switch(CharacteriseDaughters(motherLogical))
1973  {
1974  case kNormal:
1975  if ( pVoxelHeader )
1976  {
1977 #ifdef G4NEW_SAFETY
1978  G4double safetyTwo = fpVoxelSafety->ComputeSafety(localPoint,
1979  *motherPhysical, pMaxLength);
1980  newSafety= safetyTwo; // Faster and best available
1981 #else
1982  G4double safetyOldVoxel;
1983  safetyOldVoxel =
1984  fvoxelNav.ComputeSafety(localPoint,fHistory,pMaxLength);
1985  newSafety= safetyOldVoxel;
1986 #endif
1987  }
1988  else
1989  {
1990  newSafety=fnormalNav.ComputeSafety(localPoint,fHistory,pMaxLength);
1991  }
1992  break;
1993  case kParameterised:
1994  if( GetDaughtersRegularStructureId(motherLogical) != 1 )
1995  {
1996  newSafety=fparamNav.ComputeSafety(localPoint,fHistory,pMaxLength);
1997  }
1998  else // Regular structure
1999  {
2000  newSafety=fregularNav.ComputeSafety(localPoint,fHistory,pMaxLength);
2001  }
2002  break;
2003  case kReplica:
2004  G4Exception("G4ITNavigator2::ComputeSafety()", "GeomNav0001",
2005  FatalException, "Not applicable for replicated volumes.");
2006  break;
2007  }
2008  }
2009  else
2010  {
2011  newSafety = freplicaNav.ComputeSafety(pGlobalpoint, localPoint,
2012  fHistory, pMaxLength);
2013  }
2014 
2015  if (keepState)
2016  {
2017  *fpNavigatorState = *savedState;
2018  delete savedState;
2019  // RestoreSavedState();
2020  // This now overwrites the values of the Safety 'sphere' (correction)
2021  }
2022 
2023  // Remember last safety origin & value
2024  //
2025  // We overwrite the Safety 'sphere' - keeping old behaviour
2026  fPreviousSftOrigin = pGlobalpoint;
2027  fPreviousSafety = newSafety;
2028  }
2029 
2030 #ifdef G4DEBUG_NAVIGATION
2031  if( fVerbose > 1 )
2032  {
2033  G4cout << " ---- Exiting ComputeSafety " << G4endl;
2034  if( fVerbose > 2 ) { PrintState(); }
2035  G4cout << " Returned value of Safety = " << newSafety << G4endl;
2036  }
2037  G4cout.precision(oldcoutPrec);
2038 #endif
2039 
2040  return newSafety;
2041 }
2042 
2043 
2044 // ********************************************************************
2045 // RecheckDistanceToCurrentBoundary
2046 //
2047 // Trial method for checking potential displacement for MS
2048 // Check position aDisplacedGlobalpoint, to see whether it is in the
2049 // current volume (mother).
2050 // If in mother, check distance to boundary along aNewDirection.
2051 // If in entering daughter, check distance back to boundary.
2052 // NOTE:
2053 // Can be called only after ComputeStep() is called - before ReLocation
2054 // Deals only with current volume (and potentially entered)
2055 // Caveats
2056 // First VERSION: Does not consider daughter volumes if it remained in mother
2057 // neither for checking whether it has exited current (mother) volume,
2058 // nor for determining distance to exit this (mother) volume.
2059 // ********************************************************************
2060 //
2062  const G4ThreeVector &aDisplacedGlobalPoint,
2063  const G4ThreeVector &aNewDirection,
2064  const G4double ProposedMove,
2065  G4double *prDistance,
2066  G4double *prNewSafety) const
2067 {
2068  G4ThreeVector localPosition = ComputeLocalPoint(aDisplacedGlobalPoint);
2069  G4ThreeVector localDirection = ComputeLocalAxis(aNewDirection);
2070  // G4double Step = kInfinity;
2071 
2072  G4bool validExitNormal;
2073  G4ThreeVector exitNormal;
2074  // Check against mother solid
2075  G4VPhysicalVolume *motherPhysical = fHistory.GetTopVolume();
2076  G4LogicalVolume *motherLogical = motherPhysical->GetLogicalVolume();
2077 
2078 #ifdef CHECK_ORDER_OF_METHODS
2080  {
2081  G4Exception("G4Navigator::RecheckDistanceToCurrentBoundary()",
2082  "GeomNav0001", FatalException,
2083  "Method must be called after ComputeStep(), before call to LocateMethod.");
2084  }
2085 #endif
2086 
2087  EInside locatedDaug; // = kUndefined;
2088  G4double daughterStep= DBL_MAX;
2089  G4double daughterSafety= DBL_MAX;
2090 
2091  if( fEnteredDaughter )
2092  {
2093  if( motherLogical->CharacteriseDaughters() ==kReplica ) { return false; }
2094 
2095  // Track arrived at boundary of a daughter volume at
2096  // the last call of ComputeStep().
2097  // In case the proposed displaced point is inside this daughter,
2098  // it must backtrack at least to the entry point.
2099  // NOTE: No check is made against other daughter volumes. It is
2100  // assumed that the proposed displacement is small enough that
2101  // this is not needed.
2102 
2103  // Must check boundary of current daughter
2104  //
2105  G4VPhysicalVolume *candPhysical= fBlockedPhysicalVolume;
2106  G4LogicalVolume *candLogical= candPhysical->GetLogicalVolume();
2107  G4VSolid *candSolid= candLogical->GetSolid();
2108 
2109  G4AffineTransform nextLevelTrf(candPhysical->GetRotation(),
2110  candPhysical->GetTranslation());
2111 
2112  G4ThreeVector dgPosition= nextLevelTrf.TransformPoint(localPosition);
2113  G4ThreeVector dgDirection= nextLevelTrf.TransformAxis(localDirection);
2114  locatedDaug = candSolid->Inside(dgPosition);
2115 
2116  if( locatedDaug == kInside )
2117  {
2118  // Reverse direction - and find first exit. ( Is it valid?)
2119  // Must backtrack
2120  G4double distanceBackOut =
2121  candSolid->DistanceToOut(dgPosition,
2122  - dgDirection, // Reverse direction
2123  true, &validExitNormal, &exitNormal);
2124  daughterStep= - distanceBackOut;
2125  // No check is made whether the particle could have arrived at
2126  // at this location without encountering another volume or
2127  // a different psurface of the current volume
2128  if( prNewSafety )
2129  {
2130  daughterSafety= candSolid->DistanceToOut(dgPosition);
2131  }
2132  }
2133  else
2134  {
2135  if( locatedDaug == kOutside )
2136  {
2137  // See whether it still intersects it
2138  //
2139  daughterStep= candSolid->DistanceToIn(dgPosition,
2140  dgDirection);
2141  if( prNewSafety )
2142  {
2143  daughterSafety= candSolid->DistanceToIn(dgPosition);
2144  }
2145  }
2146  else
2147  {
2148  // The point remains on the surface of candidate solid
2149  //
2150  daughterStep= 0.0;
2151  daughterSafety= 0.0;
2152  }
2153  }
2154 
2155  // If trial point is in daughter (or on its surface) we have the
2156  // answer, the rest is not relevant
2157  //
2158  if( locatedDaug != kOutside )
2159  {
2160  *prDistance= daughterStep;
2161  if( prNewSafety ) { *prNewSafety= daughterSafety; }
2162  return true;
2163  }
2164  // If ever extended, so that some type of mother cut daughter,
2165  // this would change
2166  }
2167 
2168  G4VSolid *motherSolid= motherLogical->GetSolid();
2169 
2170  G4double motherStep= DBL_MAX, motherSafety= DBL_MAX;
2171 
2172  // Check distance to boundary of mother
2173  //
2174  if ( fHistory.GetTopVolumeType()!=kReplica )
2175  {
2176  EInside locatedMoth = motherSolid->Inside(localPosition);
2177 
2178  if( locatedMoth == kInside )
2179  {
2180  motherSafety= motherSolid->DistanceToOut(localPosition);
2181  if( ProposedMove >= motherSafety )
2182  {
2183  motherStep= motherSolid->DistanceToOut(localPosition,
2184  localDirection,
2185  true, &validExitNormal, &exitNormal);
2186  }
2187  else
2188  {
2189  motherStep= ProposedMove;
2190  }
2191  }
2192  else if( locatedMoth == kOutside)
2193  {
2194  motherSafety= motherSolid->DistanceToIn(localPosition);
2195  if( ProposedMove >= motherSafety )
2196  {
2197  motherStep= - motherSolid->DistanceToIn(localPosition,
2198  -localDirection);
2199  }
2200  }
2201  else
2202  {
2203  motherSafety= 0.0;
2204  *prDistance= 0.0; // On surface - no move // motherStep;
2205  if( prNewSafety ) { *prNewSafety= motherSafety; }
2206  return false;
2207  }
2208  }
2209  else
2210  {
2211  return false;
2212  }
2213 
2214  *prDistance= std::min( motherStep, daughterStep );
2215  if( prNewSafety )
2216  {
2217  *prNewSafety= std::min( motherSafety, daughterSafety );
2218  }
2219  return true;
2220 }
2221 
2222 
2223 // ********************************************************************
2224 // CreateTouchableHistoryHandle
2225 // ********************************************************************
2226 //
2228 {
2231 }
2232 
2233 // ********************************************************************
2234 // PrintState
2235 // ********************************************************************
2236 //
2238 {
2240  G4int oldcoutPrec = G4cout.precision(4);
2241  if( fVerbose >= 4 )
2242  {
2243  G4cout << "The current state of G4Navigator is: " << G4endl;
2244  G4cout << " ValidExitNormal= " << fValidExitNormal // << G4endl
2245  << " ExitNormal = " << fExitNormal // << G4endl
2246  << " Exiting = " << fExiting // << G4endl
2247  << " Entering = " << fEntering // << G4endl
2248  << " BlockedPhysicalVolume= " ;
2249  if (fBlockedPhysicalVolume==0)
2250  G4cout << "None";
2251  else
2252  G4cout << fBlockedPhysicalVolume->GetName();
2253  G4cout << G4endl
2254  << " BlockedReplicaNo = " << fBlockedReplicaNo // << G4endl
2255  << " LastStepWasZero = " << fLastStepWasZero // << G4endl
2256  << G4endl;
2257  }
2258  if( ( 1 < fVerbose) && (fVerbose < 4) )
2259  {
2260  G4cout << G4endl; // Make sure to line up
2261  G4cout << std::setw(30) << " ExitNormal " << " "
2262  << std::setw( 5) << " Valid " << " "
2263  << std::setw( 9) << " Exiting " << " "
2264  << std::setw( 9) << " Entering" << " "
2265  << std::setw(15) << " Blocked:Volume " << " "
2266  << std::setw( 9) << " ReplicaNo" << " "
2267  << std::setw( 8) << " LastStepZero " << " "
2268  << G4endl;
2269  G4cout << "( " << std::setw(7) << fExitNormal.x()
2270  << ", " << std::setw(7) << fExitNormal.y()
2271  << ", " << std::setw(7) << fExitNormal.z() << " ) "
2272  << std::setw( 5) << fValidExitNormal << " "
2273  << std::setw( 9) << fExiting << " "
2274  << std::setw( 9) << fEntering << " ";
2275  if ( fBlockedPhysicalVolume==0 )
2276  G4cout << std::setw(15) << "None";
2277  else
2278  G4cout << std::setw(15)<< fBlockedPhysicalVolume->GetName();
2279  G4cout << std::setw( 9) << fBlockedReplicaNo << " "
2280  << std::setw( 8) << fLastStepWasZero << " "
2281  << G4endl;
2282  }
2283  if( fVerbose > 2 )
2284  {
2285  G4cout.precision(8);
2286  G4cout << " Current Localpoint = " << fLastLocatedPointLocal << G4endl;
2287  G4cout << " PreviousSftOrigin = " << fPreviousSftOrigin << G4endl;
2288  G4cout << " PreviousSafety = " << fPreviousSafety << G4endl;
2289  }
2290  G4cout.precision(oldcoutPrec);
2291 }
2292 
2293 // ********************************************************************
2294 // ComputeStepLog
2295 // ********************************************************************
2296 //
2298  G4double moveLenSq) const
2299 {
2301  // The following checks only make sense if the move is larger
2302  // than the tolerance.
2303 
2304  static const G4double fAccuracyForWarning = kCarTolerance,
2305  fAccuracyForException = 1000*kCarTolerance;
2306 
2307  G4ThreeVector OriginalGlobalpoint = fHistory.GetTopTransform().Inverse().
2308  TransformPoint(fLastLocatedPointLocal);
2309 
2310  G4double shiftOriginSafSq = (fPreviousSftOrigin-pGlobalpoint).mag2();
2311 
2312  // Check that the starting point of this step is
2313  // within the isotropic safety sphere of the last point
2314  // to a accuracy/precision given by fAccuracyForWarning.
2315  // If so give warning.
2316  // If it fails by more than fAccuracyForException exit with error.
2317  //
2318  if( shiftOriginSafSq >= sqr(fPreviousSafety) )
2319  {
2320  G4double shiftOrigin = std::sqrt(shiftOriginSafSq);
2321  G4double diffShiftSaf = shiftOrigin - fPreviousSafety;
2322 
2323  if( diffShiftSaf > fAccuracyForWarning )
2324  {
2325  G4int oldcoutPrec= G4cout.precision(8);
2326  G4int oldcerrPrec= G4cerr.precision(10);
2327  std::ostringstream message, suggestion;
2328  message << "Accuracy error or slightly inaccurate position shift."
2329  << G4endl
2330  << " The Step's starting point has moved "
2331  << std::sqrt(moveLenSq)/mm << " mm " << G4endl
2332  << " since the last call to a Locate method." << G4endl
2333  << " This has resulted in moving "
2334  << shiftOrigin/mm << " mm "
2335  << " from the last point at which the safety "
2336  << " was calculated " << G4endl
2337  << " which is more than the computed safety= "
2338  << fPreviousSafety/mm << " mm at that point." << G4endl
2339  << " This difference is "
2340  << diffShiftSaf/mm << " mm." << G4endl
2341  << " The tolerated accuracy is "
2342  << fAccuracyForException/mm << " mm.";
2343 
2344  suggestion << " ";
2345  static G4ThreadLocal G4int warnNow = 0;
2346  if( ((++warnNow % 100) == 1) )
2347  {
2348  message << G4endl
2349  << " This problem can be due to either " << G4endl
2350  << " - a process that has proposed a displacement"
2351  << " larger than the current safety , or" << G4endl
2352  << " - inaccuracy in the computation of the safety";
2353  suggestion << "We suggest that you " << G4endl
2354  << " - find i) what particle is being tracked, and "
2355  << " ii) through what part of your geometry " << G4endl
2356  << " for example by re-running this event with "
2357  << G4endl
2358  << " /tracking/verbose 1 " << G4endl
2359  << " - check which processes you declare for"
2360  << " this particle (and look at non-standard ones)"
2361  << G4endl
2362  << " - in case, create a detailed logfile"
2363  << " of this event using:" << G4endl
2364  << " /tracking/verbose 6 ";
2365  }
2366  G4Exception("G4ITNavigator2::ComputeStep()",
2367  "GeomNav1002", JustWarning,
2368  message, G4String(suggestion.str()));
2369  G4cout.precision(oldcoutPrec);
2370  G4cerr.precision(oldcerrPrec);
2371  }
2372 #ifdef G4DEBUG_NAVIGATION
2373  else
2374  {
2375  G4cerr << "WARNING - G4ITNavigator2::ComputeStep()" << G4endl
2376  << " The Step's starting point has moved "
2377  << std::sqrt(moveLenSq) << "," << G4endl
2378  << " which has taken it to the limit of"
2379  << " the current safety. " << G4endl;
2380  }
2381 #endif
2382  }
2383  G4double safetyPlus = fPreviousSafety + fAccuracyForException;
2384  if ( shiftOriginSafSq > sqr(safetyPlus) )
2385  {
2386  std::ostringstream message;
2387  message << "May lead to a crash or unreliable results." << G4endl
2388  << " Position has shifted considerably without"
2389  << " notifying the navigator !" << G4endl
2390  << " Tolerated safety: " << safetyPlus << G4endl
2391  << " Computed shift : " << shiftOriginSafSq;
2392  G4Exception("G4ITNavigator2::ComputeStep()", "GeomNav1002",
2393  JustWarning, message);
2394  }
2395 }
2396 
2397 // ********************************************************************
2398 // Operator <<
2399 // ********************************************************************
2400 //
2401 std::ostream& operator << (std::ostream &os,const G4ITNavigator2 &n)
2402 {
2403  // Old version did only the following:
2404  // os << "Current History: " << G4endl << n.fHistory;
2405  // Old behaviour is recovered for fVerbose = 0
2406 
2407  // Adapted from G4ITNavigator2::PrintState() const
2408 
2409  G4int oldcoutPrec = os.precision(4);
2410  if( n.fVerbose >= 4 )
2411  {
2412  os << "The current state of G4ITNavigator2 is: " << G4endl;
2413  os << " ValidExitNormal= " << n.fValidExitNormal << G4endl
2414  << " ExitNormal = " << n.fExitNormal << G4endl
2415  << " Exiting = " << n.fExiting << G4endl
2416  << " Entering = " << n.fEntering << G4endl
2417  << " BlockedPhysicalVolume= " ;
2418  if (n.fBlockedPhysicalVolume==0)
2419  os << "None";
2420  else
2421  os << n.fBlockedPhysicalVolume->GetName();
2422  os << G4endl
2423  << " BlockedReplicaNo = " << n.fBlockedReplicaNo << G4endl
2424  << " LastStepWasZero = " << n.fLastStepWasZero << G4endl
2425  << G4endl;
2426  }
2427  if( ( 1 < n.fVerbose) && (n.fVerbose < 4) )
2428  {
2429  os << G4endl; // Make sure to line up
2430  os << std::setw(30) << " ExitNormal " << " "
2431  << std::setw( 5) << " Valid " << " "
2432  << std::setw( 9) << " Exiting " << " "
2433  << std::setw( 9) << " Entering" << " "
2434  << std::setw(15) << " Blocked:Volume " << " "
2435  << std::setw( 9) << " ReplicaNo" << " "
2436  << std::setw( 8) << " LastStepZero " << " "
2437  << G4endl;
2438  os << "( " << std::setw(7) << n.fExitNormal.x()
2439  << ", " << std::setw(7) << n.fExitNormal.y()
2440  << ", " << std::setw(7) << n.fExitNormal.z() << " ) "
2441  << std::setw( 5) << n.fValidExitNormal << " "
2442  << std::setw( 9) << n.fExiting << " "
2443  << std::setw( 9) << n.fEntering << " ";
2444  if ( n.fBlockedPhysicalVolume==0 )
2445  { os << std::setw(15) << "None"; }
2446  else
2447  { os << std::setw(15)<< n.fBlockedPhysicalVolume->GetName(); }
2448  os << std::setw( 9) << n.fBlockedReplicaNo << " "
2449  << std::setw( 8) << n.fLastStepWasZero << " "
2450  << G4endl;
2451  }
2452  if( n.fVerbose > 2 )
2453  {
2454  os.precision(8);
2455  os << " Current Localpoint = " << n.fLastLocatedPointLocal << G4endl;
2456  os << " PreviousSftOrigin = " << n.fPreviousSftOrigin << G4endl;
2457  os << " PreviousSafety = " << n.fPreviousSafety << G4endl;
2458  }
2459  if( n.fVerbose > 3 || n.fVerbose == 0 )
2460  {
2461  os << "Current History: " << G4endl << n.fHistory;
2462  }
2463 
2464  os.precision(oldcoutPrec);
2465  return os;
2466 }
2467 
2468 //------------------------------------------------------------------------------
2469 
2470 EInside
2472 {
2473  const G4AffineTransform& transform = GetGlobalToLocalTransform();
2474  G4ThreeVector localPoint(transform.TransformPoint(globalPoint));
2475 
2476  G4VSolid* solid = fHistory.GetTopVolume()->GetLogicalVolume()->GetSolid();
2477 
2478  return solid->Inside(localPoint);
2479 }
2480 
2481 //------------------------------------------------------------------------------
2482 
2484 {
2485  const G4AffineTransform& local2Global = GetLocalToGlobalTransform();
2486  G4VSolid* solid = fHistory.GetTopVolume()->GetLogicalVolume()->GetSolid();
2487 
2488  G4VoxelLimits voxelLimits; // Defaults to "infinite" limits.
2489  G4double vmin, vmax;
2490  G4AffineTransform dummy;
2491  std::vector<std::vector<double> > fExtend;
2492 
2493  solid->CalculateExtent(kXAxis,voxelLimits,dummy,vmin,vmax);
2494  fExtend[kXAxis][BoundingBox::kMin] = vmin;
2495  fExtend[kXAxis][BoundingBox::kMax] = vmax;
2496 
2497  solid->CalculateExtent(kYAxis,voxelLimits,dummy,vmin,vmax);
2498  fExtend[kYAxis][BoundingBox::kMin] = vmin;
2499  fExtend[kYAxis][BoundingBox::kMax] = vmax;
2500 
2501  solid->CalculateExtent(kZAxis,voxelLimits,dummy,vmin,vmax);
2502  fExtend[kZAxis][BoundingBox::kMin] = vmin;
2503  fExtend[kZAxis][BoundingBox::kMax] = vmax;
2504 
2505  G4ThreeVector rndmPos;
2506 
2507  while(1)
2508  {
2509  for(size_t i = 0 ; i < 3 ; ++i)
2510  {
2511  double min = fExtend[i][BoundingBox::kMin];
2512  double max = fExtend[i][BoundingBox::kMax];
2513  rndmPos[i] = G4UniformRand()*(max-min)+min;
2514  }
2515 
2516  if(solid->Inside(rndmPos) == kInside) break;
2517  }
2518 
2519  _rndmPoint = local2Global.TransformPoint(rndmPos);
2520 }
2521 
G4String GetName() const
#define fBlockedPhysicalVolume
G4SmartVoxelHeader * GetVoxelHeader() const
#define fValidExitNormal
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const =0
G4ITNavigatorState_Lock2 * GetNavigatorState()
virtual G4bool LevelLocate(G4NavigationHistory &history, const G4VPhysicalVolume *blockedVol, const G4int blockedNum, const G4ThreeVector &globalPoint, const G4ThreeVector *globalDirection, const G4bool pLocatedOnEdge, G4ThreeVector &localPoint)
#define fLastStepEndPointLocal
G4double ComputeStep(const G4ThreeVector &globalPoint, const G4ThreeVector &globalDirection, const G4double currentProposedStepLength, G4double &newSafety, G4NavigationHistory &history, G4bool &validExitNormal, G4ThreeVector &exitNormal, G4bool &exiting, G4bool &entering, G4VPhysicalVolume *(*pBlockedPhysical), G4int &blockedReplicaNo)
const G4ThreeVector & GetTranslation() const
static const G4double kInfinity
Definition: geomdefs.hh:42
virtual G4bool IsNested() const
virtual G4double ComputeStep(const G4ThreeVector &pGlobalPoint, const G4ThreeVector &pDirection, const G4double pCurrentProposedStepLength, G4double &pNewSafety)
#define fExiting
virtual void SetupHierarchy()
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
void UpdateMaterial(G4Material *pMaterial)
CLHEP::Hep3Vector G4ThreeVector
CLHEP::HepRotation G4RotationMatrix
#define fLocatedOnEdge
virtual void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
Definition: G4VSolid.cc:137
virtual G4VSolid * ComputeSolid(const G4int, G4VPhysicalVolume *)
void CheckNavigatorState() const
G4VoxelSafety * fpVoxelSafety
G4double GetSurfaceTolerance() const
void ResetNavigatorState()
void ResetStackAndState()
std::ostream & operator<<(std::ostream &os, const G4ITNavigator2 &n)
G4double ComputeSafety(const G4ThreeVector &localPoint, const G4VPhysicalVolume &currentPhysical, G4double maxLength=DBL_MAX)
G4TouchableHistory * CreateTouchableHistory() const
G4VPhysicalVolume * fTopPhysical
G4double CheckNextStep(const G4ThreeVector &pGlobalPoint, const G4ThreeVector &pDirection, const G4double pCurrentProposedStepLength, G4double &pNewSafety)
#define fLastTriedStepComputation
void SetSolid(G4VSolid *pSolid)
G4NavigatorState * fpNavigatorState
G4SmartVoxelNode * ParamVoxelLocate(G4SmartVoxelHeader *pHead, const G4ThreeVector &localPoint)
G4double ComputeStep(const G4ThreeVector &globalPoint, const G4ThreeVector &globalDirection, const G4ThreeVector &localPoint, const G4ThreeVector &localDirection, const G4double currentProposedStepLength, G4double &newSafety, G4NavigationHistory &history, G4bool &validExitNormal, G4bool &calculatedExitNormal, G4ThreeVector &exitNormal, G4bool &exiting, G4bool &entering, G4VPhysicalVolume *(*pBlockedPhysical), G4int &blockedReplicaNo)
virtual G4GeometryType GetEntityType() const =0
#define fWasLimitedByGeometry
#define G4ThreadLocal
Definition: tls.hh:89
G4VPhysicalVolume * NewNavigatorStateAndLocate(const G4ThreeVector &p, const G4ThreeVector &direction)
#define fExitNormalGlobalFrame
virtual ~G4ITNavigator2()
G4ReferenceCountedHandle< G4TouchableHistory > G4TouchableHistoryHandle
G4bool LevelLocate(G4NavigationHistory &history, const G4VPhysicalVolume *blockedVol, const G4int blockedNum, const G4ThreeVector &globalPoint, const G4ThreeVector *globalDirection, const G4bool pLocatedOnEdge, G4ThreeVector &localPoint)
void ComputeTransformation(const G4int replicaNo, G4VPhysicalVolume *pVol, G4ThreeVector &point) const
int G4int
Definition: G4Types.hh:78
G4int fAbandonThreshold_NoZeroSteps
#define fExitNormal
#define fPushed
#define fLastStepWasZero
G4double ComputeSafety(const G4ThreeVector &globalpoint, const G4NavigationHistory &history, const G4double pMaxLength=DBL_MAX)
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:77
G4ThreeVector ComputeLocalAxis(const G4ThreeVector &pVec) const
void SetNormalNavigation(G4NormalNavigation *fnormnav)
void GetRandomInCurrentVolume(G4ThreeVector &rndmPoint) const
G4double ComputeSafety(const G4ThreeVector &localPoint, const G4NavigationHistory &history, const G4double pProposedMaxLength=DBL_MAX)
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
const G4String & GetName() const
virtual EInside Inside(const G4ThreeVector &p) const =0
#define fEnteredDaughter
G4double ComputeStep(const G4ThreeVector &globalPoint, const G4ThreeVector &globalDirection, const G4double currentProposedStepLength, G4double &newSafety, G4NavigationHistory &history, G4bool &validExitNormal, G4ThreeVector &exitNormal, G4bool &exiting, G4bool &entering, G4VPhysicalVolume *(*pBlockedPhysical), G4int &blockedReplicaNo)
virtual void LocateGlobalPointWithinVolume(const G4ThreeVector &position)
#define fEntering
bool G4bool
Definition: G4Types.hh:79
virtual G4int GetRegularStructureId() const =0
virtual G4VPhysicalVolume * LocateGlobalPointAndSetup(const G4ThreeVector &point, const G4ThreeVector *direction=0, const G4bool pRelativeSearch=true, const G4bool ignoreDirection=true)
virtual G4double ComputeSafety(const G4ThreeVector &globalpoint, const G4NavigationHistory &history, const G4double pMaxLength=DBL_MAX)
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const =0
virtual G4VPVParameterisation * GetParameterisation() const =0
#define fChangedGrandMotherRefFrame
G4RegularNavigation fregularNav
G4NormalNavigation fnormalNav
G4bool LevelLocate(G4NavigationHistory &history, const G4VPhysicalVolume *blockedVol, const G4int blockedNum, const G4ThreeVector &globalPoint, const G4ThreeVector *globalDirection, const G4bool pLocatedOnEdge, G4ThreeVector &localPoint)
EVolume CharacteriseDaughters() const
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
#define fCalculatedExitNormal
EVolume VolumeType(const G4VPhysicalVolume *pVol) const
G4ReplicaNavigation freplicaNav
G4double kCarTolerance
#define G4DEBUG_NAVIGATION
const G4int n
virtual G4ThreeVector GetGlobalExitNormal(const G4ThreeVector &point, G4bool *valid)
#define CheckNavigatorStateIsValid()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define fStepEndPoint
G4ParameterisedNavigation fparamNav
#define fHistory
virtual G4bool RecheckDistanceToCurrentBoundary(const G4ThreeVector &pGlobalPoint, const G4ThreeVector &pDirection, const G4double CurrentProposedStepLength, G4double *prDistance, G4double *prNewSafety=0) const
#define fGrandMotherExitNormal
static const double perMillion
Definition: G4SIunits.hh:331
#define fBlockedReplicaNo
EInside InsideCurrentVolume(const G4ThreeVector &globalPoint) const
G4ThreeVector TransformPoint(const G4ThreeVector &vec) const
G4int fActionThreshold_NoZeroSteps
void PrintState() const
G4double ComputeStep(const G4ThreeVector &localPoint, const G4ThreeVector &localDirection, const G4double currentProposedStepLength, G4double &newSafety, G4NavigationHistory &history, G4bool &validExitNormal, G4ThreeVector &exitNormal, G4bool &exiting, G4bool &entering, G4VPhysicalVolume *(*pBlockedPhysical), G4int &blockedReplicaNo)
G4LogicalVolume * GetLogicalVolume() const
G4ThreeVector ComputeLocalPoint(const G4ThreeVector &rGlobPoint) const
G4VoxelNavigation fvoxelNav
#define fPreviousSftOrigin
T max(const T t1, const T t2)
brief Return the largest of the two arguments
EInside
Definition: geomdefs.hh:58
virtual G4ThreeVector GetLocalExitNormal(G4bool *valid)
virtual void ComputeTransformation(const G4int, G4VPhysicalVolume *) const =0
#define fPreviousSafety
G4ThreeVector TransformAxis(const G4ThreeVector &axis) const
void SetNavigatorState(G4ITNavigatorState_Lock2 *)
virtual G4int GetCopyNo() const =0
void ComputeStepLog(const G4ThreeVector &pGlobalpoint, G4double moveLenSq) const
const G4AffineTransform & GetGlobalToLocalTransform() const
G4int GetDaughtersRegularStructureId(const G4LogicalVolume *pLog) const
const G4RotationMatrix * GetRotation() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
Definition: Step.hh:41
G4SmartVoxelNode * VoxelLocate(G4SmartVoxelHeader *pHead, const G4ThreeVector &localPoint)
#define fExitedMother
virtual G4double ComputeSafety(const G4ThreeVector &globalpoint, const G4double pProposedMaxLength=DBL_MAX, const G4bool keepState=true)
G4int GetVerboseLevel() const
G4bool EnteredDaughterVolume() const
G4double ComputeSafety(const G4ThreeVector &localPoint, const G4NavigationHistory &history, const G4double pProposedMaxLength=DBL_MAX)
#define G4endl
Definition: G4ios.hh:61
#define fLocatedOutsideWorld
virtual G4ThreeVector GetLocalExitNormalAndCheck(const G4ThreeVector &point, G4bool *valid)
const G4NavigationHistory * GetHistory() const
T sqr(const T &x)
Definition: templates.hh:145
virtual void ResetState()
EVolume
Definition: geomdefs.hh:68
G4int MoveUpHistory(G4int num_levels=1)
const G4String & GetName() const
virtual G4VPhysicalVolume * ResetHierarchyAndLocate(const G4ThreeVector &point, const G4ThreeVector &direction, const G4TouchableHistory &h)
G4double ComputeSafety(const G4ThreeVector &globalPoint, const G4ThreeVector &localPoint, G4NavigationHistory &history, const G4double pProposedMaxLength=DBL_MAX)
G4bool LevelLocate(G4NavigationHistory &history, const G4VPhysicalVolume *blockedVol, const G4int blockedNum, const G4ThreeVector &globalPoint, const G4ThreeVector *globalDirection, const G4bool pLocatedOnEdge, G4ThreeVector &localPoint)
#define fLastLocatedPointLocal
double G4double
Definition: G4Types.hh:76
G4bool LevelLocate(G4NavigationHistory &history, const G4VPhysicalVolume *blockedVol, const G4int blockedNum, const G4ThreeVector &globalPoint, const G4ThreeVector *globalDirection, const G4bool pLocatedOnEdge, G4ThreeVector &localPoint)
void NewNavigatorState()
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const =0
virtual G4double ComputeStep(const G4ThreeVector &globalPoint, const G4ThreeVector &globalDirection, const G4double currentProposedStepLength, G4double &newSafety, G4NavigationHistory &history, G4bool &validExitNormal, G4ThreeVector &exitNormal, G4bool &exiting, G4bool &entering, G4VPhysicalVolume *(*pBlockedPhysical), G4int &blockedReplicaNo)
void SetVerboseLevel(G4int level)
EInside BackLocate(G4NavigationHistory &history, const G4ThreeVector &globalPoint, G4ThreeVector &localPoint, const G4bool &exiting, G4bool &notKnownInside) const
#define DBL_MAX
Definition: templates.hh:83
const G4AffineTransform GetLocalToGlobalTransform() const
static const double mm
Definition: G4SIunits.hh:114
virtual G4bool CheckOverlaps(G4int res=1000, G4double tol=0., G4bool verbose=true, G4int errMax=1)
static G4GeometryTolerance * GetInstance()
G4VSolid * GetSolid() const
#define fNumberZeroSteps
virtual G4TouchableHistoryHandle CreateTouchableHistoryHandle() const
G4GLOB_DLL std::ostream G4cerr
EVolume CharacteriseDaughters(const G4LogicalVolume *pLog) const
G4AffineTransform GetMotherToDaughterTransform(G4VPhysicalVolume *dVolume, G4int dReplicaNo, EVolume dVolumeType)