Geant4  10.02
G4ITPathFinder.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 // $Id: G4ITPathFinder.hh 69058 2013-04-17 09:08:25Z gcosmo $
27 //
28 // class G4ITPathFinder
29 //
30 // Class description:
31 //
32 // G4ITPathFinder is a duplicated version of G4ITPathFinder
33 //
34 // This class directs the lock-stepped propagation of a track in the
35 // 'mass' and other parallel geometries. It ensures that tracking
36 // in a magnetic field sees these parallel geometries at each trial step,
37 // and that the earliest boundary limits the step.
38 //
39 // For the movement in field, it relies on the class G4PropagatorInField
40 //
41 // History:
42 // -------
43 // 7.10.05 John Apostolakis, Draft design
44 // 26.04.06 John Apostolakis, Revised design and first implementation
45 // ---------------------------------------------------------------------------
46 
47 #ifndef G4ITPATHFINDER_HH
48 #define G4ITPATHFINDER_HH 1
49 
50 #include <vector>
51 #include "G4Types.hh"
52 
53 #include "G4FieldTrack.hh"
54 
56 class G4ITNavigator;
57 
58 #include "G4TouchableHandle.hh"
59 #include "G4FieldTrack.hh"
60 #include "G4ITMultiNavigator.hh"
61 #include "G4TrackState.hh"
62 
64 class G4ITPathFinder;
65 
66 // Global state (retained during stepping for one track)
67 // State changed in a step computation
68 template<>
69 class G4TrackState<G4ITPathFinder> : public G4TrackStateBase<G4ITPathFinder>
70 {
71  friend class G4ITPathFinder;
72 
73 protected:
74  G4bool fNewTrack; // Flag a new track (ensure first step)
75 
76  ELimited fLimitedStep[G4ITNavigator::fMaxNav];
77  G4bool fLimitTruth[G4ITNavigator::fMaxNav];
78  G4double fCurrentStepSize[G4ITNavigator::fMaxNav];
79  G4int fNoGeometriesLimiting; // How many processes contribute to limit
80 
81  G4ThreeVector fPreSafetyLocation; // last initial position for which safety evaluated
82  G4double fPreSafetyMinValue; // /\ corresponding value of full safety
83  G4double fPreSafetyValues[ G4ITNavigator::fMaxNav ]; // Safeties for the above point
84  // This part of the state can be retained for severall calls --> CARE
85 
86  G4ThreeVector fPreStepLocation; // point where last ComputeStep called
87  G4double fMinSafety_PreStepPt; // /\ corresponding value of full safety
88  G4double fCurrentPreStepSafety[ G4ITNavigator::fMaxNav ]; // Safeties for the above point
89  // This changes at each step,
90  // so it can differ when steps inside min-safety are made
91 
92  G4bool fPreStepCenterRenewed; // Whether PreSafety coincides with PreStep point
93 
94  G4double fMinStep; // As reported by Navigators -- can be kInfinity
95  G4double fTrueMinStep; // Corrected in case >= proposed
96 
97  // State after calling 'locate'
98 
99  G4VPhysicalVolume* fLocatedVolume[G4ITNavigator::fMaxNav];
101 
102  // State after calling 'ComputeStep' (others member variables will be affected)
103  G4FieldTrack fEndState; // Point, velocity, ... at proposed step end
104  G4bool fFieldExertedForce; // In current proposed step
105 
106  G4bool fRelocatedPoint; // Signals that point was or is being moved
107  // from the position of the last location
108  // or the endpoint resulting from ComputeStep
109  // -- invalidates fEndState
110 
111  // State for 'ComputeSafety' and related methods
112  G4ThreeVector fSafetyLocation; // point where ComputeSafety is called
113  G4double fMinSafety_atSafLocation; // /\ corresponding value of safety
114  G4double fNewSafetyComputed[ G4ITNavigator::fMaxNav ]; // Safeties for last ComputeSafety
115 
116  // State for Step numbers
118 
119 public:
121 
124  fEndState( G4ThreeVector(), G4ThreeVector(), 0., 0., 0., 0., 0.),
125  fFieldExertedForce(false),
126  fRelocatedPoint(true),
127  fLastStepNo(-1), fCurrentStepNo(-1) {
128 
129  G4ThreeVector Big3Vector( kInfinity, kInfinity, kInfinity );
130  fLastLocatedPosition= Big3Vector;
131  fSafetyLocation= Big3Vector;
132  fPreSafetyLocation= Big3Vector;
133  fPreStepLocation= Big3Vector;
134 
135  fPreSafetyMinValue= -1.0;
136  fMinSafety_PreStepPt= -1.0;
137  fMinSafety_atSafLocation= -1.0;
138  fMinStep= -1.0;
139  fTrueMinStep= -1.0;
140  fPreStepCenterRenewed= false;
141  fNewTrack= false;
142  fNoGeometriesLimiting= 0;
143 
144  for( G4int num=0; num< G4ITNavigator::fMaxNav; ++num )
145  {
146  fLimitTruth[num] = false;
147  fLimitedStep[num] = kUndefLimited;
148  fCurrentStepSize[num] = -1.0;
149  fLocatedVolume[num] = 0;
150  fPreSafetyValues[num]= -1.0;
151  fCurrentPreStepSafety[num] = -1.0;
152  fNewSafetyComputed[num]= -1.0;
153  }
154  }
155 };
156 
157 class G4ITPathFinder : public G4TrackStateDependent<G4ITPathFinder>
158 {
159 
160 public: // with description
161 
162  static G4ITPathFinder* GetInstance();
163  //
164  // Retrieve singleton instance
165 
166  G4double ComputeStep( const G4FieldTrack &pFieldTrack,
167  G4double pCurrentProposedStepLength,
168  G4int navigatorId, // Identifies the geometry
169  G4int stepNo, // See next step/check
170  G4double &pNewSafety, // Only for this geometry
171  ELimited &limitedStep,
172  G4FieldTrack &EndState,
173  G4VPhysicalVolume* currentVolume );
174  //
175  // Compute the next geometric Step -- Curved or linear
176  // If it is called with a larger 'stepNo' it will execute a new step;
177  // if 'stepNo' is same as last call, then the results for
178  // the geometry with Id. number 'navigatorId' will be returned.
179 
180  void Locate( const G4ThreeVector& position,
181  const G4ThreeVector& direction,
182  G4bool relativeSearch=true);
183  //
184  // Make primary relocation of global point in all navigators,
185  // and update them.
186 
187  void ReLocate( const G4ThreeVector& position );
188  //
189  // Make secondary relocation of global point (within safety only)
190  // in all navigators, and update them.
191 
192  void PrepareNewTrack( const G4ThreeVector& position,
193  const G4ThreeVector& direction,
194  G4VPhysicalVolume* massStartVol=0);
195  //
196  // Check and cache set of active navigators.
197 
199  inline G4VPhysicalVolume* GetLocatedVolume( G4int navId ) const;
200 
201  // -----------------------------------------------------------------
202 
203  inline G4bool IsParticleLooping() const;
204 
205  inline G4double GetCurrentSafety() const;
206  // Minimum value of safety after last ComputeStep
207  inline G4double GetMinimumStep() const;
208  // Get the minimum step size from the last ComputeStep call
209  // - in case full step is taken, this is kInfinity
210  inline unsigned int GetNumberGeometriesLimitingStep() const;
211 
212  G4double ComputeSafety( const G4ThreeVector& globalPoint);
213  // Recompute safety for the relevant point the endpoint of the last step!!
214  // Maintain vector of individual safety values (for next method)
215 
216  G4double ObtainSafety( G4int navId, G4ThreeVector& globalCenterPoint );
217  // Obtain safety for navigator/geometry navId for last point 'computed'
218  // --> last point for which ComputeSafety was called
219  // Returns the point (center) for which this safety is valid
220 
221  void EnableParallelNavigation( G4bool enableChoice=true );
222  //
223  // Must call it to ensure that G4ITNavigator is prepared,
224  // especially for curved tracks. If true it switches PropagatorInField
225  // to use MultiNavigator. Must call it with false to undo (=PiF use
226  // Navigator for tracking!)
227 
228  inline G4int SetVerboseLevel(G4int lev=-1);
229 
230 public: // with description
231 
232  inline G4int GetMaxLoopCount() const;
233  inline void SetMaxLoopCount( G4int new_max );
234  //
235  // A maximum for the number of steps that a (looping) particle can take.
236 
237 public: // without description
238 
239  inline void MovePoint();
240  //
241  // Signal that location will be moved -- internal use primarily
242 
243  // To provide best compatibility between Coupled and Old Transportation
244  // the next two methods are provided:
245  G4double LastPreSafety( G4int navId, G4ThreeVector& globalCenterPoint, G4double& minSafety );
246  // Obtain last safety needed in ComputeStep (for geometry navId)
247  // --> last point at which ComputeStep recalculated safety
248  // Returns the point (center) for which this safety is valid
249  // and also the minimum safety over all navigators (ie full)
250 
252  // Tell G4ITNavigator to copy PostStep Safety to PreSafety (for use at next step)
253 
255  // Convert ELimited to string
256 
257 protected: // without description
258 
259  G4double DoNextLinearStep( const G4FieldTrack &FieldTrack,
260  G4double proposedStepLength);
261 
262  G4double DoNextCurvedStep( const G4FieldTrack &FieldTrack,
263  G4double proposedStepLength,
264  G4VPhysicalVolume* pCurrentPhysVolume);
265 
266  void WhichLimited();
267  void PrintLimited();
268  //
269  // Print key details out - for debugging
270 
271  // void ClearState();
272  //
273  // Clear all the State of this class and its current associates
274 
276  //
277  // Whether use safety to discard unneccesary calls to navigator
278 
279  void ReportMove( const G4ThreeVector& OldV, const G4ThreeVector& NewV, const G4String& Quantity ) const;
280  // Helper method to report movement (likely of initial point)
281 
282 protected:
283 
284  G4ITPathFinder(); // Singleton
285  ~G4ITPathFinder();
286 
287  inline G4ITNavigator* GetNavigator(G4int n) const;
288 
289 private:
290 
291  // ----------------------------------------------------------------------
292  // DATA Members
293  // ----------------------------------------------------------------------
294 
296  //
297  // Object that enables G4PropagatorInField to see many geometries
298 
300 
301  G4ITNavigator* fpNavigator[G4ITNavigator::fMaxNav];
302 
303  G4int fVerboseLevel; // For debuging purposes
304 
305  G4ITTransportationManager* fpTransportManager; // Cache for frequent use
306  // G4PropagatorInField* fpFieldPropagator;
307 
309 
311 
312 };
313 
314 // ********************************************************************
315 // Inline methods.
316 // ********************************************************************
317 
319 {
320  G4VPhysicalVolume* vol=0;
321  if( (navId < G4ITNavigator::fMaxNav) && (navId >=0) ) { vol= fpTrackState->fLocatedVolume[navId]; }
322  return vol;
323 }
324 
326 {
327  G4int old= fVerboseLevel; fVerboseLevel= newLevel; return old;
328 }
329 
331 {
332  return fpTrackState->fMinStep;
333 }
334 
336 {
337  unsigned int noGeometries=fpTrackState->fNoGeometriesLimiting;
338  return noGeometries;
339 }
340 
342 {
343  return fpTrackState->fMinSafety_PreStepPt;
344 }
345 
347 {
348  fpTrackState->fRelocatedPoint= true;
349 }
350 
351 inline G4ITNavigator* G4ITPathFinder::GetNavigator(G4int n) const
352 {
353  if( (n>fNoActiveNavigators)||(n<0)) { n=0; }
354  return fpNavigator[n];
355 }
356 
357 inline G4double G4ITPathFinder::ObtainSafety( G4int navId, G4ThreeVector& globalCenterPoint )
358 {
359  globalCenterPoint= fpTrackState->fSafetyLocation;
360  // navId = std::min( navId, fMaxNav-1 );
361  return fpTrackState->fNewSafetyComputed[ navId ];
362 }
363 
365  G4ThreeVector& globalCenterPoint,
367 {
368  globalCenterPoint= fpTrackState->fPreSafetyLocation;
369  minSafety= fpTrackState->fPreSafetyMinValue;
370  // navId = std::min( navId, fMaxNav-1 );
371  return fpTrackState->fPreSafetyValues[ navId ];
372 }
373 #endif
G4ITNavigator * GetNavigator(G4int n) const
void ReLocate(const G4ThreeVector &position)
G4ITMultiNavigator * fpMultiNavigator
static G4ITPathFinder * GetInstance()
#define fCurrentPreStepSafety
G4double ObtainSafety(G4int navId, G4ThreeVector &globalCenterPoint)
G4double GetMinimumStep() const
G4double LastPreSafety(G4int navId, G4ThreeVector &globalCenterPoint, G4double &minSafety)
static const G4double kInfinity
Definition: geomdefs.hh:42
CLHEP::Hep3Vector G4ThreeVector
G4bool IsParticleLooping() const
#define fPreSafetyValues
void EnableParallelNavigation(G4bool enableChoice=true)
void Locate(const G4ThreeVector &position, const G4ThreeVector &direction, G4bool relativeSearch=true)
G4TouchableHandle CreateTouchableHandle(G4int navId) const
ELimited
void PrepareNewTrack(const G4ThreeVector &position, const G4ThreeVector &direction, G4VPhysicalVolume *massStartVol=0)
#define fLocatedVolume
#define G4ThreadLocal
Definition: tls.hh:89
G4double kCarTolerance
G4double GetCurrentSafety() const
int G4int
Definition: G4Types.hh:78
#define fFieldExertedForce
G4VPhysicalVolume * GetLocatedVolume(G4int navId) const
G4double DoNextLinearStep(const G4FieldTrack &FieldTrack, G4double proposedStepLength)
#define position
Definition: xmlparse.cc:622
#define fLimitTruth
static G4ThreadLocal G4ITPathFinder * fpPathFinder
G4bool UseSafetyForOptimization(G4bool)
bool G4bool
Definition: G4Types.hh:79
void ReportMove(const G4ThreeVector &OldV, const G4ThreeVector &NewV, const G4String &Quantity) const
G4int GetMaxLoopCount() const
G4double ComputeSafety(const G4ThreeVector &globalPoint)
const G4int n
static const G4double minSafety
G4double DoNextCurvedStep(const G4FieldTrack &FieldTrack, G4double proposedStepLength, G4VPhysicalVolume *pCurrentPhysVolume)
#define fCurrentStepNo
#define fLimitedStep
G4String & LimitedString(ELimited lim)
#define fCurrentStepSize
void SetMaxLoopCount(G4int new_max)
unsigned int GetNumberGeometriesLimitingStep() const
#define fLastStepNo
#define fEndState
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
double G4double
Definition: G4Types.hh:76
#define fRelocatedPoint
G4int SetVerboseLevel(G4int lev=-1)
void PushPostSafetyToPreSafety()
G4ITNavigator * fpNavigator[G4ITNavigator::fMaxNav]
G4ITTransportationManager * fpTransportManager
#define fNewSafetyComputed