source: trunk/source/geometry/navigation/src/G4MultiNavigator.cc @ 1347

Last change on this file since 1347 was 1347, checked in by garnier, 13 years ago

geant4 tag 9.4
File size: 21.1 KB
Line 
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: G4MultiNavigator.cc,v 1.11 2010/09/06 09:49:15 gcosmo Exp $
28// GEANT4 tag $ Name:  $
29//
30// class G4PathFinder Implementation
31//
32// Author:  John Apostolakis, November 2006
33// --------------------------------------------------------------------
34
35#include "G4MultiNavigator.hh"
36
37class G4FieldManager;
38
39#include "G4Navigator.hh"
40#include "G4PropagatorInField.hh"
41#include "G4TransportationManager.hh"
42
43#include <iomanip>
44
45// ********************************************************************
46// Constructor
47// ********************************************************************
48//
49G4MultiNavigator::G4MultiNavigator() 
50  : G4Navigator(), fLastMassWorld(0)
51{
52  fNoActiveNavigators= 0; 
53  G4ThreeVector Big3Vector( kInfinity, kInfinity, kInfinity ); 
54  fLastLocatedPosition = Big3Vector;
55  fSafetyLocation  = Big3Vector;
56  fPreStepLocation = Big3Vector;
57
58  fMinSafety_PreStepPt=  -1.0; 
59  fMinSafety_atSafLocation= -1.0; 
60  fMinSafety= -kInfinity; 
61  fTrueMinStep= fMinStep= -kInfinity; 
62
63  for(register int num=0; num< fMaxNav; ++num )
64  {
65    fpNavigator[num] =  0;   
66    fLimitTruth[num] = false;
67    fLimitedStep[num] = kUndefLimited;
68    fCurrentStepSize[num] = fNewSafety[num] = -1.0; 
69    fLocatedVolume[num] = 0; 
70  }
71
72  pTransportManager= G4TransportationManager::GetTransportationManager();
73
74  G4Navigator* massNav= pTransportManager->GetNavigatorForTracking();
75  if( massNav )
76  { 
77    G4VPhysicalVolume* pWorld= massNav->GetWorldVolume(); 
78    if( pWorld )
79    { 
80      SetWorldVolume( pWorld ); 
81      fLastMassWorld = pWorld; 
82    }
83  }
84}
85
86G4MultiNavigator::~G4MultiNavigator() 
87{
88}
89
90G4double G4MultiNavigator::ComputeStep(const G4ThreeVector &pGlobalPoint,
91                                       const G4ThreeVector &pDirection,
92                                       const G4double       proposedStepLength,
93                                             G4double      &pNewSafety)
94{
95  G4double safety= 0.0, step=0.0;
96  G4double minSafety= kInfinity, minStep= kInfinity;
97
98#ifdef G4DEBUG_NAVIGATION
99  if( fVerbose > 2 )
100  {
101    G4cout << " G4MultiNavigator::ComputeStep : entered " << G4endl;
102    G4cout << "   Input position= " << pGlobalPoint
103           << "   direction= "      << pDirection         << G4endl;
104    G4cout << "   Requested step= " << proposedStepLength << G4endl;
105  }
106#endif
107
108  std::vector<G4Navigator*>::iterator pNavigatorIter;
109
110  pNavigatorIter= pTransportManager-> GetActiveNavigatorsIterator();
111
112  G4ThreeVector initialPosition = pGlobalPoint;
113  G4ThreeVector initialDirection= pDirection;
114
115  for( register int num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num )
116  {
117     safety= kInfinity;
118
119     step= (*pNavigatorIter)->ComputeStep( initialPosition, 
120                                           initialDirection,
121                                           proposedStepLength,
122                                           safety ); 
123     if( safety < minSafety ){ minSafety = safety; } 
124     if( step < minStep )    { minStep= step; } 
125
126     fCurrentStepSize[num] = step; 
127     fNewSafety[num]= safety; 
128      // This is currently the safety from the last sub-step
129
130#ifdef G4DEBUG_NAVIGATION
131     if( fVerbose > 2 )
132     {
133       G4cout << "G4MultiNavigator::ComputeStep : Navigator ["
134              << num << "] -- step size " << step
135              << " safety= " << safety << G4endl;
136     }
137#endif
138  } 
139
140  // Save safety value, related position
141  //
142  fPreStepLocation     = initialPosition; 
143  fMinSafety_PreStepPt = minSafety;
144  fMinStep = minStep; 
145
146  if( fMinStep == kInfinity )
147  {
148     fTrueMinStep = proposedStepLength;   //  Use this below for endpoint !!
149  }
150  else
151  {
152     fTrueMinStep = minStep;
153  }
154
155#ifdef G4DEBUG_NAVIGATION
156  if( fVerbose > 1 )
157  {
158    G4ThreeVector endPosition = initialPosition+fTrueMinStep*initialDirection;
159
160    G4int oldPrec = G4cout.precision(8); 
161    G4cout << "G4MultiNavigator::ComputeStep : "
162           << " initialPosition = " << initialPosition
163           << " and endPosition = " << endPosition << G4endl;
164    G4cout.precision( oldPrec );
165  }
166#endif
167
168  pNewSafety = minSafety; 
169
170  this->WhichLimited(); 
171
172#ifdef G4DEBUG_NAVIGATION
173  if( fVerbose > 2 )
174  {
175    G4cout << " G4MultiNavigator::ComputeStep : exits returning "
176           << minStep << G4endl;
177  }
178#endif
179
180  return minStep;  // must return kInfinity if do not limit step
181}
182
183// ----------------------------------------------------------------------
184
185G4double
186G4MultiNavigator::ObtainFinalStep( G4int     navigatorId, 
187                                   G4double &pNewSafety,  // for this geometry
188                                   G4double &minStep,
189                                   ELimited &limitedStep) 
190{
191  if( navigatorId > fNoActiveNavigators )
192  { 
193     G4cerr << "ERROR - G4MultiNavigator::ObtainFinalStep()"
194            << "        Navigator Id = " << navigatorId
195            << "        No Active = " << fNoActiveNavigators << " . "
196            << G4endl;
197     G4Exception("G4MultiNavigator::ObtainFinalStep()", "InvalidSetup",
198                 FatalException, "Bad Navigator Id" ); 
199  }
200
201  // Prepare the information to return
202  //
203  pNewSafety  = fNewSafety[ navigatorId ]; 
204  limitedStep = fLimitedStep[ navigatorId ];
205  minStep= fMinStep; 
206
207#ifdef G4DEBUG_NAVIGATION
208  if( fVerbose > 1 )
209  { 
210     G4cout << " G4MultiNavigator::ComputeStep returns "
211            << fCurrentStepSize[ navigatorId ]
212            << " for Navigator " << navigatorId
213            << " Limited step = " << limitedStep
214            << " Safety(mm) = " << pNewSafety / mm << G4endl; 
215  }
216#endif
217
218  return fCurrentStepSize[ navigatorId ];
219}
220
221// ----------------------------------------------------------------------
222
223void G4MultiNavigator::PrepareNewTrack( const G4ThreeVector position, 
224                                        const G4ThreeVector direction )
225{
226#ifdef G4DEBUG_NAVIGATION
227  if( fVerbose > 1 )
228  {
229    G4cout << " Entered G4MultiNavigator::PrepareNewTrack() " << G4endl;
230  }
231#endif
232
233  G4MultiNavigator::PrepareNavigators(); 
234
235  LocateGlobalPointAndSetup( position, &direction, false, false );   
236  //
237  // The first location for each Navigator must be non-relative
238  // or else call ResetStackAndState() for each Navigator
239  // Use direction to get correct side of boundary (ignore dir= false)
240}
241
242// ----------------------------------------------------------------------
243
244void G4MultiNavigator::PrepareNavigators()
245{
246  // Key purposes:
247  //   - Check and cache set of active navigators
248  //   - Reset state for new track
249
250#ifdef G4DEBUG_NAVIGATION
251  if( fVerbose > 1 )
252  {
253    G4cout << " Entered G4MultiNavigator::PrepareNavigators() " << G4endl;
254  }
255#endif
256
257  // Message the transportation-manager to find active navigators
258
259  std::vector<G4Navigator*>::iterator pNavigatorIter; 
260  fNoActiveNavigators=  pTransportManager-> GetNoActiveNavigators();
261
262  if( fNoActiveNavigators > fMaxNav )
263  {
264    G4cerr << "ERROR - G4MultiNavigator::PrepareNavigators()"
265           << "        Too many active Navigators (worlds): "
266           << fNoActiveNavigators << G4endl
267           << "        which is more than the number allowed: "
268           << fMaxNav << " !" << G4endl;
269    G4Exception("G4MultiNavigator::PrepareNavigators()", "InvalidSetup", 
270                FatalException, "Too many active Navigators / worlds !"); 
271  }
272
273  pNavigatorIter= pTransportManager-> GetActiveNavigatorsIterator();
274  for( register int num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num )
275  {
276     fpNavigator[num] =  *pNavigatorIter;   
277     fLimitTruth[num] = false;
278     fLimitedStep[num] = kDoNot;
279     fCurrentStepSize[num] = 0.0; 
280     fLocatedVolume[num] = 0; 
281  }
282  fWasLimitedByGeometry = false; 
283
284  // Check the world volume of the mass navigator
285  // in case a call to SetWorldVolume() changed it
286
287  G4VPhysicalVolume* massWorld = GetWorldVolume();
288
289  if( (massWorld != fLastMassWorld) && (massWorld!=0) )
290  { 
291     // Pass along change to Mass Navigator
292     fpNavigator[0] -> SetWorldVolume( massWorld );
293
294#ifdef G4DEBUG_NAVIGATION
295     if( fVerbose > 0 )
296     { 
297       G4cout << " G4MultiNavigator::PrepareNavigators() changed world volume " 
298              << " for mass geometry to " << massWorld->GetName() << G4endl; 
299     }
300#endif
301
302     fLastMassWorld = massWorld;
303  }
304}
305
306// ----------------------------------------------------------------------
307
308G4VPhysicalVolume* 
309G4MultiNavigator::LocateGlobalPointAndSetup(const G4ThreeVector& position,
310                                            const G4ThreeVector* pDirection,
311                                            const G4bool pRelativeSearch,
312                                            const G4bool ignoreDirection )
313{
314  // Locate the point in each geometry
315
316  G4ThreeVector direction(0.0, 0.0, 0.0);
317  G4bool relative = pRelativeSearch; 
318  std::vector<G4Navigator*>::iterator pNavIter
319    = pTransportManager->GetActiveNavigatorsIterator(); 
320
321  if( pDirection ) { direction = *pDirection; }
322
323#ifdef G4DEBUG_NAVIGATION
324  if( fVerbose > 2 )
325  {
326    G4cout << " Entered G4MultiNavigator::LocateGlobalPointAndSetup() "
327           << G4endl;
328    G4cout << "   Locating at position: " << position
329           << ", with direction: " << direction << G4endl
330           << "   Relative: " << relative
331           << ", ignore direction: " << ignoreDirection << G4endl;
332    G4cout << "   Number of active navigators: " << fNoActiveNavigators
333           << G4endl;
334  }
335#endif
336
337  for ( register int num=0; num< fNoActiveNavigators ; ++pNavIter,++num )
338  {
339     if( fWasLimitedByGeometry && fLimitTruth[num] )
340     { 
341        (*pNavIter)->SetGeometricallyLimitedStep(); 
342     }
343
344     G4VPhysicalVolume *pLocated
345       = (*pNavIter)->LocateGlobalPointAndSetup( position, &direction,
346                                                 relative, ignoreDirection );   
347     // Set the state related to the location
348     //
349     fLocatedVolume[num] = pLocated; 
350
351     // Clear state related to the step
352     //
353     fLimitedStep[num]   = kDoNot; 
354     fCurrentStepSize[num] = 0.0;     
355     fLimitTruth[ num ] = false;   // Always clear on locating (see Navigator)
356   
357#ifdef G4DEBUG_NAVIGATION
358     if( fVerbose > 2 )
359     {
360       G4cout << " Located in world: " << num << ", at: " << position << G4endl
361              << " Used geomLimStp: " << fLimitTruth[num]
362              << ", found in volume: " << pLocated << G4endl; 
363       G4cout << " Name = '" ;       
364       if( pLocated )
365       { 
366         G4cout << pLocated->GetName() << "'"; 
367         G4cout << " - CopyNo= " << pLocated->GetCopyNo(); 
368       }
369       else
370       { 
371         G4cout <<  "Null'   Id: Not-Set "; 
372       }
373       G4cout << G4endl; 
374     }
375#endif
376  }
377
378  fWasLimitedByGeometry = false;   // Clear on locating
379  G4VPhysicalVolume* volMassLocated= fLocatedVolume[0]; 
380
381  return volMassLocated;
382}
383
384// ----------------------------------------------------------------------
385
386void
387G4MultiNavigator::LocateGlobalPointWithinVolume(const G4ThreeVector& position)
388{
389  // Relocate the point in each geometry
390
391  std::vector<G4Navigator*>::iterator pNavIter
392    = pTransportManager->GetActiveNavigatorsIterator(); 
393
394#ifdef G4DEBUG_NAVIGATION
395  if( fVerbose > 2 )
396  {
397    G4cout << " Entered G4MultiNavigator::ReLocate() " << G4endl
398           << "  Re-locating at position: " << position  << G4endl; 
399  }
400#endif
401
402  for ( register int num=0; num< fNoActiveNavigators ; ++pNavIter,++num )
403  {
404     //  ... none limited the step
405
406     (*pNavIter)->LocateGlobalPointWithinVolume( position ); 
407
408     // Clear state related to the step
409     //
410     fLimitedStep[num]     = kDoNot; 
411     fCurrentStepSize[num] = 0.0;     
412
413     fLimitTruth[ num ] = false;   // Always clear on locating (see Navigator)
414  }
415  fWasLimitedByGeometry = false;   // Clear on locating
416  fLastLocatedPosition  = position; 
417}
418
419// ----------------------------------------------------------------------
420
421G4double G4MultiNavigator::ComputeSafety( const G4ThreeVector& position,
422                                          const G4double       maxDistance,
423                                          const G4bool         state)
424{
425    // Recompute safety for the relevant point
426
427    G4double minSafety = kInfinity, safety = kInfinity;
428 
429    std::vector<G4Navigator*>::iterator pNavigatorIter;
430    pNavigatorIter= pTransportManager-> GetActiveNavigatorsIterator();
431
432    for( register int num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num )
433    {
434       safety = (*pNavigatorIter)->ComputeSafety( position, maxDistance, state);
435       if( safety < minSafety ) { minSafety = safety; } 
436    } 
437
438    fSafetyLocation = position;
439    fMinSafety_atSafLocation = minSafety;
440
441#ifdef G4DEBUG_NAVIGATION
442    if( fVerbose > 1 )
443    { 
444      G4cout << " G4MultiNavigator::ComputeSafety - returns: " 
445             << minSafety << ", at location: " << position << G4endl;
446    }
447#endif
448    return minSafety; 
449}
450
451// -----------------------------------------------------------------------
452
453G4TouchableHistoryHandle
454G4MultiNavigator::CreateTouchableHistoryHandle() const
455{
456  G4Exception( "G4MultiNavigator::CreateTouchableHistoryHandle()", 
457               "215-TouchableFromWrongNavigator", FatalException, 
458               "Getting a touchable from G4MultiNavigator is not defined."); 
459
460  G4TouchableHistory* touchHist;
461  touchHist= fpNavigator[0] -> CreateTouchableHistory(); 
462
463  G4VPhysicalVolume* locatedVolume= fLocatedVolume[0]; 
464  if( locatedVolume == 0 )
465  {
466    // Workaround to ensure that the touchable is fixed !! // TODO: fix
467    //
468    touchHist->UpdateYourself( locatedVolume, touchHist->GetHistory() );
469  }
470   
471  return G4TouchableHistoryHandle(touchHist); 
472}
473
474// -----------------------------------------------------------------------
475
476void G4MultiNavigator::WhichLimited()
477{
478  // Flag which processes limited the step
479
480  G4int last=-1; 
481  const G4int IdTransport= 0;  // Id of Mass Navigator !!
482  G4int noLimited=0; 
483  ELimited shared= kSharedOther; 
484
485#ifdef G4DEBUG_NAVIGATION
486  if( fVerbose > 2 )
487  {
488    G4cout << " Entered G4MultiNavigator::WhichLimited() " << G4endl;
489  }
490#endif
491
492  // Assume that [IdTransport] is Mass / Transport
493  //
494  G4bool transportLimited = (fCurrentStepSize[IdTransport] == fMinStep)
495                         && ( fMinStep!= kInfinity); 
496  if( transportLimited )
497  { 
498     shared= kSharedTransport;
499  }
500
501  for ( register int num= 0; num < fNoActiveNavigators; num++ )
502  { 
503    G4bool limitedStep;
504
505    G4double step= fCurrentStepSize[num]; 
506
507    limitedStep = ( step == fMinStep ) && ( step != kInfinity); 
508   
509    fLimitTruth[ num ] = limitedStep; 
510    if( limitedStep )
511    {
512      noLimited++; 
513      fLimitedStep[num] = shared;
514      last= num; 
515    }
516    else
517    {
518      fLimitedStep[num] = kDoNot;
519    }
520  }
521  if( (last > -1) && (noLimited == 1 ) )
522  {
523    fLimitedStep[ last ] = kUnique; 
524  }
525}
526
527// -----------------------------------------------------------------------
528
529void
530G4MultiNavigator::PrintLimited()
531{
532  // Report results -- for checking   
533
534  static G4String StrDoNot("DoNot"), StrUnique("Unique"),
535                  StrUndefined("Undefined"),
536                  StrSharedTransport("SharedTransport"),
537                  StrSharedOther("SharedOther");
538  G4cout << "### G4MultiNavigator::PrintLimited() reports: " << G4endl;
539  G4cout << "    Minimum step (true): " << fTrueMinStep
540         << ", reported min: " << fMinStep << G4endl; 
541
542#ifdef G4DEBUG_NAVIGATION
543  if(fVerbose>=2)
544  {
545    G4cout << std::setw(5) << " NavId"  << " "
546           << std::setw(12) << " step-size " << " "
547           << std::setw(12) << " raw-size "  << " "
548           << std::setw(12) << " pre-safety " << " " 
549           << std::setw(15) << " Limited / flag"  << " "
550           << std::setw(15) << "  World "  << " "
551           << G4endl; 
552  }
553#endif
554
555  for ( register int num= 0; num < fNoActiveNavigators; num++ )
556  { 
557    G4double rawStep = fCurrentStepSize[num]; 
558    G4double stepLen = fCurrentStepSize[num]; 
559    if( stepLen > fTrueMinStep )
560    { 
561      stepLen = fTrueMinStep;     // did not limit (went as far as asked)
562    }
563    G4int oldPrec= G4cout.precision(9); 
564
565    G4cout << std::setw(5) << num  << " "
566           << std::setw(12) << stepLen << " "
567           << std::setw(12) << rawStep << " "
568           << std::setw(12) << fNewSafety[num] << " "
569           << std::setw(5) << (fLimitTruth[num] ? "YES" : " NO") << " ";
570    G4String limitedStr;
571    switch ( fLimitedStep[num] )
572    {
573      case kDoNot          : limitedStr= StrDoNot; break;
574      case kUnique         : limitedStr = StrUnique; break; 
575      case kSharedTransport: limitedStr= StrSharedTransport; break; 
576      case kSharedOther    : limitedStr = StrSharedOther; break;
577      default              : limitedStr = StrUndefined; break;
578    }
579    G4cout << " " << std::setw(15) << limitedStr << " "; 
580    G4cout.precision(oldPrec); 
581
582    G4Navigator *pNav= fpNavigator[ num ]; 
583    G4String  WorldName( "Not-Set" ); 
584    if (pNav)
585    {
586       G4VPhysicalVolume *pWorld= pNav->GetWorldVolume(); 
587       if( pWorld )
588       {
589           WorldName = pWorld->GetName(); 
590       }
591    }
592    G4cout << " " << WorldName ; 
593    G4cout << G4endl;
594  }
595}
596 
597
598// -----------------------------------------------------------------------
599
600void G4MultiNavigator::ResetState()
601{
602   fWasLimitedByGeometry= false; 
603
604   G4Exception("G4MultiNavigator::ResetState()", "217-NotImplemented", 
605               FatalException, 
606               "Cannot call ResetState() for navigators of G4MultiNavigator.");
607   
608   std::vector<G4Navigator*>::iterator pNavigatorIter;
609   pNavigatorIter= pTransportManager-> GetActiveNavigatorsIterator();
610   for( register int num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num )
611   {
612       //  (*pNavigatorIter)->ResetState();  // KEEP THIS comment !!!
613   } 
614}
615
616// -----------------------------------------------------------------------
617
618void G4MultiNavigator::SetupHierarchy()
619{
620  G4Exception( "G4MultiNavigator::SetupHierarchy()", 
621               "217-NotImplemented", FatalException, 
622           "Cannot call SetupHierarchy() for navigators of G4MultiNavigator."); 
623}
624
625// -----------------------------------------------------------------------
626
627void G4MultiNavigator::CheckMassWorld()
628{
629   G4VPhysicalVolume* navTrackWorld=
630     pTransportManager->GetNavigatorForTracking()->GetWorldVolume();
631
632   if( navTrackWorld != fLastMassWorld )
633   { 
634      G4Exception( "G4MultiNavigator::CheckMassWorld()",
635                   "220-InvalidSetup", FatalException, 
636                   "Mass world pointer has been changed." ); 
637   }
638}
639
640// -----------------------------------------------------------------------
641
642G4VPhysicalVolume*
643G4MultiNavigator::ResetHierarchyAndLocate(const G4ThreeVector &point,
644                                          const G4ThreeVector &direction,
645                                          const G4TouchableHistory &MassHistory)
646{
647   // Reset geometry for all -- and use the touchable for the mass history
648
649   G4VPhysicalVolume* massVolume=0; 
650   G4Navigator* pMassNavigator= fpNavigator[0]; 
651
652   if( pMassNavigator )
653   {
654      massVolume= pMassNavigator->ResetHierarchyAndLocate( point, direction,
655                                                           MassHistory); 
656   }
657   else
658   {
659      G4Exception("G4MultiNavigator::ResetHierarchyAndLocate()",
660                  "218-TooEarlyToReset", FatalException,
661                  "Cannot reset hierarchy before navigators are initialised.");
662   }
663
664   std::vector<G4Navigator*>::iterator pNavIter= 
665       pTransportManager->GetActiveNavigatorsIterator(); 
666
667   for ( register int num=0; num< fNoActiveNavigators ; ++pNavIter,++num )
668   {
669      G4bool relativeSearch, ignoreDirection; 
670
671      (*pNavIter)-> LocateGlobalPointAndSetup( point, 
672                                               &direction, 
673                                               relativeSearch=false,
674                                               ignoreDirection=false);
675   }
676   return massVolume; 
677}
Note: See TracBrowser for help on using the repository browser.