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