Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/geometry/navigation/include/G4VIntersectionLocator.hh

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 // Class G4VIntersectionLocator 
 27 //
 28 // class description:
 29 //
 30 // Base class for the calculation of the intersection point with a boundary 
 31 // when PropagationInField is used.
 32 // Gives possibility to choose the method of intersection; concrete locators
 33 // implemented are: G4SimpleLocator, G4MultiLevelLocator, G4BrentLocator.
 34 //
 35 // Key Method: EstimateIntersectionPoint()
 36 
 37 // 27.10.08 - John Apostolakis, Tatiana Nikitina: Design and implementation 
 38 // ---------------------------------------------------------------------------
 39 #ifndef G4VINTERSECTIONLOCATOR_HH
 40 #define G4VINTERSECTIONLOCATOR_HH
 41 
 42 #include "G4Types.hh" 
 43 #include "G4ThreeVector.hh"
 44 #include "G4FieldTrack.hh"
 45 
 46 #include "G4Navigator.hh"
 47 #include "G4ChordFinder.hh"
 48 
 49 class G4VIntersectionLocator
 50  {
 51    public:  // with description 
 52  
 53      G4VIntersectionLocator(G4Navigator *theNavigator);
 54        // Constructor
 55      virtual ~G4VIntersectionLocator();
 56        // Default destructor
 57      
 58      virtual G4bool EstimateIntersectionPoint( 
 59          const G4FieldTrack&       curveStartPointTangent,  // A
 60          const G4FieldTrack&       curveEndPointTangent,    // B
 61          const G4ThreeVector&      trialPoint,              // E
 62                G4FieldTrack&       intersectPointTangent,   // Output
 63                G4bool&             recalculatedEndPoint,    // Out
 64                G4double&           fPreviousSafety,         // In/Out
 65                G4ThreeVector&      fPreviousSftOrigin) = 0; // In/Out   
 66        // If such an intersection exists, this function calculates the
 67        // intersection point of the true path of the particle with the surface
 68        // of the current volume (or of one of its daughters). 
 69        // Should use lateral displacement as measure of convergence
 70        // NOTE: changes the safety!
 71 
 72      void printStatus( const G4FieldTrack& startFT,
 73                        const G4FieldTrack& currentFT, 
 74                              G4double      requestStep, 
 75                              G4double      safety,
 76                              G4int         stepNum);
 77        // Print Method, useful mostly for debugging
 78 
 79      inline G4bool IntersectChord( const G4ThreeVector&  StartPointA,
 80                                    const G4ThreeVector&  EndPointB,
 81                                    G4double&      NewSafety,
 82                                    G4double&      PreviousSafety,    // In/Out
 83                                    G4ThreeVector& PreviousSftOrigin, // In/Out
 84                                    G4double&      LinearStepLength,
 85                                    G4ThreeVector& IntersectionPoint,
 86                                    G4bool*        calledNavigator = nullptr );
 87        // Intersect the chord from StartPointA to EndPointB and return
 88        // whether an intersection occurred. NOTE: changes the Safety!
 89 
 90      inline void SetEpsilonStepFor( G4double EpsilonStep );
 91      inline void SetDeltaIntersectionFor( G4double deltaIntersection );
 92      inline void SetNavigatorFor( G4Navigator* fNavigator );
 93      inline void SetChordFinderFor(G4ChordFinder* fCFinder );
 94        // These parameters must be set at each step, in case they were changed
 95        // Note: This simple approach ensures that all scenarios are considered. 
 96        //       Future refinement may identify which are invariant during a 
 97        //       track, run or event
 98 
 99      inline void  SetVerboseFor(G4int fVerbose);
100      inline G4int GetVerboseFor();
101        // Controling verbosity enables checking of the locating of intersections
102 
103   public:  // without description
104 
105     // Additional inline Set/Get methods for parameters, dependent objects
106 
107     inline G4double       GetDeltaIntersectionFor();
108     inline G4double       GetEpsilonStepFor();
109     inline G4Navigator*   GetNavigatorFor();
110     inline G4ChordFinder* GetChordFinderFor();
111 
112     inline void SetSafetyParametersFor(G4bool UseSafety );
113 
114     inline void AddAdjustementOfFoundIntersection(G4bool UseCorrection);
115     inline G4bool GetAdjustementOfFoundIntersection();
116       // Methods to be made Obsolete - replaced by methods below
117     inline void AdjustIntersections(G4bool UseCorrection); 
118     inline G4bool AreIntersectionsAdjusted(){ return fUseNormalCorrection; }  
119       // Change adjustment flag  ( New Interface ) 
120 
121     static void printStatus( const G4FieldTrack& startFT,
122                              const G4FieldTrack& currentFT, 
123                                    G4double      requestStep, 
124                                    G4double      safety,
125                                    G4int         stepNum,
126                                    std::ostream& oss,
127                                    G4int         verboseLevel );
128       // Print Method for any ostream - e.g. cerr -- and for G4Exception
129 
130     inline void SetCheckMode( G4bool value ) { fCheckMode = value; }
131     inline G4bool GetCheckMode()             { return fCheckMode; }
132 
133   protected:  // with description
134 
135     G4FieldTrack ReEstimateEndpoint( const G4FieldTrack& CurrentStateA,  
136                                      const G4FieldTrack& EstimtdEndStateB,
137                                            G4double linearDistSq, // not used
138                                            G4double curveDist );  // not used
139       // Return new estimate for state after curveDist starting from
140       // CurrentStateA, to replace EstimtdEndStateB, and report displacement
141       // (if field is compiled verbose)
142 
143     G4bool CheckAndReEstimateEndpoint( const G4FieldTrack& CurrentStartA,  
144                                        const G4FieldTrack& EstimatedEndB,
145                                              G4FieldTrack& RevisedEndPoint,
146                                              G4int &       errorCode);
147       // Check whether EndB is too far from StartA to be reached 
148       // and if, re-estimate new value for EndB (return in RevisedEndPoint)
149       // Report error if  EndB is before StartA (in curve length)
150       // In that case return errorCode = 2.
151 
152     G4ThreeVector GetSurfaceNormal(const G4ThreeVector& CurrentInt_Point,
153                                          G4bool& validNormal);
154       // Position *must* be the intersection point from last call
155       // to G4Navigator's ComputeStep  (via IntersectChord )
156       // Will try to use cached (last) value in Navigator for speed, 
157       // if it was kept and valid.
158       // Value returned is in global coordinates.
159       // It does NOT guarantee to obtain Normal. This can happen eg if:
160       //  - the "Intersection" Point is not on a surface, potentially due to
161       //  - inaccuracies in the transformations used, or
162       //  - issues with the Solid.
163 
164     G4ThreeVector GetGlobalSurfaceNormal(const G4ThreeVector& CurrentE_Point,
165                                                G4bool& validNormal);
166       // Return the SurfaceNormal of Intersecting Solid in global coordinates
167       // Costlier then GetSurfaceNormal
168 
169     G4bool AdjustmentOfFoundIntersection(const G4ThreeVector& A,
170                                          const G4ThreeVector& CurrentE_Point, 
171                                          const G4ThreeVector& CurrentF_Point,
172                                          const G4ThreeVector& MomentumDir,
173                                          const G4bool         IntersectAF, 
174                                                G4ThreeVector& IntersectionPoint,
175                                                G4double&      NewSafety,
176                                                G4double&      fPrevSafety,
177                                                G4ThreeVector& fPrevSftOrigin );
178       // Optional method for adjustment of located intersection point
179       // using the surface-normal
180   
181     void ReportTrialStep( G4int step_no, 
182                           const G4ThreeVector& ChordAB_v,
183                           const G4ThreeVector& ChordEF_v,
184                           const G4ThreeVector& NewMomentumDir,
185                           const G4ThreeVector& NormalAtEntry,
186                           G4bool validNormal   );
187       // Print a three-line report on the current "sub-step", ie trial
188       // intersection
189 
190     G4bool LocateGlobalPointWithinVolumeAndCheck( const G4ThreeVector& pos );
191       // Locate point using navigator - updates state of Navigator.
192       // By default, it assumes that the point is inside the current volume,
193       // and returns true.
194       // In check mode, checks whether the point is *inside* the volume.
195       //   If it is inside, it returns true.
196       //   If not, issues a warning and returns false.
197 
198     void LocateGlobalPointWithinVolumeCheckAndReport( const G4ThreeVector& pos,
199                                             const G4String& CodeLocationInfo,
200                                                   G4int     CheckMode );
201       // As above, but report information about code location.
202       // If CheckMode > 1, report extra information.
203 
204   protected:  // without description
205 
206     // Auxiliary methods -- to report issues
207 
208     void ReportReversedPoints( std::ostringstream& ossMsg,
209                                const G4FieldTrack& StartPointVel, 
210                                const G4FieldTrack& EndPointVel,
211                                      G4double NewSafety, G4double epsStep,
212                                const G4FieldTrack& CurrentA_PointVelocity,
213                                const G4FieldTrack& CurrentB_PointVelocity,
214                                const G4FieldTrack& SubStart_PointVelocity,
215                                const G4ThreeVector& CurrentE_Point,
216                                const G4FieldTrack& ApproxIntersecPointV,
217                                G4int sbstp_no, G4int sbstp_no_p, G4int depth );
218       // Build error messsage (in ossMsg) to report that point 'B' has
219       // gone past 'A'
220 
221     void ReportProgress( std::ostream& oss,
222                          const G4FieldTrack& StartPointVel, 
223                          const G4FieldTrack& EndPointVel,
224                                G4int         substep_no, 
225                          const G4FieldTrack& A_PtVel,    // G4double safetyA
226                          const G4FieldTrack& B_PtVel,  
227                                G4double      safetyLast,
228                                G4int         depth= -1 );
229       // Report the current status / progress in finding the first intersection
230 
231      void ReportImmediateHit( const char*          MethodName, 
232                               const G4ThreeVector& StartPosition, 
233                               const G4ThreeVector& TrialPoint, 
234                                     G4double       tolerance,
235                                  unsigned long int numCalls );
236       // Report case: trial point is 'close' to start, within tolerance
237     
238   private:  // no description
239 
240     G4ThreeVector GetLocalSurfaceNormal(const G4ThreeVector& CurrentE_Point,
241                                               G4bool& validNormal);
242       // Return the SurfaceNormal of Intersecting Solid  in local coordinates
243 
244     G4ThreeVector GetLastSurfaceNormal( const G4ThreeVector& intersectPoint,
245                                               G4bool& validNormal) const; 
246       // Position *must* be the intersection point from last call
247       // to G4Navigator's ComputeStep  (via IntersectChord )
248       // Temporary - will use the same method in the Navigator
249 
250   protected:
251 
252     G4double kCarTolerance;                  // Constant
253 
254     G4int    fVerboseLevel = 0;              // For debugging
255     G4bool   fUseNormalCorrection = false;   // Configuration parameter
256     G4bool   fCheckMode = false;
257     G4bool   fiUseSafety = false;    // Whether to use safety for 'fast steps'
258    
259     G4Navigator* fiNavigator;
260 
261     G4ChordFinder* fiChordFinder = nullptr;  // Overridden at each step
262     G4double fiEpsilonStep = -1.0;           // Overridden at each step
263     G4double fiDeltaIntersection = -1.0;     // Overridden at each step
264       // Parameters set at each physical step by calling method 
265       // by G4PropagatorInField
266 
267     G4Navigator *fHelpingNavigator;
268       // Helper for location
269 
270     G4TouchableHistory *fpTouchable = nullptr;
271       // Touchable history hook
272 };
273 
274 #include "G4VIntersectionLocator.icc"
275 
276 #endif
277