Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // Class G4PropagatorInField 26 // Class G4PropagatorInField 27 // 27 // 28 // class description: 28 // class description: 29 // 29 // 30 // This class performs the navigation/propagat 30 // This class performs the navigation/propagation of a particle/track 31 // in a magnetic field. The field is in genera 31 // in a magnetic field. The field is in general non-uniform. 32 // For the calculation of the path, it relies 32 // For the calculation of the path, it relies on the class G4ChordFinder. 33 33 34 // History: 34 // History: 35 // ------- 35 // ------- 36 // 25.10.96 John Apostolakis, design and impl 36 // 25.10.96 John Apostolakis, design and implementation 37 // 25.03.97 John Apostolakis, adaptation for 37 // 25.03.97 John Apostolakis, adaptation for G4Transportation and cleanup 38 // 8.11.02 John Apostolakis, changes to enab 38 // 8.11.02 John Apostolakis, changes to enable use of safety in intersecting 39 // ------------------------------------------- 39 // --------------------------------------------------------------------------- 40 #ifndef G4PropagatorInField_hh 40 #ifndef G4PropagatorInField_hh 41 #define G4PropagatorInField_hh 1 41 #define G4PropagatorInField_hh 1 42 42 43 #include "G4Types.hh" 43 #include "G4Types.hh" 44 44 45 #include <vector> 45 #include <vector> 46 46 47 #include "G4FieldTrack.hh" 47 #include "G4FieldTrack.hh" 48 #include "G4FieldManager.hh" 48 #include "G4FieldManager.hh" 49 #include "G4VIntersectionLocator.hh" 49 #include "G4VIntersectionLocator.hh" 50 50 51 class G4ChordFinder; 51 class G4ChordFinder; 52 52 53 class G4Navigator; 53 class G4Navigator; 54 class G4VPhysicalVolume; 54 class G4VPhysicalVolume; 55 class G4VCurvedTrajectoryFilter; 55 class G4VCurvedTrajectoryFilter; 56 56 57 class G4PropagatorInField 57 class G4PropagatorInField 58 { 58 { 59 59 60 public: // with description 60 public: // with description 61 61 62 G4PropagatorInField( G4Navigator* theNaviga 62 G4PropagatorInField( G4Navigator* theNavigator, 63 G4FieldManager* detect 63 G4FieldManager* detectorFieldMgr, 64 G4VIntersectionLocator 64 G4VIntersectionLocator* vLocator = nullptr ); 65 ~G4PropagatorInField(); 65 ~G4PropagatorInField(); 66 66 67 G4double ComputeStep( G4FieldTrack& pFieldT 67 G4double ComputeStep( G4FieldTrack& pFieldTrack, 68 G4double pCurrentProp 68 G4double pCurrentProposedStepLength, 69 G4double& pNewSafety, 69 G4double& pNewSafety, 70 G4VPhysicalVolume* pP 70 G4VPhysicalVolume* pPhysVol = nullptr, 71 G4bool canRelaxDeltaC 71 G4bool canRelaxDeltaChord = false); 72 // Compute the next geometric Step 72 // Compute the next geometric Step 73 73 74 inline G4ThreeVector EndPosition() const; 74 inline G4ThreeVector EndPosition() const; 75 inline G4ThreeVector EndMomentumDir() const 75 inline G4ThreeVector EndMomentumDir() const; 76 inline G4bool IsParticleLooping() co 76 inline G4bool IsParticleLooping() const; 77 // Return the state after the Step 77 // Return the state after the Step 78 78 79 inline G4double GetEpsilonStep() const; 79 inline G4double GetEpsilonStep() const; 80 // Relative accuracy for current Step (Ca 80 // Relative accuracy for current Step (Calc.) 81 inline void SetEpsilonStep(G4double new 81 inline void SetEpsilonStep(G4double newEps); 82 // The ratio DeltaOneStep()/h_current_ste 82 // The ratio DeltaOneStep()/h_current_step 83 83 84 G4FieldManager* FindAndSetFieldManager(G4VP 84 G4FieldManager* FindAndSetFieldManager(G4VPhysicalVolume* pCurrentPhysVol); 85 // Set (and return) the correct field man 85 // Set (and return) the correct field manager (global or local), 86 // if it exists. 86 // if it exists. 87 // Should be called before ComputeStep is 87 // Should be called before ComputeStep is called; 88 // Currently, ComputeStep will call it, i 88 // Currently, ComputeStep will call it, if it has not been called. 89 89 90 inline G4ChordFinder* GetChordFinder(); 90 inline G4ChordFinder* GetChordFinder(); 91 91 92 G4int SetVerboseLevel( G4int verbose 92 G4int SetVerboseLevel( G4int verbose ); 93 inline G4int GetVerboseLevel() const; 93 inline G4int GetVerboseLevel() const; 94 inline G4int Verbose() const; 94 inline G4int Verbose() const; 95 inline void CheckMode(G4bool mode); 95 inline void CheckMode(G4bool mode); 96 96 97 inline void SetVerboseTrace( G4bool enabl 97 inline void SetVerboseTrace( G4bool enable ); 98 inline G4bool GetVerboseTrace(); 98 inline G4bool GetVerboseTrace(); 99 // Tracing key parts of Compute Step 99 // Tracing key parts of Compute Step 100 100 101 inline G4int GetMaxLoopCount() const; 101 inline G4int GetMaxLoopCount() const; 102 inline void SetMaxLoopCount( G4int new_max 102 inline void SetMaxLoopCount( G4int new_max ); 103 // A maximum for the number of substeps t 103 // A maximum for the number of substeps that a particle can take. 104 // Above this number it is signaled as 104 // Above this number it is signaled as 'looping'. 105 105 106 void printStatus( const G4FieldTrack& 106 void printStatus( const G4FieldTrack& startFT, 107 const G4FieldTrack& 107 const G4FieldTrack& currentFT, 108 G4double 108 G4double requestStep, 109 G4double 109 G4double safety, 110 G4int 110 G4int step, 111 G4VPhysicalVolume* 111 G4VPhysicalVolume* startVolume); 112 // Print Method - useful mostly for debug 112 // Print Method - useful mostly for debugging. 113 113 114 inline G4FieldTrack GetEndState() const; 114 inline G4FieldTrack GetEndState() const; 115 115 116 inline G4double GetMinimumEpsilonStep() con 116 inline G4double GetMinimumEpsilonStep() const; // Min for relative accuracy 117 inline void SetMinimumEpsilonStep( G4do 117 inline void SetMinimumEpsilonStep( G4double newEpsMin ); // of any step 118 inline G4double GetMaximumEpsilonStep() con 118 inline G4double GetMaximumEpsilonStep() const; 119 inline void SetMaximumEpsilonStep( G4do 119 inline void SetMaximumEpsilonStep( G4double newEpsMax ); 120 // The 4 above methods are now obsolescen 120 // The 4 above methods are now obsolescent but *for now* will work 121 // They are being replaced by same-name m 121 // They are being replaced by same-name methods in G4FieldManager, 122 // allowing the specialisation in differe 122 // allowing the specialisation in different volumes. 123 // Their new behaviour is to change the v 123 // Their new behaviour is to change the values for the global field 124 // manager 124 // manager 125 125 126 void SetLargestAcceptableStep( G4double << 126 inline void SetLargestAcceptableStep( G4double newBigDist ); 127 G4double GetLargestAcceptableStep(); << 127 inline G4double GetLargestAcceptableStep(); 128 void ResetLargestAcceptableStep(); << 129 // Obtain / change the size of the larges << 130 // Reset method uses the world volume's << 131 << 132 G4double GetMaxStepSizeMultiplier(); << 133 void SetMaxStepSizeMultiplier(G4double << 134 // Control extra Multiplier parameter for << 135 G4double GetMinBigDistance(); << 136 void SetMinBigDistance(G4double val); << 137 // Control minimum 'directional' distance << 138 128 139 void SetTrajectoryFilter(G4VCurvedTrajector 129 void SetTrajectoryFilter(G4VCurvedTrajectoryFilter* filter); 140 // Set the filter that examines & stores 130 // Set the filter that examines & stores 'intermediate' 141 // curved trajectory points. Currently o 131 // curved trajectory points. Currently only position is stored. 142 132 143 std::vector<G4ThreeVector>* GimmeTrajectory 133 std::vector<G4ThreeVector>* GimmeTrajectoryVectorAndForgetIt() const; 144 // Access the points which have passed by 134 // Access the points which have passed by the filter. 145 // Responsibility for deleting the points 135 // Responsibility for deleting the points lies with the client. 146 // This method MUST BE called exactly ONC 136 // This method MUST BE called exactly ONCE per step. 147 137 148 void ClearPropagatorState(); 138 void ClearPropagatorState(); 149 // Clear all the State of this class and 139 // Clear all the State of this class and its current associates 150 // --> the current field manager & chord 140 // --> the current field manager & chord finder will also be called 151 141 152 inline void SetDetectorFieldManager( G4Fiel 142 inline void SetDetectorFieldManager( G4FieldManager* newGlobalFieldManager ); 153 // Update this (dangerous) state -- for t 143 // Update this (dangerous) state -- for the time being 154 144 155 inline void SetUseSafetyForOptimization( 145 inline void SetUseSafetyForOptimization( G4bool ); 156 inline G4bool GetUseSafetyForOptimization() 146 inline G4bool GetUseSafetyForOptimization(); 157 // Toggle & view parameter for using safe 147 // Toggle & view parameter for using safety to discard 158 // unneccesary calls to navigator (thus ' 148 // unneccesary calls to navigator (thus 'optimising' performance) 159 inline G4bool IntersectChord( const G4Three 149 inline G4bool IntersectChord( const G4ThreeVector& StartPointA, 160 const G4Three 150 const G4ThreeVector& EndPointB, 161 G4doubl 151 G4double& NewSafety, 162 G4doubl 152 G4double& LinearStepLength, 163 G4Three 153 G4ThreeVector& IntersectionPoint); 164 // Intersect the chord from StartPointA t 154 // Intersect the chord from StartPointA to EndPointB 165 // and return whether an intersection occ 155 // and return whether an intersection occurred 166 // NOTE: Safety is changed! 156 // NOTE: Safety is changed! 167 157 168 inline G4bool IsFirstStepInVolume(); 158 inline G4bool IsFirstStepInVolume(); 169 inline G4bool IsLastStepInVolume(); 159 inline G4bool IsLastStepInVolume(); 170 inline void PrepareNewTrack(); 160 inline void PrepareNewTrack(); 171 161 172 inline G4VIntersectionLocator* GetIntersect 162 inline G4VIntersectionLocator* GetIntersectionLocator(); 173 inline void SetIntersectionLocator(G4VInter 163 inline void SetIntersectionLocator(G4VIntersectionLocator* pLocator ); 174 // Change or get the object which calcula 164 // Change or get the object which calculates the exact 175 // intersection point with the next bound 165 // intersection point with the next boundary 176 166 177 inline G4int GetIterationsToIncreaseChordDi 167 inline G4int GetIterationsToIncreaseChordDistance() const; 178 inline void SetIterationsToIncreaseChordDi 168 inline void SetIterationsToIncreaseChordDistance(G4int numIters); 179 // Control the parameter which enables th 169 // Control the parameter which enables the temporary 'relaxation' 180 // which ensures that chord segments ar 170 // which ensures that chord segments are short enough so that 181 // their sagitta is small than delta-ch 171 // their sagitta is small than delta-chord parameter. 182 // The Set method increases the value of 172 // The Set method increases the value of delta-chord temporarily, 183 // doubling it once the number of itera 173 // doubling it once the number of iterations substeps reach 184 // value of 'IncreaseChordDistanceThres 174 // value of 'IncreaseChordDistanceThreshold'. It is also doubled 185 // again every time the iteration count 175 // again every time the iteration count reaches a multiple of this 186 // value. 176 // value. 187 // Note: delta-chord is reset to its orig 177 // Note: delta-chord is reset to its original value at the end of 188 // each call to ComputeStep. 178 // each call to ComputeStep. 189 179 190 public: // without description 180 public: // without description 191 181 192 inline G4double GetDeltaIntersection() cons 182 inline G4double GetDeltaIntersection() const; 193 inline G4double GetDeltaOneStep() const; 183 inline G4double GetDeltaOneStep() const; 194 184 195 inline G4FieldManager* GetCurrentFieldManag 185 inline G4FieldManager* GetCurrentFieldManager(); 196 inline G4EquationOfMotion* GetCurrentEquati 186 inline G4EquationOfMotion* GetCurrentEquationOfMotion(); 197 // Auxiliary methods - their results can 187 // Auxiliary methods - their results can/will change during propagation 198 188 199 inline void SetNavigatorForPropagating(G4Na 189 inline void SetNavigatorForPropagating(G4Navigator* SimpleOrMultiNavigator); 200 inline G4Navigator* GetNavigatorForPropagat 190 inline G4Navigator* GetNavigatorForPropagating(); 201 191 202 inline void SetThresholdNoZeroStep( G4int n 192 inline void SetThresholdNoZeroStep( G4int noAct, 203 G4int n 193 G4int noHarsh, 204 G4int n 194 G4int noAbandon ); 205 inline G4int GetThresholdNoZeroSteps( G4int 195 inline G4int GetThresholdNoZeroSteps( G4int i ); 206 196 207 inline G4double GetZeroStepThreshold(); 197 inline G4double GetZeroStepThreshold(); 208 inline void SetZeroStepThreshold( G4dou 198 inline void SetZeroStepThreshold( G4double newLength ); 209 199 210 void RefreshIntersectionLocator(); 200 void RefreshIntersectionLocator(); 211 // Update the Locator with parameters fro 201 // Update the Locator with parameters from this class 212 // and from current field manager 202 // and from current field manager 213 203 214 protected: // without description 204 protected: // without description 215 205 216 void PrintStepLengthDiagnostic( G4double 206 void PrintStepLengthDiagnostic( G4double currentProposedStepLength, 217 G4double 207 G4double decreaseFactor, 218 G4double 208 G4double stepTrial, 219 const G4FieldTrac 209 const G4FieldTrack& aFieldTrack); 220 210 221 void ReportLoopingParticle( G4int count, G 211 void ReportLoopingParticle( G4int count, G4double StepTaken, 222 G4double stepRe 212 G4double stepRequest, const char* methodName, 223 const G4ThreeVe << 213 G4ThreeVector momentumVec, 224 G4VPhysicalVolu 214 G4VPhysicalVolume* physVol); 225 void ReportStuckParticle(G4int noZeroSteps, 215 void ReportStuckParticle(G4int noZeroSteps, G4double proposedStep, 226 G4double lastTried 216 G4double lastTriedStep, G4VPhysicalVolume* physVol); 227 217 228 private: 218 private: 229 219 230 // ---------------------------------------- 220 // ---------------------------------------------------------------------- 231 // DATA Members 221 // DATA Members 232 // ---------------------------------------- 222 // ---------------------------------------------------------------------- 233 223 234 // ======================================= 224 // ================================================================== 235 // INVARIANTS - Must not change during tra 225 // INVARIANTS - Must not change during tracking 236 226 237 // ** PARAMETERS ----------- 227 // ** PARAMETERS ----------- 238 G4int fMax_loop_count = 1000; 228 G4int fMax_loop_count = 1000; 239 // Limit for the number of sub-steps take 229 // Limit for the number of sub-steps taken in one call to ComputeStep 240 G4int fIncreaseChordDistanceThreshold = 100 230 G4int fIncreaseChordDistanceThreshold = 100; 241 G4bool fUseSafetyForOptimisation = true; 231 G4bool fUseSafetyForOptimisation = true; 242 // (false) is less sensitive to incorrect 232 // (false) is less sensitive to incorrect safety 243 233 244 // Thresholds for identifying "abnormal" c 234 // Thresholds for identifying "abnormal" cases - which cause looping 245 // 235 // 246 G4int fActionThreshold_NoZeroSteps = 2; 236 G4int fActionThreshold_NoZeroSteps = 2; // Threshold # - above it act 247 G4int fSevereActionThreshold_NoZeroSteps = 237 G4int fSevereActionThreshold_NoZeroSteps = 10; // Threshold # to act harshly 248 G4int fAbandonThreshold_NoZeroSteps = 50; 238 G4int fAbandonThreshold_NoZeroSteps = 50; // Threshold # to abandon 249 G4double fZeroStepThreshold = 0.0; 239 G4double fZeroStepThreshold = 0.0; 250 // Threshold *length* for counting of tin 240 // Threshold *length* for counting of tiny or 'zero' steps 251 << 241 252 // Parameters related to handling of very l << 253 // occur typically in large volumes with << 254 G4double fLargestAcceptableStep; 242 G4double fLargestAcceptableStep; 255 // Maximum size of a step - for optimizat 243 // Maximum size of a step - for optimization (and to avoid problems) 256 G4double fMaxStepSizeMultiplier = 3; << 257 // Multiplier for directional exit distan << 258 G4double fMinBigDistance= 100. ; // * CLHEP << 259 // Minimum distance added to directional << 260 // ** End of PARAMETERS ----- 244 // ** End of PARAMETERS ----- 261 245 262 G4double kCarTolerance; 246 G4double kCarTolerance; 263 // Geometrical tolerance defining surfa 247 // Geometrical tolerance defining surface thickness 264 248 265 G4bool fAllocatedLocator; 249 G4bool fAllocatedLocator; // Book-keeping 266 250 267 // --------------------------------------- 251 // -------------------------------------------------------- 268 // ** Dependent Objects - to which work is 252 // ** Dependent Objects - to which work is delegated 269 253 270 G4FieldManager* fDetectorFieldMgr; 254 G4FieldManager* fDetectorFieldMgr; 271 // The Field Manager of the whole Dete 255 // The Field Manager of the whole Detector. (default) 272 256 273 G4VIntersectionLocator* fIntersectionLocato 257 G4VIntersectionLocator* fIntersectionLocator; 274 // Refines candidate intersection 258 // Refines candidate intersection 275 259 276 G4VCurvedTrajectoryFilter* fpTrajectoryFilt 260 G4VCurvedTrajectoryFilter* fpTrajectoryFilter = nullptr; 277 // The filter encapsulates the algorithm 261 // The filter encapsulates the algorithm which selects which 278 // intermediate points should be stored i 262 // intermediate points should be stored in a trajectory. 279 // When it is NULL, no intermediate point 263 // When it is NULL, no intermediate points will be stored. 280 // Else PIF::ComputeStep must submit (all 264 // Else PIF::ComputeStep must submit (all) intermediate 281 // points it calculates, to this filter. 265 // points it calculates, to this filter. (jacek 04/11/2002) 282 266 283 G4Navigator* fNavigator; 267 G4Navigator* fNavigator; 284 // Set externally - only by tracking / ru 268 // Set externally - only by tracking / run manager 285 // 269 // 286 // ** End of Dependent Objects ----------- 270 // ** End of Dependent Objects ---------------------------- 287 271 288 // End of INVARIANTS 272 // End of INVARIANTS 289 // ======================================= 273 // ================================================================== 290 274 291 // STATE information 275 // STATE information 292 // ----------------- 276 // ----------------- 293 G4FieldManager* fCurrentFieldMgr; 277 G4FieldManager* fCurrentFieldMgr; 294 // The Field Manager of the current volu 278 // The Field Manager of the current volume (may be the global) 295 G4bool fSetFieldMgr = false; // Has it bee 279 G4bool fSetFieldMgr = false; // Has it been set for the current step? 296 280 297 // Parameters of current step 281 // Parameters of current step 298 G4double fEpsilonStep; // Relati 282 G4double fEpsilonStep; // Relative accuracy of current Step 299 G4FieldTrack End_PointAndTangent; // End po 283 G4FieldTrack End_PointAndTangent; // End point storage 300 G4bool fParticleIsLooping = false; 284 G4bool fParticleIsLooping = false; 301 G4int fNoZeroStep = 0; // Count 285 G4int fNoZeroStep = 0; // Count of zero Steps 302 286 303 // State used for Optimisation 287 // State used for Optimisation 304 G4double fFull_CurveLen_of_LastAttempt = -1 288 G4double fFull_CurveLen_of_LastAttempt = -1; 305 G4double fLast_ProposedStepLength = -1; 289 G4double fLast_ProposedStepLength = -1; 306 // Previous step information -- for use i 290 // Previous step information -- for use in adjust step size 307 G4ThreeVector fPreviousSftOrigin; 291 G4ThreeVector fPreviousSftOrigin; 308 G4double fPreviousSafety = 0.0; 292 G4double fPreviousSafety = 0.0; 309 // Last safety origin & value: for optimi 293 // Last safety origin & value: for optimisation 310 294 311 G4int fVerboseLevel = 0; 295 G4int fVerboseLevel = 0; 312 G4bool fVerbTracePiF = false; 296 G4bool fVerbTracePiF = false; 313 G4bool fCheck = false; 297 G4bool fCheck = false; 314 // For debugging purposes 298 // For debugging purposes 315 299 316 G4bool fFirstStepInVolume = true; 300 G4bool fFirstStepInVolume = true; 317 G4bool fLastStepInVolume = true; 301 G4bool fLastStepInVolume = true; 318 G4bool fNewTrack = true; 302 G4bool fNewTrack = true; 319 }; 303 }; 320 304 321 // Inline methods 305 // Inline methods 322 // 306 // 323 #include "G4PropagatorInField.icc" 307 #include "G4PropagatorInField.icc" 324 308 325 #endif 309 #endif 326 310