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 G4SafetyHelper 27 // 28 // Class description: 29 // 30 // This class is a helper for physics processes which require 31 // knowledge of the safety, and the step size for the 'mass' geometry 32 33 // First version: J.Apostolakis, July 5th, 2006 34 // -------------------------------------------------------------------- 35 #ifndef G4SAFETYHELPER_HH 36 #define G4SAFETYHELPER_HH 1 37 38 #include <vector> 39 40 #include "G4Types.hh" 41 #include "G4ThreeVector.hh" 42 #include "G4Navigator.hh" 43 44 class G4PathFinder; 45 46 class G4SafetyHelper 47 { 48 public: // with description 49 50 G4SafetyHelper(); 51 ~G4SafetyHelper(); 52 // Constructor and destructor 53 54 G4double CheckNextStep( const G4ThreeVector& position, 55 const G4ThreeVector& direction, 56 const G4double currentMaxStep, 57 G4double& newSafety ); 58 // Return linear step for mass geometry 59 60 G4double ComputeSafety( const G4ThreeVector& pGlobalPoint, 61 G4double maxRadius = DBL_MAX ); 62 // Return safety for all geometries. 63 // 64 // The 2nd argument is the radius of your interest (e.g. maximum 65 // displacement). Giving this you can reduce the average computational 66 // cost. If the second argument is not given, this is the real 67 // isotropic safety 68 69 void Locate(const G4ThreeVector& pGlobalPoint, 70 const G4ThreeVector& direction); 71 // Locate the point for all geometries 72 73 void ReLocateWithinVolume(const G4ThreeVector& pGlobalPoint ); 74 // Relocate the point in the volume of interest 75 76 inline void EnableParallelNavigation(G4bool parallel); 77 // To have parallel worlds considered, must be true. 78 // Alternative is to use single (mass) Navigator directly 79 80 void InitialiseNavigator(); 81 // Check for new navigator for tracking, and reinitialise pointer 82 83 inline G4int SetVerboseLevel( G4int lev ); 84 inline G4VPhysicalVolume* GetWorldVolume(); 85 inline void SetCurrentSafety(G4double val, const G4ThreeVector& pos); 86 87 public: // without description 88 89 void InitialiseHelper(); 90 91 private: 92 93 G4PathFinder* fpPathFinder = nullptr; 94 G4Navigator* fpMassNavigator = nullptr; 95 96 G4bool fUseParallelGeometries = false; 97 // Flag whether to use PathFinder or single (mass) Navigator directly 98 // By default, one geometry only 99 G4bool fFirstCall = true; 100 // Flag of first call 101 G4int fVerbose = 0; 102 // Whether to print warning in case of move outside safety 103 104 // State used during tracking -- for optimisation 105 106 G4ThreeVector fLastSafetyPosition; 107 G4double fLastSafety = 0.0; 108 109 // const G4double fRecomputeFactor = 0.0; 110 // parameter for further optimisation: 111 // if ( move < fact*safety ) do fast recomputation of safety 112 113 // End State (tracking) 114 }; 115 116 // -------------------------------------------------------------------- 117 // Inline definitions 118 // -------------------------------------------------------------------- 119 120 inline G4int G4SafetyHelper::SetVerboseLevel( G4int lev ) 121 { 122 G4int oldlv = fVerbose; 123 fVerbose = lev; 124 return oldlv; 125 } 126 127 inline 128 void G4SafetyHelper::EnableParallelNavigation(G4bool parallel) 129 { 130 fUseParallelGeometries = parallel; 131 } 132 133 inline 134 G4VPhysicalVolume* G4SafetyHelper::GetWorldVolume() 135 { 136 return fpMassNavigator->GetWorldVolume(); 137 } 138 139 inline 140 void G4SafetyHelper::SetCurrentSafety(G4double val, const G4ThreeVector& pos) 141 { 142 fLastSafety = val; 143 fLastSafetyPosition = pos; 144 } 145 146 #endif 147