0032882: Modeling Data - Extrema curve/curve cannot find all solutions (OCCT 7.2...

Extrema/Extrema_GenExtCC.gxx - estimation of Lipchitz constant is improved
Extrema_GlobOptFuncCC.cxx - function value is changed

LocOpe/LocOpe_WiresOnShape.cxx - small correction to fix regression

lowalgos/extcc/bug32882 - new test case is added

some test were updated according new behavior of extrema algo

diff --git a/src/Extrema/Extrema_GenExtCC.gxx b/src/Extrema/Extrema_GenExtCC.gxx
index 10e497d2e28d474d0735781c585bd1de5b059bea..86719def542764f8c79b196f4592bbd55b68e127 100644 (file)
--- a/src/Extrema/Extrema_GenExtCC.gxx
+++ b/src/Extrema/Extrema_GenExtCC.gxx
@@ -26,6 +26,7 @@
 #include <Precision.hxx>
 #include <NCollection_Vector.hxx>
 #include <NCollection_CellFilter.hxx>
+#include <GCPnts_AbscissaPoint.hxx>
 
 // Comparator, used in std::sort.
 static Standard_Boolean comp(const gp_XY& theA,
@@ -47,6 +48,57 @@ static Standard_Boolean comp(const gp_XY& theA,
   return Standard_False;
 }
 
+static void ChangeIntervals(Handle(TColStd_HArray1OfReal)& theInts, const Standard_Integer theNbInts)
+{
+  Standard_Integer aNbInts = theInts->Length() - 1;
+  Standard_Integer aNbAdd = theNbInts - aNbInts;
+  Handle(TColStd_HArray1OfReal) aNewInts = new TColStd_HArray1OfReal(1, theNbInts + 1);
+  Standard_Integer aNbLast = theInts->Length();
+  Standard_Integer i;
+  if (aNbInts == 1)
+  {
+    aNewInts->SetValue(1, theInts->First());
+    aNewInts->SetValue(theNbInts + 1, theInts->Last());
+    Standard_Real dt = (theInts->Last() - theInts->First()) / theNbInts;
+    Standard_Real t = theInts->First() + dt;
+    for (i = 2; i <= theNbInts; ++i, t += dt)
+    {
+      aNewInts->SetValue(i, t);
+    }
+    theInts = aNewInts;
+    return;
+  }
+  //
+  for (i = 1; i <= aNbLast; ++i)
+  {
+    aNewInts->SetValue(i, theInts->Value(i));
+  }
+  //
+  while (aNbAdd > 0)
+  {
+    Standard_Real anLIntMax = -1.;
+    Standard_Integer aMaxInd = -1;
+    for (i = 1; i < aNbLast; ++i)
+    {
+      Standard_Real anL = aNewInts->Value(i + 1) - aNewInts->Value(i);
+      if (anL > anLIntMax)
+      {
+        anLIntMax = anL;
+        aMaxInd = i;
+      }
+    }
+
+    Standard_Real t = (aNewInts->Value(aMaxInd + 1) + aNewInts->Value(aMaxInd)) / 2.;
+    for (i = aNbLast; i > aMaxInd; --i)
+    {
+      aNewInts->SetValue(i + 1, aNewInts->Value(i));
+    }
+    aNbLast++;
+    aNbAdd--;
+    aNewInts->SetValue(aMaxInd + 1, t);
+  }
+  theInts = aNewInts;
+}
 class Extrema_CCPointsInspector : public NCollection_CellFilter_InspectorXY
 {
 public:
@@ -111,7 +163,6 @@ static Standard_Real ProjPOnC(const Pnt& theP,
       if (aD < aDist)
         aDist = aD;
     }
-    aDist = sqrt(aDist);
   }
   return aDist;
 }
@@ -229,14 +280,79 @@ void Extrema_GenExtCC::Perform()
     aNbInter[1] = C2.NbIntervals(aContinuity);
   }
 
-  TColStd_Array1OfReal anIntervals1(1, aNbInter[0] + 1);
-  TColStd_Array1OfReal anIntervals2(1, aNbInter[1] + 1);
-  C1.Intervals(anIntervals1, aContinuity);
-  C2.Intervals(anIntervals2, aContinuity);
+  Standard_Real anL[2];
+  Standard_Integer indmax = -1, indmin = -1;
+  const Standard_Real mult = 20.;
+  if (!(Precision::IsInfinite(C1.FirstParameter()) || Precision::IsInfinite(C1.LastParameter()) ||
+    Precision::IsInfinite(C2.FirstParameter()) || Precision::IsInfinite(C2.LastParameter())))
+  {
+    anL[0] = GCPnts_AbscissaPoint::Length(C1);
+    anL[1] = GCPnts_AbscissaPoint::Length(C2);
+    if (anL[0] / aNbInter[0] > mult * anL[1] / aNbInter[1])
+    {
+      indmax = 0;
+      indmin = 1;
+    }
+    else if (anL[1] / aNbInter[1] > mult * anL[0] / aNbInter[0])
+    {
+      indmax = 1;
+      indmin = 0;
+    }
+  }
+  Standard_Integer aNbIntOpt = 0;
+  if (indmax >= 0)
+  {
+    aNbIntOpt = RealToInt(anL[indmax] * aNbInter[indmin] / anL[indmin] / (mult / 4.)) + 1;
+    if (aNbIntOpt > 100 || aNbIntOpt < aNbInter[indmax])
+    {
+      indmax = -1;
+    }
+    else
+    {
+      if (aNbIntOpt * aNbInter[indmin] > 100)
+      {
+        aNbIntOpt = 100 / aNbInter[indmin];
+        if (aNbIntOpt < aNbInter[indmax])
+        {
+          indmax = -1;
+        }
+      }
+    }
+  }
+
+  Handle(TColStd_HArray1OfReal) anIntervals1 = new TColStd_HArray1OfReal(1, aNbInter[0] + 1);
+  Handle(TColStd_HArray1OfReal) anIntervals2 = new TColStd_HArray1OfReal(1, aNbInter[1] + 1);
+  C1.Intervals(anIntervals1->ChangeArray1(), aContinuity);
+  C2.Intervals(anIntervals2->ChangeArray1(), aContinuity);
+  if (indmax >= 0)
+  {
+    if (indmax == 0)
+    {
+      //Change anIntervals1
+      ChangeIntervals(anIntervals1, aNbIntOpt);
+      aNbInter[0] = anIntervals1->Length() - 1;
+    }
+    else
+    {
+      //Change anIntervals2;
+      ChangeIntervals(anIntervals2, aNbIntOpt);
+      aNbInter[1] = anIntervals2->Length() - 1;
+    }
+  }
+  if (C1.IsClosed() && aNbInter[0] == 1)
+  {
+    ChangeIntervals(anIntervals1, 3);
+    aNbInter[0] = anIntervals1->Length() - 1;
+  }
+  if (C2.IsClosed() && aNbInter[1] == 1)
+  {
+    ChangeIntervals(anIntervals2, 3);
+    aNbInter[1] = anIntervals2->Length() - 1;
+  }
 
   // Lipchitz constant computation.
   const Standard_Real aMaxLC = 10000.;
-  Standard_Real aLC = 9.0; // Default value.
+  Standard_Real aLC = 100.0; // Default value.
   const Standard_Real aMaxDer1 = 1.0 / C1.Resolution(1.0);
   const Standard_Real aMaxDer2 = 1.0 / C2.Resolution(1.0);
   Standard_Real aMaxDer = Max(aMaxDer1, aMaxDer2) * Sqrt(2.0);
@@ -276,6 +392,43 @@ void Extrema_GenExtCC::Perform()
   }
 
   Extrema_GlobOptFuncCCC2 aFunc (C1, C2);
+  if (aLC < aMaxLC || aMaxDer > aMaxLC)
+  {
+    //Estimation of Lipschitz constant by gradient of optimization function 
+    //using sampling in parameter space.
+    math_Vector aT(1, 2), aG(1, 2);
+    Standard_Real aF, aMaxG = 0.;
+    Standard_Real t1, t2, dt1, dt2;
+    Standard_Integer n1 = 21, n2 = 21, i1, i2;
+    dt1 = (C1.LastParameter() - C1.FirstParameter()) / (n1 - 1);
+    dt2 = (C2.LastParameter() - C2.FirstParameter()) / (n2 - 1);
+    for (i1 = 1, t1 = C1.FirstParameter(); i1 <= n1; ++i1, t1 += dt1)
+    {
+      aT(1) = t1;
+      for (i2 = 1, t2 = C2.FirstParameter(); i2 <= n2; ++i2, t2 += dt2)
+      {
+        aT(2) = t2;
+        aFunc.Values(aT, aF, aG);
+        Standard_Real aMod = aG(1)*aG(1) + aG(2)*aG(2);
+        aMaxG = Max(aMaxG, aMod);
+      }
+    }
+    aMaxG = Sqrt(aMaxG);
+    if (aMaxG > aMaxDer)
+    {
+      aLC = Min(aMaxG, aMaxLC);
+      isConstLockedFlag = Standard_True;
+    }
+    if (aMaxG > 100. * aMaxLC)
+    {
+      aLC = 100. * aMaxLC;
+      isConstLockedFlag = Standard_True;
+    }
+    else if (aMaxG < 0.1 * aMaxDer)
+    {
+      isConstLockedFlag = Standard_True;
+    }
+  }
   math_GlobOptMin aFinder(&aFunc, myLowBorder, myUppBorder, aLC);
   aFinder.SetLipConstState(isConstLockedFlag);
   aFinder.SetContinuity(aContinuity == GeomAbs_C2 ? 2 : 1);
@@ -286,8 +439,8 @@ void Extrema_GenExtCC::Perform()
   aFinder.SetFunctionalMinimalValue(0.0); // Best distance cannot be lower than 0.0.
 
   // Size computed to have cell index inside of int32 value.
-  const Standard_Real aCellSize = Max(Max(anIntervals1.Last() - anIntervals1.First(),
-                                      anIntervals2.Last() - anIntervals2.First())
+  const Standard_Real aCellSize = Max(Max(anIntervals1->Last() - anIntervals1->First(),
+                                      anIntervals2->Last() - anIntervals2->First())
                                       * Precision::PConfusion() / (2.0 * Sqrt(2.0)),
                                       Precision::PConfusion());
   Extrema_CCPointsInspector anInspector(aCellSize);
@@ -303,10 +456,10 @@ void Extrema_GenExtCC::Perform()
   {
     for(j = 1; j <= aNbInter[1]; j++)
     {
-      aFirstBorderInterval(1) = anIntervals1(i);
-      aFirstBorderInterval(2) = anIntervals2(j); 
-      aSecondBorderInterval(1) = anIntervals1(i + 1);
-      aSecondBorderInterval(2) = anIntervals2(j + 1);
+      aFirstBorderInterval(1) = anIntervals1->Value(i);
+      aFirstBorderInterval(2) = anIntervals2->Value(j);
+      aSecondBorderInterval(1) = anIntervals1->Value(i + 1);
+      aSecondBorderInterval(2) = anIntervals2->Value(j + 1);
 
       aFinder.SetLocalParams(aFirstBorderInterval, aSecondBorderInterval);
       aFinder.Perform(GetSingleSolutionFlag());
@@ -515,6 +668,7 @@ Standard_Boolean Extrema_GenExtCC::IsDone() const
 //=======================================================================
 Standard_Boolean Extrema_GenExtCC::IsParallel() const 
 {
+  if (!IsDone()) throw StdFail_NotDone();
   return myParallel; 
 }
 
@@ -524,8 +678,7 @@ Standard_Boolean Extrema_GenExtCC::IsParallel() const
 //=======================================================================
 Standard_Integer Extrema_GenExtCC::NbExt() const
 {
-  StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::NbExt()")
-
+  if (!IsDone()) throw StdFail_NotDone();
   return myPoints1.Length();
 }
 
@@ -535,8 +688,10 @@ Standard_Integer Extrema_GenExtCC::NbExt() const
 //=======================================================================
 Standard_Real Extrema_GenExtCC::SquareDistance(const Standard_Integer N) const
 {
-  StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::SquareDistance()")
-  Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::SquareDistance()")
+  if (N < 1 || N > NbExt())
+  {
+    throw Standard_OutOfRange();
+  }
 
   return Tool1::Value(*((Curve1*)myC[0]), myPoints1(N)).SquareDistance(Tool2::Value(*((Curve2*)myC[1]), myPoints2(N)));
 }
@@ -549,8 +704,10 @@ void Extrema_GenExtCC::Points(const Standard_Integer N,
                               POnC& P1,
                               POnC& P2) const
 {
-  StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::Points()")
-  Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::Points()")
+  if (N < 1 || N > NbExt())
+  {
+    throw Standard_OutOfRange();
+  }
 
   P1.SetValues(myPoints1(N), Tool1::Value(*((Curve1*)myC[0]), myPoints1(N)));
   P2.SetValues(myPoints2(N), Tool2::Value(*((Curve2*)myC[1]), myPoints2(N)));

diff --git a/src/Extrema/Extrema_GlobOptFuncCC.cxx b/src/Extrema/Extrema_GlobOptFuncCC.cxx
index 5246d3c05a6369e8ff76a9e854d45ae70d2b9ac5..1c0f1f07cc7767976d6c34c079b1866649999a63 100644 (file)
--- a/src/Extrema/Extrema_GlobOptFuncCC.cxx
+++ b/src/Extrema/Extrema_GlobOptFuncCC.cxx
@@ -45,7 +45,7 @@ static Standard_Boolean _Value(const Adaptor3d_Curve& C1,
     return Standard_False;
   }
 
-  F = C2.Value(v).Distance(C1.Value(u));
+  F = C2.Value(v).SquareDistance(C1.Value(u));
   return Standard_True;
 }
 
@@ -66,7 +66,7 @@ static Standard_Boolean _Value(const Adaptor2d_Curve2d& C1,
     return Standard_False;
   }
 
-  F = C2.Value(v).Distance(C1.Value(u));
+  F = C2.Value(v).SquareDistance(C1.Value(u));
   return Standard_True;
 }
 
@@ -91,13 +91,14 @@ static Standard_Boolean _Gradient(const Adaptor3d_Curve& C1,
 
   C1.D1(X(1), C1D0, C1D1);
   C2.D1(X(2), C2D0, C2D1);
-
+  
   G(1) = - (C2D0.X() - C1D0.X()) * C1D1.X() 
          - (C2D0.Y() - C1D0.Y()) * C1D1.Y() 
          - (C2D0.Z() - C1D0.Z()) * C1D1.Z();
   G(2) =   (C2D0.X() - C1D0.X()) * C2D1.X() 
          + (C2D0.Y() - C1D0.Y()) * C2D1.Y() 
          + (C2D0.Z() - C1D0.Z()) * C2D1.Z();
+  G *= 2.;
   return Standard_True;
 }
 
@@ -123,8 +124,11 @@ static Standard_Boolean _Gradient(const Adaptor2d_Curve2d& C1,
 
   G(1) = - (C2D0.X() - C1D0.X()) * C1D1.X() 
          - (C2D0.Y() - C1D0.Y()) * C1D1.Y();
+
   G(2) =   (C2D0.X() - C1D0.X()) * C2D1.X() 
          + (C2D0.Y() - C1D0.Y()) * C2D1.Y();
+  G *= 2.;
+
   return Standard_True;
 }
 
@@ -168,6 +172,7 @@ static Standard_Boolean _Hessian (const Adaptor3d_Curve& C1,
            + (C2D0.X() - C1D0.X()) * C2D2.X() 
            + (C2D0.Y() - C1D0.Y()) * C2D2.Y() 
            + (C2D0.Z() - C1D0.Z()) * C2D2.Z();
+  H *= 2.;
   return Standard_True;
 }
 
@@ -206,10 +211,11 @@ static Standard_Boolean _Hessian (const Adaptor2d_Curve2d& C1,
            + C2D1.Y() * C2D1.Y() 
            + (C2D0.X() - C1D0.X()) * C2D2.X() 
            + (C2D0.Y() - C1D0.Y()) * C2D2.Y();
+  H *= 2.;
   return Standard_True;
 }
 
-// C0
+//C0
 
 //=======================================================================
 //function : Extrema_GlobOptFuncCCC0
@@ -405,6 +411,5 @@ Standard_Boolean Extrema_GlobOptFuncCCC2::Values(const math_Vector& X,Standard_R
   else
     isHessianComputed = _Hessian(*myC1_2d, *myC2_2d, X, H);
 
-
   return (Value(X, F) && Gradient(X, G) && isHessianComputed);
 }

diff --git a/src/IntTools/IntTools_FaceFace.cxx b/src/IntTools/IntTools_FaceFace.cxx
index 8726f0ad567ee237f99ac07a488b33cbd57ff20b..6527e5715f25b4721b9fe8438409eb7b94ca3d2b 100644 (file)
--- a/src/IntTools/IntTools_FaceFace.cxx
+++ b/src/IntTools/IntTools_FaceFace.cxx
@@ -1188,7 +1188,7 @@ void IntTools_FaceFace::MakeCurve(const Standard_Integer Index,
                              bAvoidLineConstructor, 
                              myTol,
                              aSeqOfL, 
-                             aReachedTol,
+                             aReachedTol, // obsolete parameter
                              myContext);
       //
       aNbSeqOfL=aSeqOfL.Length();
@@ -1196,7 +1196,7 @@ void IntTools_FaceFace::MakeCurve(const Standard_Integer Index,
       Standard_Real aTolC = 0.;
       if (bIsDecomposited) {
         nbiter=aNbSeqOfL;
-        aTolC = aReachedTol;
+        aTolC = Precision::Confusion();
       }
       else {
         nbiter=1;

diff --git a/src/IntTools/IntTools_WLineTool.cxx b/src/IntTools/IntTools_WLineTool.cxx
index f36c78833bff2c4cf07fdb9d86f11d04615a5339..b9a964f220c51a7e9851ca054ac8eb9a1ac1639e 100644 (file)
--- a/src/IntTools/IntTools_WLineTool.cxx
+++ b/src/IntTools/IntTools_WLineTool.cxx
@@ -229,203 +229,6 @@ Standard_Boolean IntTools_WLineTool::NotUseSurfacesForApprox(const TopoDS_Face&
 
 /////////////////////// DecompositionOfWLine ////////////////////////////
 
-//=======================================================================
-//function : CheckTangentZonesExist
-//purpose  : static subfunction in ComputeTangentZones
-//=======================================================================
-static
-Standard_Boolean CheckTangentZonesExist(const Handle(GeomAdaptor_HSurface)& theSurface1,
-                                        const Handle(GeomAdaptor_HSurface)& theSurface2)
-{
-  if ( ( theSurface1->GetType() != GeomAbs_Torus ) ||
-       ( theSurface2->GetType() != GeomAbs_Torus ) )
-    return Standard_False;
-
-  gp_Torus aTor1 = theSurface1->Torus();
-  gp_Torus aTor2 = theSurface2->Torus();
-
-  if ( aTor1.Location().Distance( aTor2.Location() ) > Precision::Confusion() )
-    return Standard_False;
-
-  if ( ( fabs( aTor1.MajorRadius() - aTor2.MajorRadius() ) > Precision::Confusion() ) ||
-       ( fabs( aTor1.MinorRadius() - aTor2.MinorRadius() ) > Precision::Confusion() ) )
-    return Standard_False;
-
-  if ( ( aTor1.MajorRadius() < aTor1.MinorRadius() ) ||
-       ( aTor2.MajorRadius() < aTor2.MinorRadius() ) )
-    return Standard_False;
-
-  return Standard_True;
-}
-
-
-//=======================================================================
-//function : ComputeTangentZones
-//purpose  : static subfunction in DecompositionOfWLine
-//=======================================================================
-static
-Standard_Integer ComputeTangentZones( const Handle(GeomAdaptor_HSurface)& theSurface1,
-                                     const Handle(GeomAdaptor_HSurface)&  theSurface2,
-                                     const TopoDS_Face&                   theFace1,
-                                     const TopoDS_Face&                   theFace2,
-                                     Handle(TColgp_HArray1OfPnt2d)&       theResultOnS1,
-                                     Handle(TColgp_HArray1OfPnt2d)&       theResultOnS2,
-                                     Handle(TColStd_HArray1OfReal)&       theResultRadius,
-                                     const Handle(IntTools_Context)& aContext)
-{
-  Standard_Integer aResult = 0;
-  if ( !CheckTangentZonesExist( theSurface1, theSurface2 ) )
-    return aResult;
-
-
-  TColgp_SequenceOfPnt2d aSeqResultS1, aSeqResultS2;
-  TColStd_SequenceOfReal aSeqResultRad;
-
-  gp_Torus aTor1 = theSurface1->Torus();
-  gp_Torus aTor2 = theSurface2->Torus();
-
-  gp_Ax2 anax1( aTor1.Location(), aTor1.Axis().Direction() );
-  gp_Ax2 anax2( aTor2.Location(), aTor2.Axis().Direction() );
-  Standard_Integer j = 0;
-
-  for ( j = 0; j < 2; j++ ) {
-    Standard_Real aCoef = ( j == 0 ) ? -1 : 1;
-    Standard_Real aRadius1 = fabs(aTor1.MajorRadius() + aCoef * aTor1.MinorRadius());
-    Standard_Real aRadius2 = fabs(aTor2.MajorRadius() + aCoef * aTor2.MinorRadius());
-
-    gp_Circ aCircle1( anax1, aRadius1 );
-    gp_Circ aCircle2( anax2, aRadius2 );
-
-    // roughly compute radius of tangent zone for perpendicular case
-    Standard_Real aCriteria = Precision::Confusion() * 0.5;
-
-    Standard_Real aT1 = aCriteria;
-    Standard_Real aT2 = aCriteria;
-    if ( j == 0 ) {
-      // internal tangency
-      Standard_Real aR = ( aRadius1 > aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
-      //aT1 = aCriteria * aCriteria + aR * aR - ( aR - aCriteria ) * ( aR - aCriteria );
-      aT1 = 2. * aR * aCriteria;
-      aT2 = aT1;
-    }
-    else {
-      // external tangency
-      Standard_Real aRb = ( aRadius1 > aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
-      Standard_Real aRm = ( aRadius1 < aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
-      Standard_Real aDelta = aRb - aCriteria;
-      aDelta *= aDelta;
-      aDelta -= aRm * aRm;
-      aDelta /= 2. * (aRb - aRm);
-      aDelta -= 0.5 * (aRb - aRm);
-      
-      aT1 = 2. * aRm * (aRm - aDelta);
-      aT2 = aT1;
-    }
-    aCriteria = ( aT1 > aT2) ? aT1 : aT2;
-    if ( aCriteria > 0 )
-      aCriteria = sqrt( aCriteria );
-
-    if ( aCriteria > 0.5 * aTor1.MinorRadius() ) {
-      // too big zone -> drop to minimum
-      aCriteria = Precision::Confusion();
-    }
-
-    GeomAdaptor_Curve aC1( new Geom_Circle(aCircle1) );
-    GeomAdaptor_Curve aC2( new Geom_Circle(aCircle2) );
-    Extrema_ExtCC anExtrema(aC1, aC2, 0, 2. * M_PI, 0, 2. * M_PI, 
-                            Precision::PConfusion(), Precision::PConfusion());
-            
-    if ( anExtrema.IsDone() ) {
-
-      Standard_Integer i = 0;
-      for ( i = 1; i <= anExtrema.NbExt(); i++ ) {
-        if ( anExtrema.SquareDistance(i) > aCriteria * aCriteria )
-          continue;
-
-        Extrema_POnCurv P1, P2;
-        anExtrema.Points( i, P1, P2 );
-
-        Standard_Boolean bFoundResult = Standard_True;
-        gp_Pnt2d pr1, pr2;
-
-        Standard_Integer surfit = 0;
-        for ( surfit = 0; surfit < 2; surfit++ ) {
-          GeomAPI_ProjectPointOnSurf& aProjector = 
-            (surfit == 0) ? aContext->ProjPS(theFace1) : aContext->ProjPS(theFace2);
-
-          gp_Pnt aP3d = (surfit == 0) ? P1.Value() : P2.Value();
-          aProjector.Perform(aP3d);
-
-          if(!aProjector.IsDone())
-            bFoundResult = Standard_False;
-          else {
-            if(aProjector.LowerDistance() > aCriteria) {
-              bFoundResult = Standard_False;
-            }
-            else {
-              Standard_Real foundU = 0, foundV = 0;
-              aProjector.LowerDistanceParameters(foundU, foundV);
-              if ( surfit == 0 )
-                pr1 = gp_Pnt2d( foundU, foundV );
-              else
-                pr2 = gp_Pnt2d( foundU, foundV );
-            }
-          }
-        }
-        if ( bFoundResult ) {
-          aSeqResultS1.Append( pr1 );
-          aSeqResultS2.Append( pr2 );
-          aSeqResultRad.Append( aCriteria );
-
-          // torus is u and v periodic
-          const Standard_Real twoPI = M_PI + M_PI;
-          Standard_Real arr1tmp[2] = {pr1.X(), pr1.Y()};
-          Standard_Real arr2tmp[2] = {pr2.X(), pr2.Y()};
-
-          // iteration on period bounds
-          for ( Standard_Integer k1 = 0; k1 < 2; k1++ ) {
-            Standard_Real aBound = ( k1 == 0 ) ? 0 : twoPI;
-            Standard_Real aShift = ( k1 == 0 ) ? twoPI : -twoPI;
-
-            // iteration on surfaces
-            for ( Standard_Integer k2 = 0; k2 < 2; k2++ ) {
-              Standard_Real* arr1 = ( k2 == 0 ) ? arr1tmp : arr2tmp;
-              Standard_Real* arr2 = ( k2 != 0 ) ? arr1tmp : arr2tmp;
-              TColgp_SequenceOfPnt2d& aSeqS1 = ( k2 == 0 ) ? aSeqResultS1 : aSeqResultS2; 
-              TColgp_SequenceOfPnt2d& aSeqS2 = ( k2 != 0 ) ? aSeqResultS1 : aSeqResultS2; 
-
-              if (fabs(arr1[0] - aBound) < Precision::PConfusion()) {
-                aSeqS1.Append( gp_Pnt2d( arr1[0] + aShift, arr1[1] ) );
-                aSeqS2.Append( gp_Pnt2d( arr2[0], arr2[1] ) );
-                aSeqResultRad.Append( aCriteria );
-              }
-              if (fabs(arr1[1] - aBound) < Precision::PConfusion()) {
-                aSeqS1.Append( gp_Pnt2d( arr1[0], arr1[1] + aShift) );
-                aSeqS2.Append( gp_Pnt2d( arr2[0], arr2[1] ) );
-                aSeqResultRad.Append( aCriteria );
-              }
-            }
-          } //
-        }
-      }
-    }
-  }
-  aResult = aSeqResultRad.Length();
-
-  if ( aResult > 0 ) {
-    theResultOnS1 = new TColgp_HArray1OfPnt2d( 1, aResult );
-    theResultOnS2 = new TColgp_HArray1OfPnt2d( 1, aResult );
-    theResultRadius = new TColStd_HArray1OfReal( 1, aResult );
-
-    for ( Standard_Integer i = 1 ; i <= aResult; i++ ) {
-      theResultOnS1->SetValue( i, aSeqResultS1.Value(i) );
-      theResultOnS2->SetValue( i, aSeqResultS2.Value(i) );
-      theResultRadius->SetValue( i, aSeqResultRad.Value(i) );
-    }
-  }
-  return aResult;
-}
-
 //=======================================================================
 //function : IsPointOnBoundary
 //purpose  : static subfunction in DecompositionOfWLine
@@ -457,27 +260,6 @@ Standard_Boolean IsPointOnBoundary(const Standard_Real theParameter,
   return bRet;
 }
 
-//=======================================================================
-//function : IsInsideTanZone
-//purpose  : Check if point is inside a radial tangent zone.
-//           static subfunction in DecompositionOfWLine and FindPoint
-//=======================================================================
-static
-Standard_Boolean IsInsideTanZone(const gp_Pnt2d&     thePoint,
-                                 const gp_Pnt2d&     theTanZoneCenter,
-                                 const Standard_Real theZoneRadius,
-                                 Handle(GeomAdaptor_HSurface) theGASurface)
-{
-  Standard_Real aUResolution = theGASurface->UResolution( theZoneRadius );
-  Standard_Real aVResolution = theGASurface->VResolution( theZoneRadius );
-  Standard_Real aRadiusSQR = ( aUResolution < aVResolution ) ? aUResolution : aVResolution;
-  aRadiusSQR *= aRadiusSQR;
-  if ( thePoint.SquareDistance( theTanZoneCenter ) <= aRadiusSQR )
-    return Standard_True;
-
-  return Standard_False;
-}
-
 //=======================================================================
 //function : AdjustByNeighbour
 //purpose  : static subfunction in DecompositionOfWLine
@@ -652,72 +434,6 @@ Standard_Boolean FindPoint(const gp_Pnt2d&     theFirstPoint,
   return Standard_False;
 }
 
-//=======================================================================
-//function : FindPoint
-//purpose  : Find point on the boundary of radial tangent zone
-//           static subfunction in DecompositionOfWLine
-//=======================================================================
-static
-Standard_Boolean FindPoint(const gp_Pnt2d&     theFirstPoint,
-                           const gp_Pnt2d&     theLastPoint,
-                           const Standard_Real theUmin, 
-                           const Standard_Real theUmax,
-                           const Standard_Real theVmin,
-                           const Standard_Real theVmax,
-                           const gp_Pnt2d&     theTanZoneCenter,
-                           const Standard_Real theZoneRadius,
-                           Handle(GeomAdaptor_HSurface) theGASurface,
-                           gp_Pnt2d&           theNewPoint) {
-  theNewPoint = theLastPoint;
-
-  if ( !IsInsideTanZone( theLastPoint, theTanZoneCenter, theZoneRadius, theGASurface) )
-    return Standard_False;
-
-  Standard_Real aUResolution = theGASurface->UResolution( theZoneRadius );
-  Standard_Real aVResolution = theGASurface->VResolution( theZoneRadius );
-
-  Standard_Real aRadius = ( aUResolution < aVResolution ) ? aUResolution : aVResolution;
-  gp_Ax22d anAxis( theTanZoneCenter, gp_Dir2d(1, 0), gp_Dir2d(0, 1) );
-  gp_Circ2d aCircle( anAxis, aRadius );
-  
-  //
-  gp_Vec2d aDir( theLastPoint.XY() - theFirstPoint.XY() );
-  Standard_Real aLength = aDir.Magnitude();
-  if ( aLength <= gp::Resolution() )
-    return Standard_False;
-  gp_Lin2d aLine( theFirstPoint, aDir );
-
-  //
-  Handle(Geom2d_Line) aCLine = new Geom2d_Line( aLine );
-  Handle(Geom2d_TrimmedCurve) aC1 = new Geom2d_TrimmedCurve( aCLine, 0, aLength );
-  Handle(Geom2d_Circle) aC2 = new Geom2d_Circle( aCircle );
-
-  Standard_Real aTol = aRadius * 0.001;
-  aTol = ( aTol < Precision::PConfusion() ) ? Precision::PConfusion() : aTol;
-
-  Geom2dAPI_InterCurveCurve anIntersector;
-  anIntersector.Init( aC1, aC2, aTol );
-
-  if ( anIntersector.NbPoints() == 0 )
-    return Standard_False;
-
-  Standard_Boolean aFound = Standard_False;
-  Standard_Real aMinDist = aLength * aLength;
-  Standard_Integer i = 0;
-  for ( i = 1; i <= anIntersector.NbPoints(); i++ ) {
-    gp_Pnt2d aPInt = anIntersector.Point( i );
-    if ( aPInt.SquareDistance( theFirstPoint ) < aMinDist ) {
-      if ( ( aPInt.X() >= theUmin ) && ( aPInt.X() <= theUmax ) &&
-           ( aPInt.Y() >= theVmin ) && ( aPInt.Y() <= theVmax ) ) {
-        theNewPoint = aPInt;
-        aFound = Standard_True;
-      }
-    }
-  }
-
-  return aFound;
-}
-
 //=======================================================================
 //function : DecompositionOfWLine
 //purpose  : 
@@ -732,7 +448,7 @@ Standard_Boolean IntTools_WLineTool::
                        const Standard_Boolean                         theAvoidLConstructor,
                        const Standard_Real                            theTol,
                        IntPatch_SequenceOfLine&                       theNewLines,
-                       Standard_Real&                                 theReachedTol3d,
+                       Standard_Real&                                 /*theReachedTol3d*/,
                        const Handle(IntTools_Context)& aContext) 
 {
   Standard_Boolean bRet, bAvoidLineConstructor;
@@ -758,13 +474,7 @@ Standard_Boolean IntTools_WLineTool::
   TColStd_Array1OfListOfInteger anArrayOfLines(1, aNbPnts); 
   TColStd_Array1OfInteger       anArrayOfLineType(1, aNbPnts);
   TColStd_ListOfInteger aListOfPointIndex;
-  
-  Handle(TColgp_HArray1OfPnt2d) aTanZoneS1;
-  Handle(TColgp_HArray1OfPnt2d) aTanZoneS2;
-  Handle(TColStd_HArray1OfReal) aTanZoneRadius;
-  Standard_Integer aNbZone = ComputeTangentZones( theSurface1, theSurface2, theFace1, theFace2,
-                                                 aTanZoneS1, aTanZoneS2, aTanZoneRadius, aContext);
-  
+    
   //
   nblines=0;
   aTol=Precision::Confusion();
@@ -835,24 +545,6 @@ Standard_Boolean IntTools_WLineTool::
           bIsCurrentPointOnBoundary = Standard_True;
           break;
         }
-        else {
-          // check if a point belong to a tangent zone. Begin
-          Standard_Integer zIt = 0;
-          for ( zIt = 1; zIt <= aNbZone; zIt++ ) {
-            gp_Pnt2d aPZone = (i == 0) ? aTanZoneS1->Value(zIt) : aTanZoneS2->Value(zIt);
-            Standard_Real aZoneRadius = aTanZoneRadius->Value(zIt);
-
-            if ( IsInsideTanZone(gp_Pnt2d( U, V ), aPZone, aZoneRadius, aGASurface ) ) {
-              // set boundary flag to split the curve by a tangent zone
-              bIsPointOnBoundary = Standard_True;
-              bIsCurrentPointOnBoundary = Standard_True;
-              if ( theReachedTol3d < aZoneRadius ) {
-                theReachedTol3d = aZoneRadius;
-              }
-              break;
-            }
-          }
-        }
       }//for(j = 0; j < 2; j++) {
 
       if(bIsCurrentPointOnBoundary){
@@ -931,7 +623,7 @@ Standard_Boolean IntTools_WLineTool::
         Standard_Integer nbboundaries = 0;
 
         Standard_Boolean bIsNearBoundary = Standard_False;
-        Standard_Integer aZoneIndex = 0;
+        //Standard_Integer aZoneIndex = 0;
         Standard_Integer bIsUBoundary = Standard_False; // use if nbboundaries == 1
         Standard_Integer bIsFirstBoundary = Standard_False; // use if nbboundaries == 1
         
@@ -981,32 +673,6 @@ Standard_Boolean IntTools_WLineTool::
           }
         }
 
-        // check if a point belong to a tangent zone. Begin
-        for ( Standard_Integer zIt = 1; zIt <= aNbZone; zIt++ ) {
-          gp_Pnt2d aPZone = (surfit == 0) ? aTanZoneS1->Value(zIt) : aTanZoneS2->Value(zIt);
-          Standard_Real aZoneRadius = aTanZoneRadius->Value(zIt);
-
-          Standard_Integer aneighbourpointindex1 = (j == 0) ? iFirst : iLast;
-          const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex1);
-          Standard_Real nU1, nV1;
-            
-          if(surfit == 0)
-            aNeighbourPoint.ParametersOnS1(nU1, nV1);
-          else
-            aNeighbourPoint.ParametersOnS2(nU1, nV1);
-          gp_Pnt2d ap1(nU1, nV1);
-          gp_Pnt2d ap2 = AdjustByNeighbour( ap1, gp_Pnt2d( U, V ), aGASurface );
-
-
-          if ( IsInsideTanZone( ap2, aPZone, aZoneRadius, aGASurface ) ) {
-            aZoneIndex = zIt;
-            bIsNearBoundary = Standard_True;
-            if ( theReachedTol3d < aZoneRadius ) {
-              theReachedTol3d = aZoneRadius;
-            }
-          }
-        }
-        // check if a point belong to a tangent zone. End
         Standard_Boolean bComputeLineEnd = Standard_False;
 
         if(nbboundaries == 2) {
@@ -1145,20 +811,7 @@ Standard_Boolean IntTools_WLineTool::
             gp_Pnt2d ap1(nU1, nV1);
             gp_Pnt2d ap2;
 
-
-            if ( aZoneIndex ) {
-              // exclude point from a tangent zone
-              anewpoint = AdjustByNeighbour( ap1, anewpoint, aGASurface );
-              gp_Pnt2d aPZone = (surfit == 0) ? aTanZoneS1->Value(aZoneIndex) : aTanZoneS2->Value(aZoneIndex);
-              Standard_Real aZoneRadius = aTanZoneRadius->Value(aZoneIndex);
-
-              if ( FindPoint(ap1, anewpoint, umin, umax, vmin, vmax, 
-                             aPZone, aZoneRadius, aGASurface, ap2) ) {
-                anewpoint = ap2;
-                found = Standard_True;
-              }
-            }
-            else if ( aGASurface->IsUPeriodic() || aGASurface->IsVPeriodic() ) {
+            if ( aGASurface->IsUPeriodic() || aGASurface->IsVPeriodic() ) {
               // re-compute point near boundary if shifted on a period
               ap2 = AdjustByNeighbour( ap1, anewpoint, aGASurface );
 

diff --git a/src/LocOpe/LocOpe_WiresOnShape.cxx b/src/LocOpe/LocOpe_WiresOnShape.cxx
index 76fcc8191fb45ae2065f34ac1d48062493d13cf4..450ef7dea9bbeb80967f70bf04987eb4988da504 100644 (file)
--- a/src/LocOpe/LocOpe_WiresOnShape.cxx
+++ b/src/LocOpe/LocOpe_WiresOnShape.cxx
@@ -1300,6 +1300,10 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge,
   Standard_Real aTolV[2];
   aTolV[0] =BRep_Tool::Tolerance(theVertices[0]);
   aTolV[1] =BRep_Tool::Tolerance(theVertices[1]);
+  Standard_Real ext = 16.; // = 4 * 4 - to avoid creating microedges, area around vertices is increased
+                           // up to 4 vertex tolerance. Such approach is usual for other topological
+                           // algorithms, for example, Boolean Operations.
+  Standard_Real aTolVExt[2] = { ext * aTolV[0] * aTolV[0], ext * aTolV[1] * aTolV[1] };
 
   BRepAdaptor_Curve2d thePCurve(theEdge, theFace);
   Bnd_Box2d theBox;
@@ -1309,6 +1313,9 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge,
   Standard_Real  aFpar, aLpar;
   const Handle(Geom_Curve)& theCurve = BRep_Tool::Curve(theEdge, thePar[0], thePar[1]);
   GeomAdaptor_Curve theGAcurve(theCurve, thePar[0], thePar[1]);
+  Standard_Real aTolV2d[2] = { theGAcurve.Resolution(aTolV[0]), theGAcurve.Resolution(aTolV[1]) };
+  aTolV2d[0] = Max(aTolV2d[0], Precision::PConfusion());
+  aTolV2d[1] = Max(aTolV2d[1], Precision::PConfusion());
   Standard_Real aDistMax = Precision::Confusion() * Precision::Confusion();
   TopExp_Explorer Explo(theFace, TopAbs_EDGE);
   for (; Explo.More(); Explo.Next())
@@ -1335,6 +1342,23 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge,
       isOverlapped = Standard_True;
       return;
     }
+    // Check extremity distances
+    Standard_Real dists[4];
+    gp_Pnt aP11, aP12, aP21, aP22;
+    anExtrema.TrimmedSquareDistances(dists[0], dists[1], dists[2], dists[3],
+      aP11, aP12, aP21, aP22);
+    for (i = 0; i < 4; ++i)
+    {
+      if (i < 2)
+        j = 0;
+      else
+        j = 1;
+      if (dists[i] < aTolVExt[j] / ext)
+      {
+        return;
+      }
+    }
+
     for (i = 1; i <= aNbExt; i++)
     {
       Standard_Real aDist = anExtrema.SquareDistance(i);
@@ -1347,7 +1371,7 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge,
       Standard_Real anIntPar = aPOnC2.Parameter();
       for (j = 0; j < 2; j++) //try to find intersection on an extremity of "theEdge"
       {
-        if (Abs(theIntPar - thePar[j]) <= Precision::PConfusion())
+        if (Abs(theIntPar - thePar[j]) <= aTolV2d[j])
           break;
       }
       //intersection found in the middle of the edge
@@ -1356,10 +1380,10 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge,
         gp_Pnt aPoint = aCurve->Value(anIntPar);
         gp_Pnt aPointInt = theCurve->Value(theIntPar);
     
-        if (aPointInt.SquareDistance(thePnt[0]) > aTolV[0] * aTolV[0] &&
-          aPointInt.SquareDistance(thePnt[1]) > aTolV[1] * aTolV[1] &&
-          aPoint.SquareDistance(thePnt[0]) > aTolV[0] * aTolV[0] &&
-          aPoint.SquareDistance(thePnt[1]) > aTolV[1] * aTolV[1])
+        if (aPointInt.SquareDistance(thePnt[0]) > aTolVExt[0] &&
+          aPointInt.SquareDistance(thePnt[1]) > aTolVExt[1]  &&
+          aPoint.SquareDistance(thePnt[0]) > aTolVExt[0]  &&
+          aPoint.SquareDistance(thePnt[1]) > aTolVExt[1])
         {
           SplitPars.Append(theIntPar);
           if( aDist > aDistMax)

diff --git a/tests/bugs/modalg_5/bug23706_10 b/tests/bugs/modalg_5/bug23706_10
index d80c11f0e4aeb0eae688a9fd0c67487c412f667c..6a5d69e19df0d5783004fbfcd50a34cf2ae2600a 100755 (executable)
--- a/tests/bugs/modalg_5/bug23706_10
+++ b/tests/bugs/modalg_5/bug23706_10
@@ -13,7 +13,7 @@ set info [extrema r3 r4]
 
 if {[regexp "ext_1" $info]} {
   set dist [lindex [length ext_1] end]
-  if { $dist > 4.0e-13 } {
+  if { $dist > 5.0e-11 } {
     puts "Error: Extrema distance is too big"
   }
 } else {

diff --git a/tests/bugs/modalg_5/bug23706_11 b/tests/bugs/modalg_5/bug23706_11
index b1afdfe5a0393a3c4719feaa019a99c2015318f0..9df8d5c5f4beecd77cc51334c38ebdb4cc0cedd2 100755 (executable)
--- a/tests/bugs/modalg_5/bug23706_11
+++ b/tests/bugs/modalg_5/bug23706_11
@@ -10,8 +10,11 @@ bsplinecurve r1 2 5 1 3 2 1 3 1 4 1 5 3 2 5 3 1 3 7 3 1 4 8 3 1 4 8 3 1 5 9 3 1
 bsplinecurve r2 2 5 2 3 2.5 1 3 1 3.5 1 4 3 -1 2 3 1 1 11 3 1 3 9 3 1 3 9 3 1 3 9 3 1 5 7 3 1 7 4 3 1
 set info [extrema r1 r2]
 
-if { [llength $info] != 1 } {
-    puts "Error : Extrema is wrong"
+if {[regexp "ext_1" $info]} {
+  set dist [lindex [length ext_1] end]
+  if { $dist > 1.0e-8 } {
+    puts "Error: Extrema distance is too big"
+  }
 } else {
-    puts "OK: Extrema is valid"
+  puts "Error: Extrema is not found"
 }

diff --git a/tests/bugs/modalg_6/bug27665 b/tests/bugs/modalg_6/bug27665
index 85fe8294f20a42f186a7156127bdae464f27090d..7fa5fc9be794f5950270dea92a2b18c92acdc7ab 100644 (file)
--- a/tests/bugs/modalg_6/bug27665
+++ b/tests/bugs/modalg_6/bug27665
@@ -5,7 +5,7 @@ puts ""
 #################################################
 # BrepExrtrema_DistShapeShape bad performance on OCCT 6.7.0
 #################################################
-cpulimit 100
+cpulimit 500
 restore [locate_data_file bug27665_wircmpd.brep] w
 explode w
 

diff --git a/tests/lowalgos/extcc/bug29858_03 b/tests/lowalgos/extcc/bug29858_03
new file mode 100644 (file)
index 0000000..51702e2
--- /dev/null
+++ b/tests/lowalgos/extcc/bug29858_03
@@ -0,0 +1,29 @@
+puts "========"
+puts "OCC29858"
+puts "========"
+puts ""
+#################################################
+# Regression in GeomAPI_ExtremaCurveCurve
+#################################################
+
+# Read input
+restore [locate_data_file bug29858_03_e1.brep] e1
+restore [locate_data_file bug29858_03_e2.brep] e2
+
+# Extract geometry from topology
+mkcurve c1 e1
+mkcurve c2 e2
+
+# Run extrema
+set info [extrema c1 c2]
+
+# Check result
+if {[regexp "ext_1" $info]} {
+  set dist [lindex [length ext_1] end]
+  if { $dist > 1.0e-10 } {
+    puts "Error: Extrema distance is too big"
+  }
+} else {
+  puts "Error: Extrema is not found"
+}
+

diff --git a/tests/lowalgos/extcc/bug32882 b/tests/lowalgos/extcc/bug32882
new file mode 100644 (file)
index 0000000..9e91bbe
--- /dev/null
+++ b/tests/lowalgos/extcc/bug32882
@@ -0,0 +1,36 @@
+puts "========"
+puts "OCC32882: Modeling Data - Extrema curve/curve cannot find all solutions"
+puts "========"
+puts ""
+
+
+# Read input
+restore [locate_data_file bug32882.brep] cc
+explode cc
+
+# Extract geometry from topology
+mkcurve c1 cc_1
+mkcurve c2 cc_2
+mkcurve c3 cc_3
+
+# Run extrema c1/c3
+set info [extrema c1 c3]
+
+# Check number of solution
+if { [llength $info] != 3 } {
+   puts "Error: Invalid extrema number in extrema c1-c3 output"
+}
+
+# Check result
+checklength ext_1 -l 2.929642751e-14 -eps .01
+checklength ext_2 -l 3.480934286e-14 -eps .01
+checklength ext_3 -l 3.177643716e-14 -eps .01
+
+# Run extrema c3/c2
+set info [extrema c3 c2]
+# Check number of solutions
+if { [llength $info] != 2 } {
+   puts "Error: Invalid extrema number in extrema c3-c2 output"
+}
+checklength ext_1 -l 5.684341886e-14 -eps .01
+checklength ext_2 -l 2.929642751e-14 -eps .01