| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 // class G4MultiNavigator Implementation
27 //
28 // Author: John Apostolakis, November 2006
29 // -------------------------------------------
30
31 #include <iomanip>
32
33 #include "G4MultiNavigator.hh"
34
35 class G4FieldManager;
36
37 #include "G4SystemOfUnits.hh"
38 #include "G4Navigator.hh"
39 #include "G4PropagatorInField.hh"
40 #include "G4TransportationManager.hh"
41
42 // -------------------------------------------
43
44 G4MultiNavigator::G4MultiNavigator()
45 {
46 G4ThreeVector Big3Vector( kInfinity, kInfini
47 fLastLocatedPosition = Big3Vector;
48 fSafetyLocation = Big3Vector;
49 fPreStepLocation = Big3Vector;
50
51 for(auto num=0; num< fMaxNav; ++num )
52 {
53 fpNavigator[num] = nullptr;
54 fLimitTruth[num] = false;
55 fLimitedStep[num] = kUndefLimited;
56 fCurrentStepSize[num] = fNewSafety[num] =
57 fLocatedVolume[num] = nullptr;
58 }
59
60 pTransportManager= G4TransportationManager::
61
62 G4Navigator* massNav= pTransportManager->Get
63 if( massNav != nullptr )
64 {
65 G4VPhysicalVolume* pWorld= massNav->GetWor
66 if( pWorld != nullptr )
67 {
68 SetWorldVolume( pWorld );
69 fLastMassWorld = pWorld;
70 }
71 }
72 }
73
74 // -------------------------------------------
75
76 G4MultiNavigator::~G4MultiNavigator() = defaul
77
78 // -------------------------------------------
79
80 G4double G4MultiNavigator::ComputeStep(const G
81 const G
82 const G
83 G
84 {
85 G4double safety= 0.0, step=0.0;
86 G4double minSafety= kInfinity, minStep= kInf
87
88 fNoLimitingStep= -1;
89 fIdNavLimiting= -1; // Reset for new ste
90
91 #ifdef G4DEBUG_NAVIGATION
92 if( fVerbose > 2 )
93 {
94 G4cout << " G4MultiNavigator::ComputeStep
95 G4cout << " Input position= " << pGlobal
96 << " direction= " << pDirect
97 G4cout << " Requested step= " << propose
98 }
99 #endif
100
101 std::vector<G4Navigator*>::iterator pNavigat
102
103 pNavigatorIter= pTransportManager-> GetActiv
104
105 G4ThreeVector initialPosition = pGlobalPoint
106 G4ThreeVector initialDirection= pDirection;
107
108 for( auto num=0; num< fNoActiveNavigators; +
109 {
110 safety= kInfinity;
111
112 step= (*pNavigatorIter)->ComputeStep( ini
113 ini
114 pro
115 saf
116 if( safety < minSafety ){ minSafety = saf
117 if( step < minStep ) { minStep= step;
118
119 fCurrentStepSize[num] = step;
120 fNewSafety[num]= safety;
121 // This is currently the safety from the
122
123 #ifdef G4DEBUG_NAVIGATION
124 if( fVerbose > 2 )
125 {
126 G4cout << "G4MultiNavigator::ComputeSte
127 << num << "] -- step size " << s
128 << " safety= " << safety << G4en
129 }
130 #endif
131 }
132
133 // Save safety value, related position
134 //
135 fPreStepLocation = initialPosition;
136 fMinSafety_PreStepPt = minSafety;
137 fMinStep = minStep;
138
139 if( fMinStep == kInfinity )
140 {
141 fTrueMinStep = proposedStepLength; //
142 }
143 else
144 {
145 fTrueMinStep = minStep;
146 }
147
148 #ifdef G4DEBUG_NAVIGATION
149 if( fVerbose > 1 )
150 {
151 G4ThreeVector endPosition = initialPositio
152
153 G4int oldPrec = G4cout.precision(8);
154 G4cout << "G4MultiNavigator::ComputeStep :
155 << " initialPosition = " << initial
156 << " and endPosition = " << endPosi
157 G4cout.precision( oldPrec );
158 }
159 #endif
160
161 pNewSafety = minSafety;
162
163 this->WhichLimited();
164
165 #ifdef G4DEBUG_NAVIGATION
166 if( fVerbose > 2 )
167 {
168 G4cout << " G4MultiNavigator::ComputeStep
169 << minStep << G4endl;
170 }
171 #endif
172
173 return minStep; // must return kInfinity if
174 }
175
176 // -------------------------------------------
177
178 G4double
179 G4MultiNavigator::ObtainFinalStep( G4int n
180 G4double& p
181 G4double& m
182 ELimited& l
183 {
184 if( navigatorId > fNoActiveNavigators )
185 {
186 std::ostringstream message;
187 message << "Bad Navigator Id!" << G4endl
188 << " Navigator Id = " << n
189 << " No Active = " << fNoA
190 G4Exception("G4MultiNavigator::ObtainFina
191 FatalException, message);
192 }
193
194 // Prepare the information to return
195 //
196 pNewSafety = fNewSafety[ navigatorId ];
197 limitedStep = fLimitedStep[ navigatorId ];
198 minStep= fMinStep;
199
200 #ifdef G4DEBUG_NAVIGATION
201 if( fVerbose > 1 )
202 {
203 G4cout << " G4MultiNavigator::ComputeStep
204 << fCurrentStepSize[ navigatorId ]
205 << " for Navigator " << navigatorI
206 << " Limited step = " << limitedSt
207 << " Safety(mm) = " << pNewSafety
208 }
209 #endif
210
211 return fCurrentStepSize[ navigatorId ];
212 }
213
214 // -------------------------------------------
215
216 void G4MultiNavigator::PrepareNewTrack( const
217 const
218 {
219 #ifdef G4DEBUG_NAVIGATION
220 if( fVerbose > 1 )
221 {
222 G4cout << " Entered G4MultiNavigator::Prep
223 }
224 #endif
225
226 G4MultiNavigator::PrepareNavigators();
227
228 LocateGlobalPointAndSetup( position, &direct
229 //
230 // The first location for each Navigator mus
231 // or else call ResetStackAndState() for eac
232 // Use direction to get correct side of boun
233 }
234
235 // -------------------------------------------
236
237 void G4MultiNavigator::PrepareNavigators()
238 {
239 // Key purposes:
240 // - Check and cache set of active navigat
241 // - Reset state for new track
242
243 #ifdef G4DEBUG_NAVIGATION
244 if( fVerbose > 1 )
245 {
246 G4cout << " Entered G4MultiNavigator::Prep
247 }
248 #endif
249
250 // Message the transportation-manager to fin
251
252 std::vector<G4Navigator*>::const_iterator pN
253 fNoActiveNavigators = (G4int)pTransportManag
254
255 if( fNoActiveNavigators > fMaxNav )
256 {
257 std::ostringstream message;
258 message << "Too many active Navigators / w
259 << " Active Navigators (wor
260 << fNoActiveNavigators << G4endl
261 << " which is more than the
262 << fMaxNav << " !";
263 G4Exception("G4MultiNavigator::PrepareNavi
264 FatalException, message);
265 }
266
267 pNavigatorIter= pTransportManager->GetActive
268 for( auto num=0; num< fNoActiveNavigators; +
269 {
270 fpNavigator[num] = *pNavigatorIter;
271 fLimitTruth[num] = false;
272 fLimitedStep[num] = kDoNot;
273 fCurrentStepSize[num] = 0.0;
274 fLocatedVolume[num] = nullptr;
275 }
276 fWasLimitedByGeometry = false;
277
278 // Check the world volume of the mass naviga
279 // in case a call to SetWorldVolume() change
280
281 G4VPhysicalVolume* massWorld = GetWorldVolum
282
283 if( (massWorld != fLastMassWorld) && (massWo
284 {
285 // Pass along change to Mass Navigator
286 fpNavigator[0] -> SetWorldVolume( massWor
287
288 #ifdef G4DEBUG_NAVIGATION
289 if( fVerbose > 0 )
290 {
291 G4cout << " G4MultiNavigator::PrepareNa
292 << " for mass geometry to " << m
293 }
294 #endif
295
296 fLastMassWorld = massWorld;
297 }
298 }
299
300 // -------------------------------------------
301
302 G4VPhysicalVolume*
303 G4MultiNavigator::LocateGlobalPointAndSetup(co
304 co
305 co
306 co
307 {
308 // Locate the point in each geometry
309
310 G4ThreeVector direction(0.0, 0.0, 0.0);
311 G4bool relative = pRelativeSearch;
312 auto pNavIter = pTransportManager->GetActive
313
314 if( pDirection != nullptr ) { direction = *p
315
316 #ifdef G4DEBUG_NAVIGATION
317 if( fVerbose > 2 )
318 {
319 G4cout << " Entered G4MultiNavigator::Loca
320 << G4endl;
321 G4cout << " Locating at position: " << p
322 << ", with direction: " << directio
323 << " Relative: " << relative
324 << ", ignore direction: " << ignore
325 G4cout << " Number of active navigators:
326 << G4endl;
327 }
328 #endif
329
330 for ( auto num=0; num< fNoActiveNavigators ;
331 {
332 if( fWasLimitedByGeometry && fLimitTruth[
333 {
334 (*pNavIter)->SetGeometricallyLimitedSt
335 }
336
337 G4VPhysicalVolume *pLocated
338 = (*pNavIter)->LocateGlobalPointAndSetu
339
340 // Set the state related to the location
341 //
342 fLocatedVolume[num] = pLocated;
343
344 // Clear state related to the step
345 //
346 fLimitedStep[num] = kDoNot;
347 fCurrentStepSize[num] = 0.0;
348 fLimitTruth[ num ] = false; // Always c
349
350 #ifdef G4DEBUG_NAVIGATION
351 if( fVerbose > 2 )
352 {
353 G4cout << " Located in world: " << num
354 << " Used geomLimStp: " << fLimi
355 << ", found in volume: " << pLoc
356 G4cout << " Name = '" ;
357 if( pLocated )
358 {
359 G4cout << pLocated->GetName() << "'";
360 G4cout << " - CopyNo= " << pLocated->
361 }
362 else
363 {
364 G4cout << "Null' Id: Not-Set ";
365 }
366 G4cout << G4endl;
367 }
368 #endif
369 }
370
371 fWasLimitedByGeometry = false; // Clear on
372 G4VPhysicalVolume* volMassLocated = fLocated
373
374 return volMassLocated;
375 }
376
377 // -------------------------------------------
378
379 void
380 G4MultiNavigator::LocateGlobalPointWithinVolum
381 {
382 // Relocate the point in each geometry
383
384 auto pNavIter = pTransportManager->GetActive
385
386 #ifdef G4DEBUG_NAVIGATION
387 if( fVerbose > 2 )
388 {
389 G4cout << " Entered G4MultiNavigator::ReLo
390 << " Re-locating at position: " <<
391 }
392 #endif
393
394 for ( auto num=0; num< fNoActiveNavigators ;
395 {
396 // ... none limited the step
397
398 (*pNavIter)->LocateGlobalPointWithinVolum
399
400 // Clear state related to the step
401 //
402 fLimitedStep[num] = kDoNot;
403 fCurrentStepSize[num] = 0.0;
404
405 fLimitTruth[ num ] = false; // Always c
406 }
407 fWasLimitedByGeometry = false; // Clear on
408 fLastLocatedPosition = position;
409 }
410
411 // -------------------------------------------
412
413 G4double G4MultiNavigator::ComputeSafety( cons
414 cons
415 cons
416 {
417 // Recompute safety for the relevant point
418
419 G4double minSafety = kInfinity, safety = k
420
421 std::vector<G4Navigator*>::iterator pNavig
422 pNavigatorIter= pTransportManager-> GetAct
423
424 for( auto num=0; num< fNoActiveNavigators;
425 {
426 safety = (*pNavigatorIter)->ComputeSafe
427 if( safety < minSafety ) { minSafety =
428 }
429
430 fSafetyLocation = position;
431 fMinSafety_atSafLocation = minSafety;
432
433 #ifdef G4DEBUG_NAVIGATION
434 if( fVerbose > 1 )
435 {
436 G4cout << " G4MultiNavigator::ComputeSaf
437 << minSafety << ", at location: "
438 }
439 #endif
440 return minSafety;
441 }
442
443 // -------------------------------------------
444
445 G4TouchableHandle G4MultiNavigator::CreateTouc
446 {
447 G4Exception( "G4MultiNavigator::CreateToucha
448 "GeomNav0001", FatalException,
449 "Getting a touchable from G4Mul
450
451 G4TouchableHistory* touchHist;
452 touchHist= fpNavigator[0] -> CreateTouchable
453
454 G4VPhysicalVolume* locatedVolume= fLocatedVo
455 if( locatedVolume == nullptr )
456 {
457 // Workaround to ensure that the touchable
458 //
459 touchHist->UpdateYourself( locatedVolume,
460 }
461
462 return {touchHist};
463 }
464
465 // -------------------------------------------
466
467 void G4MultiNavigator::WhichLimited()
468 {
469 // Flag which processes limited the step
470
471 G4int last=-1;
472 const G4int IdTransport= 0; // Id of Mass N
473 G4int noLimited=0;
474 ELimited shared= kSharedOther;
475
476 #ifdef G4DEBUG_NAVIGATION
477 if( fVerbose > 2 )
478 {
479 G4cout << " Entered G4MultiNavigator::Whic
480 }
481 #endif
482
483 // Assume that [IdTransport] is Mass / Trans
484 //
485 G4bool transportLimited = (fCurrentStepSize[
486 && ( fMinStep!= kInfi
487 if( transportLimited )
488 {
489 shared = kSharedTransport;
490 }
491
492 for ( auto num = 0; num < fNoActiveNavigator
493 {
494 G4bool limitedStep;
495
496 G4double step = fCurrentStepSize[num];
497
498 limitedStep = ( step == fMinStep ) && ( st
499
500 fLimitTruth[ num ] = limitedStep;
501 if( limitedStep )
502 {
503 ++noLimited;
504 fLimitedStep[num] = shared;
505 last= num;
506 }
507 else
508 {
509 fLimitedStep[num] = kDoNot;
510 }
511 }
512 if( (last > -1) && (noLimited == 1 ) )
513 {
514 fLimitedStep[ last ] = kUnique;
515 fIdNavLimiting = last;
516 }
517
518 fNoLimitingStep = noLimited;
519
520 return;
521 }
522
523 // -------------------------------------------
524
525 void
526 G4MultiNavigator::PrintLimited()
527 {
528 // Report results -- for checking
529
530 static const G4String StrDoNot("DoNot"), Str
531 StrUndefined("Undefined"),
532 StrSharedTransport("SharedTr
533 StrSharedOther("SharedOther"
534 G4cout << "### G4MultiNavigator::PrintLimite
535 G4cout << " Minimum step (true): " << fTr
536 << ", reported min: " << fMinStep <<
537
538 #ifdef G4DEBUG_NAVIGATION
539 if(fVerbose>=2)
540 {
541 G4cout << std::setw(5) << " NavId" << " "
542 << std::setw(12) << " step-size " <
543 << std::setw(12) << " raw-size " <
544 << std::setw(12) << " pre-safety "
545 << std::setw(15) << " Limited / fla
546 << std::setw(15) << " World " <<
547 << G4endl;
548 }
549 #endif
550
551 for ( auto num = 0; num < fNoActiveNavigator
552 {
553 G4double rawStep = fCurrentStepSize[num];
554 G4double stepLen = fCurrentStepSize[num];
555 if( stepLen > fTrueMinStep )
556 {
557 stepLen = fTrueMinStep; // did not l
558 }
559 G4long oldPrec = G4cout.precision(9);
560
561 G4cout << std::setw(5) << num << " "
562 << std::setw(12) << stepLen << " "
563 << std::setw(12) << rawStep << " "
564 << std::setw(12) << fNewSafety[num]
565 << std::setw(5) << (fLimitTruth[num
566 G4String limitedStr;
567 switch ( fLimitedStep[num] )
568 {
569 case kDoNot : limitedStr = StrD
570 case kUnique : limitedStr = StrU
571 case kSharedTransport: limitedStr = StrS
572 case kSharedOther : limitedStr = StrS
573 default : limitedStr = StrU
574 }
575 G4cout << " " << std::setw(15) << limitedS
576 G4cout.precision(oldPrec);
577
578 G4Navigator *pNav = fpNavigator[ num ];
579 G4String WorldName( "Not-Set" );
580 if (pNav != nullptr)
581 {
582 G4VPhysicalVolume *pWorld = pNav->GetWo
583 if( pWorld != nullptr )
584 {
585 WorldName = pWorld->GetName();
586 }
587 }
588 G4cout << " " << WorldName ;
589 G4cout << G4endl;
590 }
591 }
592
593
594 // -------------------------------------------
595
596 void G4MultiNavigator::ResetState()
597 {
598 fWasLimitedByGeometry = false;
599
600 G4Exception("G4MultiNavigator::ResetState()
601 FatalException,
602 "Cannot reset state for navigat
603 /*
604 std::vector<G4Navigator*>::iterator pNaviga
605 pNavigatorIter = pTransportManager->GetActi
606 for( auto num = 0; num< fNoActiveNavigators
607 {
608 // (*pNavigatorIter)->ResetState(); /
609 }
610 */
611 }
612
613 // -------------------------------------------
614
615 void G4MultiNavigator::SetupHierarchy()
616 {
617 G4Exception( "G4MultiNavigator::SetupHierarc
618 "GeomNav0001", FatalException,
619 "Cannot setup hierarchy for nav
620 }
621
622 // -------------------------------------------
623
624 void G4MultiNavigator::CheckMassWorld()
625 {
626 G4VPhysicalVolume* navTrackWorld =
627 pTransportManager->GetNavigatorForTrackin
628
629 if( navTrackWorld != fLastMassWorld )
630 {
631 G4Exception( "G4MultiNavigator::CheckMas
632 "GeomNav0003", FatalExcepti
633 "Mass world pointer has bee
634 }
635 }
636
637 // -------------------------------------------
638
639 G4VPhysicalVolume*
640 G4MultiNavigator::ResetHierarchyAndLocate(cons
641 cons
642 cons
643 {
644 // Reset geometry for all -- and use the to
645
646 G4VPhysicalVolume* massVolume = nullptr;
647 G4Navigator* pMassNavigator = fpNavigator[0
648
649 if( pMassNavigator != nullptr )
650 {
651 massVolume= pMassNavigator->ResetHierarc
652
653 }
654 else
655 {
656 G4Exception("G4MultiNavigator::ResetHier
657 "GeomNav0002", FatalExceptio
658 "Cannot reset hierarchy befo
659 }
660
661 auto pNavIter= pTransportManager->GetActive
662
663 for ( auto num = 0; num < fNoActiveNavigato
664 {
665 G4bool relativeSearch, ignoreDirection;
666
667 (*pNavIter)-> LocateGlobalPointAndSetup(
668
669
670
671 }
672 return massVolume;
673 }
674
675 // -------------------------------------------
676
677 G4ThreeVector
678 G4MultiNavigator::GetGlobalExitNormal(const G4
679 G4bool*
680 {
681 G4ThreeVector normalGlobalCrd(0.0, 0.0, 0.0)
682 G4bool isObtained = false;
683 // These default values will be used if fNoL
684 G4int firstNavigatorId = -1;
685 G4bool oneObtained = false;
686
687 if( fNoLimitingStep == 1 )
688 {
689 // Only message the Navigator which limite
690 normalGlobalCrd = fpNavigator[ fIdNavLimit
691 ->GetGlobalExitNormal( arg
692 *argpObtained = isObtained;
693 }
694 else
695 {
696 if( fNoLimitingStep > 1 )
697 {
698 auto pNavIter= pTransportManager->GetAct
699
700 for ( auto num = 0; num < fNoActiveNavig
701 {
702 G4ThreeVector oneNormal;
703 if( fLimitTruth[ num ] ) // Did this
704 {
705 G4ThreeVector newNormal =
706 (*pNavIter)->GetGlobalExitNormal(
707 if( oneObtained )
708 {
709 // Keep first one - only if it is
710 if( !isObtained && (newNormal.mag2
711 {
712 normalGlobalCrd = newNormal;
713 isObtained = oneObtained;
714 firstNavigatorId = num;
715 }
716 else
717 {
718 // Check for clash
719 G4double dotNewPrevious = newNor
720 G4double productMagSq = normalGl
721 if( productMagSq > 0.0 )
722 {
723 G4double productMag = std::sqr
724 dotNewPrevious /= productMag;
725 if( dotNewPrevious < (1 - perT
726 {
727 *argpObtained = false;
728
729 if( fVerbose > 2 ) // dotN
730 {
731 std::ostringstream message
732 message << "Clash of Norma
733 << G4endl
734 << " Previo
735 << firstNavigatorI
736 << " Curren
737 << num << G4endl;
738 message << " Dot product
739 << dotNewPrevious
740 message << " Normal
741 << normalGlobalCrd
742 message << " Normal
743 G4Exception("G4MultiNaviga
744 "GeomNav0002",
745 }
746 }
747 else
748 {
749 // Close agreement - Do not
750 }
751 }
752 }
753 }
754 }
755 } // end for over the Navigators
756
757 // Report if no Normal was obtained
758 if( !oneObtained )
759 {
760 std::ostringstream message;
761 message << "No Normal obtained despite
762 << " candidate Navigators limi
763 G4Exception("G4MultiNavigator::GetGlob
764 JustWarning, message);
765 }
766
767 } // end if ( fNoLimiting > 1 )
768 } // end else
769
770 *argpObtained = isObtained;
771 return normalGlobalCrd;
772 }
773
774 // -------------------------------------------
775
776 G4ThreeVector
777 G4MultiNavigator::GetLocalExitNormal(G4bool* a
778 {
779 // If it is the mass navigator, then expect
780 G4ThreeVector normalGlobalCrd(0.0, 0.0, 0.0)
781 G4bool isObtained = false;
782 // These default values will be used if fNoL
783
784 if( fNoLimitingStep == 1 )
785 {
786 // Only message the Navigator which limite
787 normalGlobalCrd = fpNavigator[ fIdNavLimit
788 ->GetLocalExitNormal( &isO
789 *argpObtained = isObtained;
790
791 static G4ThreadLocal G4int numberWarnings
792 G4int noWarningsStart = 10, noModuloWarnin
793 ++numberWarnings;
794 if( (numberWarnings < noWarningsStart )
795 || (numberWarnings%noModuloWarnings == 0)
796 {
797 std::ostringstream message;
798 message << "Cannot obtain normal in loca
799 << "coordinate systems." << G4en
800 G4Exception("G4MultiNavigator::GetGlobal
801 JustWarning, message);
802 }
803 }
804 else
805 {
806 if( fNoLimitingStep > 1 )
807 {
808 std::ostringstream message;
809 message << "Cannot obtain normal in loca
810 << "coordinate systems." << G4en
811 G4Exception("G4MultiNavigator::GetGlobal
812 FatalException, message);
813 }
814 }
815
816 *argpObtained = isObtained;
817 return normalGlobalCrd;
818 }
819
820 // -------------------------------------------
821
822 G4ThreeVector
823 G4MultiNavigator::GetLocalExitNormalAndCheck(c
824
825 {
826 return G4MultiNavigator::GetLocalExitNormal(
827 }
828