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