From: JGV <>
Date: Tue, 26 Jul 2011 17:25:28 +0000 (+0000)
Subject: OCC22610 The algorithm GeomAPI_ProjectPointOnSurf produces wrong results
X-Git-Tag: V6_5_2~63
X-Git-Url: http://git.dev.opencascade.org/gitweb/?p=occt.git;a=commitdiff_plain;h=38f33510d6b3eee07add7d304b3d263e53870ecd

OCC22610 The algorithm GeomAPI_ProjectPointOnSurf produces wrong results
---

diff --git a/src/Extrema/Extrema_GenExtPS.cxx b/src/Extrema/Extrema_GenExtPS.cxx
index 0cf3008b47..b76cf7882b 100755
--- a/src/Extrema/Extrema_GenExtPS.cxx
+++ b/src/Extrema/Extrema_GenExtPS.cxx
@@ -419,71 +419,82 @@ void Extrema_GenExtPS::FindSolution(const gp_Pnt& P, const math_Vector& UV, cons
   if (myF.HasDegIso())
     aNbMaxIter = 150;
 
+  gp_Pnt PStart = myS->Value(UV(1), UV(2));
+  Standard_Real DistStart = P.SquareDistance(PStart);
+  Standard_Real DistSol = DistStart;
+  
   math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
-  if(S.IsDone()) {
-    root = S.Root();
-    myF.Value(root, errors);
-    if(f == Extrema_ExtFlag_MIN) 
+  Standard_Boolean ToResolveOnSubgrid = Standard_False;
+  if (f == Extrema_ExtFlag_MIN)
+  {
+    if(S.IsDone())
     {
-      if(Abs(errors(1)) > eps || Abs(errors(2)) > eps) {
+      root = S.Root();
+      myF.Value(root, errors);
+      gp_Pnt PSol = myS->Value(root(1), root(2));
+      DistSol = P.SquareDistance(PSol);
+      if(Abs(errors(1)) > eps || Abs(errors(2)) > eps || DistStart < DistSol)
         //try to improve solution on subgrid of sample points
-        gp_Pnt PSol = myS->Value(root(1), root(2));
-        Standard_Real DistSol = P.SquareDistance(PSol);
-
-        Standard_Real u1 = Max(UV(1) - PasU, myumin), u2 = Min(UV(1) + PasU, myusup);
-        Standard_Real v1 = Max(UV(2) - PasV, myvmin), v2 = Min(UV(2) + PasV, myvsup);
-
-        if(u2 - u1 < 2.*PasU) {
-          if(Abs(u1 - myumin) < 1.e-9) {
-            u2 = u1 + 2.*PasU;
-            u2 = Min(u2, myusup);
-          }
-          if(Abs(u2 - myusup) < 1.e-9) {
-            u1 = u2 - 2.*PasU;
-            u1 = Max(u1, myumin);
-          }
+        ToResolveOnSubgrid = Standard_True;
+    }
+    else
+      ToResolveOnSubgrid = Standard_True;
+    
+    if (ToResolveOnSubgrid)
+    {
+      Standard_Real u1 = Max(UV(1) - PasU, myumin), u2 = Min(UV(1) + PasU, myusup);
+      Standard_Real v1 = Max(UV(2) - PasV, myvmin), v2 = Min(UV(2) + PasV, myvsup);
+      
+      if(u2 - u1 < 2.*PasU) {
+        if(Abs(u1 - myumin) < 1.e-9) {
+          u2 = u1 + 2.*PasU;
+          u2 = Min(u2, myusup);
         }
-
-        if(v2 - v1 < 2.*PasV) {
-          if(Abs(v1 - myvmin) < 1.e-9) {
-            v2 = v1 + 2.*PasV;
-            v2 = Min(v2, myvsup);
-          }
-          if(Abs(v2 - myvsup) < 1.e-9) {
-            v1 = v2 - 2.*PasV;
-            v1 = Max(v1, myvmin);
-          }
+        if(Abs(u2 - myusup) < 1.e-9) {
+          u1 = u2 - 2.*PasU;
+          u1 = Max(u1, myumin);
         }
-
-        Standard_Real du = (u2 - u1)/(nbsubsample-1);
-        Standard_Real dv = (v2 - v1)/(nbsubsample-1);
-        Standard_Real u, v;
-        Standard_Real dist;
-
-        Standard_Boolean NewSolution = Standard_False;
-        Standard_Integer Nu, Nv;
-        for (Nu = 1, u = u1; Nu < nbsubsample; Nu++, u += du) {
-          for (Nv = 1, v = v1; Nv < nbsubsample; Nv++, v += dv) {
-            gp_Pnt Puv = myS->Value(u, v);
-            dist = P.SquareDistance(Puv);
-
-            if(dist < DistSol) {
-              UV(1) = u;
-              UV(2) = v;
-              NewSolution = Standard_True;
-              DistSol = dist;
-            }
-          }
+      }
+      
+      if(v2 - v1 < 2.*PasV) {
+        if(Abs(v1 - myvmin) < 1.e-9) {
+          v2 = v1 + 2.*PasV;
+          v2 = Min(v2, myvsup);
         }
-
-        if(NewSolution) {
-          //try to precise
-          math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
+        if(Abs(v2 - myvsup) < 1.e-9) {
+          v1 = v2 - 2.*PasV;
+          v1 = Max(v1, myvmin);
         }
-
       }
-    }
-  }
+      
+      Standard_Real du = (u2 - u1)/(nbsubsample-1);
+      Standard_Real dv = (v2 - v1)/(nbsubsample-1);
+      Standard_Real u, v;
+      Standard_Real dist;
+      
+      Standard_Boolean NewSolution = Standard_False;
+      Standard_Integer Nu, Nv;
+      for (Nu = 1, u = u1; Nu < nbsubsample; Nu++, u += du) {
+        for (Nv = 1, v = v1; Nv < nbsubsample; Nv++, v += dv) {
+          gp_Pnt Puv = myS->Value(u, v);
+          dist = P.SquareDistance(Puv);
+          
+          if(dist < DistSol) {
+            UV(1) = u;
+            UV(2) = v;
+            NewSolution = Standard_True;
+            DistSol = dist;
+          }
+        }
+      }
+      
+      if(NewSolution) {
+        //try to precise
+        math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
+      }
+    } //end of if (ToResolveOnSubgrid)
+  } //end of if (f == Extrema_ExtFlag_MIN)
+  
   myDone = Standard_True;
 }
 
diff --git a/src/GeometryTest/FILES b/src/GeometryTest/FILES
index 9b93feba37..d79bef44e3 100755
--- a/src/GeometryTest/FILES
+++ b/src/GeometryTest/FILES
@@ -6,5 +6,4 @@ GeometryTest_API2dCommands.cxx
 GeometryTest_APICommands.cxx
 GeometryTest_ContinuityCommands.cxx
 GeometryTest_PolyCommands.cxx
-
-
+GeometryTest_TestProjCommands.cxx
diff --git a/src/GeometryTest/GeometryTest.cdl b/src/GeometryTest/GeometryTest.cdl
index 41990bd233..baadcce4b0 100755
--- a/src/GeometryTest/GeometryTest.cdl
+++ b/src/GeometryTest/GeometryTest.cdl
@@ -55,4 +55,7 @@ is
 	---Purpose: defines     command  to    test  the    polyhedral
 	--          triangulations and the polygons from the Poly package.
     
+    TestProjCommands(I : in out Interpretor from Draw);
+	---Purpose: defines commands to test projection of geometric objects
+    
 end GeometryTest;
diff --git a/src/GeometryTest/GeometryTest.cxx b/src/GeometryTest/GeometryTest.cxx
index 314069aa25..448803af38 100755
--- a/src/GeometryTest/GeometryTest.cxx
+++ b/src/GeometryTest/GeometryTest.cxx
@@ -24,6 +24,7 @@ void GeometryTest::AllCommands(Draw_Interpretor& theCommands)
 //  GeometryTest::API2dCommands(theCommands);
   GeometryTest::APICommands(theCommands);
   GeometryTest::ContinuityCommands(theCommands);
+  GeometryTest::TestProjCommands(theCommands);
   // define the TCL variable Draw_GEOMETRY
   //char* com = "set Draw_GEOMETRY 1";
   //theCommands.Eval(com);
diff --git a/src/GeometryTest/GeometryTest_TestProjCommands.cxx b/src/GeometryTest/GeometryTest_TestProjCommands.cxx
new file mode 100644
index 0000000000..fc3466afa1
--- /dev/null
+++ b/src/GeometryTest/GeometryTest_TestProjCommands.cxx
@@ -0,0 +1,112 @@
+// File:      GeometryTest_TestProjCommands.cxx
+// Created:   30.06.11 17:21:02
+// Author:    jgv@ROLEX
+// Copyright: Open CASCADE 2011
+
+#include <GeometryTest.hxx>
+#include <Geom_Curve.hxx>
+#include <Geom_Surface.hxx>
+#include <GeomAPI_ProjectPointOnSurf.hxx>
+
+#include <DBRep.hxx>
+#include <DBRep_DrawableShape.hxx>
+#include <Draw_Interpretor.hxx>
+#include <Draw_Appli.hxx>
+#include <DrawTrSurf.hxx>
+#include <Draw_Marker3D.hxx>
+
+#include <stdio.h>
+#ifdef WNT
+Standard_IMPORT Draw_Viewer dout;
+#endif
+
+//=======================================================================
+//function : xdistcs
+//purpose  : 
+//=======================================================================
+static Standard_Integer xdistcs(Draw_Interpretor& , Standard_Integer n, const char** a)
+{
+  if (n<5) {
+    cout<<" Use xdistcs c s t1 t2 nbp"<<endl;
+    return 0;
+  }
+  //
+  Standard_Boolean bRet;
+  Standard_Integer i, aNbP, iSize;
+  Standard_Real aTol, aD, aT, aT1, aT2, dT;
+  gp_Pnt aP;
+  Handle(Geom_Curve) aC;
+  Handle(Geom_Surface) aS;
+  GeomAPI_ProjectPointOnSurf aPPS;
+  Handle(Draw_Marker3D) aMr;
+  Draw_Color aColor(Draw_rouge);
+  //
+  aTol=1.e-7;
+  //
+  aC=DrawTrSurf::GetCurve(a[1]);
+  if (aC.IsNull()) {
+    cout<<a[1]<<" is null curve"<<endl;
+    return 0;
+  }
+  //
+  aS=DrawTrSurf::GetSurface(a[2]);
+  if (aS.IsNull()) {
+    cout<<a[2]<<" is null"<<endl;
+    return 0;
+  }
+  //
+  aT1=atof(a[3]);
+  aT2=atof(a[4]);
+  //
+  aNbP=10;
+  if (n>5) {
+    aNbP=atoi(a[5]);
+  }
+  //
+  iSize=3;
+  //
+  dT=(aT2-aT1)/(aNbP-1);
+  for (i=0; i<aNbP; ++i) {
+    aT=aT1+i*dT;
+    if (i==aNbP-1) {
+      aT=aT2;
+    }
+    //
+    aC->D0(aT, aP);
+    aPPS.Init(aP, aS, aTol);
+    bRet=aPPS.IsDone();
+    if (!bRet) {
+      cout<<" GeomAPI_ProjectPointOnSurf failed"<<endl;
+      return 0;
+    }
+    //
+    aD=aPPS.LowerDistance();
+    printf(" T=%lg\tD=%lg\n", aT, aD);
+    //
+    aMr=new Draw_Marker3D(aP, Draw_Plus, aColor, iSize);
+    dout << aMr;
+  }
+  //
+  return 0;
+}
+
+//=======================================================================
+//function : TestProjCommands
+//purpose  : 
+//=======================================================================
+
+void  GeometryTest::TestProjCommands(Draw_Interpretor& theCommands)
+{
+  
+  static Standard_Boolean loaded = Standard_False;
+  if (loaded) return;
+  loaded = Standard_True;
+  
+  DrawTrSurf::BasicCommands(theCommands);
+  
+  const char* g;
+  
+  g = "Testing of projection (geometric objects)";
+  
+  theCommands.Add("xdistcs", "xdistcs c s t1 t2 nbp", __FILE__, xdistcs, g);
+}