Geant4 Cross Reference |
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 // 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, kInfinity, kInfinity ); 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] = -1.0; 57 fLocatedVolume[num] = nullptr; 58 } 59 60 pTransportManager= G4TransportationManager::GetTransportationManager(); 61 62 G4Navigator* massNav= pTransportManager->GetNavigatorForTracking(); 63 if( massNav != nullptr ) 64 { 65 G4VPhysicalVolume* pWorld= massNav->GetWorldVolume(); 66 if( pWorld != nullptr ) 67 { 68 SetWorldVolume( pWorld ); 69 fLastMassWorld = pWorld; 70 } 71 } 72 } 73 74 // ----------------------------------------------------------------------- 75 76 G4MultiNavigator::~G4MultiNavigator() = default; 77 78 // ----------------------------------------------------------------------- 79 80 G4double G4MultiNavigator::ComputeStep(const G4ThreeVector& pGlobalPoint, 81 const G4ThreeVector& pDirection, 82 const G4double proposedStepLength, 83 G4double& pNewSafety) 84 { 85 G4double safety= 0.0, step=0.0; 86 G4double minSafety= kInfinity, minStep= kInfinity; 87 88 fNoLimitingStep= -1; 89 fIdNavLimiting= -1; // Reset for new step 90 91 #ifdef G4DEBUG_NAVIGATION 92 if( fVerbose > 2 ) 93 { 94 G4cout << " G4MultiNavigator::ComputeStep : entered " << G4endl; 95 G4cout << " Input position= " << pGlobalPoint 96 << " direction= " << pDirection << G4endl; 97 G4cout << " Requested step= " << proposedStepLength << G4endl; 98 } 99 #endif 100 101 std::vector<G4Navigator*>::iterator pNavigatorIter; 102 103 pNavigatorIter= pTransportManager-> GetActiveNavigatorsIterator(); 104 105 G4ThreeVector initialPosition = pGlobalPoint; 106 G4ThreeVector initialDirection= pDirection; 107 108 for( auto num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num ) 109 { 110 safety= kInfinity; 111 112 step= (*pNavigatorIter)->ComputeStep( initialPosition, 113 initialDirection, 114 proposedStepLength, 115 safety ); 116 if( safety < minSafety ){ minSafety = safety; } 117 if( step < minStep ) { minStep= step; } 118 119 fCurrentStepSize[num] = step; 120 fNewSafety[num]= safety; 121 // This is currently the safety from the last sub-step 122 123 #ifdef G4DEBUG_NAVIGATION 124 if( fVerbose > 2 ) 125 { 126 G4cout << "G4MultiNavigator::ComputeStep : Navigator [" 127 << num << "] -- step size " << step 128 << " safety= " << safety << G4endl; 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; // Use this below for endpoint !! 142 } 143 else 144 { 145 fTrueMinStep = minStep; 146 } 147 148 #ifdef G4DEBUG_NAVIGATION 149 if( fVerbose > 1 ) 150 { 151 G4ThreeVector endPosition = initialPosition+fTrueMinStep*initialDirection; 152 153 G4int oldPrec = G4cout.precision(8); 154 G4cout << "G4MultiNavigator::ComputeStep : " 155 << " initialPosition = " << initialPosition 156 << " and endPosition = " << endPosition << G4endl; 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 : exits returning " 169 << minStep << G4endl; 170 } 171 #endif 172 173 return minStep; // must return kInfinity if do not limit step 174 } 175 176 // ---------------------------------------------------------------------- 177 178 G4double 179 G4MultiNavigator::ObtainFinalStep( G4int navigatorId, 180 G4double& pNewSafety, // for this geometry 181 G4double& minStep, 182 ELimited& limitedStep) 183 { 184 if( navigatorId > fNoActiveNavigators ) 185 { 186 std::ostringstream message; 187 message << "Bad Navigator Id!" << G4endl 188 << " Navigator Id = " << navigatorId 189 << " No Active = " << fNoActiveNavigators << "."; 190 G4Exception("G4MultiNavigator::ObtainFinalStep()", "GeomNav0002", 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 returns " 204 << fCurrentStepSize[ navigatorId ] 205 << " for Navigator " << navigatorId 206 << " Limited step = " << limitedStep 207 << " Safety(mm) = " << pNewSafety / mm << G4endl; 208 } 209 #endif 210 211 return fCurrentStepSize[ navigatorId ]; 212 } 213 214 // ---------------------------------------------------------------------- 215 216 void G4MultiNavigator::PrepareNewTrack( const G4ThreeVector& position, 217 const G4ThreeVector direction ) 218 { 219 #ifdef G4DEBUG_NAVIGATION 220 if( fVerbose > 1 ) 221 { 222 G4cout << " Entered G4MultiNavigator::PrepareNewTrack() " << G4endl; 223 } 224 #endif 225 226 G4MultiNavigator::PrepareNavigators(); 227 228 LocateGlobalPointAndSetup( position, &direction, false, false ); 229 // 230 // The first location for each Navigator must be non-relative 231 // or else call ResetStackAndState() for each Navigator 232 // Use direction to get correct side of boundary (ignore dir= false) 233 } 234 235 // ---------------------------------------------------------------------- 236 237 void G4MultiNavigator::PrepareNavigators() 238 { 239 // Key purposes: 240 // - Check and cache set of active navigators 241 // - Reset state for new track 242 243 #ifdef G4DEBUG_NAVIGATION 244 if( fVerbose > 1 ) 245 { 246 G4cout << " Entered G4MultiNavigator::PrepareNavigators() " << G4endl; 247 } 248 #endif 249 250 // Message the transportation-manager to find active navigators 251 252 std::vector<G4Navigator*>::const_iterator pNavigatorIter; 253 fNoActiveNavigators = (G4int)pTransportManager-> GetNoActiveNavigators(); 254 255 if( fNoActiveNavigators > fMaxNav ) 256 { 257 std::ostringstream message; 258 message << "Too many active Navigators / worlds !" << G4endl 259 << " Active Navigators (worlds): " 260 << fNoActiveNavigators << G4endl 261 << " which is more than the number allowed: " 262 << fMaxNav << " !"; 263 G4Exception("G4MultiNavigator::PrepareNavigators()", "GeomNav0002", 264 FatalException, message); 265 } 266 267 pNavigatorIter= pTransportManager->GetActiveNavigatorsIterator(); 268 for( auto num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num ) 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 navigator 279 // in case a call to SetWorldVolume() changed it 280 281 G4VPhysicalVolume* massWorld = GetWorldVolume(); 282 283 if( (massWorld != fLastMassWorld) && (massWorld!=nullptr) ) 284 { 285 // Pass along change to Mass Navigator 286 fpNavigator[0] -> SetWorldVolume( massWorld ); 287 288 #ifdef G4DEBUG_NAVIGATION 289 if( fVerbose > 0 ) 290 { 291 G4cout << " G4MultiNavigator::PrepareNavigators() changed world volume " 292 << " for mass geometry to " << massWorld->GetName() << G4endl; 293 } 294 #endif 295 296 fLastMassWorld = massWorld; 297 } 298 } 299 300 // ---------------------------------------------------------------------- 301 302 G4VPhysicalVolume* 303 G4MultiNavigator::LocateGlobalPointAndSetup(const G4ThreeVector& position, 304 const G4ThreeVector* pDirection, 305 const G4bool pRelativeSearch, 306 const G4bool ignoreDirection ) 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->GetActiveNavigatorsIterator(); 313 314 if( pDirection != nullptr ) { direction = *pDirection; } 315 316 #ifdef G4DEBUG_NAVIGATION 317 if( fVerbose > 2 ) 318 { 319 G4cout << " Entered G4MultiNavigator::LocateGlobalPointAndSetup() " 320 << G4endl; 321 G4cout << " Locating at position: " << position 322 << ", with direction: " << direction << G4endl 323 << " Relative: " << relative 324 << ", ignore direction: " << ignoreDirection << G4endl; 325 G4cout << " Number of active navigators: " << fNoActiveNavigators 326 << G4endl; 327 } 328 #endif 329 330 for ( auto num=0; num< fNoActiveNavigators ; ++pNavIter,++num ) 331 { 332 if( fWasLimitedByGeometry && fLimitTruth[num] ) 333 { 334 (*pNavIter)->SetGeometricallyLimitedStep(); 335 } 336 337 G4VPhysicalVolume *pLocated 338 = (*pNavIter)->LocateGlobalPointAndSetup( position, &direction, 339 relative, ignoreDirection ); 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 clear on locating (see Navigator) 349 350 #ifdef G4DEBUG_NAVIGATION 351 if( fVerbose > 2 ) 352 { 353 G4cout << " Located in world: " << num << ", at: " << position << G4endl 354 << " Used geomLimStp: " << fLimitTruth[num] 355 << ", found in volume: " << pLocated << G4endl; 356 G4cout << " Name = '" ; 357 if( pLocated ) 358 { 359 G4cout << pLocated->GetName() << "'"; 360 G4cout << " - CopyNo= " << pLocated->GetCopyNo(); 361 } 362 else 363 { 364 G4cout << "Null' Id: Not-Set "; 365 } 366 G4cout << G4endl; 367 } 368 #endif 369 } 370 371 fWasLimitedByGeometry = false; // Clear on locating 372 G4VPhysicalVolume* volMassLocated = fLocatedVolume[0]; 373 374 return volMassLocated; 375 } 376 377 // ---------------------------------------------------------------------- 378 379 void 380 G4MultiNavigator::LocateGlobalPointWithinVolume(const G4ThreeVector& position) 381 { 382 // Relocate the point in each geometry 383 384 auto pNavIter = pTransportManager->GetActiveNavigatorsIterator(); 385 386 #ifdef G4DEBUG_NAVIGATION 387 if( fVerbose > 2 ) 388 { 389 G4cout << " Entered G4MultiNavigator::ReLocate() " << G4endl 390 << " Re-locating at position: " << position << G4endl; 391 } 392 #endif 393 394 for ( auto num=0; num< fNoActiveNavigators ; ++pNavIter,++num ) 395 { 396 // ... none limited the step 397 398 (*pNavIter)->LocateGlobalPointWithinVolume( position ); 399 400 // Clear state related to the step 401 // 402 fLimitedStep[num] = kDoNot; 403 fCurrentStepSize[num] = 0.0; 404 405 fLimitTruth[ num ] = false; // Always clear on locating (see Navigator) 406 } 407 fWasLimitedByGeometry = false; // Clear on locating 408 fLastLocatedPosition = position; 409 } 410 411 // ---------------------------------------------------------------------- 412 413 G4double G4MultiNavigator::ComputeSafety( const G4ThreeVector& position, 414 const G4double maxDistance, 415 const G4bool state) 416 { 417 // Recompute safety for the relevant point 418 419 G4double minSafety = kInfinity, safety = kInfinity; 420 421 std::vector<G4Navigator*>::iterator pNavigatorIter; 422 pNavigatorIter= pTransportManager-> GetActiveNavigatorsIterator(); 423 424 for( auto num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num ) 425 { 426 safety = (*pNavigatorIter)->ComputeSafety( position, maxDistance, state); 427 if( safety < minSafety ) { minSafety = safety; } 428 } 429 430 fSafetyLocation = position; 431 fMinSafety_atSafLocation = minSafety; 432 433 #ifdef G4DEBUG_NAVIGATION 434 if( fVerbose > 1 ) 435 { 436 G4cout << " G4MultiNavigator::ComputeSafety - returns: " 437 << minSafety << ", at location: " << position << G4endl; 438 } 439 #endif 440 return minSafety; 441 } 442 443 // ----------------------------------------------------------------------- 444 445 G4TouchableHandle G4MultiNavigator::CreateTouchableHistoryHandle() const 446 { 447 G4Exception( "G4MultiNavigator::CreateTouchableHistoryHandle()", 448 "GeomNav0001", FatalException, 449 "Getting a touchable from G4MultiNavigator is not defined."); 450 451 G4TouchableHistory* touchHist; 452 touchHist= fpNavigator[0] -> CreateTouchableHistory(); 453 454 G4VPhysicalVolume* locatedVolume= fLocatedVolume[0]; 455 if( locatedVolume == nullptr ) 456 { 457 // Workaround to ensure that the touchable is fixed !! // TODO: fix 458 // 459 touchHist->UpdateYourself( locatedVolume, touchHist->GetHistory() ); 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 Navigator !! 473 G4int noLimited=0; 474 ELimited shared= kSharedOther; 475 476 #ifdef G4DEBUG_NAVIGATION 477 if( fVerbose > 2 ) 478 { 479 G4cout << " Entered G4MultiNavigator::WhichLimited() " << G4endl; 480 } 481 #endif 482 483 // Assume that [IdTransport] is Mass / Transport 484 // 485 G4bool transportLimited = (fCurrentStepSize[IdTransport] == fMinStep) 486 && ( fMinStep!= kInfinity); 487 if( transportLimited ) 488 { 489 shared = kSharedTransport; 490 } 491 492 for ( auto num = 0; num < fNoActiveNavigators; ++num ) 493 { 494 G4bool limitedStep; 495 496 G4double step = fCurrentStepSize[num]; 497 498 limitedStep = ( step == fMinStep ) && ( step != kInfinity); 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"), StrUnique("Unique"), 531 StrUndefined("Undefined"), 532 StrSharedTransport("SharedTransport"), 533 StrSharedOther("SharedOther"); 534 G4cout << "### G4MultiNavigator::PrintLimited() reports: " << G4endl; 535 G4cout << " Minimum step (true): " << fTrueMinStep 536 << ", reported min: " << fMinStep << G4endl; 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 / flag" << " " 546 << std::setw(15) << " World " << " " 547 << G4endl; 548 } 549 #endif 550 551 for ( auto num = 0; num < fNoActiveNavigators; ++num ) 552 { 553 G4double rawStep = fCurrentStepSize[num]; 554 G4double stepLen = fCurrentStepSize[num]; 555 if( stepLen > fTrueMinStep ) 556 { 557 stepLen = fTrueMinStep; // did not limit (went as far as asked) 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] ? "YES" : " NO") << " "; 566 G4String limitedStr; 567 switch ( fLimitedStep[num] ) 568 { 569 case kDoNot : limitedStr = StrDoNot; break; 570 case kUnique : limitedStr = StrUnique; break; 571 case kSharedTransport: limitedStr = StrSharedTransport; break; 572 case kSharedOther : limitedStr = StrSharedOther; break; 573 default : limitedStr = StrUndefined; break; 574 } 575 G4cout << " " << std::setw(15) << limitedStr << " "; 576 G4cout.precision(oldPrec); 577 578 G4Navigator *pNav = fpNavigator[ num ]; 579 G4String WorldName( "Not-Set" ); 580 if (pNav != nullptr) 581 { 582 G4VPhysicalVolume *pWorld = pNav->GetWorldVolume(); 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()", "GeomNav0001", 601 FatalException, 602 "Cannot reset state for navigators of G4MultiNavigator."); 603 /* 604 std::vector<G4Navigator*>::iterator pNavigatorIter; 605 pNavigatorIter = pTransportManager->GetActiveNavigatorsIterator(); 606 for( auto num = 0; num< fNoActiveNavigators; ++pNavigatorIter,++num ) 607 { 608 // (*pNavigatorIter)->ResetState(); // KEEP THIS comment !!! 609 } 610 */ 611 } 612 613 // ----------------------------------------------------------------------- 614 615 void G4MultiNavigator::SetupHierarchy() 616 { 617 G4Exception( "G4MultiNavigator::SetupHierarchy()", 618 "GeomNav0001", FatalException, 619 "Cannot setup hierarchy for navigators of G4MultiNavigator."); 620 } 621 622 // ----------------------------------------------------------------------- 623 624 void G4MultiNavigator::CheckMassWorld() 625 { 626 G4VPhysicalVolume* navTrackWorld = 627 pTransportManager->GetNavigatorForTracking()->GetWorldVolume(); 628 629 if( navTrackWorld != fLastMassWorld ) 630 { 631 G4Exception( "G4MultiNavigator::CheckMassWorld()", 632 "GeomNav0003", FatalException, 633 "Mass world pointer has been changed." ); 634 } 635 } 636 637 // ----------------------------------------------------------------------- 638 639 G4VPhysicalVolume* 640 G4MultiNavigator::ResetHierarchyAndLocate(const G4ThreeVector& point, 641 const G4ThreeVector& direction, 642 const G4TouchableHistory& MassHistory) 643 { 644 // Reset geometry for all -- and use the touchable for the mass history 645 646 G4VPhysicalVolume* massVolume = nullptr; 647 G4Navigator* pMassNavigator = fpNavigator[0]; 648 649 if( pMassNavigator != nullptr ) 650 { 651 massVolume= pMassNavigator->ResetHierarchyAndLocate( point, direction, 652 MassHistory); 653 } 654 else 655 { 656 G4Exception("G4MultiNavigator::ResetHierarchyAndLocate()", 657 "GeomNav0002", FatalException, 658 "Cannot reset hierarchy before navigators are initialised."); 659 } 660 661 auto pNavIter= pTransportManager->GetActiveNavigatorsIterator(); 662 663 for ( auto num = 0; num < fNoActiveNavigators ; ++pNavIter,++num ) 664 { 665 G4bool relativeSearch, ignoreDirection; 666 667 (*pNavIter)-> LocateGlobalPointAndSetup( point, 668 &direction, 669 relativeSearch=false, 670 ignoreDirection=false); 671 } 672 return massVolume; 673 } 674 675 // ----------------------------------------------------------------------- 676 677 G4ThreeVector 678 G4MultiNavigator::GetGlobalExitNormal(const G4ThreeVector& argPoint, 679 G4bool* argpObtained) // obtained valid 680 { 681 G4ThreeVector normalGlobalCrd(0.0, 0.0, 0.0); 682 G4bool isObtained = false; 683 // These default values will be used if fNoLimitingStep== 0 684 G4int firstNavigatorId = -1; 685 G4bool oneObtained = false; 686 687 if( fNoLimitingStep == 1 ) 688 { 689 // Only message the Navigator which limited the step! 690 normalGlobalCrd = fpNavigator[ fIdNavLimiting ] 691 ->GetGlobalExitNormal( argPoint, &isObtained ); 692 *argpObtained = isObtained; 693 } 694 else 695 { 696 if( fNoLimitingStep > 1 ) 697 { 698 auto pNavIter= pTransportManager->GetActiveNavigatorsIterator(); 699 700 for ( auto num = 0; num < fNoActiveNavigators ; ++pNavIter, ++num ) 701 { 702 G4ThreeVector oneNormal; 703 if( fLimitTruth[ num ] ) // Did this geometry limit the step ? 704 { 705 G4ThreeVector newNormal = 706 (*pNavIter)->GetGlobalExitNormal( argPoint, &oneObtained ); 707 if( oneObtained ) 708 { 709 // Keep first one - only if it is valid (ie not null) 710 if( !isObtained && (newNormal.mag2() != 0.0) ) 711 { 712 normalGlobalCrd = newNormal; 713 isObtained = oneObtained; 714 firstNavigatorId = num; 715 } 716 else 717 { 718 // Check for clash 719 G4double dotNewPrevious = newNormal.dot( normalGlobalCrd ); 720 G4double productMagSq = normalGlobalCrd.mag2()*newNormal.mag2(); 721 if( productMagSq > 0.0 ) 722 { 723 G4double productMag = std::sqrt( productMagSq ); 724 dotNewPrevious /= productMag; // Normalise 725 if( dotNewPrevious < (1 - perThousand) ) 726 { 727 *argpObtained = false; 728 729 if( fVerbose > 2 ) // dotNewPrevious <= 0.0 ) 730 { 731 std::ostringstream message; 732 message << "Clash of Normal from different Navigators!" 733 << G4endl 734 << " Previous Navigator Id = " 735 << firstNavigatorId << G4endl 736 << " Current Navigator Id = " 737 << num << G4endl; 738 message << " Dot product of 2 normals = " 739 << dotNewPrevious << G4endl; 740 message << " Normal (previous) = " 741 << normalGlobalCrd << G4endl; 742 message << " Normal (current) = " << newNormal << G4endl; 743 G4Exception("G4MultiNavigator::GetGlobalExitNormal()", 744 "GeomNav0002", JustWarning, message); 745 } 746 } 747 else 748 { 749 // Close agreement - Do not change 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 having " << fNoLimitingStep 762 << " candidate Navigators limiting the step!" << G4endl; 763 G4Exception("G4MultiNavigator::GetGlobalExitNormal()", "GeomNav0002", 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* argpObtained) 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 fNoLimitingStep== 0 783 784 if( fNoLimitingStep == 1 ) 785 { 786 // Only message the Navigator which limited the step! 787 normalGlobalCrd = fpNavigator[ fIdNavLimiting ] 788 ->GetLocalExitNormal( &isObtained ); 789 *argpObtained = isObtained; 790 791 static G4ThreadLocal G4int numberWarnings = 0; 792 G4int noWarningsStart = 10, noModuloWarnings = 100; 793 ++numberWarnings; 794 if( (numberWarnings < noWarningsStart ) 795 || (numberWarnings%noModuloWarnings == 0) ) 796 { 797 std::ostringstream message; 798 message << "Cannot obtain normal in local coordinates of two or more " 799 << "coordinate systems." << G4endl; 800 G4Exception("G4MultiNavigator::GetGlobalExitNormal()", "GeomNav0002", 801 JustWarning, message); 802 } 803 } 804 else 805 { 806 if( fNoLimitingStep > 1 ) 807 { 808 std::ostringstream message; 809 message << "Cannot obtain normal in local coordinates of two or more " 810 << "coordinate systems." << G4endl; 811 G4Exception("G4MultiNavigator::GetGlobalExitNormal()", "GeomNav0002", 812 FatalException, message); 813 } 814 } 815 816 *argpObtained = isObtained; 817 return normalGlobalCrd; 818 } 819 820 // ----------------------------------------------------------------------- 821 822 G4ThreeVector 823 G4MultiNavigator::GetLocalExitNormalAndCheck(const G4ThreeVector&, // point, 824 G4bool* obtained) 825 { 826 return G4MultiNavigator::GetLocalExitNormal( obtained ); 827 } 828