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