0028211: Modeling Algorithms - Boolean fuse operation produces incorrect result
authormsv <msv@opencascade.com>
Fri, 29 Dec 2017 14:44:42 +0000 (17:44 +0300)
committerapn <apn@opencascade.com>
Fri, 12 Jan 2018 11:54:59 +0000 (14:54 +0300)
Correct procedure of initialization of BRepTopAdaptor_FClass2d and IntTools_FClass2d classifiers so as to produce more tight polygon in the case of self-intersections on very thin faces.

The idea is concluded in checking the condition:
defl < 2 * S / P, where S - is the surface area computed on produced polygon, P - its perimeter, defl - deflection computed on it.
If the condition is not true the polygon is discretized again using QuasiUniformDeflection tool.

src/BRepTopAdaptor/BRepTopAdaptor_FClass2d.cxx
src/IntTools/IntTools_FClass2d.cxx
tests/bugs/modalg_7/bug28211_1 [new file with mode: 0644]
tests/bugs/modalg_7/bug28211_2 [new file with mode: 0644]
tests/de/step_5/A1

index 527693d..1ab9716 100644 (file)
@@ -30,6 +30,7 @@
 #include <ElCLib.hxx>
 #include <Geom2dInt_Geom2dCurveTool.hxx>
 #include <GeomAbs_SurfaceType.hxx>
+#include <GCPnts_QuasiUniformDeflection.hxx>
 #include <gp_Pnt.hxx>
 #include <gp_Pnt2d.hxx>
 #include <Precision.hxx>
@@ -298,7 +299,7 @@ BRepTopAdaptor_FClass2d::BRepTopAdaptor_FClass2d(const TopoDS_Face& aFace,const
              PClass(im1)=SeqPnt2d.Value(im1);
              PClass(nbpnts)=SeqPnt2d.Value(nbpnts);
 
-
+              Standard_Real aPer = 0.;
 //           for(Standard_Integer ii=1; ii<nbpnts; ii++,im0++,im1++,im2++)
              for(Standard_Integer ii=1; ii<nbpnts; ii++,im0++,im1++)
                { 
@@ -310,10 +311,86 @@ BRepTopAdaptor_FClass2d::BRepTopAdaptor_FClass2d(const TopoDS_Face& aFace,const
 //               Standard_Real N = A.Magnitude() * B.Magnitude();
 
                  square += (PClass(im0).X()-PClass(im1).X())*(PClass(im0).Y()+PClass(im1).Y())*.5; 
+                  aPer += (PClass(im0).XY() - PClass(im1).XY()).Modulus();
 
 //               if(N>1e-16){ Standard_Real a=A.Angle(B); angle+=a; }
                }
 
+              Standard_Real anExpThick = Max(2. * Abs(square) / aPer, 1e-7);
+              Standard_Real aDefl = Max(FlecheU, FlecheV);
+              Standard_Real aDiscrDefl = Min(aDefl*0.1, anExpThick * 10.);
+              while (aDefl > anExpThick && aDiscrDefl > 1e-7)
+              {
+                // Deflection of the polygon is too much for this ratio of area and perimeter,
+                // and this might lead to self-intersections.
+                // Discretize the wire more tightly to eliminate the error.
+                firstpoint = 1;
+                SeqPnt2d.Clear();
+                FlecheU = 0.0;
+                FlecheV = 0.0;
+                for (WireExplorer.Init(TopoDS::Wire(FaceExplorer.Current()), Face);
+                  WireExplorer.More(); WireExplorer.Next())
+                {
+                  edge = WireExplorer.Current();
+                  Or = edge.Orientation();
+                  if (Or == TopAbs_FORWARD || Or == TopAbs_REVERSED)
+                  {
+                    Standard_Real pfbid, plbid;
+                    BRep_Tool::Range(edge, Face, pfbid, plbid);
+                    if (Abs(plbid - pfbid) < 1.e-9) continue;
+                    BRepAdaptor_Curve2d C(edge, Face);
+                    GCPnts_QuasiUniformDeflection aDiscr(C, aDiscrDefl);
+                    if (!aDiscr.IsDone())
+                      break;
+                    Standard_Integer nbp = aDiscr.NbPoints();
+                    Standard_Integer iStep = 1, i = 1, iEnd = nbp + 1;
+                    if (Or == TopAbs_REVERSED)
+                    {
+                      iStep = -1;
+                      i = nbp;
+                      iEnd = 0;
+                    }
+                    if (firstpoint == 2)
+                      i += iStep;
+                    for (; i != iEnd; i += iStep)
+                    {
+                      gp_Pnt2d aP2d = C.Value(aDiscr.Parameter(i));
+                      SeqPnt2d.Append(aP2d);
+                    }
+                    if (nbp > 2)
+                    {
+                      Standard_Integer ii = SeqPnt2d.Length();
+                      gp_Lin2d Lin(SeqPnt2d(ii - 2), gp_Dir2d(gp_Vec2d(SeqPnt2d(ii - 2), SeqPnt2d(ii))));
+                      Standard_Real ul = ElCLib::Parameter(Lin, SeqPnt2d(ii - 1));
+                      gp_Pnt2d Pp = ElCLib::Value(ul, Lin);
+                      Standard_Real dU = Abs(Pp.X() - SeqPnt2d(ii - 1).X());
+                      Standard_Real dV = Abs(Pp.Y() - SeqPnt2d(ii - 1).Y());
+                      if (dU > FlecheU) FlecheU = dU;
+                      if (dV > FlecheV) FlecheV = dV;
+                    }
+                    firstpoint = 2;
+                  }
+                }
+                nbpnts = SeqPnt2d.Length();
+                PClass.Resize(1, nbpnts, Standard_False);
+                im1 = nbpnts - 1;
+                im0 = 1;
+                PClass(im1) = SeqPnt2d.Value(im1);
+                PClass(nbpnts) = SeqPnt2d.Value(nbpnts);
+                square = 0.;
+                aPer = 0.;
+                for (Standard_Integer ii = 1; ii<nbpnts; ii++, im0++, im1++)
+                {
+                  if (im1 >= nbpnts) im1 = 1;
+                  PClass(ii) = SeqPnt2d.Value(ii);
+                  square += (PClass(im0).X() - PClass(im1).X())*(PClass(im0).Y() + PClass(im1).Y())*.5;
+                  aPer += (PClass(im0).XY() - PClass(im1).XY()).Modulus();
+                }
+
+                anExpThick = Max(2. * Abs(square) / aPer, 1e-7);
+                aDefl = Max(FlecheU, FlecheV);
+                aDiscrDefl = Min(aDiscrDefl * 0.1, anExpThick * 10.);
+              }
       
              //-- FlecheU*=10.0;
              //-- FlecheV*=10.0;
index 6e1edf5..13863dd 100644 (file)
@@ -26,6 +26,7 @@
 #include <Geom2dInt_Geom2dCurveTool.hxx>
 #include <GeomAbs_SurfaceType.hxx>
 #include <GeomInt.hxx>
+#include <GCPnts_QuasiUniformDeflection.hxx>
 #include <gp_Pnt.hxx>
 #include <gp_Pnt2d.hxx>
 #include <IntTools_FClass2d.hxx>
 #include <TopoDS_Wire.hxx>
 
 #include <stdio.h>
+
+//#define DEBUG_PCLASS_POLYGON
+#ifdef DEBUG_PCLASS_POLYGON
+#include <DrawTrSurf.hxx>
+#include <Geom2d_BSplineCurve.hxx>
+#endif
+
 //=======================================================================
 //function : IntTools_FClass2d:IntTools:_FClass2d
 //purpose  : 
@@ -386,7 +394,6 @@ void IntTools_FClass2d::Init(const TopoDS_Face& aFace,
         Standard_Integer im2=nbpnts-2;
         Standard_Integer im1=nbpnts-1;
         Standard_Integer im0=1;
-        Standard_Integer ii;
         Standard_Real    angle = 0.0;
         Standard_Real aX0, aY0, aX1, aY1, aS;
         //
@@ -397,7 +404,9 @@ void IntTools_FClass2d::Init(const TopoDS_Face& aFace,
         PClass(im2)=SeqPnt2d.Value(im2);
         PClass(im1)=SeqPnt2d.Value(im1);
         PClass(nbpnts)=SeqPnt2d.Value(nbpnts);
-        for(ii=1; ii<nbpnts; ii++,im0++,im1++,im2++) { 
+        Standard_Real aPer = 0.;
+        for (Standard_Integer ii = 1; ii<nbpnts; ii++, im0++, im1++, im2++)
+        {
           if(im2>=nbpnts) im2=1;
           if(im1>=nbpnts) im1=1;
           PClass(ii)=SeqPnt2d.Value(ii);
@@ -408,6 +417,7 @@ void IntTools_FClass2d::Init(const TopoDS_Face& aFace,
           aP2D0.Coord(aX0, aY0);
           aP2D1.Coord(aX1, aY1);
           aS=aS+(aY0+aY1)*(aX1-aX0); 
+          aPer += aP2D1.Distance(aP2D0);
 
           gp_Vec2d A(PClass(im2),PClass(im1));
           gp_Vec2d B(PClass(im1),PClass(im0));
@@ -455,6 +465,106 @@ void IntTools_FClass2d::Init(const TopoDS_Face& aFace,
         if (!iFlag) {
           angle = 0.; 
         }
+#ifdef DEBUG_PCLASS_POLYGON
+        TColStd_Array1OfReal aKnots(1, nbpnts);
+        TColStd_Array1OfInteger aMults(1, nbpnts);
+        for (int i = 1; i <= nbpnts; i++)
+        {
+          aKnots(i) = i;
+          aMults(i) = 1;
+        }
+        aMults(1) = aMults(nbpnts) = 2;
+        Handle(Geom2d_BSplineCurve) aPol = new Geom2d_BSplineCurve(PClass, aKnots, aMults, 1);
+        DrawTrSurf::Set("pol", aPol);
+#endif
+
+        Standard_Real anExpThick = Max(2. * Abs(aS) / aPer, 1e-7);
+        Standard_Real aDefl = Max(FlecheU, FlecheV);
+        Standard_Real aDiscrDefl = Min(aDefl*0.1, anExpThick * 10.);
+        while (aDefl > anExpThick && aDiscrDefl > 1e-7)
+        {
+          // Deflection of the polygon is too much for this ratio of area and perimeter,
+          // and this might lead to self-intersections.
+          // Discretize the wire more tightly to eliminate the error.
+          firstpoint = 1;
+          SeqPnt2d.Clear();
+          FlecheU = 0.0;
+          FlecheV = 0.0;
+          for (aWExp.Init(TopoDS::Wire(aExpF.Current()), Face);
+            aWExp.More(); aWExp.Next())
+          {
+            edge = aWExp.Current();
+            Or = edge.Orientation();
+            if (Or == TopAbs_FORWARD || Or == TopAbs_REVERSED)
+            {
+              BRep_Tool::Range(edge, Face, pfbid, plbid);
+              if (Abs(plbid - pfbid) < 1.e-9) continue;
+              BRepAdaptor_Curve2d C(edge, Face);
+              GCPnts_QuasiUniformDeflection aDiscr(C, aDiscrDefl);
+              if (!aDiscr.IsDone())
+                break;
+              Standard_Integer nbp = aDiscr.NbPoints();
+              Standard_Integer iStep = 1, i = 1, iEnd = nbp + 1;
+              if (Or == TopAbs_REVERSED)
+              {
+                iStep = -1;
+                i = nbp;
+                iEnd = 0;
+              }
+              if (firstpoint == 2)
+                i += iStep;
+              for (; i != iEnd; i += iStep)
+              {
+                gp_Pnt2d aP2d = C.Value(aDiscr.Parameter(i));
+                SeqPnt2d.Append(aP2d);
+              }
+              if (nbp > 2)
+              {
+                Standard_Integer ii = SeqPnt2d.Length();
+                gp_Lin2d Lin(SeqPnt2d(ii - 2), gp_Dir2d(gp_Vec2d(SeqPnt2d(ii - 2), SeqPnt2d(ii))));
+                Standard_Real ul = ElCLib::Parameter(Lin, SeqPnt2d(ii - 1));
+                gp_Pnt2d Pp = ElCLib::Value(ul, Lin);
+                Standard_Real dU = Abs(Pp.X() - SeqPnt2d(ii - 1).X());
+                Standard_Real dV = Abs(Pp.Y() - SeqPnt2d(ii - 1).Y());
+                if (dU > FlecheU) FlecheU = dU;
+                if (dV > FlecheV) FlecheV = dV;
+              }
+              firstpoint = 2;
+            }
+          }
+          nbpnts = SeqPnt2d.Length();
+          PClass.Resize(1, nbpnts, Standard_False);
+          im1 = nbpnts - 1;
+          im0 = 1;
+          PClass(im1) = SeqPnt2d.Value(im1);
+          PClass(nbpnts) = SeqPnt2d.Value(nbpnts);
+          aS = 0.;
+          aPer = 0.;
+          for (Standard_Integer ii = 1; ii<nbpnts; ii++, im0++, im1++)
+          {
+            if (im1 >= nbpnts) im1 = 1;
+            PClass(ii) = SeqPnt2d.Value(ii);
+            aS += (PClass(im1).X() - PClass(im0).X())*(PClass(im0).Y() + PClass(im1).Y())*.5;
+            aPer += (PClass(im0).XY() - PClass(im1).XY()).Modulus();
+          }
+#ifdef DEBUG_PCLASS_POLYGON
+          TColStd_Array1OfReal aKnots(1, nbpnts);
+          TColStd_Array1OfInteger aMults(1, nbpnts);
+          for (int i = 1; i <= nbpnts; i++)
+          {
+            aKnots(i) = i;
+            aMults(i) = 1;
+          }
+          aMults(1) = aMults(nbpnts) = 2;
+          Handle(Geom2d_BSplineCurve) aPol = new Geom2d_BSplineCurve(PClass, aKnots, aMults, 1);
+          DrawTrSurf::Set("pol1", aPol);
+#endif
+
+          anExpThick = Max(2. * Abs(aS) / aPer, 1e-7);
+          aDefl = Max(FlecheU, FlecheV);
+          aDiscrDefl = Min(aDiscrDefl * 0.1, anExpThick * 10.);
+        }
+
         if(aS>0.){
           myIsHole=Standard_False;
         }
diff --git a/tests/bugs/modalg_7/bug28211_1 b/tests/bugs/modalg_7/bug28211_1
new file mode 100644 (file)
index 0000000..bba45fb
--- /dev/null
@@ -0,0 +1,25 @@
+puts "========"
+puts "OCC28211"
+puts "========"
+puts ""
+#################################################
+# Modeling Algorithms - Boolean fuse operation produces incorrect result
+#################################################
+
+restore [locate_data_file bug28211_DDJ_BLD1_B_Pos3.brep] a
+
+explode a so
+bclearobjects
+bcleartools
+baddobjects a_1
+baddtools a_2
+bfillds
+
+# FUSE
+bbop result 1
+
+checkshape result
+checknbshapes result -solid 1 -shell 1
+checkprops result -v 13953.1 -s 4096.56
+
+checkview -display result -2d -path ${imagedir}/${test_image}.png
\ No newline at end of file
diff --git a/tests/bugs/modalg_7/bug28211_2 b/tests/bugs/modalg_7/bug28211_2
new file mode 100644 (file)
index 0000000..4fd6012
--- /dev/null
@@ -0,0 +1,26 @@
+puts "========"
+puts "OCC28211"
+puts "========"
+puts ""
+#################################################
+# Modeling Algorithms - Boolean fuse operation produces incorrect result
+#################################################
+
+restore [locate_data_file bug28211_MHX_VKT_WS_Pos3_source.brep] s
+restore [locate_data_file bug28211_MHX_VKT_WS_Pos3_tool.brep] t
+
+explode s f
+explode t f
+
+bnondestructive 1
+
+# Before the fix, exception in bsection algo was thrown.
+bsection r s_91 t_4
+
+bcut result s t
+
+checkshape result
+checknbshapes result -solid 4 -shell 4 -face 105 -wire 105
+checkprops result -v 226838 -s 36317.7
+
+checkview -display result -2d -path ${imagedir}/${test_image}.png
\ No newline at end of file
index a2df6c0..8d6a9c4 100755 (executable)
@@ -1,14 +1,16 @@
 # !!!! This file is generated automatically, do not edit manually! See end script
 
+puts "TODO CR28805 ALL: STATSHAPE : Faulty " 
+
 set LinuxDiff 3
 set filename Z8INV5.stp
 
 set ref_data {
 DATA        : Faulties = 0  ( 0 )  Warnings = 0  ( 0 )  Summary  = 0  ( 0 )
 TPSTAT      : Faulties = 0  ( 0 )  Warnings = 119  ( 620 )  Summary  = 119  ( 620 )
-CHECKSHAPE  : Wires    = 15  ( 16 )  Faces    = 17  ( 17 )  Shells   = 1  ( 1 )   Solids   = 0 ( 0 )
-NBSHAPES    : Solid    = 22  ( 22 )  Shell    = 24  ( 24 )  Face     = 1519  ( 1519 )   Summary  = 11224  ( 11212 )
-STATSHAPE   : Solid    = 22  ( 22 )  Shell    = 24  ( 24 )  Face     = 1519  ( 1519 )   FreeWire = 0  ( 0 )   FreeEdge  = 0 ( 0 )   SharedEdge = 4794  ( 4787 )
+CHECKSHAPE  : Wires    = 15  ( 16 )  Faces    = 17  ( 18 )  Shells   = 1  ( 1 )   Solids   = 0 ( 0 )
+NBSHAPES    : Solid    = 22  ( 22 )  Shell    = 25  ( 24 )  Face     = 1520  ( 1519 )   Summary  = 11227  ( 11212 )
+STATSHAPE   : Solid    = 22  ( 22 )  Shell    = 25  ( 24 )  Face     = 1520  ( 1519 )   FreeWire = 0  ( 0 )   FreeEdge  = 0 ( 0 )   SharedEdge = 4794  ( 4787 )
 TOLERANCE   : MaxTol   =    7.499684301  (    7.499684301 )  AvgTol   =   0.03452487556  (   0.03544461059 )
 LABELS      : N0Labels = 25  ( 25 )  N1Labels = 23  ( 23 )  N2Labels = 0  ( 0 )   TotalLabels = 48  ( 48 )   NameLabels = 48  ( 48 )   ColorLabels = 0  ( 0 )   LayerLabels = 0  ( 0 )
 PROPS       : Centroid = 0  ( 0 )  Volume   = 0  ( 0 )  Area     = 0  ( 0 )