0022775: [CEA 1091] Add an option on GetNonBlocks to retrieve quadrangular faces...

[modules/geom.git] / src / BlockFix / BlockFix_CheckTool.cxx

diff --git a/src/BlockFix/BlockFix_CheckTool.cxx b/src/BlockFix/BlockFix_CheckTool.cxx
index 32c0714be12bc8ddfe723351794ca404b496ab0a..36e527550605b9afac35c33dc370a3fa6b2eef37 100644 (file)
--- a/src/BlockFix/BlockFix_CheckTool.cxx
+++ b/src/BlockFix/BlockFix_CheckTool.cxx
@@ -28,6 +28,11 @@
 
 #include <BRep_Tool.hxx>
 
+#include <Geom_Curve.hxx>
+
+#include <gp_Pnt.hxx>
+#include <gp_Vec.hxx>
+
 #include <TopExp.hxx>
 #include <TopExp_Explorer.hxx>
 
@@ -36,8 +41,10 @@
 #include <TopoDS_Face.hxx>
 #include <TopoDS_Solid.hxx>
 #include <TopoDS_Shape.hxx>
+#include <TopoDS_Vertex.hxx>
 
 #include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
+#include <TopTools_IndexedMapOfShape.hxx>
 #include <TopTools_MapOfShape.hxx>
 #include <TopTools_ListOfShape.hxx>
 #include <TopTools_ListIteratorOfListOfShape.hxx>
@@ -49,6 +56,7 @@
 BlockFix_CheckTool::BlockFix_CheckTool( )
 {
   myHasCheck = Standard_False;
+  myAngTolerance = -1.;
   myPossibleBlocks.Clear();
 }
 
@@ -63,6 +71,17 @@ void BlockFix_CheckTool::SetShape(const TopoDS_Shape& aShape)
   myPossibleBlocks.Clear();
 }
 
+//=======================================================================
+//function : SetAngTolerance
+//purpose  :
+//=======================================================================
+void BlockFix_CheckTool::SetAngTolerance(const Standard_Real theTolerance)
+{
+  myHasCheck     = Standard_False;
+  myAngTolerance = theTolerance;
+  myPossibleBlocks.Clear();
+}
+
 //=======================================================================
 //function : Perform
 //purpose  :
@@ -159,7 +178,6 @@ void BlockFix_CheckTool::Perform()
     if (nbe < 12)
       IsBlock = Standard_False;
     if (nbe > 12) {
-      IsBlock = Standard_False;
       // check edges unification
       // creating map of edge faces
       TopTools_IndexedDataMapOfShapeListOfShape aMapEdgeFaces;
@@ -194,9 +212,19 @@ void BlockFix_CheckTool::Perform()
         Standard_Integer i = 1;
         for (; i <= aMapFacesEdges.Extent(); i++) {
           const TopTools_ListOfShape& ListEdges = aMapFacesEdges.FindFromIndex(i);
-          if (ListEdges.Extent() > 1) break;
+          if (ListEdges.Extent() > 1) {
+            if (myAngTolerance < 0.) {
+              break;
+            }
+
+            // Check if edges have C1 continuity.
+            if (!isC1(ListEdges)) {
+              break;
+            }
+          }
         }
         if (i <= aMapFacesEdges.Extent()) {
+          IsBlock = Standard_False;
           MayBeUE = Standard_True;
           break;
         }
@@ -265,3 +293,84 @@ void BlockFix_CheckTool::DumpCheckResult(Standard_OStream& S) const
     S<<"             number of impossible blocks = "<<nbtmp<<endl;
   }
 }
+
+//=======================================================================
+//function : isC1
+//purpose  :
+//=======================================================================
+Standard_Boolean BlockFix_CheckTool::isC1
+          (const TopTools_ListOfShape &theEdges) const
+{
+  // Fill the map vertex - list of ancestor edges
+  TopTools_IndexedDataMapOfShapeListOfShape aMapVE;
+  TopTools_ListIteratorOfListOfShape        anIter(theEdges);
+  TopTools_MapOfShape                       aMapFence;
+  Standard_Integer                          i;
+  Standard_Integer                          aNbVtx;
+
+  for (; anIter.More(); anIter.Next()) {
+    TopTools_IndexedMapOfShape  aMapVtx;
+    const TopoDS_Shape         &anEdge = anIter.Value();
+
+    if (aMapFence.Add(anEdge)) {
+      TopExp::MapShapes(anEdge, TopAbs_VERTEX, aMapVtx);
+      aNbVtx = aMapVtx.Extent();
+
+      for (i = 1; i <= aNbVtx; ++i) {
+        const TopoDS_Shape &aVtx = aMapVtx.FindKey(i);
+
+        if (!aMapVE.Contains(aVtx)) {
+          aMapVE.Add(aVtx, TopTools_ListOfShape());
+        }
+
+        aMapVE.ChangeFromKey(aVtx).Append(anEdge);
+      }
+    }
+  }
+
+  // Check C1 continuity.
+  Standard_Integer aNbEnds = 0;
+
+  for (i = 1, aNbVtx = aMapVE.Extent(); i <= aNbVtx; ++i) {
+    const TopTools_ListOfShape &anEdges  = aMapVE.FindFromIndex(i);
+    Standard_Integer            aNbEdges = anEdges.Extent();
+
+    if (aNbEdges == 1) {
+      ++aNbEnds;
+    } else if (aNbEdges == 2) {
+      TopoDS_Vertex      aCommonVtx = TopoDS::Vertex(aMapVE.FindKey(i));
+      TopoDS_Edge        anEdge1    = TopoDS::Edge(anEdges.First());
+      TopoDS_Edge        anEdge2    = TopoDS::Edge(anEdges.Last());
+      Standard_Real      aParam1    = BRep_Tool::Parameter(aCommonVtx, anEdge1);
+      Standard_Real      aParam2    = BRep_Tool::Parameter(aCommonVtx, anEdge2);
+      Standard_Real      aPar[2];
+      Handle(Geom_Curve) aCurve1    =
+        BRep_Tool::Curve(anEdge1, aPar[0], aPar[1]);
+      Handle(Geom_Curve) aCurve2    =
+        BRep_Tool::Curve(anEdge2, aPar[0], aPar[1]);
+      gp_Pnt             aPnt;
+      gp_Vec             aVec1;
+      gp_Vec             aVec2;
+
+      aCurve1->D1(aParam1, aPnt, aVec1);
+      aCurve2->D1(aParam2, aPnt, aVec2);
+
+      if (anEdge1.Orientation() != anEdge2.Orientation()) {
+        // Orientations are different. One vector should be reversed.
+        aVec1.Reverse();
+      }
+
+      const Standard_Real anAngle = aVec1.Angle(aVec2);
+
+      if (anAngle > myAngTolerance) {
+        // There is no C1 continuity.
+        break;
+      }
+    } else {
+      // Non-manifold case.
+      break;
+    }
+  }
+
+  return (i > aNbVtx && aNbEnds == 2);
+}