include/Geant4/G4VIntersectionLocator.hh

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 // G4VIntersectionLocator
0027 //
0028 // Class description:
0029 //
0030 // Base class for the calculation of the intersection point with a boundary
0031 // when PropagationInField is used.
0032 // Gives possibility to choose the method of intersection; concrete locators
0033 // implemented are: G4SimpleLocator, G4MultiLevelLocator, G4BrentLocator.
0034 // Key Method: EstimateIntersectionPoint().
0035
0036 // Authors: John Apostolakis, Tatiana Nikitina (CERN), 27 October 2008
0037 // --------------------------------------------------------------------
0038 #ifndef G4VINTERSECTIONLOCATOR_HH
0039 #define G4VINTERSECTIONLOCATOR_HH 1
0040
0041 #include "G4Types.hh"
0042 #include "G4ThreeVector.hh"
0043 #include "G4FieldTrack.hh"
0044
0045 #include "G4Navigator.hh"
0046 #include "G4ChordFinder.hh"
0047
0048 /**
0049  * @brief G4VIntersectionLocator is a base class for the calculation of the
0050  * intersection point with a boundary when PropagationInField is used.
0051  * It gives the possibility to choose the method of intersection.
0052  */
0053
0054 class G4VIntersectionLocator
0055 {
0056   public:
0057
0058     /**
0059      * Constructor and virtual Destructor.
0060      */
0061     G4VIntersectionLocator(G4Navigator *theNavigator);
0062     virtual ~G4VIntersectionLocator();
0063
0064     /**
0065      * If such an intersection exists, this method calculates the intersection
0066      * point of the true path of the particle with the surface of the current
0067      * volume (or of one of its daughters).
0068      * Should use lateral displacement as measure of convergence.
0069      *  @note Changes the safety!
0070      *  @param[in] curveStartPointTangent Start point tangent track.
0071      *  @param[in] curveEndPointTangent End point tangent track.
0072      *  @param[in] trialPoint Trial point.
0073      *  @param[out] intersectPointTangent Intersection point tangent track.
0074      *  @param[out] recalculatedEndPoint Flagging if end point was recomputed.
0075      *  @param[in,out] fPreviousSafety Previous safety distance.
0076      *  @param[in,out] fPreviousSftOrigin Previous safety point origin.
0077      *  @returns Whether intersection exists or not.
0078      */
0079     virtual G4bool EstimateIntersectionPoint(
0080                      const G4FieldTrack&  curveStartPointTangent,  // A
0081                      const G4FieldTrack&  curveEndPointTangent,    // B
0082                      const G4ThreeVector& trialPoint,              // E
0083                            G4FieldTrack&  intersectPointTangent,   // Output
0084                            G4bool&        recalculatedEndPoint,    // Out
0085                            G4double&      fPreviousSafety,         // In/Out
0086                            G4ThreeVector& fPreviousSftOrigin) = 0; // In/Out
0087
0088     /**
0089      * Intersects the chord from StartPointA to EndPointB and returns
0090      * whether an intersection occurred.
0091      *  @note Changes the safety!
0092      *  @param[in] StartPointA Chord starting point.
0093      *  @param[in] EndPointB Chord end point.
0094      *  @param[out] NewSafety New calculated safety distance.
0095      *  @param[in,out] PreviousSafety Previous safety distance.
0096      *  @param[in,out] PreviousSftOrigin Previous safety point origin.
0097      *  @param[out] LinearStepLength Linear chord length.
0098      *  @param[out] IntersectionPoint Intersection point.
0099      *  @param[in,out] calledNavigator Pointer to flag indicating if the
0100      *                 navigator has been called or not.
0101      *  @returns Whether intersection exists or not.
0102      */
0103     inline G4bool IntersectChord( const G4ThreeVector&  StartPointA,
0104                                   const G4ThreeVector&  EndPointB,
0105                                   G4double&      NewSafety,
0106                                   G4double&      PreviousSafety,    // In/Out
0107                                   G4ThreeVector& PreviousSftOrigin, // In/Out
0108                                   G4double&      LinearStepLength,
0109                                   G4ThreeVector& IntersectionPoint,
0110                                   G4bool*        calledNavigator = nullptr );
0111
0112     /**
0113      * Setters for parameters which must be set at each step, in case they are
0114      * changed.
0115      *  @note This simple approach ensures that all scenarios are considered.
0116      *        Future refinement may identify which are invariant during a
0117      *        track, run or event.
0118      */
0119     inline void SetEpsilonStepFor( G4double EpsilonStep );
0120     inline void SetDeltaIntersectionFor( G4double deltaIntersection );
0121     inline void SetNavigatorFor( G4Navigator* fNavigator );
0122     inline void SetChordFinderFor(G4ChordFinder* fCFinder );
0123
0124     /**
0125      * Verbosity control.
0126      * Controlling verbosity enables checking of the locating of intersections.
0127      */
0128     inline void  SetVerboseFor(G4int fVerbose);
0129     inline G4int GetVerboseFor();
0130
0131     /**
0132      * Additional inline Get/Set methods for parameters, dependent objects.
0133      */
0134     inline G4double       GetDeltaIntersectionFor();
0135     inline G4double       GetEpsilonStepFor();
0136     inline G4Navigator*   GetNavigatorFor();
0137     inline G4ChordFinder* GetChordFinderFor();
0138     inline void SetSafetyParametersFor(G4bool UseSafety );
0139
0140     /**
0141      * Adjustment flag accessor/modifier.
0142      */
0143     inline void AdjustIntersections(G4bool UseCorrection);
0144     inline G4bool AreIntersectionsAdjusted(){ return fUseNormalCorrection; }
0145
0146     /**
0147      * Adjustment flag accessor/modifier.
0148      *  @deprecated Replaced by methods above. To be removed in future releases.
0149      */
0150     inline void AddAdjustementOfFoundIntersection(G4bool UseCorrection);
0151     inline G4bool GetAdjustementOfFoundIntersection();
0152
0153     /**
0154      * Dumps status of propagator to any ostream.
0155      */
0156     static void printStatus( const G4FieldTrack& startFT,
0157                              const G4FieldTrack& currentFT,
0158                                    G4double      requestStep,
0159                                    G4double      safety,
0160                                    G4int         stepNum,
0161                                    std::ostream& oss,
0162                                    G4int         verboseLevel );
0163
0164     /**
0165      * Dumps status of propagator to cout, useful mostly for debugging.
0166      */
0167     void printStatus( const G4FieldTrack& startFT,
0168                       const G4FieldTrack& currentFT,
0169                             G4double      requestStep,
0170                             G4double      safety,
0171                             G4int         stepNum);
0172
0173     /**
0174      * Sets/gets check mode.
0175      * When enabled, uses additional verifications and stricter condictions
0176      * for ensuring correctness. Effective only when G4VERBOSE is enabled.
0177      */
0178     inline void SetCheckMode( G4bool value ) { fCheckMode = value; }
0179     inline G4bool GetCheckMode()             { return fCheckMode; }
0180
0181   protected:
0182
0183     /**
0184      * Returns new estimate for state after curveDist starting from
0185      * CurrentStateA, to replace EstimtdEndStateB, and reports displacement
0186      * (if field is compiled verbose).
0187      *  @param[in] CurrentStateA Start point tangent track.
0188      *  @param[in] EstimtdEndStateB Estimated end point tangent track.
0189      *  @param[in] linearDistSq Not used.
0190      *  @param[in] curveDist Not used.
0191      *  @returns New estimate for state.
0192      */
0193     G4FieldTrack ReEstimateEndpoint( const G4FieldTrack& CurrentStateA,
0194                                      const G4FieldTrack& EstimtdEndStateB,
0195                                            G4double linearDistSq, // not used
0196                                            G4double curveDist );  // not used
0197
0198     /**
0199      * Checks whether EndB is too far from StartA to be reached and if,
0200      * re-estimates new value for EndB (return in RevisedEndPoint).
0201      * Reports error if EndB is before StartA (in curve length)
0202      * In that case return errorCode = 2.
0203      *  @param[in] CurrentStartA Start point tangent track.
0204      *  @param[in] EstimatedEndB Estimated end point tangent track.
0205      *  @param[in] RevisedEndPoint Revised end point tangent track.
0206      *  @param[in] errorCode Error code (0=OK, 1=coincident points, 2=error).
0207      *  @returns true if end point has been revised.
0208      */
0209     G4bool CheckAndReEstimateEndpoint( const G4FieldTrack& CurrentStartA,
0210                                        const G4FieldTrack& EstimatedEndB,
0211                                              G4FieldTrack& RevisedEndPoint,
0212                                              G4int&        errorCode);
0213
0214     /**
0215      * Returns the surface normal. Position *must* be the intersection point
0216      * from last call to G4Navigator's ComputeStep() (via IntersectChord).
0217      * It tries to use cached (last) value in Navigator for speed, if it was
0218      * kept and valid. The value returned is in global coordinates.
0219      *  @note It does NOT guarantee to obtain Normal. This can happen e.g.
0220      *        if the "Intersection" Point is not on a surface, potentially
0221      *        due to either inaccuracies in the transformations used, or
0222      *        issues with the Solid.
0223      *  @param[in] CurrentInt_Point Current point.
0224      *  @param[in,out] validNormal Flagging if normal is a valid vector.
0225      *  @returns The surface normal vector in local coordinates.
0226      */
0227     G4ThreeVector GetSurfaceNormal(const G4ThreeVector& CurrentInt_Point,
0228                                          G4bool& validNormal);
0229
0230     /**
0231      * Returns the surface normal of the Intersecting Solid in global
0232      * coordinates.
0233      *  @note This method is costlier then GetSurfaceNormal().
0234      *  @param[in] CurrentInt_Point Current point.
0235      *  @param[in,out] validNormal Flagging if normal is a valid vector.
0236      *  @returns The surface normal vector in global coordinates.
0237      */
0238     G4ThreeVector GetGlobalSurfaceNormal(const G4ThreeVector& CurrentE_Point,
0239                                                G4bool& validNormal);
0240     /**
0241      * Optional method for adjustment of located intersection point using
0242      * the surface-normal.
0243      *  @param[in] A Chord starting point.
0244      *  @param[in] CurrentE_Point E Chord point.
0245      *  @param[in] CurrentF_Point F Chord point.
0246      *  @param[in] MomentumDir Momentum direction.
0247      *  @param[in] IntersectAF First part intersecting?
0248      *  @param[in,out] IntersectionPoint Intersection point tangent track.
0249      *  @param[in,out] NewSafety New safety distance.
0250      *  @param[in,out] fPrevSafety Previous safety distance.
0251      *  @param[in,out] fPrevSftOrigin Previous safety point origin.
0252      *  @returns Whether intersection exists or not.
0253      */
0254     G4bool AdjustmentOfFoundIntersection(const G4ThreeVector& A,
0255                                          const G4ThreeVector& CurrentE_Point,
0256                                          const G4ThreeVector& CurrentF_Point,
0257                                          const G4ThreeVector& MomentumDir,
0258                                          const G4bool         IntersectAF,
0259                                                G4ThreeVector& IntersectionPoint,
0260                                                G4double&      NewSafety,
0261                                                G4double&      fPrevSafety,
0262                                                G4ThreeVector& fPrevSftOrigin );
0263
0264     /**
0265      * Prints a three-line report on the current "sub-step",
0266      * i.e. trial intersection.
0267      *  @param[in] step_no Step number.
0268      *  @param[in] ChordAB_v AB chord.
0269      *  @param[in] ChordEF_v EF chord.
0270      *  @param[in] NewMomentumDir Momentum direction.
0271      *  @param[in] NormalAtEntry Normal vector at entry.
0272      *  @param[in] validNormal Validity flag for normal vector at E.
0273      */
0274     void ReportTrialStep( G4int step_no,
0275                           const G4ThreeVector& ChordAB_v,
0276                           const G4ThreeVector& ChordEF_v,
0277                           const G4ThreeVector& NewMomentumDir,
0278                           const G4ThreeVector& NormalAtEntry,
0279                           G4bool validNormal   );
0280
0281     /**
0282      * Locates a point using the navigator and updates the state of Navigator.
0283      * By default, it assumes that the point is inside the current volume,
0284      * and returns true.
0285      * In check mode, it checks whether the point is *inside* the volume.
0286      * If it is inside, it returns true.
0287      * If not, issues a warning and returns false.
0288      *  @param[in] pos The point to locate.
0289      *  @returns If a point is inside the volume or not.
0290      */
0291     G4bool LocateGlobalPointWithinVolumeAndCheck( const G4ThreeVector& pos );
0292
0293     /**
0294      * Locates a point using the navigator and updates the state of Navigator,
0295      * but report information about code location.
0296      * If CheckMode > 1, report extra information.
0297      *  @param[in] pos The point to locate.
0298      *  @param[in] CodeLocationInfo String for code location info.
0299      *  @param[in] CheckMode Not used.
0300      */
0301     void LocateGlobalPointWithinVolumeCheckAndReport( const G4ThreeVector& pos,
0302                                             const G4String& CodeLocationInfo,
0303                                                   G4int     CheckMode );
0304
0305     // Auxiliary methods -- to report issues
0306
0307     /**
0308      * Builds error message (in ossMsg) to report that point 'B' has
0309      * gone past 'A'.
0310      */
0311     void ReportReversedPoints( std::ostringstream& ossMsg,
0312                                const G4FieldTrack& StartPointVel,
0313                                const G4FieldTrack& EndPointVel,
0314                                      G4double NewSafety, G4double epsStep,
0315                                const G4FieldTrack& CurrentA_PointVelocity,
0316                                const G4FieldTrack& CurrentB_PointVelocity,
0317                                const G4FieldTrack& SubStart_PointVelocity,
0318                                const G4ThreeVector& CurrentE_Point,
0319                                const G4FieldTrack& ApproxIntersecPointV,
0320                                G4int sbstp_no, G4int sbstp_no_p, G4int depth );
0321
0322     /**
0323      * Reports the current status / progress in finding the first intersection.
0324      */
0325     void ReportProgress( std::ostream& oss,
0326                          const G4FieldTrack& StartPointVel,
0327                          const G4FieldTrack& EndPointVel,
0328                                G4int         substep_no,
0329                          const G4FieldTrack& A_PtVel,    // G4double safetyA
0330                          const G4FieldTrack& B_PtVel,
0331                                G4double      safetyLast,
0332                                G4int         depth= -1 );
0333
0334     /**
0335      * Report case: trial point is 'close' to start, within tolerance.
0336      */
0337      void ReportImmediateHit( const char*          MethodName,
0338                               const G4ThreeVector& StartPosition,
0339                               const G4ThreeVector& TrialPoint,
0340                                     G4double       tolerance,
0341                                  unsigned long int numCalls );
0342
0343   private:
0344
0345     /**
0346      * Returns the SurfaceNormal of the Intersecting Solid in local coordinates.
0347      */
0348     G4ThreeVector GetLocalSurfaceNormal(const G4ThreeVector& CurrentE_Point,
0349                                               G4bool& validNormal);
0350
0351     /**
0352      * Position *must* be the intersection point from last call
0353      * to G4Navigator's ComputeStep  (via IntersectChord).
0354      */
0355     G4ThreeVector GetLastSurfaceNormal( const G4ThreeVector& intersectPoint,
0356                                               G4bool& validNormal) const;
0357
0358   protected:
0359
0360     G4double kCarTolerance;                  // Constant
0361
0362     G4int    fVerboseLevel = 0;              // For debugging
0363     G4bool   fUseNormalCorrection = false;   // Configuration parameter
0364     G4bool   fCheckMode = false;
0365     G4bool   fiUseSafety = false;    // Whether to use safety for 'fast steps'
0366
0367     G4Navigator* fiNavigator;
0368
0369     /**
0370      * Parameters set at each physical step by G4PropagatorInField.
0371      */
0372     G4ChordFinder* fiChordFinder = nullptr;  // Overridden at each step
0373     G4double fiEpsilonStep = -1.0;           // Overridden at each step
0374     G4double fiDeltaIntersection = -1.0;     // Overridden at each step
0375
0376     /** Helper for location. */
0377     G4Navigator *fHelpingNavigator;
0378
0379     /** Touchable history hook. */
0380     G4TouchableHistory *fpTouchable = nullptr;
0381 };
0382
0383 #include "G4VIntersectionLocator.icc"
0384
0385 #endif