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