0024154: Wrong result of CUT operation
[occt.git] / src / BOPTools / BOPTools_AlgoTools.cxx
index 070c5ad..c9e24b8 100644 (file)
 #include <gp_Pnt.hxx>
 #include <gp_XYZ.hxx>
 #include <gp_Pnt2d.hxx>
+#include <gp_Cylinder.hxx>
+#include <gp_Cone.hxx>
+#include <gp_Sphere.hxx>
+#include <gp_Torus.hxx>
 //
 #include <Geom2d_Curve.hxx>
 #include <Geom_Surface.hxx>
+#include <Geom_Plane.hxx>
 #include <Geom_TrimmedCurve.hxx>
 #include <Geom_Curve.hxx>
 #include <GeomAPI_ProjectPointOnSurf.hxx>
 #include <BRep_Tool.hxx>
 #include <BRepLib.hxx>
 #include <BRepAdaptor_Curve2d.hxx>
+#include <BRepAdaptor_Surface.hxx>
 #include <BRepClass3d_SolidClassifier.hxx>
 #include <TopExp.hxx>
 #include <TopExp_Explorer.hxx>
 //
 #include <IntTools_Tools.hxx>
 //
+#include <BOPTools.hxx>
+#include <BOPTools_CoupleOfShape.hxx>
+#include <BOPTools_ListOfCoupleOfShape.hxx>
 #include <BOPTools_AlgoTools2D.hxx>
 #include <BOPTools_AlgoTools3D.hxx>
 //
 #include <BOPCol_IndexedMapOfShape.hxx>
 #include <BOPCol_MapOfShape.hxx>
 //
-#include <BOPTools.hxx>
-#include <BOPTools_CoupleOfShape.hxx>
-#include <BOPTools_ListOfCoupleOfShape.hxx>
-#include <Geom_SurfaceOfLinearExtrusion.hxx>
-#include <GeomAdaptor_Surface.hxx>
-#include <gp_Cylinder.hxx>
-#include <Geom_CylindricalSurface.hxx>
-#include <gp_Lin.hxx>
 #include <BOPInt_ShrunkRange.hxx>
+//
 
 static
   Standard_Real AngleWithRef(const gp_Dir& theD1,
@@ -85,16 +87,33 @@ static
   TopAbs_Orientation Orientation(const TopoDS_Edge& anE,
                                  const TopoDS_Face& aF);
 
-
+static
+  void GetFaceDir(const TopoDS_Edge& aE,
+                  const TopoDS_Face& aF,
+                  const gp_Pnt& aP,
+                  const Standard_Real aT,
+                  const gp_Dir& aDTgt,
+                  gp_Dir& aDN,
+                  gp_Dir& aDB,
+                  Handle(BOPInt_Context)& theContext,
+                  GeomAPI_ProjectPointOnSurf& aProjPL);
+static
+  Standard_Boolean FindPointInFace(const TopoDS_Edge& aE,
+                                   const TopoDS_Face& aF,
+                                   const gp_Pnt& aP,
+                                   gp_Dir& aDB,
+                                   gp_Pnt& aPOut,
+                                   Handle(BOPInt_Context)& theContext,
+                                   GeomAPI_ProjectPointOnSurf& aProjPL);
 
 //=======================================================================
 // function: MakeConnexityBlocks
 // purpose: 
 //=======================================================================
-  void BOPTools_AlgoTools::MakeConnexityBlocks (const TopoDS_Shape& theS,
-                                            const TopAbs_ShapeEnum theType1,
-                                            const TopAbs_ShapeEnum theType2,
-                                            BOPCol_ListOfShape& theLCB)
+void BOPTools_AlgoTools::MakeConnexityBlocks (const TopoDS_Shape& theS,
+                                              const TopAbs_ShapeEnum theType1,
+                                              const TopAbs_ShapeEnum theType2,
+                                              BOPCol_ListOfShape& theLCB)
 {
   Standard_Integer  aNbF, aNbAdd, aNbAdd1, i;
   BRep_Builder aBB;
@@ -121,7 +140,7 @@ static
     aMAdd.Clear();
     aMAdd.Add(aF1);
     //
-    while(1) {
+    for(;;) {
       aMAdd1.Clear();
       //
       aNbAdd = aMAdd.Extent();
@@ -177,7 +196,7 @@ static
 // function: OrientFacesOnShell
 // purpose: 
 //=======================================================================
-  void BOPTools_AlgoTools::OrientFacesOnShell (TopoDS_Shape& aShell)
+void BOPTools_AlgoTools::OrientFacesOnShell (TopoDS_Shape& aShell)
 {
   Standard_Boolean bIsProcessed1, bIsProcessed2;
   Standard_Integer i, aNbE, aNbF, j;
@@ -305,8 +324,8 @@ static
 //function : Orientation
 //purpose  :
 //=======================================================================
-  TopAbs_Orientation Orientation(const TopoDS_Edge& anE,
-                                 const TopoDS_Face& aF)
+TopAbs_Orientation Orientation(const TopoDS_Edge& anE,
+                               const TopoDS_Face& aF)
 {
   TopAbs_Orientation anOr=TopAbs_INTERNAL;
 
@@ -329,10 +348,10 @@ static
 // function: MakeConnexityBlock.
 // purpose: 
 //=======================================================================
-  void BOPTools_AlgoTools::MakeConnexityBlock (BOPCol_ListOfShape& theLFIn,
-                                           BOPCol_IndexedMapOfShape& theMEAvoid,
-                                           BOPCol_ListOfShape& theLCB,
-                                           const Handle(NCollection_BaseAllocator)& theAllocator)
+void BOPTools_AlgoTools::MakeConnexityBlock (BOPCol_ListOfShape& theLFIn,
+                                             BOPCol_IndexedMapOfShape& theMEAvoid,
+                                             BOPCol_ListOfShape& theLCB,
+                                             const Handle(NCollection_BaseAllocator)& theAllocator)
 {
   Standard_Integer  aNbF, aNbAdd1, aNbAdd, i;
   TopExp_Explorer aExp;
@@ -355,7 +374,7 @@ static
   const TopoDS_Shape& aF1=theLFIn.First();
   aMAdd.Add(aF1);
   //
-  while(1) {
+  for(;;) {
     aMAdd1.Clear();
     aNbAdd = aMAdd.Extent();
     for (i=1; i<=aNbAdd; ++i) {
@@ -409,10 +428,10 @@ static
 // function:  ComputeStateByOnePoint
 // purpose: 
 //=======================================================================
-  TopAbs_State BOPTools_AlgoTools::ComputeStateByOnePoint(const TopoDS_Shape& theS,
-                                                      const TopoDS_Solid& theRef,
-                                                      const Standard_Real theTol,
-                                                      Handle(BOPInt_Context)& theContext)
+TopAbs_State BOPTools_AlgoTools::ComputeStateByOnePoint(const TopoDS_Shape& theS,
+                                                        const TopoDS_Solid& theRef,
+                                                        const Standard_Real theTol,
+                                                        Handle(BOPInt_Context)& theContext)
 {
   TopAbs_State aState;
   TopAbs_ShapeEnum aType;
@@ -434,11 +453,11 @@ static
 // function:  ComputeState
 // purpose: 
 //=======================================================================
-  TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Face& theF,
-                                            const TopoDS_Solid& theRef,
-                                            const Standard_Real theTol,
-                                            BOPCol_IndexedMapOfShape& theBounds,
-                                            Handle(BOPInt_Context)& theContext)
+TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Face& theF,
+                                              const TopoDS_Solid& theRef,
+                                              const Standard_Real theTol,
+                                              BOPCol_IndexedMapOfShape& theBounds,
+                                              Handle(BOPInt_Context)& theContext)
 {
   TopAbs_State aState;
   TopExp_Explorer aExp; 
@@ -476,10 +495,10 @@ static
 // function:  ComputeState
 // purpose: 
 //=======================================================================
-  TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Vertex& theV,
-                                            const TopoDS_Solid& theRef,
-                                            const Standard_Real theTol,
-                                            Handle(BOPInt_Context)& theContext)
+TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Vertex& theV,
+                                              const TopoDS_Solid& theRef,
+                                              const Standard_Real theTol,
+                                              Handle(BOPInt_Context)& theContext)
 {
   TopAbs_State aState;
   gp_Pnt aP3D; 
@@ -492,10 +511,10 @@ static
 // function:  ComputeState
 // purpose: 
 //=======================================================================
-  TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Edge& theE,
-                                            const TopoDS_Solid& theRef,
-                                            const Standard_Real theTol,
-                                            Handle(BOPInt_Context)& theContext)
+TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Edge& theE,
+                                              const TopoDS_Solid& theRef,
+                                              const Standard_Real theTol,
+                                              Handle(BOPInt_Context)& theContext)
 {
   Standard_Real aT1, aT2, aT = 0.;
   TopAbs_State aState;
@@ -542,10 +561,10 @@ static
 // function:  ComputeState
 // purpose: 
 //=======================================================================
-  TopAbs_State BOPTools_AlgoTools::ComputeState(const gp_Pnt& theP,
-                                            const TopoDS_Solid& theRef,
-                                            const Standard_Real theTol,
-                                            Handle(BOPInt_Context)& theContext)
+TopAbs_State BOPTools_AlgoTools::ComputeState(const gp_Pnt& theP,
+                                              const TopoDS_Solid& theRef,
+                                              const Standard_Real theTol,
+                                              Handle(BOPInt_Context)& theContext)
 {
   TopAbs_State aState;
   //
@@ -560,22 +579,33 @@ static
 //function : IsInternalFace
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsInternalFace(const TopoDS_Face& theFace,
-                                                      const TopoDS_Solid& theSolid,
-                                                      BOPCol_IndexedDataMapOfShapeListOfShape& theMEF,
-                                                      const Standard_Real theTol,
-                                                      Handle(BOPInt_Context)& theContext)
+Standard_Integer BOPTools_AlgoTools::IsInternalFace
+  (const TopoDS_Face& theFace,
+   const TopoDS_Solid& theSolid,
+   BOPCol_IndexedDataMapOfShapeListOfShape& theMEF,
+   const Standard_Real theTol,
+   Handle(BOPInt_Context)& theContext)
 {
-  Standard_Boolean bRet, bDegenerated;
-  Standard_Integer aNbF;
+  Standard_Boolean bDegenerated;
+  Standard_Integer aNbF, iRet, iFound;
   TopAbs_Orientation aOr;
-  TopoDS_Edge aEL;
+  TopoDS_Edge aE1;
   TopExp_Explorer aExp;
   BOPCol_ListIteratorOfListOfShape aItF;
   //
-  bRet=Standard_False;
+  // For all invoked functions: [::IsInternalFace(...)] 
+  // the returned value iRet means:
+  // iRet=0;  - state is not IN
+  // iRet=1;  - state is IN
+  // iRet=2;  - state can not be found by the method of angles
+  //
+  // For this function the returned value iRet means:
+  // iRet=0;  - state is not IN
+  // iRet=1;  - state is IN
   //
+  iRet=0; 
   // 1 Try to find an edge from theFace in theMEF
+  iFound=0;
   aExp.Init(theFace, TopAbs_EDGE);
   for(; aExp.More(); aExp.Next()) {
     const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aExp.Current()));
@@ -583,6 +613,8 @@ static
       continue;
     }
     //
+    ++iFound;
+    //
     aOr=aE.Orientation();
     if (aOr==TopAbs_INTERNAL) {
       continue;
@@ -595,64 +627,85 @@ static
     BOPCol_ListOfShape& aLF=theMEF.ChangeFromKey(aE);
     aNbF=aLF.Extent();
     if (!aNbF) {
-      return bRet; // it can not be so
+      return iRet; // it can not be so
     }
+    //
     else if (aNbF==1) {
       // aE is internal edge on aLF.First()
       const TopoDS_Face& aF1=(*(TopoDS_Face*)(&aLF.First()));
-      bRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aF1, aF1, theContext);
-      return bRet;
+      BOPTools_AlgoTools::GetEdgeOnFace(aE, aF1, aE1);
+      if (aE1.Orientation()!=TopAbs_INTERNAL) {
+        iRet=2; 
+        break;
+      }
+      //
+      iRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aF1, aF1, theContext);
+      break;
     }
+    //
     else if (aNbF==2) {
       const TopoDS_Face& aF1=(*(TopoDS_Face*)(&aLF.First()));
       const TopoDS_Face& aF2=(*(TopoDS_Face*)(&aLF.Last()));
       //
       if (aF2.IsSame(aF1) && BRep_Tool::IsClosed(aE, aF1)) {
         // treat as it was for 1 face
-        bRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aF1, aF2, theContext);
-        return bRet;
+        iRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aF1, aF2, theContext);
+        break;
       }
     }
+    //
     if (aNbF%2) {
-      return bRet; // it can not be so
+      iRet=0;
+      return iRet; // it can not be so
     }
     else { // aNbF=2,4,6,8,...
-      bRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aLF, theContext);
-      return bRet;
+      iRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aLF, theContext);
+      break;
     }
   }//for(; aExp.More(); aExp.Next()) {
   //
+  if (!iFound) {
+    // the face has no shared edges with the solid
+    iRet=2;
+  }
+  //
+  if (iRet!=2) {
+    return iRet;
+  }
+  //
   //========================================
   // 2. Classify face using classifier
   //
   TopAbs_State aState;
   BOPCol_IndexedMapOfShape aBounds;
   //
+  BOPTools::MapShapes(theSolid, TopAbs_EDGE, aBounds);
+  //
   aState=BOPTools_AlgoTools::ComputeState(theFace, theSolid, theTol, aBounds, theContext);
-  bRet=(aState==TopAbs_IN);
   //
-  return bRet;
+  iRet=(aState==TopAbs_IN)? 1 : 0;
+  //
+  return iRet;
 }
 //=======================================================================
 //function : IsInternalFace
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsInternalFace(const TopoDS_Face& theFace,
-                                                      const TopoDS_Edge& theEdge,
-                                                      BOPCol_ListOfShape& theLF,
-                                                      Handle(BOPInt_Context)& theContext)
+Standard_Integer BOPTools_AlgoTools::IsInternalFace(const TopoDS_Face& theFace,
+                                                    const TopoDS_Edge& theEdge,
+                                                    BOPCol_ListOfShape& theLF,
+                                                    Handle(BOPInt_Context)& theContext)
 {
-  Standard_Boolean bRet;
-  Standard_Boolean aNbF;
+  Standard_Integer aNbF, iRet;
   //
-  bRet=Standard_False;
+  iRet=0;
   //
   aNbF=theLF.Extent();
   if (aNbF==2) {
     const TopoDS_Face& aF1=(*(TopoDS_Face*)(&theLF.First()));
     const TopoDS_Face& aF2=(*(TopoDS_Face*)(&theLF.Last()));
-    bRet=BOPTools_AlgoTools::IsInternalFace(theFace, theEdge, aF1, aF2, theContext);
-    return bRet;
+    iRet=BOPTools_AlgoTools::IsInternalFace(theFace, theEdge, aF1, aF2, theContext);
+    return iRet;
   }
   //
   else {
@@ -667,25 +720,27 @@ static
       //
       const TopoDS_Face& aF1=(*(TopoDS_Face*)(&aCSFF.Shape1()));
       const TopoDS_Face& aF2=(*(TopoDS_Face*)(&aCSFF.Shape2()));
-      bRet=BOPTools_AlgoTools::IsInternalFace(theFace, theEdge, aF1, aF2, theContext);
-      if (bRet) {
-        return bRet;
+      iRet=BOPTools_AlgoTools::IsInternalFace(theFace, theEdge, aF1, aF2, theContext);
+      if (iRet) {
+        return iRet;
       }
     }
   }
-  return bRet;
+  return iRet;
 }
 //=======================================================================
 //function : IsInternalFace
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsInternalFace(const TopoDS_Face& theFace,
-                                                      const TopoDS_Edge& theEdge,
-                                                      const TopoDS_Face& theFace1,
-                                                      const TopoDS_Face& theFace2,
-                                                      Handle(BOPInt_Context)& theContext)
+Standard_Integer BOPTools_AlgoTools::IsInternalFace
+  (const TopoDS_Face& theFace,
+   const TopoDS_Edge& theEdge,
+   const TopoDS_Face& theFace1,
+   const TopoDS_Face& theFace2,
+   Handle(BOPInt_Context)& theContext)
 {
   Standard_Boolean bRet;
+  Standard_Integer iRet;
   TopoDS_Edge aE1, aE2;
   TopoDS_Face aFOff;
   BOPTools_ListOfCoupleOfShape theLCSOff;
@@ -714,91 +769,103 @@ static
   aCS2.SetShape2(theFace2);
   theLCSOff.Append(aCS2);
   //
-  GetFaceOff(aE1, theFace1, theLCSOff, aFOff, theContext);
+  bRet=GetFaceOff(aE1, theFace1, theLCSOff, aFOff, theContext);
   //
-  bRet = theFace.IsEqual(aFOff);
-  return bRet;
+  iRet=0; // theFace is not internal
+  if (theFace.IsEqual(aFOff)) {
+    // theFace is internal
+    iRet=1;  
+    if (!bRet) {
+      // theFace seems to be internal  
+      iRet=2;
+    }
+  }
+  return iRet;
 }
 //=======================================================================
 //function : GetFaceOff
 //purpose  : 
 //=======================================================================
-  void BOPTools_AlgoTools::GetFaceOff(const TopoDS_Edge& theE1,
-                                      const TopoDS_Face& theF1,
-                                      BOPTools_ListOfCoupleOfShape& theLCSOff,
-                                      TopoDS_Face& theFOff,
-                                      Handle(BOPInt_Context)& theContext)
+Standard_Boolean BOPTools_AlgoTools::GetFaceOff
+  (const TopoDS_Edge& theE1,
+   const TopoDS_Face& theF1,
+   BOPTools_ListOfCoupleOfShape& theLCSOff,
+   TopoDS_Face& theFOff,
+   Handle(BOPInt_Context)& theContext)
 {
+  Standard_Boolean bRet;
   Standard_Real aT, aT1, aT2, aAngle, aTwoPI, aAngleMin;
+  Standard_Real aUmin, aUsup, aVmin, aVsup;
   gp_Pnt aPn1, aPn2, aPx;
+  gp_Dir aDN1, aDN2, aDBF, aDBF2, aDTF;
   gp_Vec aVTgt;
-  gp_Dir aDN1, aDN2;
+  TopAbs_Orientation aOr;
   Handle(Geom_Curve)aC3D;
+  Handle(Geom_Plane) aPL;
   BOPTools_ListIteratorOfListOfCoupleOfShape aIt;
+  GeomAPI_ProjectPointOnSurf aProjPL;
   //
   aAngleMin=100.;
   aTwoPI=M_PI+M_PI;
   aC3D =BRep_Tool::Curve(theE1, aT1, aT2);
-  //BRep_Tool::Range(theE1, aT1, aT2);
   aT=BOPTools_AlgoTools2D::IntermediatePoint(aT1, aT2);
   aC3D->D0(aT, aPx);
-  // Ref
+  //
   BOPTools_AlgoTools2D::EdgeTangent(theE1, aT, aVTgt);
-  gp_Dir aDTtgt(aVTgt);
-  aDTtgt.Reverse();
-  Handle(Geom_Plane) aPL;
-  aPL = new Geom_Plane(aPx, aDTtgt);
-  // N1
-  BOPTools_AlgoTools3D::GetApproxNormalToFaceOnEdge(theE1, theF1, aT, aPn1, aDN1, theContext);
-  //
-  gp_Pnt aPFx, aPF2x;
-  CorrectPoint(aPn1, aPL, aPFx);
-  gp_Vec aVBF (aPx, aPFx );
-  gp_Dir aDTF;
-  gp_Dir aDBF (aVBF);
+  gp_Dir aDTgt(aVTgt), aDTgt2;
+  aOr = theE1.Orientation();
+  //
+  aPL = new Geom_Plane(aPx, aDTgt);
+  aPL->Bounds(aUmin, aUsup, aVmin, aVsup);
+  aProjPL.Init(aPL, aUmin, aUsup, aVmin, aVsup);
+  //
+  GetFaceDir(theE1, theF1, aPx, aT, aDTgt, aDN1, aDBF, theContext, aProjPL);
+  //
   aDTF=aDN1^aDBF;
   //
+  bRet=Standard_True;
   aIt.Initialize(theLCSOff);
   for (; aIt.More(); aIt.Next()) {
     const BOPTools_CoupleOfShape& aCS=aIt.Value();
     const TopoDS_Edge& aE2=(*(TopoDS_Edge*)(&aCS.Shape1()));
     const TopoDS_Face& aF2=(*(TopoDS_Face*)(&aCS.Shape2()));
     //
-    /*if (aF2==theF1) {
-      aAngle=M_PI;
-    }
-    else if (aF2.IsSame(theF1)) {
-      aAngle=aTwoPI;
-    }
-    else {*/
-    if (!theE1.IsEqual(aE2) || 
-        !GetProjectPoint(theF1, aPn1, aF2, aPn2, aDN2, theContext)) {
-      BOPTools_AlgoTools3D::GetApproxNormalToFaceOnEdge (aE2, aF2, aT, aPn2, aDN2, theContext);
-    }
-    CorrectPoint(aPn2, aPL, aPF2x);
-    gp_Vec aVBF2(aPx, aPF2x);
-    gp_Dir aDBF2(aVBF2);
+    aDTgt2 = (aE2.Orientation()==aOr) ? aDTgt : aDTgt.Reversed();
+    GetFaceDir(aE2, aF2, aPx, aT, aDTgt2, aDN2, aDBF2, theContext, aProjPL);
     //Angle
     aAngle=AngleWithRef(aDBF, aDBF2, aDTF);
     //
     if(aAngle<0.) {
       aAngle=aTwoPI+aAngle;
     }
-    //}
-  
+    //
+    if (aAngle<Precision::Angular()) {
+      if (aF2==theF1) {
+        aAngle=M_PI;
+      }
+      else if (aF2.IsSame(theF1)) {
+        aAngle=aTwoPI;
+      }
+    }
+    //
     if (aAngle<aAngleMin){
       aAngleMin=aAngle;
       theFOff=aF2;
     }
+    else if (aAngle==aAngleMin) {
+      // the minimal angle can not be found
+      bRet=Standard_False;  
+    }
   }
+  return bRet;
 }
 //=======================================================================
 //function : GetEdgeOff
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::GetEdgeOff(const TopoDS_Edge& theE1,
-                                              const TopoDS_Face& theF2,
-                                              TopoDS_Edge& theE2)
+Standard_Boolean BOPTools_AlgoTools::GetEdgeOff(const TopoDS_Edge& theE1,
+                                                const TopoDS_Face& theF2,
+                                                TopoDS_Edge& theE2)
 {
   Standard_Boolean bFound;
   TopAbs_Orientation aOr1, aOr1C, aOr2;
@@ -827,12 +894,12 @@ static
 //function : AreFacesSameDomain
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::AreFacesSameDomain(const TopoDS_Face& theF1,
-                                                          const TopoDS_Face& theF2,
-                                                          Handle(BOPInt_Context)& theContext)
+Standard_Boolean BOPTools_AlgoTools::AreFacesSameDomain(const TopoDS_Face& theF1,
+                                                        const TopoDS_Face& theF2,
+                                                        Handle(BOPInt_Context)& theContext)
 {
   Standard_Boolean bFlag;
-  Standard_Integer iP;
+  Standard_Integer iErr;
   Standard_Real aTolF1, aTolF2, aTol;
   gp_Pnt2d aP2D;
   gp_Pnt aP;
@@ -849,27 +916,22 @@ static
   //
   aTolF1=BRep_Tool::Tolerance(aF1);
   // 1
-  iP=0;
   aExp.Init(aF1, TopAbs_EDGE);
   for (; aExp.More(); aExp.Next()) {
     aE1=(*(TopoDS_Edge*)(&aExp.Current()));
     if (!BRep_Tool::Degenerated(aE1)) {
-      iP=1;
-      //break;
       Standard_Real aTolE = BRep_Tool::Tolerance(aE1);
       aTolF1 = (aTolE > aTolF1) ? aTolE : aTolF1;
     }
   }
-  if (!iP) {
-    return bFlag;
-  }
-  //
   // 2
   aTolF2=BRep_Tool::Tolerance(aF2);
   aTol=aTolF1+aTolF2;
   //
-  BOPTools_AlgoTools3D::PointNearEdge(aE1, aF1, aP2D, aP, theContext);
-  bFlag=theContext->IsValidPointForFace(aP, aF2, aTol);
+  iErr = BOPTools_AlgoTools3D::PointInFace(aF1, aP, aP2D, theContext);
+  if (!iErr) {
+    bFlag=theContext->IsValidPointForFace(aP, aF2, aTol);
+  }
   //
   return bFlag;
 }
@@ -878,9 +940,9 @@ static
 //function : CheckSameGeom
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::CheckSameGeom(const TopoDS_Face& theF1,
-                                                 const TopoDS_Face& theF2,
-                                                 Handle(BOPInt_Context)& theContext)
+Standard_Boolean BOPTools_AlgoTools::CheckSameGeom(const TopoDS_Face& theF1,
+                                                   const TopoDS_Face& theF2,
+                                                   Handle(BOPInt_Context)& theContext)
 {
   Standard_Boolean bRet;
   Standard_Real aTolF1, aTolF2, aTol;
@@ -907,8 +969,8 @@ static
 // function: Sense
 // purpose: 
 //=======================================================================
-  Standard_Integer BOPTools_AlgoTools::Sense (const TopoDS_Face& theF1,
-                                          const TopoDS_Face& theF2)
+Standard_Integer BOPTools_AlgoTools::Sense (const TopoDS_Face& theF1,
+                                            const TopoDS_Face& theF2)
 {
   Standard_Integer iSense=0;
   gp_Dir aDNF1, aDNF2;
@@ -953,9 +1015,9 @@ static
 // function: IsSplitToReverse
 // purpose: 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Shape& theSp,
-                                                    const TopoDS_Shape& theSr,
-                                                    Handle(BOPInt_Context)& theContext)
+Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Shape& theSp,
+                                                      const TopoDS_Shape& theSr,
+                                                      Handle(BOPInt_Context)& theContext)
 {
   Standard_Boolean bRet;
   TopAbs_ShapeEnum aType;
@@ -987,9 +1049,9 @@ static
 //function :IsSplitToReverse
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Face& theFSp,
-                                                    const TopoDS_Face& theFSr,
-                                                    Handle(BOPInt_Context)& theContext)
+Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Face& theFSp,
+                                                      const TopoDS_Face& theFSr,
+                                                      Handle(BOPInt_Context)& theContext)
 {
   Standard_Boolean bRet, bFound, bInFace;
   Standard_Real aT1, aT2, aT, aU, aV, aScPr;
@@ -1076,9 +1138,9 @@ static
 //function :IsSplitToReverse
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Edge& aEF1,
-                                                    const TopoDS_Edge& aEF2,
-                                                    Handle(BOPInt_Context)& theContext)
+Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Edge& aEF1,
+                                                      const TopoDS_Edge& aEF2,
+                                                      Handle(BOPInt_Context)& theContext)
 {
   Standard_Boolean bRet, bIsDegenerated;
   //
@@ -1127,8 +1189,8 @@ static
 //function : IsHole
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsHole(const TopoDS_Shape& aW,
-                                          const TopoDS_Shape& aFace)
+Standard_Boolean BOPTools_AlgoTools::IsHole(const TopoDS_Shape& aW,
+                                            const TopoDS_Shape& aFace)
 {
   Standard_Boolean bIsHole;
   Standard_Integer i, aNbS;
@@ -1200,8 +1262,8 @@ static
 // function: MakeContainer
 // purpose: 
 //=======================================================================
-  void BOPTools_AlgoTools::MakeContainer(const TopAbs_ShapeEnum theType,
-                                     TopoDS_Shape& theC)
+void BOPTools_AlgoTools::MakeContainer(const TopAbs_ShapeEnum theType,
+                                       TopoDS_Shape& theC)
 {
   BRep_Builder aBB;
   //
@@ -1250,12 +1312,12 @@ static
 // function: MakePCurve
 // purpose: 
 //=======================================================================
-  void  BOPTools_AlgoTools::MakePCurve(const TopoDS_Edge& aE,
-                                       const TopoDS_Face& aF1,
-                                       const TopoDS_Face& aF2,
-                                       const IntTools_Curve& aIC,
-                                       const Standard_Boolean bPC1,
-                                      const Standard_Boolean bPC2)
+void BOPTools_AlgoTools::MakePCurve(const TopoDS_Edge& aE,
+                                    const TopoDS_Face& aF1,
+                                    const TopoDS_Face& aF2,
+                                    const IntTools_Curve& aIC,
+                                    const Standard_Boolean bPC1,
+                                    const Standard_Boolean bPC2)
 
 {
   Standard_Integer i;
@@ -1311,13 +1373,13 @@ static
 // function: MakeEdge
 // purpose: 
 //=======================================================================
-  void  BOPTools_AlgoTools::MakeEdge(const IntTools_Curve& theIC,
-                                 const TopoDS_Vertex& theV1,
-                                 const Standard_Real theT1,
-                                 const TopoDS_Vertex& theV2,
-                                 const Standard_Real theT2,
-                                 const Standard_Real theTolR3D,
-                                 TopoDS_Edge& theE)
+void BOPTools_AlgoTools::MakeEdge(const IntTools_Curve& theIC,
+                                  const TopoDS_Vertex& theV1,
+                                  const Standard_Real theT1,
+                                  const TopoDS_Vertex& theV2,
+                                  const Standard_Real theT2,
+                                  const Standard_Real theTolR3D,
+                                  TopoDS_Edge& theE)
 {
   Standard_Real aTolV;
   BRep_Builder aBB;
@@ -1340,9 +1402,9 @@ static
 // function: ComputeVV
 // purpose: 
 //=======================================================================
-  Standard_Integer BOPTools_AlgoTools::ComputeVV(const TopoDS_Vertex& aV1, 
-                                             const gp_Pnt& aP2,
-                                             const Standard_Real aTolP2)
+Standard_Integer BOPTools_AlgoTools::ComputeVV(const TopoDS_Vertex& aV1, 
+                                               const gp_Pnt& aP2,
+                                               const Standard_Real aTolP2)
 {
   Standard_Real aTolV1, aTolSum, aTolSum2, aD2;
   gp_Pnt aP1;
@@ -1364,8 +1426,8 @@ static
 // function: ComputeVV
 // purpose: 
 //=======================================================================
-  Standard_Integer BOPTools_AlgoTools::ComputeVV(const TopoDS_Vertex& aV1, 
-                                             const TopoDS_Vertex& aV2)
+Standard_Integer BOPTools_AlgoTools::ComputeVV(const TopoDS_Vertex& aV1, 
+                                               const TopoDS_Vertex& aV2)
 {
   Standard_Real aTolV1, aTolV2, aTolSum, aTolSum2, aD2;
   gp_Pnt aP1, aP2;
@@ -1388,8 +1450,8 @@ static
 // function: MakeVertex
 // purpose : 
 //=======================================================================
-  void BOPTools_AlgoTools::MakeVertex(BOPCol_ListOfShape& aLV,
-                                  TopoDS_Vertex& aVnew)
+void BOPTools_AlgoTools::MakeVertex(BOPCol_ListOfShape& aLV,
+                                    TopoDS_Vertex& aVnew)
 {
   Standard_Integer aNb;
   Standard_Real aTi, aDi, aDmax;
@@ -1432,9 +1494,9 @@ static
 //function : GetEdgeOnFace
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::GetEdgeOnFace(const TopoDS_Edge& theE1,
-                                                 const TopoDS_Face& theF2,
-                                                 TopoDS_Edge& theE2)
+Standard_Boolean BOPTools_AlgoTools::GetEdgeOnFace(const TopoDS_Edge& theE1,
+                                                   const TopoDS_Face& theF2,
+                                                   TopoDS_Edge& theE2)
 {
   Standard_Boolean bFound;
   TopoDS_Iterator aItF, aItW;
@@ -1467,7 +1529,7 @@ Standard_Boolean FindFacePairs (const TopoDS_Edge& theE,
 {
   Standard_Boolean bFound;
   Standard_Integer i, aNbCEF;
-  TopAbs_Orientation aOr, aOrC;
+  TopAbs_Orientation aOr, aOrC = TopAbs_FORWARD;
   BOPCol_MapOfShape aMFP;
   TopoDS_Face aF1, aF2;
   TopoDS_Edge aEL, aE1;
@@ -1601,30 +1663,14 @@ Standard_Real fsqrt(Standard_Real val)
   return (double)u.val;
 }
 
-//=======================================================================
-//function : CorrectPoint
-//purpose  : 
-//=======================================================================
-  void BOPTools_AlgoTools::CorrectPoint(const gp_Pnt& thePnt,
-                                    const Handle(Geom_Plane)& thePL,
-                                    gp_Pnt& theNewPnt)
-{
-  theNewPnt = thePnt;
-  GeomAPI_ProjectPointOnSurf aProjector;
-  aProjector.Init(thePnt, thePL);
-  if (aProjector.NbPoints()) {
-    theNewPnt = aProjector.NearestPoint();
-  }
-}
-
 //=======================================================================
 // function: IsBlockInOnFace
 // purpose: 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsBlockInOnFace (const IntTools_Range& aShrR,
-                                                        const TopoDS_Face& aF,
-                                                        const TopoDS_Edge& aE1,
-                                                        Handle(BOPInt_Context)& aContext)
+Standard_Boolean BOPTools_AlgoTools::IsBlockInOnFace (const IntTools_Range& aShrR,
+                                                      const TopoDS_Face& aF,
+                                                      const TopoDS_Edge& aE1,
+                                                      Handle(BOPInt_Context)& aContext)
 {
   Standard_Boolean bFlag;
   Standard_Real f1, l1, ULD, VLD;
@@ -1707,79 +1753,12 @@ Standard_Real fsqrt(Standard_Real val)
   return bFlag;
 }
 
-//=======================================================================
-//function : GetProjectPoint
-//purpose  : 
-//=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::GetProjectPoint(const TopoDS_Face& aF,
-                                                       const gp_Pnt& aPF, 
-                                                       const TopoDS_Face& aF1,
-                                                       gp_Pnt& aPF1, 
-                                                       gp_Dir& aDNF1,
-                                                       Handle(BOPInt_Context)& aContext)
-{
-  Standard_Boolean bRet = Standard_False;
-  Handle(Geom_Surface) aS, aS1;
-  GeomAbs_SurfaceType aTS, aTS1;
-  Handle(Geom_CylindricalSurface) aCS, aCS1;
-  Standard_Real aR, aR1, dR, aU, aV, aDMin;
-  //
-  aS = BRep_Tool::Surface(aF);
-  aS1 = BRep_Tool::Surface(aF1);
-  GeomAdaptor_Surface aGAS(aS), aGAS1(aS1);
-  //
-  aTS = aGAS.GetType();
-  aTS1 = aGAS1.GetType();
-  //
-  if (!(aTS == GeomAbs_Cylinder && aTS1 == GeomAbs_Cylinder)) {
-    return bRet;
-  }
-  //
-  aCS = Handle(Geom_CylindricalSurface)::DownCast(aS);
-  aCS1 = Handle(Geom_CylindricalSurface)::DownCast(aS1);
-  //
-  if (aCS.IsNull() || aCS1.IsNull()) {
-    return bRet;
-  }
-  //
-  if (!aCS->Axis().IsParallel(aCS1->Axis(), Precision::Angular())) {
-    return bRet;
-  }
-  //
-  aR = aCS->Radius();
-  aR1 = aCS1->Radius();
-  dR = fabs(aR - aR1);
-  //
-  if (dR < Precision::Angular()) {
-    return bRet;
-  }
-    
-  gp_Lin aL(aCS->Axis()), aL1(aCS1->Axis());    
-  if (Abs(aL.Distance(aL1) - dR) < Precision::Confusion()) {
-    GeomAPI_ProjectPointOnSurf& aProjector=aContext->ProjPS(aF1);
-    aProjector.Perform(aPF);
-    if (!aProjector.IsDone()) {
-      return bRet;
-    }
-    aDMin = aProjector.LowerDistance();
-    if (aDMin > dR) {
-      return bRet;
-    }
-    aProjector.LowerDistanceParameters(aU, aV);
-    Handle(Geom_Surface) aS = BRep_Tool::Surface(aF1);
-    BOPTools_AlgoTools3D::GetNormalToSurface (aS, aU, aV, aDNF1);
-    aS->D0(aU, aV, aPF1);
-    bRet = !bRet;
-  }
-  return bRet;
-}
-
 //=======================================================================
 //function : IsMicroEdge
 //purpose  : 
 //=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::IsMicroEdge(const TopoDS_Edge& aE,
-                                                   const Handle(BOPInt_Context)& aCtx) 
+Standard_Boolean BOPTools_AlgoTools::IsMicroEdge(const TopoDS_Edge& aE,
+                                                 const Handle(BOPInt_Context)& aCtx) 
 {
   Standard_Boolean bRet;
   Standard_Integer iErr;
@@ -1812,77 +1791,116 @@ Standard_Real fsqrt(Standard_Real val)
   return bRet;
 }
 
-/*
 //=======================================================================
-//function : AreFacesSameDomain
-//purpose  : 
-//=======================================================================
-  Standard_Boolean BOPTools_AlgoTools::AreFacesSameDomain(const TopoDS_Face& theF1,
-                                                      const TopoDS_Face& theF2,
-                                                      Handle(BOPInt_Context)& theContext)
-  {
-  Standard_Boolean bFlag;
-  
-  Standard_Integer aNbE1, aNbE2;
-  Standard_Real aTolF1, aTolF2, aTol;
-  gp_Pnt2d aP2D;
-  gp_Pnt aP;
-  TopoDS_Face aF1, aF2;
-  TopExp_Explorer aExp;
-  BOPCol_MapOfShape aME1, aME2;
-  BOPCol_MapIteratorOfMapOfShape aIt;
-  //
-  bFlag=Standard_False;
-  //
-  aF1=theF1;
-  aF1.Orientation(TopAbs_FORWARD);
-  aF2=theF2;
-  aF2.Orientation(TopAbs_FORWARD);
-  //
-  // 1
-  aExp.Init(aF1, TopAbs_EDGE);
-  for (; aExp.More(); aExp.Next()) {
-    const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aExp.Current()));
-    if (!BRep_Tool::Degenerated(aE)) {
-      aME1.Add(aE);
-    }
-  }
-  //
-  aExp.Init(aF2, TopAbs_EDGE);
-  for (; aExp.More(); aExp.Next()) {
-    const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aExp.Current()));
-    if (!BRep_Tool::Degenerated(aE)) {
-      if (!aME1.Contains(aE)) {
-      return bFlag;
-      }
-      aME2.Add(aE);
-    }
+//function : GetFaceDir
+//purpose  : Get binormal direction for the face in the point aP
+//=======================================================================
+void GetFaceDir(const TopoDS_Edge& aE,
+                const TopoDS_Face& aF,
+                const gp_Pnt& aP,
+                const Standard_Real aT,
+                const gp_Dir& aDTgt,
+                gp_Dir& aDN,
+                gp_Dir& aDB,
+                Handle(BOPInt_Context)& theContext,
+                GeomAPI_ProjectPointOnSurf& aProjPL)
+{
+  BOPTools_AlgoTools3D::GetNormalToFaceOnEdge(aE, aF, aT, aDN);
+  if (aF.Orientation()==TopAbs_REVERSED){
+    aDN.Reverse();
   }
+  aDB = aDN^aDTgt;
   //
-  aNbE1=aME1.Extent();
-  aNbE2=aME2.Extent();
-  //
-  if(!aNbE1 || !aNbE2){
-    return bFlag;
+  gp_Pnt aPx;
+  if (!FindPointInFace(aE, aF, aP, aDB, aPx, theContext, aProjPL)) {
+    BOPTools_AlgoTools3D::GetApproxNormalToFaceOnEdge(aE, aF, aT, aPx, aDN, theContext);
+    aProjPL.Perform(aPx);
+    aPx = aProjPL.NearestPoint();
+    gp_Vec aVec(aP, aPx);
+    aDB.SetXYZ(aVec.XYZ());
   }
-  //
-  if(aNbE1!=aNbE2) {
-    return bFlag;
+}
+
+//=======================================================================
+//function : FindPointInFace
+//purpose  : Find a point in the face in direction of <aDB>
+//=======================================================================
+Standard_Boolean FindPointInFace(const TopoDS_Edge& aE,
+                                 const TopoDS_Face& aF,
+                                 const gp_Pnt& aP,
+                                 gp_Dir& aDB,
+                                 gp_Pnt& aPOut,
+                                 Handle(BOPInt_Context)& theContext,
+                                 GeomAPI_ProjectPointOnSurf& aProjPL) 
+{
+  Standard_Integer aNbItMax;
+  Standard_Real aDt, aDtMin, aTolE, aTolF, aDist;
+  Standard_Boolean bRet;
+  gp_Pnt aP1;
+  BRepAdaptor_Surface aBAS;
+  //
+  bRet = Standard_False;
+  aTolE = BRep_Tool::Tolerance(aE);
+  aTolF = BRep_Tool::Tolerance(aF);
+  aDt = 2*(aTolE+aTolF);
+  aBAS.Initialize(aF, Standard_False);
+  //
+  aDtMin=5.e-4;
+  if (aDt < aDtMin) {
+    Standard_Real aR;
+    GeomAbs_SurfaceType aSType=aBAS.GetType();
+    switch (aSType) {
+    case GeomAbs_Cylinder:
+      aR = aBAS.Cylinder().Radius();
+      break;
+    case GeomAbs_Cone:
+      aR = aBAS.Cone().RefRadius();
+      break;
+    case GeomAbs_Sphere:
+      aR = aBAS.Sphere().Radius();
+      break;
+    case GeomAbs_Torus:
+      aR = aBAS.Torus().MinorRadius();
+      break;
+    case GeomAbs_SurfaceOfRevolution:
+      aR=1.;
+      break;
+    default:
+      aR=0.001;
+    }
+    if (aR < 0.01) {
+      aDtMin=5.e-6;
+    }
+    else if (aR > 100.) {
+      aDtMin*=10;
+    }
+    if (aDt < aDtMin) {
+      aDt=aDtMin;
+    }
   }
   //
-  // 2
-  aTolF1=BRep_Tool::Tolerance(aF1);
-  aTolF2=BRep_Tool::Tolerance(aF2);
-  aTol=aTolF1+aTolF2;
+  GeomAPI_ProjectPointOnSurf& aProj=theContext->ProjPS(aF);
+  aNbItMax = 15;
   //
-  aIt.Initialize(aME1);
-  for (; aIt.More(); aIt.Next()) {
-    const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aIt.Key()));
-    BOPTools_AlgoTools3D::PointNearEdge(aE, aF1, aP2D, aP);
-    bFlag=theContext->IsValidPointForFace(aP, aF2, aTol);
-    break;
-  }
+  do {
+    aP1.SetCoord(aP.X()+aDt*aDB.X(),
+                 aP.Y()+aDt*aDB.Y(),
+                 aP.Z()+aDt*aDB.Z());
+    //
+    aProj.Perform(aP1);
+    if (!aProj.IsDone()) {
+      return bRet;
+    }
+    aPOut = aProj.NearestPoint();
+    aDist = aProj.LowerDistance();
+    //
+    aProjPL.Perform(aPOut);
+    aPOut = aProjPL.NearestPoint();
+    //
+    gp_Vec aV(aP, aPOut);
+    aDB.SetXYZ(aV.XYZ());
+  } while (aDist>Precision::Angular() && --aNbItMax);
   //
-  return bFlag;
+  bRet = aDist < Precision::Angular();
+  return bRet;
 }
-*/