0028599: Replacement of old Boolean operations with new ones in BRepProj_Projection...
[occt.git] / src / Extrema / Extrema_ExtCS.cxx
index 3bc81ae..797e001 100644 (file)
 
 //  Modified by skv - Thu Jul  7 12:29:34 2005 OCC9134
 
-#include <Extrema_ExtCS.ixx>
+#include <Adaptor3d_Curve.hxx>
+#include <Adaptor3d_Surface.hxx>
+#include <Bnd_Box.hxx>
+#include <BndLib_AddSurface.hxx>
+#include <ElCLib.hxx>
+#include <ElSLib.hxx>
+#include <Extrema_ExtCS.hxx>
+#include <Extrema_ExtPElC.hxx>
+#include <Extrema_ExtPElS.hxx>
+#include <Extrema_ExtPS.hxx>
 #include <Extrema_GenExtCS.hxx>
-#include <StdFail_NotDone.hxx>
-#include <Standard_NotImplemented.hxx>
-#include <StdFail_InfiniteSolutions.hxx>
-#include <Precision.hxx>
+#include <Extrema_POnCurv.hxx>
+#include <Extrema_POnSurf.hxx>
 #include <GeomAbs_CurveType.hxx>
-#include <gp_Pnt.hxx>
-#include <gp_Pln.hxx>
-#include <gp_Cylinder.hxx>
 #include <gp_Cone.hxx>
+#include <gp_Cylinder.hxx>
+#include <gp_Lin.hxx>
+#include <gp_Pln.hxx>
+#include <gp_Pnt.hxx>
 #include <gp_Sphere.hxx>
 #include <gp_Torus.hxx>
-#include <gp_Lin.hxx>
-
-#include <ElCLib.hxx>
-#include <Extrema_ExtPElC.hxx>
-#include <Bnd_Box.hxx>
-#include <BndLib_AddSurface.hxx>
-#include <Extrema_ExtPElS.hxx>
+#include <Precision.hxx>
+#include <Standard_NotImplemented.hxx>
+#include <Standard_OutOfRange.hxx>
+#include <Standard_TypeMismatch.hxx>
+#include <StdFail_InfiniteSolutions.hxx>
+#include <StdFail_NotDone.hxx>
 #include <TColStd_Array1OfReal.hxx>
-#include <Extrema_ExtPS.hxx>
-#include <ElSLib.hxx>
 
 Extrema_ExtCS::Extrema_ExtCS() 
 {
@@ -108,6 +113,9 @@ void Extrema_ExtCS::Perform(const Adaptor3d_Curve& C,
   NbT = NbU = NbV = 10;
   GeomAbs_CurveType myCtype  = C.GetType();
 
+  myDone = Standard_False;
+  // Try analytic computation of extrema
+  Standard_Boolean isComputeAnalytic = Standard_True;
 
   switch(myCtype) {
 
@@ -124,6 +132,7 @@ void Extrema_ExtCS::Perform(const Adaptor3d_Curve& C,
       case GeomAbs_Plane:
         myExtElCS.Perform(C.Line(), myS->Plane());
         if (myExtElCS.IsParallel())   break;
+        Standard_FALLTHROUGH
 
       case GeomAbs_Torus:
       case GeomAbs_Cone:
@@ -187,7 +196,10 @@ void Extrema_ExtCS::Perform(const Adaptor3d_Curve& C,
 
           }
 
-
+          if (myS->IsUPeriodic())
+            NbU = 13;
+          if (myS->IsVPeriodic())
+            NbV = 13;
 
           Extrema_GenExtCS Ext(C, *myS, NbT, NbU, NbV, cfirst, clast, ufirst, ulast,
             vfirst, vlast, mytolC, mytolS);
@@ -224,7 +236,13 @@ void Extrema_ExtCS::Perform(const Adaptor3d_Curve& C,
         myExtElCS.Perform(C.Circle(), myS->Plane());
         break;
       }
+      else if (myStype == GeomAbs_Sphere)
+      {
+        myExtElCS.Perform(C.Circle(), myS->Sphere());
+        break;
+      }
     }
+    Standard_FALLTHROUGH
   case GeomAbs_Hyperbola: 
     {
       if(myCtype == GeomAbs_Hyperbola && myStype == GeomAbs_Plane) {
@@ -233,157 +251,228 @@ void Extrema_ExtCS::Perform(const Adaptor3d_Curve& C,
         break;
       }
     }
+    Standard_FALLTHROUGH
   default:
     {
-      Extrema_GenExtCS Ext;
-      Ext.Initialize(*myS, NbU, NbV, mytolS);
-      if(myCtype == GeomAbs_Hyperbola) {
-        Standard_Real tmin = Max(-20., C.FirstParameter());
-        Standard_Real tmax = Min(20., C.LastParameter());
-        Ext.Perform(C, NbT, tmin, tmax, mytolC); // to avoid overflow
-      }
-      else {
-        if(myCtype == GeomAbs_Circle && NbT < 13) {
-          NbT = 13;
-        }
-        Ext.Perform(C, NbT, mytolC);
-      }
+      isComputeAnalytic = Standard_False;
+      break;
+    }
+  }
 
-      myDone = Ext.IsDone();
-      if (myDone) {
-        Standard_Integer NbExt = Ext.NbExt();
-        Standard_Real T,U,V;
-        Extrema_POnCurv PC;
-        Extrema_POnSurf PS;
-        for (i = 1; i <= NbExt; i++) {
-          PC = Ext.PointOnCurve(i);
-          PS = Ext.PointOnSurface(i);
-          T = PC.Parameter();
+  if (isComputeAnalytic)
+  {
+    if (myExtElCS.IsDone())
+    {
+      myDone = Standard_True;
+      myIsPar = myExtElCS.IsParallel();
+      if (myIsPar)
+      {
+        mySqDist.Append(myExtElCS.SquareDistance(1));
+      }
+      else
+      {
+        Standard_Integer NbExt = myExtElCS.NbExt();
+        for (i = 1; i <= NbExt; i++)
+        {
+          Extrema_POnCurv PC;
+          Extrema_POnSurf PS;
+          myExtElCS.Points(i, PC, PS);
+          Standard_Real Ucurve = PC.Parameter();
+          Standard_Real U, V;
           PS.Parameter(U, V);
-          AddSolution(C, T, U, V, PC.Value(), PS.Value(), Ext.SquareDistance(i));
+          AddSolution(C, Ucurve, U, V, PC.Value(), PS.Value(), myExtElCS.SquareDistance(i));
         }
 
-        //Add sharp points
-        Standard_Integer SolNumber = mySqDist.Length();
-        Standard_Address CopyC = (Standard_Address)&C;
-        Adaptor3d_Curve& aC = *(Adaptor3d_Curve*)CopyC;
-        Standard_Integer NbIntervals = aC.NbIntervals(GeomAbs_C1);
-        TColStd_Array1OfReal SharpPoints(1, NbIntervals+1);
-        aC.Intervals(SharpPoints, GeomAbs_C1);
-        
-        Extrema_ExtPS aProjPS;
-        aProjPS.Initialize (*myS,
-                            myS->FirstUParameter(),
-                            myS->LastUParameter(),
-                            myS->FirstVParameter(),
-                            myS->LastVParameter(),
-                            mytolS,
-                            mytolS);
-
-        for (i = 2; i < SharpPoints.Upper(); ++i)
+        if (mySqDist.Length() == 0 && NbExt > 0)
         {
-          T = SharpPoints(i);
-          gp_Pnt aPnt = C.Value(T);
-          aProjPS.Perform (aPnt);
-          if (!aProjPS.IsDone())
-            continue;
-          Standard_Integer NbProj = aProjPS.NbExt(), jmin = 0;
-          Standard_Real MinSqDist = RealLast();
-          for (j = 1; j <= NbProj; j++)
+          // Analytical extrema seem to be out of curve/surface boundaries.
+          // Try extremity points of curve.
+          gp_Pnt aPOnC[2], aPOnS[2];
+          Standard_Real aT[2] = { myucinf, myucsup }, U[2], V[2];
+          Standard_Real aDist[2] = { -1, -1 };
+          for (i = 0; i < 2; ++i)
           {
-            Standard_Real aSqDist = aProjPS.SquareDistance(j);
-            if (aSqDist < MinSqDist)
+            if (Precision::IsInfinite(aT[i]))
+              continue;
+
+            aPOnC[i] = C.Value(aT[i]);
+            switch (myStype)
             {
-              MinSqDist = aSqDist;
-              jmin = j;
+              case GeomAbs_Plane:
+              {
+                ElSLib::Parameters(myS->Plane(), aPOnC[i], U[i], V[i]);
+                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Plane());
+                break;
+              }
+              case GeomAbs_Sphere:
+              {
+                ElSLib::Parameters(myS->Sphere(), aPOnC[i], U[i], V[i]);
+                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Sphere());
+                break;
+              }
+              case GeomAbs_Cylinder:
+              {
+                ElSLib::Parameters(myS->Cylinder(), aPOnC[i], U[i], V[i]);
+                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Cylinder());
+                break;
+              }
+              case GeomAbs_Torus:
+              {
+                ElSLib::Parameters(myS->Torus(), aPOnC[i], U[i], V[i]);
+                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Torus());
+                break;
+              }
+              case GeomAbs_Cone:
+              {
+                ElSLib::Parameters(myS->Cone(), aPOnC[i], U[i], V[i]);
+                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Cone());
+                break;
+              }
+              default:
+                continue;
             }
+
+            aDist[i] = aPOnC[i].SquareDistance(aPOnS[i]);
           }
-          if (jmin != 0)
+
+          Standard_Boolean bAdd[2] = {Standard_False, Standard_False};
+
+          // Choose solution to add
+          if (aDist[0] >= 0. && aDist[1] >= 0.)
           {
-            aProjPS.Point(jmin).Parameter(U,V);
-            AddSolution(C, T, U, V,
-                        aPnt, aProjPS.Point(jmin).Value(), MinSqDist);
+            Standard_Real aDiff = aDist[0] - aDist[1];
+            // Both computed -> take only minimal
+            if (Abs(aDiff) < Precision::Confusion())
+              // Add both
+              bAdd[0] = bAdd[1] = Standard_True;
+            else if (aDiff < 0)
+              // Add first
+              bAdd[0] = Standard_True;
+            else
+              // Add second
+              bAdd[1] = Standard_True;
+          }
+          else if (aDist[0] >= 0.)
+            // Add first
+            bAdd[0] = Standard_True;
+          else if (aDist[1] >= 0.)
+            // Add second
+            bAdd[1] = Standard_True;
+
+          for (i = 0; i < 2; ++i)
+          {
+            if (bAdd[i])
+              AddSolution(C, aT[i], U[i], V[i], aPOnC[i], aPOnS[i], aDist[i]);
           }
-        }
-        //Cut sharp solutions to keep only minimum and maximum
-        Standard_Integer imin = SolNumber + 1, imax = mySqDist.Length();
-        for (i = SolNumber + 1; i <= mySqDist.Length(); i++)
-        {
-          if (mySqDist(i) < mySqDist(imin))
-            imin = i;
-          if (mySqDist(i) > mySqDist(imax))
-            imax = i;
-        }
-        if (mySqDist.Length() > SolNumber + 2)
-        {
-          Standard_Real MinSqDist = mySqDist(imin);
-          Standard_Real MaxSqDist = mySqDist(imax);
-          Extrema_POnCurv MinPC = myPOnC(imin);
-          Extrema_POnCurv MaxPC = myPOnC(imax);
-          Extrema_POnSurf MinPS = myPOnS(imin);
-          Extrema_POnSurf MaxPS = myPOnS(imax);
-
-          mySqDist.Remove(SolNumber + 1, mySqDist.Length());
-          myPOnC.Remove(SolNumber + 1, myPOnC.Length());
-          myPOnS.Remove(SolNumber + 1, myPOnS.Length());
-
-          mySqDist.Append(MinSqDist);
-          myPOnC.Append(MinPC);
-          myPOnS.Append(MinPS);
-          mySqDist.Append(MaxSqDist);
-          myPOnC.Append(MaxPC);
-          myPOnS.Append(MaxPS);
         }
       }
       return;
     }
-    break;
   }
 
-  myDone = myExtElCS.IsDone();
+  // Elementary extrema is not done, try generic solution
+  Extrema_GenExtCS Ext;
+  Ext.Initialize(*myS, NbU, NbV, mytolS);
+  if (myCtype == GeomAbs_Hyperbola) {
+    Standard_Real tmin = Max(-20., C.FirstParameter());
+    Standard_Real tmax = Min(20., C.LastParameter());
+    Ext.Perform(C, NbT, tmin, tmax, mytolC); // to avoid overflow
+  }
+  else {
+    if ((myCtype == GeomAbs_Circle       && NbT < 13) ||
+      (myCtype == GeomAbs_BSplineCurve && NbT < 13))
+    {
+      NbT = 13;
+    }
+    Ext.Perform(C, NbT, mytolC);
+  }
+
+  myDone = Ext.IsDone();
   if (myDone) {
-    myIsPar = myExtElCS.IsParallel();
-    if (myIsPar) {
-      mySqDist.Append(myExtElCS.SquareDistance(1));
+    Standard_Integer NbExt = Ext.NbExt();
+    Standard_Real T, U, V;
+    Extrema_POnCurv PC;
+    Extrema_POnSurf PS;
+    for (i = 1; i <= NbExt; i++) {
+      PC = Ext.PointOnCurve(i);
+      PS = Ext.PointOnSurface(i);
+      T = PC.Parameter();
+      PS.Parameter(U, V);
+      AddSolution(C, T, U, V, PC.Value(), PS.Value(), Ext.SquareDistance(i));
     }
-    else {
-      Standard_Integer NbExt = myExtElCS.NbExt();
-      Standard_Real U, V;
-      for (i = 1; i <= NbExt; i++) {
-        Extrema_POnCurv PC;
-        Extrema_POnSurf PS;
-        myExtElCS.Points(i, PC, PS);
-        Standard_Real Ucurve = PC.Parameter();
-        PS.Parameter(U, V);
-        AddSolution(C, Ucurve, U, V, PC.Value(), PS.Value(), myExtElCS.SquareDistance(i));
-      }
-      if(mySqDist.Length() == 0 && NbExt > 0)
+
+    //Add sharp points
+    Standard_Integer SolNumber = mySqDist.Length();
+    Standard_Address CopyC = (Standard_Address)&C;
+    Adaptor3d_Curve& aC = *(Adaptor3d_Curve*)CopyC;
+    Standard_Integer NbIntervals = aC.NbIntervals(GeomAbs_C1);
+    TColStd_Array1OfReal SharpPoints(1, NbIntervals + 1);
+    aC.Intervals(SharpPoints, GeomAbs_C1);
+
+    Extrema_ExtPS aProjPS;
+    aProjPS.Initialize(*myS,
+      myS->FirstUParameter(),
+      myS->LastUParameter(),
+      myS->FirstVParameter(),
+      myS->LastVParameter(),
+      mytolS,
+      mytolS);
+
+    for (i = 2; i < SharpPoints.Upper(); ++i)
+    {
+      T = SharpPoints(i);
+      gp_Pnt aPnt = C.Value(T);
+      aProjPS.Perform(aPnt);
+      if (!aProjPS.IsDone())
+        continue;
+      Standard_Integer NbProj = aProjPS.NbExt(), jmin = 0;
+      Standard_Real MinSqDist = RealLast();
+      for (j = 1; j <= NbProj; j++)
       {
-        //Analytical extremas seem to be out of curve/surface boundaries.
-        //For plane it is possible to add extremity points of curve
-        if(myStype == GeomAbs_Plane)
+        Standard_Real aSqDist = aProjPS.SquareDistance(j);
+        if (aSqDist < MinSqDist)
         {
-          gp_Pln aPln = myS->Plane();
-          gp_Pnt PC, PP;
-          if(!Precision::IsInfinite(myucinf))
-          {
-            PC = C.Value(myucinf);
-            ElSLib::PlaneParameters(aPln.Position(), PC, U, V);
-            PP = ElSLib::PlaneValue(U, V, aPln.Position());
-            AddSolution(C, myucinf, U, V, PC, PP, PC.SquareDistance(PP));
-          }
-          if(!Precision::IsInfinite(myucsup))
-          {
-            PC = C.Value(myucsup);
-            ElSLib::PlaneParameters(aPln.Position(), PC, U, V);
-            PP = ElSLib::PlaneValue(U, V, aPln.Position());
-            AddSolution(C, myucsup, U, V, PC, PP, PC.SquareDistance(PP));
-          }
+          MinSqDist = aSqDist;
+          jmin = j;
         }
       }
+      if (jmin != 0)
+      {
+        aProjPS.Point(jmin).Parameter(U, V);
+        AddSolution(C, T, U, V,
+          aPnt, aProjPS.Point(jmin).Value(), MinSqDist);
+      }
+    }
+    //Cut sharp solutions to keep only minimum and maximum
+    Standard_Integer imin = SolNumber + 1, imax = mySqDist.Length();
+    for (i = SolNumber + 1; i <= mySqDist.Length(); i++)
+    {
+      if (mySqDist(i) < mySqDist(imin))
+        imin = i;
+      if (mySqDist(i) > mySqDist(imax))
+        imax = i;
+    }
+    if (mySqDist.Length() > SolNumber + 2)
+    {
+      Standard_Real MinSqDist = mySqDist(imin);
+      Standard_Real MaxSqDist = mySqDist(imax);
+      Extrema_POnCurv MinPC = myPOnC(imin);
+      Extrema_POnCurv MaxPC = myPOnC(imax);
+      Extrema_POnSurf MinPS = myPOnS(imin);
+      Extrema_POnSurf MaxPS = myPOnS(imax);
+
+      mySqDist.Remove(SolNumber + 1, mySqDist.Length());
+      myPOnC.Remove(SolNumber + 1, myPOnC.Length());
+      myPOnS.Remove(SolNumber + 1, myPOnS.Length());
+
+      mySqDist.Append(MinSqDist);
+      myPOnC.Append(MinPC);
+      myPOnS.Append(MinPS);
+      mySqDist.Append(MaxSqDist);
+      myPOnC.Append(MaxPC);
+      myPOnS.Append(MaxPS);
     }
   }
-
 }
 
 
@@ -400,16 +489,16 @@ Standard_Boolean Extrema_ExtCS::IsParallel() const
 
 Standard_Real Extrema_ExtCS::SquareDistance(const Standard_Integer N) const
 {
-  if(!myDone) StdFail_NotDone::Raise();
-  if (myIsPar && N != 1) StdFail_InfiniteSolutions::Raise();
-  if ((N < 1) || (N > mySqDist.Length())) Standard_OutOfRange::Raise();
+  if(!myDone) throw StdFail_NotDone();
+  if (myIsPar && N != 1) throw StdFail_InfiniteSolutions();
+  if ((N < 1) || (N > mySqDist.Length())) throw Standard_OutOfRange();
   return mySqDist.Value(N);
 }
 
 
 Standard_Integer Extrema_ExtCS::NbExt() const
 {
-  if(!myDone) StdFail_NotDone::Raise();
+  if(!myDone) throw StdFail_NotDone();
   return myPOnC.Length();
 }
 
@@ -419,7 +508,7 @@ void Extrema_ExtCS::Points(const Standard_Integer N,
   Extrema_POnCurv&       P1,
   Extrema_POnSurf&       P2) const
 {
-  if(!myDone) StdFail_NotDone::Raise();
+  if(!myDone) throw StdFail_NotDone();
   P1 = myPOnC.Value(N);
   P2 = myPOnS.Value(N);
 }