Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-09-19 09:27:43

0001 // Created on: 1993-12-14
0002 // Created by: Christophe MARION
0003 // Copyright (c) 1993-1999 Matra Datavision
0004 // Copyright (c) 1999-2014 OPEN CASCADE SAS
0005 //
0006 // This file is part of Open CASCADE Technology software library.
0007 //
0008 // This library is free software; you can redistribute it and/or modify it under
0009 // the terms of the GNU Lesser General Public License version 2.1 as published
0010 // by the Free Software Foundation, with special exception defined in the file
0011 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
0012 // distribution for complete text of the license and disclaimer of any warranty.
0013 //
0014 // Alternatively, this file may be used under the terms of Open CASCADE
0015 // commercial license or contractual agreement.
0016 
0017 #ifndef _Extrema_GLocateExtPC_HeaderFile
0018 #define _Extrema_GLocateExtPC_HeaderFile
0019 
0020 #include <Standard.hxx>
0021 #include <Standard_DefineAlloc.hxx>
0022 
0023 #include <GeomAbs_CurveType.hxx>
0024 #include <Precision.hxx>
0025 #include <Standard_DomainError.hxx>
0026 #include <StdFail_NotDone.hxx>
0027 #include <NCollection_Array1.hxx>
0028 
0029 #include <cmath>
0030 
0031 //! Template class for locating extremum of distance between a point and a curve.
0032 //! Calculates the distance with a close point. The close point is defined by
0033 //! the parameter value U0. The function F(u)=distance(P,C(u)) has an extremum
0034 //! when g(u)=dF/du=0. The algorithm searches a zero near the close point.
0035 //!
0036 //! @tparam TheCurve    Curve type (e.g., Adaptor3d_Curve, Adaptor2d_Curve2d)
0037 //! @tparam TheCurveTool Tool for curve operations
0038 //! @tparam ThePoint    Point type (e.g., gp_Pnt, gp_Pnt2d)
0039 //! @tparam TheVector   Vector type (e.g., gp_Vec, gp_Vec2d)
0040 //! @tparam ThePOnC     Point on curve type
0041 //! @tparam TheELPC     Extended local projection curve type
0042 //! @tparam TheLocEPC   Local extremum point curve type
0043 template <typename TheCurve,
0044           typename TheCurveTool,
0045           typename ThePoint,
0046           typename TheVector,
0047           typename ThePOnC,
0048           typename TheELPC,
0049           typename TheLocEPC>
0050 class Extrema_GLocateExtPC
0051 {
0052 public:
0053   DEFINE_STANDARD_ALLOC
0054 
0055   //! Default constructor.
0056   Extrema_GLocateExtPC()
0057       : myC(nullptr),
0058         mydist2(0.0),
0059         myismin(false),
0060         myDone(false),
0061         myumin(0.0),
0062         myusup(0.0),
0063         mytol(0.0),
0064         type(GeomAbs_OtherCurve),
0065         numberext(0)
0066   {
0067   }
0068 
0069   //! Calculates the distance with a close point.
0070   //! The close point is defined by the parameter value U0.
0071   //! TolF is used to decide to stop the iterations.
0072   //! At the nth iteration, the criteria is: abs(Un - Un-1) < TolF.
0073   Extrema_GLocateExtPC(const ThePoint& theP,
0074                        const TheCurve& theC,
0075                        const double    theU0,
0076                        const double    theTolF)
0077   {
0078     Initialize(theC,
0079                TheCurveTool::FirstParameter(theC),
0080                TheCurveTool::LastParameter(theC),
0081                theTolF);
0082     Perform(theP, theU0);
0083   }
0084 
0085   //! Calculates the distance with a close point.
0086   //! The close point is defined by the parameter value U0.
0087   //! Zeros are searched between Umin and Usup.
0088   //! TolF is used to decide to stop the iterations.
0089   //! At the nth iteration, the criteria is: abs(Un - Un-1) < TolF.
0090   Extrema_GLocateExtPC(const ThePoint& theP,
0091                        const TheCurve& theC,
0092                        const double    theU0,
0093                        const double    theUmin,
0094                        const double    theUsup,
0095                        const double    theTolF)
0096   {
0097     Initialize(theC, theUmin, theUsup, theTolF);
0098     Perform(theP, theU0);
0099   }
0100 
0101   //! Sets the fields of the algorithm.
0102   void Initialize(const TheCurve& theC,
0103                   const double    theUmin,
0104                   const double    theUsup,
0105                   const double    theTolF)
0106   {
0107     myC         = const_cast<TheCurve*>(&theC);
0108     mytol       = theTolF;
0109     myumin      = theUmin;
0110     myusup      = theUsup;
0111     type        = TheCurveTool::GetType(theC);
0112     double tolu = TheCurveTool::Resolution(theC, Precision::Confusion());
0113     if ((type == GeomAbs_BSplineCurve) || (type == GeomAbs_BezierCurve)
0114         || (type == GeomAbs_OffsetCurve) || (type == GeomAbs_OtherCurve))
0115     {
0116       myLocExtPC.Initialize(theC, theUmin, theUsup, tolu);
0117     }
0118     else
0119     {
0120       myExtremPC.Initialize(theC, theUmin, theUsup, tolu);
0121     }
0122   }
0123 
0124   //! Performs the algorithm with point P and initial parameter U0.
0125   void Perform(const ThePoint& theP, const double theU0)
0126   {
0127     int    i, i1, i2, inter;
0128     double Par, valU, valU2 = RealLast(), local_u0;
0129     double myintuinf = 0, myintusup = 0;
0130     local_u0 = theU0;
0131     switch (type)
0132     {
0133       case GeomAbs_OtherCurve:
0134       case GeomAbs_OffsetCurve:
0135       case GeomAbs_BSplineCurve: {
0136         // Search for extremum is done interval by continuous C2 interval
0137         int                        n = TheCurveTool::NbIntervals(*myC, GeomAbs_C2);
0138         NCollection_Array1<double> theInter(1, n + 1);
0139         TheCurveTool::Intervals(*myC, theInter, GeomAbs_C2);
0140         //
0141         // be gentle with the caller
0142         //
0143         if (local_u0 < myumin)
0144         {
0145           local_u0 = myumin;
0146         }
0147         else if (local_u0 > myusup)
0148         {
0149           local_u0 = myusup;
0150         }
0151         // Search for interval containing U0
0152         bool found = false;
0153         inter      = 1;
0154         while (!found && inter <= n)
0155         {
0156           myintuinf = std::max(theInter(inter), myumin);
0157           myintusup = std::min(theInter(inter + 1), myusup);
0158           if ((local_u0 >= myintuinf) && (local_u0 < myintusup))
0159             found = true;
0160           inter++;
0161         }
0162 
0163         if (found)
0164           inter--; // IFV 16.06.00 - inter is increased after found!
0165 
0166         // Try on found interval
0167         myLocExtPC.Initialize(*myC, myintuinf, myintusup, mytol);
0168         myLocExtPC.Perform(theP, local_u0);
0169         myDone = myLocExtPC.IsDone();
0170         if (myDone)
0171         {
0172           mypp    = myLocExtPC.Point();
0173           myismin = myLocExtPC.IsMin();
0174           mydist2 = myLocExtPC.SquareDistance();
0175         }
0176         else
0177         {
0178           int k = 1;
0179           // Try on neighboring intervals:
0180           i1 = inter;
0181           i2 = inter;
0182           double    s1inf, s2inf, s1sup, s2sup;
0183           ThePoint  P1;
0184           TheVector V1;
0185           TheCurveTool::D1(*myC, myintuinf, P1, V1);
0186           s2inf = (TheVector(theP, P1) * V1);
0187           TheCurveTool::D1(*myC, myintusup, P1, V1);
0188           s1sup = (TheVector(theP, P1) * V1);
0189 
0190           while (!myDone && (i2 > 0) && (i1 <= n))
0191           {
0192             i1 = inter + k;
0193             i2 = inter - k;
0194             if (i1 <= n)
0195             {
0196               myintuinf = std::max(theInter(i1), myumin);
0197               myintusup = std::min(theInter(i1 + 1), myusup);
0198               if (myintuinf < myintusup)
0199               {
0200                 TheCurveTool::D1(*myC, myintuinf, P1, V1);
0201                 s2sup = (TheVector(theP, P1) * V1);
0202                 if (Precision::IsInfinite(s2sup) || Precision::IsInfinite(s1sup))
0203                 {
0204                   break;
0205                 }
0206                 if (s1sup * s2sup <= RealEpsilon())
0207                 {
0208                   // extremum:
0209                   myDone = true;
0210                   mypp.SetValues(myintuinf, P1);
0211                   myismin = (s1sup <= 0.0);
0212                   mydist2 = theP.SquareDistance(P1);
0213                   break;
0214                 }
0215 
0216                 TheCurveTool::D1(*myC, myintusup, P1, V1);
0217                 s1sup = (TheVector(theP, P1) * V1);
0218                 myLocExtPC.Initialize(*myC, myintuinf, myintusup, mytol);
0219                 myLocExtPC.Perform(theP, (myintuinf + myintusup) * 0.5);
0220                 myDone = myLocExtPC.IsDone();
0221                 if (myDone)
0222                 {
0223                   mypp    = myLocExtPC.Point();
0224                   myismin = myLocExtPC.IsMin();
0225                   mydist2 = myLocExtPC.SquareDistance();
0226                   break;
0227                 }
0228               }
0229             }
0230 
0231             if (i2 > 0)
0232             {
0233               myintuinf = std::max(theInter(i2), myumin);
0234               myintusup = std::min(theInter(i2 + 1), myusup);
0235               if (myintuinf < myintusup)
0236               {
0237                 TheCurveTool::D1(*myC, myintusup, P1, V1);
0238                 s1inf = (TheVector(theP, P1) * V1);
0239                 if (Precision::IsInfinite(s2inf) || Precision::IsInfinite(s1inf))
0240                 {
0241                   break;
0242                 }
0243                 if (s1inf * s2inf <= RealEpsilon())
0244                 {
0245                   // extremum:
0246                   myDone = true;
0247                   mypp.SetValues(myintusup, P1);
0248                   myismin = (s1inf <= 0.0);
0249                   mydist2 = theP.SquareDistance(P1);
0250                   break;
0251                 }
0252 
0253                 TheCurveTool::D1(*myC, myintuinf, P1, V1);
0254                 s2inf = (TheVector(theP, P1) * V1);
0255                 myLocExtPC.Initialize(*myC, myintuinf, myintusup, mytol);
0256                 myLocExtPC.Perform(theP, (myintuinf + myintusup) * 0.5);
0257                 myDone = myLocExtPC.IsDone();
0258 
0259                 if (myDone)
0260                 {
0261                   mypp    = myLocExtPC.Point();
0262                   myismin = myLocExtPC.IsMin();
0263                   mydist2 = myLocExtPC.SquareDistance();
0264                   break;
0265                 }
0266               }
0267             }
0268 
0269             k++;
0270           }
0271         }
0272       }
0273 
0274       break;
0275 
0276       case GeomAbs_BezierCurve: {
0277         myLocExtPC.Perform(theP, theU0);
0278         myDone = myLocExtPC.IsDone();
0279       }
0280 
0281       break;
0282       default: {
0283         myExtremPC.Perform(theP);
0284         numberext = 0;
0285         if (myExtremPC.IsDone())
0286         {
0287           for (i = 1; i <= myExtremPC.NbExt(); i++)
0288           {
0289             Par  = myExtremPC.Point(i).Parameter();
0290             valU = std::abs(Par - theU0);
0291             if (valU <= valU2)
0292             {
0293               valU2     = valU;
0294               numberext = i;
0295               myDone    = true;
0296             }
0297           }
0298         }
0299 
0300         if (numberext == 0)
0301           myDone = false;
0302 
0303         break;
0304       }
0305     }
0306   }
0307 
0308   //! Returns True if the distance is found.
0309   bool IsDone() const { return myDone; }
0310 
0311   //! Returns the value of the extremum square distance.
0312   double SquareDistance() const
0313   {
0314     if (!IsDone())
0315     {
0316       throw StdFail_NotDone();
0317     }
0318     double d = 0;
0319     if ((type == GeomAbs_BezierCurve))
0320     {
0321       d = myLocExtPC.SquareDistance();
0322     }
0323     else if (type == GeomAbs_BSplineCurve || type == GeomAbs_OffsetCurve
0324              || type == GeomAbs_OtherCurve)
0325     {
0326       d = mydist2;
0327     }
0328     else
0329     {
0330       if (numberext != 0)
0331       {
0332         d = myExtremPC.SquareDistance(numberext);
0333       }
0334     }
0335     return d;
0336   }
0337 
0338   //! Returns True if the extremum distance is a minimum.
0339   bool IsMin() const
0340   {
0341     if (!IsDone())
0342     {
0343       throw StdFail_NotDone();
0344     }
0345     bool b = false;
0346     if ((type == GeomAbs_BezierCurve))
0347     {
0348       b = myLocExtPC.IsMin();
0349     }
0350     else if (type == GeomAbs_BSplineCurve || type == GeomAbs_OffsetCurve
0351              || type == GeomAbs_OtherCurve)
0352     {
0353       b = myismin;
0354     }
0355     else
0356     {
0357       if (numberext != 0)
0358       {
0359         b = myExtremPC.IsMin(numberext);
0360       }
0361     }
0362     return b;
0363   }
0364 
0365   //! Returns the point of the extremum distance.
0366   const ThePOnC& Point() const
0367   {
0368     if (!IsDone())
0369     {
0370       throw StdFail_NotDone();
0371     }
0372     if (type == GeomAbs_BezierCurve)
0373     {
0374       return myLocExtPC.Point();
0375     }
0376     else if (type == GeomAbs_BSplineCurve || type == GeomAbs_OffsetCurve
0377              || type == GeomAbs_OtherCurve)
0378     {
0379       return mypp;
0380     }
0381     return myExtremPC.Point(numberext);
0382   }
0383 
0384 private:
0385   ThePOnC           mypp;
0386   TheCurve*         myC;
0387   double            mydist2;
0388   bool              myismin;
0389   bool              myDone;
0390   double            myumin;
0391   double            myusup;
0392   double            mytol;
0393   TheLocEPC         myLocExtPC;
0394   TheELPC           myExtremPC;
0395   GeomAbs_CurveType type;
0396   int               numberext;
0397 };
0398 
0399 #endif // _Extrema_GLocateExtPC_HeaderFile