src/BOPTools/BOPTools_AlgoTools.cxx

   1 // Created by: Peter KURNEV
   2 // Copyright (c) 2010-2012 OPEN CASCADE SAS
   3 // Copyright (c) 2007-2010 CEA/DEN, EDF R&D, OPEN CASCADE
   4 // Copyright (c) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, CEDRAT,
   5 //                         EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   6 //
   7 // The content of this file is subject to the Open CASCADE Technology Public
   8 // License Version 6.5 (the "License"). You may not use the content of this file
   9 // except in compliance with the License. Please obtain a copy of the License
  10 // at http://www.opencascade.org and read it completely before using this file.
  11 //
  12 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
  13 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
  14 //
  15 // The Original Code and all software distributed under the License is
  16 // distributed on an "AS IS" basis, without warranty of any kind, and the
  17 // Initial Developer hereby disclaims all such warranties, including without
  18 // limitation, any warranties of merchantability, fitness for a particular
  19 // purpose or non-infringement. Please see the License for the specific terms
  20 // and conditions governing the rights and limitations under the License.
  21
  22
  23 #include <BOPTools_AlgoTools.ixx>
  24 //
  25 #include <Precision.hxx>
  26 //
  27 #include <gp_Pnt.hxx>
  28 #include <gp_XYZ.hxx>
  29 #include <gp_Pnt2d.hxx>
  30 #include <gp_Cylinder.hxx>
  31 #include <gp_Cone.hxx>
  32 #include <gp_Sphere.hxx>
  33 #include <gp_Torus.hxx>
  34 //
  35 #include <Geom2d_Curve.hxx>
  36 #include <Geom_Surface.hxx>
  37 #include <Geom_Plane.hxx>
  38 #include <Geom_TrimmedCurve.hxx>
  39 #include <Geom_Curve.hxx>
  40 #include <GeomAPI_ProjectPointOnSurf.hxx>
  41 #include <Geom2dInt_Geom2dCurveTool.hxx>
  42 //
  43 #include <TopAbs_Orientation.hxx>
  44 //
  45 #include <TopoDS_Compound.hxx>
  46 #include <TopoDS_CompSolid.hxx>
  47 #include <TopoDS_Solid.hxx>
  48 #include <TopoDS_Shell.hxx>
  49 #include <TopoDS_Wire.hxx>
  50 //
  51 #include <BRep_Builder.hxx>
  52 #include <BRep_Tool.hxx>
  53 #include <BRepLib.hxx>
  54 #include <BRepAdaptor_Curve2d.hxx>
  55 #include <BRepAdaptor_Surface.hxx>
  56 #include <BRepClass3d_SolidClassifier.hxx>
  57 #include <TopExp.hxx>
  58 #include <TopExp_Explorer.hxx>
  59 //
  60 #include <IntTools_Tools.hxx>
  61 //
  62 #include <BOPTools.hxx>
  63 #include <BOPTools_CoupleOfShape.hxx>
  64 #include <BOPTools_ListOfCoupleOfShape.hxx>
  65 #include <BOPTools_AlgoTools2D.hxx>
  66 #include <BOPTools_AlgoTools3D.hxx>
  67 //
  68 #include <BOPCol_IndexedMapOfShape.hxx>
  69 #include <BOPCol_MapOfShape.hxx>
  70 //
  71 #include <BOPInt_ShrunkRange.hxx>
  72 //
  73
  74 static
  75   Standard_Real AngleWithRef(const gp_Dir& theD1,
  76                              const gp_Dir& theD2,
  77                              const gp_Dir& theDRef);
  78 static
  79   inline Standard_Real fsqrt(Standard_Real val);
  80
  81 static
  82   Standard_Boolean FindFacePairs (const TopoDS_Edge& theE,
  83                                   const BOPCol_ListOfShape& thLF,
  84                                   BOPTools_ListOfCoupleOfShape& theLCFF,
  85                                   Handle(BOPInt_Context)& theContext);
  86 static
  87   TopAbs_Orientation Orientation(const TopoDS_Edge& anE,
  88                                  const TopoDS_Face& aF);
  89
  90 static
  91   void GetFaceDir(const TopoDS_Edge& aE,
  92                   const TopoDS_Face& aF,
  93                   const gp_Pnt& aP,
  94                   const Standard_Real aT,
  95                   const gp_Dir& aDTgt,
  96                   gp_Dir& aDN,
  97                   gp_Dir& aDB,
  98                   Handle(BOPInt_Context)& theContext,
  99                   GeomAPI_ProjectPointOnSurf& aProjPL);
 100 static
 101   Standard_Boolean FindPointInFace(const TopoDS_Edge& aE,
 102                                    const TopoDS_Face& aF,
 103                                    const gp_Pnt& aP,
 104                                    gp_Dir& aDB,
 105                                    gp_Pnt& aPOut,
 106                                    Handle(BOPInt_Context)& theContext,
 107                                    GeomAPI_ProjectPointOnSurf& aProjPL);
 108
 109 //=======================================================================
 110 // function: MakeConnexityBlocks
 111 // purpose:
 112 //=======================================================================
 113 void BOPTools_AlgoTools::MakeConnexityBlocks (const TopoDS_Shape& theS,
 114                                               const TopAbs_ShapeEnum theType1,
 115                                               const TopAbs_ShapeEnum theType2,
 116                                               BOPCol_ListOfShape& theLCB)
 117 {
 118   Standard_Integer  aNbF, aNbAdd, aNbAdd1, i;
 119   BRep_Builder aBB;
 120   TopoDS_Compound aC;
 121   TopoDS_Iterator aIt;
 122   TopExp_Explorer aExp;
 123   BOPCol_MapOfShape aMP;
 124   BOPCol_IndexedMapOfShape aMCB, aMAdd, aMAdd1;
 125   BOPCol_IndexedDataMapOfShapeListOfShape aMEF;
 126   BOPCol_ListIteratorOfListOfShape aItLF;
 127   //
 128   // 1. aMEF
 129   BOPTools::MapShapesAndAncestors(theS, theType1, theType2, aMEF);
 130   //
 131   // 2. aMCB
 132   aIt.Initialize(theS);
 133   for (; aIt.More(); aIt.Next()) {
 134     const TopoDS_Shape& aF1=aIt.Value();
 135     if (aMP.Contains(aF1)) {
 136       continue;
 137     }
 138     //
 139     aMCB.Clear();
 140     aMAdd.Clear();
 141     aMAdd.Add(aF1);
 142     //
 143     for(;;) {
 144       aMAdd1.Clear();
 145       //
 146       aNbAdd = aMAdd.Extent();
 147       for (i=1; i<=aNbAdd; ++i) {
 148         const TopoDS_Shape& aF=aMAdd(i);
 149         //
 150         aExp.Init(aF, theType1);
 151         for (; aExp.More(); aExp.Next()) {
 152           const TopoDS_Shape& aE=aExp.Current();
 153           //
 154           const BOPCol_ListOfShape& aLF=aMEF.FindFromKey(aE);
 155           aItLF.Initialize(aLF);
 156           for (; aItLF.More(); aItLF.Next()) {
 157             const TopoDS_Shape& aFx=aItLF.Value();
 158             if (aFx.IsSame(aF)) {
 159               continue;
 160             }
 161             if (aMCB.Contains(aFx)) {
 162               continue;
 163             }
 164             aMAdd1.Add(aFx);
 165           }
 166         }//for (; aExp.More(); aExp.Next()){
 167         aMCB.Add(aF);
 168       }// for (i=1; i<=aNbAdd; ++i) {
 169       //
 170       aNbAdd1=aMAdd1.Extent();
 171       if (!aNbAdd1) {
 172         break;// ->make new CB from aMCB
 173       }
 174       //
 175       aMAdd.Clear();
 176       for (i=1; i<=aNbAdd1; ++i) {
 177         const TopoDS_Shape& aFAdd = aMAdd1(i);
 178         aMAdd.Add(aFAdd);
 179       }
 180     }//while(1) {
 181     //
 182     aNbF=aMCB.Extent();
 183     if (aNbF) {
 184       aBB.MakeCompound(aC);
 185       //
 186       for (i=1; i<=aNbF; ++i) {
 187         const TopoDS_Shape& aF=aMCB(i);
 188         aBB.Add(aC, aF);
 189         aMP.Add(aF);
 190       }
 191       theLCB.Append(aC);
 192     }
 193   }// for (; aIt.More(); aIt.Next())
 194 }
 195 //=======================================================================
 196 // function: OrientFacesOnShell
 197 // purpose:
 198 //=======================================================================
 199 void BOPTools_AlgoTools::OrientFacesOnShell (TopoDS_Shape& aShell)
 200 {
 201   Standard_Boolean bIsProcessed1, bIsProcessed2;
 202   Standard_Integer i, aNbE, aNbF, j;
 203   TopAbs_Orientation anOrE1, anOrE2;
 204   TopoDS_Face aF1x, aF2x;
 205   TopoDS_Shape aShellNew;
 206   BOPCol_IndexedDataMapOfShapeListOfShape aEFMap;
 207   BOPCol_IndexedMapOfShape aProcessedFaces;
 208   BRep_Builder aBB;
 209   //
 210   BOPTools_AlgoTools::MakeContainer(TopAbs_SHELL, aShellNew);
 211   //
 212   BOPTools::MapShapesAndAncestors(aShell, TopAbs_EDGE, TopAbs_FACE, aEFMap);
 213   aNbE=aEFMap.Extent();
 214   //
 215   // One seam edge  in aEFMap contains  2 equivalent faces.
 216   for (i=1; i<=aNbE; ++i) {
 217     BOPCol_ListOfShape& aLF=aEFMap.ChangeFromIndex(i);
 218     aNbF=aLF.Extent();
 219     if (aNbF>1) {
 220       BOPCol_ListOfShape aLFTmp;
 221       BOPCol_IndexedMapOfShape aFM;
 222       //
 223       BOPCol_ListIteratorOfListOfShape anIt(aLF);
 224       for (; anIt.More(); anIt.Next()) {
 225         const TopoDS_Shape& aF=anIt.Value();
 226         if (!aFM.Contains(aF)) {
 227           aFM.Add(aF);
 228           aLFTmp.Append(aF);
 229         }
 230       }
 231       aLF.Clear();
 232       aLF=aLFTmp;
 233     }
 234   }
 235   //
 236   // Do
 237   for (i=1; i<=aNbE; ++i) {
 238     const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aEFMap.FindKey(i)));
 239     if (BRep_Tool::Degenerated(aE)) {
 240       continue;
 241     }
 242     //
 243     const BOPCol_ListOfShape& aLF=aEFMap.FindFromIndex(i);
 244     aNbF=aLF.Extent();
 245     if (aNbF!=2) {
 246       continue;
 247     }
 248     //
 249     TopoDS_Face& aF1=(*(TopoDS_Face*)(&aLF.First()));
 250     TopoDS_Face& aF2=(*(TopoDS_Face*)(&aLF.Last()));
 251     //
 252     bIsProcessed1=aProcessedFaces.Contains(aF1);
 253     bIsProcessed2=aProcessedFaces.Contains(aF2);
 254     if (bIsProcessed1 && bIsProcessed2) {
 255       continue;
 256     }
 257
 258     if (!bIsProcessed1 && !bIsProcessed2) {
 259       aProcessedFaces.Add(aF1);
 260       aBB.Add(aShellNew, aF1);
 261       bIsProcessed1=!bIsProcessed1;
 262     }
 263     //
 264     aF1x=aF1;
 265     if (bIsProcessed1) {
 266       j=aProcessedFaces.FindIndex(aF1);
 267       aF1x=(*(TopoDS_Face*)(&aProcessedFaces.FindKey(j)));
 268     }
 269     //
 270     aF2x=aF2;
 271     if (bIsProcessed2) {
 272       j=aProcessedFaces.FindIndex(aF2);
 273       aF2x=(*(TopoDS_Face*)(&aProcessedFaces.FindKey(j)));
 274     }
 275     //
 276     anOrE1=Orientation(aE, aF1x);
 277     anOrE2=Orientation(aE, aF2x);
 278     //
 279     if (bIsProcessed1 && !bIsProcessed2) {
 280       if (anOrE1==anOrE2) {
 281         if (!BRep_Tool::IsClosed(aE, aF1) &&
 282             !BRep_Tool::IsClosed(aE, aF2)) {
 283           aF2.Reverse();
 284         }
 285       }
 286       aProcessedFaces.Add(aF2);
 287       aBB.Add(aShellNew, aF2);
 288     }
 289     else if (!bIsProcessed1 && bIsProcessed2) {
 290       if (anOrE1==anOrE2) {
 291         if (!BRep_Tool::IsClosed(aE, aF1) &&
 292             !BRep_Tool::IsClosed(aE, aF2)) {
 293           aF1.Reverse();
 294         }
 295       }
 296       aProcessedFaces.Add(aF1);
 297       aBB.Add(aShellNew, aF1);
 298     }
 299   }
 300   //
 301   //
 302   for (i=1; i<=aNbE; ++i) {
 303     const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aEFMap.FindKey(i)));
 304     if (BRep_Tool::Degenerated(aE)) {
 305       continue;
 306     }
 307     //
 308     const BOPCol_ListOfShape& aLF=aEFMap.FindFromIndex(i);
 309     aNbF=aLF.Extent();
 310     if (aNbF!=2) {
 311       BOPCol_ListIteratorOfListOfShape anIt(aLF);
 312       for(; anIt.More(); anIt.Next()) {
 313         const TopoDS_Face& aF=(*(TopoDS_Face*)(&anIt.Value()));
 314         if (!aProcessedFaces.Contains(aF)) {
 315           aProcessedFaces.Add(aF);
 316           aBB.Add(aShellNew, aF);
 317         }
 318       }
 319     }
 320   }
 321   aShell=aShellNew;
 322 }
 323 //=======================================================================
 324 //function : Orientation
 325 //purpose  :
 326 //=======================================================================
 327 TopAbs_Orientation Orientation(const TopoDS_Edge& anE,
 328                                const TopoDS_Face& aF)
 329 {
 330   TopAbs_Orientation anOr=TopAbs_INTERNAL;
 331
 332   TopExp_Explorer anExp;
 333   anExp.Init(aF, TopAbs_EDGE);
 334   for (; anExp.More(); anExp.Next()) {
 335     const TopoDS_Edge& anEF1=(*(TopoDS_Edge*)(&anExp.Current()));
 336     if (anEF1.IsSame(anE)) {
 337       anOr=anEF1.Orientation();
 338       break;
 339     }
 340   }
 341   return anOr;
 342 }
 343
 344 ////////////
 345
 346
 347 //=======================================================================
 348 // function: MakeConnexityBlock.
 349 // purpose:
 350 //=======================================================================
 351 void BOPTools_AlgoTools::MakeConnexityBlock (BOPCol_ListOfShape& theLFIn,
 352                                              BOPCol_IndexedMapOfShape& theMEAvoid,
 353                                              BOPCol_ListOfShape& theLCB,
 354                                              const Handle(NCollection_BaseAllocator)& theAllocator)
 355 {
 356   Standard_Integer  aNbF, aNbAdd1, aNbAdd, i;
 357   TopExp_Explorer aExp;
 358   BOPCol_ListIteratorOfListOfShape aIt;
 359   //
 360   BOPCol_IndexedMapOfShape aMCB(100, theAllocator);
 361   BOPCol_IndexedMapOfShape aMAdd(100, theAllocator);
 362   BOPCol_IndexedMapOfShape aMAdd1(100, theAllocator);
 363   BOPCol_IndexedDataMapOfShapeListOfShape aMEF(100, theAllocator);
 364   //
 365   // 1. aMEF
 366   aNbF=theLFIn.Extent();
 367   aIt.Initialize(theLFIn);
 368   for (; aIt.More(); aIt.Next()) {
 369     const TopoDS_Shape& aF=aIt.Value();
 370     BOPTools::MapShapesAndAncestors(aF, TopAbs_EDGE, TopAbs_FACE, aMEF);
 371   }
 372   //
 373   // 2. aMCB
 374   const TopoDS_Shape& aF1=theLFIn.First();
 375   aMAdd.Add(aF1);
 376   //
 377   for(;;) {
 378     aMAdd1.Clear();
 379     aNbAdd = aMAdd.Extent();
 380     for (i=1; i<=aNbAdd; ++i) {
 381       const TopoDS_Shape& aF=aMAdd(i);
 382       //
 383       //aMAdd1.Clear();
 384       aExp.Init(aF, TopAbs_EDGE);
 385       for (; aExp.More(); aExp.Next()) {
 386         const TopoDS_Shape& aE=aExp.Current();
 387         if (theMEAvoid.Contains(aE)){
 388           continue;
 389         }
 390         //
 391         const BOPCol_ListOfShape& aLF=aMEF.FindFromKey(aE);
 392         aIt.Initialize(aLF);
 393         for (; aIt.More(); aIt.Next()) {
 394           const TopoDS_Shape& aFx=aIt.Value();
 395           if (aFx.IsSame(aF)) {
 396             continue;
 397           }
 398           if (aMCB.Contains(aFx)) {
 399             continue;
 400           }
 401           aMAdd1.Add(aFx);
 402         }
 403       }//for (; aExp.More(); aExp.Next()){
 404       aMCB.Add(aF);
 405     }// for (i=1; i<=aNbAdd; ++i) {
 406     //
 407     aNbAdd1=aMAdd1.Extent();
 408     if (!aNbAdd1) {
 409       break;
 410     }
 411     //
 412     aMAdd.Clear();
 413     for (i=1; i<=aNbAdd1; ++i) {
 414       const TopoDS_Shape& aFAdd=aMAdd1(i);
 415       aMAdd.Add(aFAdd);
 416     }
 417     //
 418   }//while(1) {
 419
 420   //
 421   aNbF=aMCB.Extent();
 422   for (i=1; i<=aNbF; ++i) {
 423     const TopoDS_Shape& aF=aMCB(i);
 424     theLCB.Append(aF);
 425   }
 426 }
 427 //=======================================================================
 428 // function:  ComputeStateByOnePoint
 429 // purpose:
 430 //=======================================================================
 431 TopAbs_State BOPTools_AlgoTools::ComputeStateByOnePoint(const TopoDS_Shape& theS,
 432                                                         const TopoDS_Solid& theRef,
 433                                                         const Standard_Real theTol,
 434                                                         Handle(BOPInt_Context)& theContext)
 435 {
 436   TopAbs_State aState;
 437   TopAbs_ShapeEnum aType;
 438   //
 439   aState=TopAbs_UNKNOWN;
 440   aType=theS.ShapeType();
 441   if (aType==TopAbs_VERTEX) {
 442     const TopoDS_Vertex& aV=(*(TopoDS_Vertex*)(&theS));
 443     aState=BOPTools_AlgoTools::ComputeState(aV, theRef, theTol, theContext);
 444   }
 445   else if (aType==TopAbs_EDGE) {
 446     const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&theS));
 447     aState=BOPTools_AlgoTools::ComputeState(aE, theRef, theTol, theContext);
 448   }
 449   return aState;
 450 }
 451
 452 //=======================================================================
 453 // function:  ComputeState
 454 // purpose:
 455 //=======================================================================
 456 TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Face& theF,
 457                                               const TopoDS_Solid& theRef,
 458                                               const Standard_Real theTol,
 459                                               BOPCol_IndexedMapOfShape& theBounds,
 460                                               Handle(BOPInt_Context)& theContext)
 461 {
 462   TopAbs_State aState;
 463   TopExp_Explorer aExp;
 464   TopoDS_Edge aE1;
 465   gp_Pnt2d aP2D;
 466   gp_Pnt aP3D;
 467   //
 468   aState=TopAbs_UNKNOWN;
 469   //
 470   aExp.Init(theF, TopAbs_EDGE);
 471   for (; aExp.More(); aExp.Next()) {
 472     const TopoDS_Edge& aSE=(*(TopoDS_Edge*)(&aExp.Current()));
 473     if (BRep_Tool::Degenerated(aSE)) {
 474       continue;
 475     }
 476     //
 477     if (!theBounds.Contains(aSE)) {
 478       const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aSE));
 479       aState=BOPTools_AlgoTools::ComputeState(aE, theRef, theTol, theContext);
 480       return aState;
 481     }
 482     if (aE1.IsNull()) {
 483       aE1=(*(TopoDS_Edge*)(&aSE));
 484     }
 485   }
 486   // !!<- process edges that are all on theRef
 487   if (!aE1.IsNull()) {
 488     BOPTools_AlgoTools3D::PointNearEdge(aE1, theF, aP2D, aP3D, theContext);
 489     aState=BOPTools_AlgoTools::ComputeState(aP3D, theRef, theTol, theContext);
 490   }
 491   //
 492   return aState;
 493 }
 494 //=======================================================================
 495 // function:  ComputeState
 496 // purpose:
 497 //=======================================================================
 498 TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Vertex& theV,
 499                                               const TopoDS_Solid& theRef,
 500                                               const Standard_Real theTol,
 501                                               Handle(BOPInt_Context)& theContext)
 502 {
 503   TopAbs_State aState;
 504   gp_Pnt aP3D;
 505   //
 506   aP3D=BRep_Tool::Pnt(theV);
 507   aState=BOPTools_AlgoTools::ComputeState(aP3D, theRef, theTol, theContext);
 508   return aState;
 509 }
 510 //=======================================================================
 511 // function:  ComputeState
 512 // purpose:
 513 //=======================================================================
 514 TopAbs_State BOPTools_AlgoTools::ComputeState(const TopoDS_Edge& theE,
 515                                               const TopoDS_Solid& theRef,
 516                                               const Standard_Real theTol,
 517                                               Handle(BOPInt_Context)& theContext)
 518 {
 519   Standard_Real aT1, aT2, aT = 0.;
 520   TopAbs_State aState;
 521   Handle(Geom_Curve) aC3D;
 522   gp_Pnt aP3D;
 523   //
 524   aC3D = BRep_Tool::Curve(theE, aT1, aT2);
 525   //
 526   if(aC3D.IsNull()) {
 527     //it means that we are in degenerated edge
 528     const TopoDS_Vertex& aV = TopExp::FirstVertex(theE);
 529     if(aV.IsNull()){
 530       return TopAbs_UNKNOWN;
 531     }
 532     aP3D=BRep_Tool::Pnt(aV);
 533   }
 534   else {//usual case
 535     Standard_Boolean bF2Inf, bL2Inf;
 536     Standard_Real dT=10.;
 537     //
 538     bF2Inf = Precision::IsNegativeInfinite(aT1);
 539     bL2Inf = Precision::IsPositiveInfinite(aT2);
 540     //
 541     if (bF2Inf && !bL2Inf) {
 542       aT=aT2-dT;
 543     }
 544     else if (!bF2Inf && bL2Inf) {
 545       aT=aT1+dT;
 546     }
 547     else if (bF2Inf && bL2Inf) {
 548       aT=0.;
 549     }
 550     else {
 551       aT=IntTools_Tools::IntermediatePoint(aT1, aT2);
 552     }
 553     aC3D->D0(aT, aP3D);
 554   }
 555   //
 556   aState=BOPTools_AlgoTools::ComputeState(aP3D, theRef, theTol, theContext);
 557   //
 558   return aState;
 559 }
 560 //=======================================================================
 561 // function:  ComputeState
 562 // purpose:
 563 //=======================================================================
 564 TopAbs_State BOPTools_AlgoTools::ComputeState(const gp_Pnt& theP,
 565                                               const TopoDS_Solid& theRef,
 566                                               const Standard_Real theTol,
 567                                               Handle(BOPInt_Context)& theContext)
 568 {
 569   TopAbs_State aState;
 570   //
 571   BRepClass3d_SolidClassifier& aSC=theContext->SolidClassifier(theRef);
 572   aSC.Perform(theP, theTol);
 573   //
 574   aState=aSC.State();
 575   //
 576   return aState;
 577 }
 578 //=======================================================================
 579 //function : IsInternalFace
 580 //purpose  :
 581 //=======================================================================
 582 Standard_Integer BOPTools_AlgoTools::IsInternalFace
 583   (const TopoDS_Face& theFace,
 584    const TopoDS_Solid& theSolid,
 585    BOPCol_IndexedDataMapOfShapeListOfShape& theMEF,
 586    const Standard_Real theTol,
 587    Handle(BOPInt_Context)& theContext)
 588 {
 589   Standard_Boolean bDegenerated;
 590   Standard_Integer aNbF, iRet, iFound;
 591   TopAbs_Orientation aOr;
 592   TopoDS_Edge aE1;
 593   TopExp_Explorer aExp;
 594   BOPCol_ListIteratorOfListOfShape aItF;
 595   //
 596   // For all invoked functions: [::IsInternalFace(...)]
 597   // the returned value iRet means:
 598   // iRet=0;  - state is not IN
 599   // iRet=1;  - state is IN
 600   // iRet=2;  - state can not be found by the method of angles
 601   //
 602   // For this function the returned value iRet means:
 603   // iRet=0;  - state is not IN
 604   // iRet=1;  - state is IN
 605   //
 606   iRet=0;
 607   // 1 Try to find an edge from theFace in theMEF
 608   iFound=0;
 609   aExp.Init(theFace, TopAbs_EDGE);
 610   for(; aExp.More(); aExp.Next()) {
 611     const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aExp.Current()));
 612     if (!theMEF.Contains(aE)) {
 613       continue;
 614     }
 615     //
 616     ++iFound;
 617     //
 618     aOr=aE.Orientation();
 619     if (aOr==TopAbs_INTERNAL) {
 620       continue;
 621     }
 622     bDegenerated=BRep_Tool::Degenerated(aE);
 623     if (bDegenerated){
 624       continue;
 625     }
 626     // aE
 627     BOPCol_ListOfShape& aLF=theMEF.ChangeFromKey(aE);
 628     aNbF=aLF.Extent();
 629     if (!aNbF) {
 630       return iRet; // it can not be so
 631     }
 632     //
 633     else if (aNbF==1) {
 634       // aE is internal edge on aLF.First()
 635       const TopoDS_Face& aF1=(*(TopoDS_Face*)(&aLF.First()));
 636       BOPTools_AlgoTools::GetEdgeOnFace(aE, aF1, aE1);
 637       if (aE1.Orientation()!=TopAbs_INTERNAL) {
 638         iRet=2;
 639         break;
 640       }
 641       //
 642       iRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aF1, aF1, theContext);
 643       break;
 644     }
 645     //
 646     else if (aNbF==2) {
 647       const TopoDS_Face& aF1=(*(TopoDS_Face*)(&aLF.First()));
 648       const TopoDS_Face& aF2=(*(TopoDS_Face*)(&aLF.Last()));
 649       //
 650       if (aF2.IsSame(aF1) && BRep_Tool::IsClosed(aE, aF1)) {
 651         // treat as it was for 1 face
 652         iRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aF1, aF2, theContext);
 653         break;
 654       }
 655     }
 656     //
 657     if (aNbF%2) {
 658       iRet=0;
 659       return iRet; // it can not be so
 660     }
 661     else { // aNbF=2,4,6,8,...
 662       iRet=BOPTools_AlgoTools::IsInternalFace(theFace, aE, aLF, theContext);
 663       break;
 664     }
 665   }//for(; aExp.More(); aExp.Next()) {
 666   //
 667   if (!iFound) {
 668     // the face has no shared edges with the solid
 669     iRet=2;
 670   }
 671   //
 672   if (iRet!=2) {
 673     return iRet;
 674   }
 675   //
 676   //========================================
 677   // 2. Classify face using classifier
 678   //
 679   TopAbs_State aState;
 680   BOPCol_IndexedMapOfShape aBounds;
 681   //
 682   BOPTools::MapShapes(theSolid, TopAbs_EDGE, aBounds);
 683   //
 684   aState=BOPTools_AlgoTools::ComputeState(theFace, theSolid, theTol, aBounds, theContext);
 685   //
 686   iRet=(aState==TopAbs_IN)? 1 : 0;
 687   //
 688   return iRet;
 689 }
 690 //=======================================================================
 691 //function : IsInternalFace
 692 //purpose  :
 693 //=======================================================================
 694 Standard_Integer BOPTools_AlgoTools::IsInternalFace(const TopoDS_Face& theFace,
 695                                                     const TopoDS_Edge& theEdge,
 696                                                     BOPCol_ListOfShape& theLF,
 697                                                     Handle(BOPInt_Context)& theContext)
 698 {
 699   Standard_Integer aNbF, iRet;
 700   //
 701   iRet=0;
 702   //
 703   aNbF=theLF.Extent();
 704   if (aNbF==2) {
 705     const TopoDS_Face& aF1=(*(TopoDS_Face*)(&theLF.First()));
 706     const TopoDS_Face& aF2=(*(TopoDS_Face*)(&theLF.Last()));
 707     iRet=BOPTools_AlgoTools::IsInternalFace(theFace, theEdge, aF1, aF2, theContext);
 708     return iRet;
 709   }
 710   //
 711   else {
 712     BOPTools_ListOfCoupleOfShape aLCFF;
 713     BOPTools_ListIteratorOfListOfCoupleOfShape aIt;
 714     //
 715     FindFacePairs(theEdge, theLF, aLCFF, theContext);
 716     //
 717     aIt.Initialize(aLCFF);
 718     for (; aIt.More(); aIt.Next()) {
 719       BOPTools_CoupleOfShape& aCSFF=aIt.ChangeValue();
 720       //
 721       const TopoDS_Face& aF1=(*(TopoDS_Face*)(&aCSFF.Shape1()));
 722       const TopoDS_Face& aF2=(*(TopoDS_Face*)(&aCSFF.Shape2()));
 723       iRet=BOPTools_AlgoTools::IsInternalFace(theFace, theEdge, aF1, aF2, theContext);
 724       if (iRet) {
 725         return iRet;
 726       }
 727     }
 728   }
 729   return iRet;
 730 }
 731 //=======================================================================
 732 //function : IsInternalFace
 733 //purpose  :
 734 //=======================================================================
 735 Standard_Integer BOPTools_AlgoTools::IsInternalFace
 736   (const TopoDS_Face& theFace,
 737    const TopoDS_Edge& theEdge,
 738    const TopoDS_Face& theFace1,
 739    const TopoDS_Face& theFace2,
 740    Handle(BOPInt_Context)& theContext)
 741 {
 742   Standard_Boolean bRet;
 743   Standard_Integer iRet;
 744   TopoDS_Edge aE1, aE2;
 745   TopoDS_Face aFOff;
 746   BOPTools_ListOfCoupleOfShape theLCSOff;
 747   BOPTools_CoupleOfShape aCS1, aCS2;
 748   //
 749   BOPTools_AlgoTools::GetEdgeOnFace(theEdge, theFace1, aE1);
 750   if (aE1.Orientation()==TopAbs_INTERNAL) {
 751     aE2=aE1;
 752     aE1.Orientation(TopAbs_FORWARD);
 753     aE2.Orientation(TopAbs_REVERSED);
 754   }
 755   else if (theFace1==theFace2) {
 756     aE2=aE1;
 757     aE1.Orientation(TopAbs_FORWARD);
 758     aE2.Orientation(TopAbs_REVERSED);
 759   }
 760   else {
 761     BOPTools_AlgoTools::GetEdgeOnFace(theEdge, theFace2, aE2);
 762   }
 763   //
 764   aCS1.SetShape1(theEdge);
 765   aCS1.SetShape2(theFace);
 766   theLCSOff.Append(aCS1);
 767   //
 768   aCS2.SetShape1(aE2);
 769   aCS2.SetShape2(theFace2);
 770   theLCSOff.Append(aCS2);
 771   //
 772   bRet=GetFaceOff(aE1, theFace1, theLCSOff, aFOff, theContext);
 773   //
 774   iRet=0; // theFace is not internal
 775   if (theFace.IsEqual(aFOff)) {
 776     // theFace is internal
 777     iRet=1;
 778     if (!bRet) {
 779       // theFace seems to be internal
 780       iRet=2;
 781     }
 782   }
 783   return iRet;
 784 }
 785 //=======================================================================
 786 //function : GetFaceOff
 787 //purpose  :
 788 //=======================================================================
 789 Standard_Boolean BOPTools_AlgoTools::GetFaceOff
 790   (const TopoDS_Edge& theE1,
 791    const TopoDS_Face& theF1,
 792    BOPTools_ListOfCoupleOfShape& theLCSOff,
 793    TopoDS_Face& theFOff,
 794    Handle(BOPInt_Context)& theContext)
 795 {
 796   Standard_Boolean bRet;
 797   Standard_Real aT, aT1, aT2, aAngle, aTwoPI, aAngleMin;
 798   Standard_Real aUmin, aUsup, aVmin, aVsup;
 799   gp_Pnt aPn1, aPn2, aPx;
 800   gp_Dir aDN1, aDN2, aDBF, aDBF2, aDTF;
 801   gp_Vec aVTgt;
 802   TopAbs_Orientation aOr;
 803   Handle(Geom_Curve)aC3D;
 804   Handle(Geom_Plane) aPL;
 805   BOPTools_ListIteratorOfListOfCoupleOfShape aIt;
 806   GeomAPI_ProjectPointOnSurf aProjPL;
 807   //
 808   aAngleMin=100.;
 809   aTwoPI=M_PI+M_PI;
 810   aC3D =BRep_Tool::Curve(theE1, aT1, aT2);
 811   aT=BOPTools_AlgoTools2D::IntermediatePoint(aT1, aT2);
 812   aC3D->D0(aT, aPx);
 813   //
 814   BOPTools_AlgoTools2D::EdgeTangent(theE1, aT, aVTgt);
 815   gp_Dir aDTgt(aVTgt), aDTgt2;
 816   aOr = theE1.Orientation();
 817   //
 818   aPL = new Geom_Plane(aPx, aDTgt);
 819   aPL->Bounds(aUmin, aUsup, aVmin, aVsup);
 820   aProjPL.Init(aPL, aUmin, aUsup, aVmin, aVsup);
 821   //
 822   GetFaceDir(theE1, theF1, aPx, aT, aDTgt, aDN1, aDBF, theContext, aProjPL);
 823   //
 824   aDTF=aDN1^aDBF;
 825   //
 826   bRet=Standard_True;
 827   aIt.Initialize(theLCSOff);
 828   for (; aIt.More(); aIt.Next()) {
 829     const BOPTools_CoupleOfShape& aCS=aIt.Value();
 830     const TopoDS_Edge& aE2=(*(TopoDS_Edge*)(&aCS.Shape1()));
 831     const TopoDS_Face& aF2=(*(TopoDS_Face*)(&aCS.Shape2()));
 832     //
 833     aDTgt2 = (aE2.Orientation()==aOr) ? aDTgt : aDTgt.Reversed();
 834     GetFaceDir(aE2, aF2, aPx, aT, aDTgt2, aDN2, aDBF2, theContext, aProjPL);
 835     //Angle
 836     aAngle=AngleWithRef(aDBF, aDBF2, aDTF);
 837     //
 838     if(aAngle<0.) {
 839       aAngle=aTwoPI+aAngle;
 840     }
 841     //
 842     if (aAngle<Precision::Angular()) {
 843       if (aF2==theF1) {
 844         aAngle=M_PI;
 845       }
 846       else if (aF2.IsSame(theF1)) {
 847         aAngle=aTwoPI;
 848       }
 849     }
 850     //
 851     if (aAngle<aAngleMin){
 852       aAngleMin=aAngle;
 853       theFOff=aF2;
 854     }
 855     else if (aAngle==aAngleMin) {
 856       // the minimal angle can not be found
 857       bRet=Standard_False;
 858     }
 859   }
 860   return bRet;
 861 }
 862 //=======================================================================
 863 //function : GetEdgeOff
 864 //purpose  :
 865 //=======================================================================
 866 Standard_Boolean BOPTools_AlgoTools::GetEdgeOff(const TopoDS_Edge& theE1,
 867                                                 const TopoDS_Face& theF2,
 868                                                 TopoDS_Edge& theE2)
 869 {
 870   Standard_Boolean bFound;
 871   TopAbs_Orientation aOr1, aOr1C, aOr2;
 872   TopExp_Explorer anExp;
 873   //
 874   bFound=Standard_False;
 875   aOr1=theE1.Orientation();
 876   aOr1C=TopAbs::Reverse(aOr1);
 877   //
 878   anExp.Init(theF2, TopAbs_EDGE);
 879   for (; anExp.More(); anExp.Next()) {
 880     const TopoDS_Edge& aEF2=(*(TopoDS_Edge*)(&anExp.Current()));
 881     if (aEF2.IsSame(theE1)) {
 882       aOr2=aEF2.Orientation();
 883       if (aOr2==aOr1C) {
 884         theE2=aEF2;
 885         bFound=!bFound;
 886         return bFound;
 887       }
 888     }
 889   }
 890   return bFound;
 891 }
 892
 893 //=======================================================================
 894 //function : AreFacesSameDomain
 895 //purpose  :
 896 //=======================================================================
 897 Standard_Boolean BOPTools_AlgoTools::AreFacesSameDomain(const TopoDS_Face& theF1,
 898                                                         const TopoDS_Face& theF2,
 899                                                         Handle(BOPInt_Context)& theContext)
 900 {
 901   Standard_Boolean bFlag;
 902   Standard_Integer iErr;
 903   Standard_Real aTolF1, aTolF2, aTol;
 904   gp_Pnt2d aP2D;
 905   gp_Pnt aP;
 906   TopoDS_Face aF1, aF2;
 907   TopoDS_Edge aE1;
 908   TopExp_Explorer aExp;
 909   //
 910   bFlag=Standard_False;
 911   //
 912   aF1=theF1;
 913   aF1.Orientation(TopAbs_FORWARD);
 914   aF2=theF2;
 915   aF2.Orientation(TopAbs_FORWARD);
 916   //
 917   aTolF1=BRep_Tool::Tolerance(aF1);
 918   // 1
 919   aExp.Init(aF1, TopAbs_EDGE);
 920   for (; aExp.More(); aExp.Next()) {
 921     aE1=(*(TopoDS_Edge*)(&aExp.Current()));
 922     if (!BRep_Tool::Degenerated(aE1)) {
 923       Standard_Real aTolE = BRep_Tool::Tolerance(aE1);
 924       aTolF1 = (aTolE > aTolF1) ? aTolE : aTolF1;
 925     }
 926   }
 927   // 2
 928   aTolF2=BRep_Tool::Tolerance(aF2);
 929   aTol=aTolF1+aTolF2;
 930   //
 931   iErr = BOPTools_AlgoTools3D::PointInFace(aF1, aP, aP2D, theContext);
 932   if (!iErr) {
 933     bFlag=theContext->IsValidPointForFace(aP, aF2, aTol);
 934   }
 935   //
 936   return bFlag;
 937 }
 938
 939 //=======================================================================
 940 //function : CheckSameGeom
 941 //purpose  :
 942 //=======================================================================
 943 Standard_Boolean BOPTools_AlgoTools::CheckSameGeom(const TopoDS_Face& theF1,
 944                                                    const TopoDS_Face& theF2,
 945                                                    Handle(BOPInt_Context)& theContext)
 946 {
 947   Standard_Boolean bRet;
 948   Standard_Real aTolF1, aTolF2, aTol;
 949   gp_Pnt2d aP2D;
 950   gp_Pnt aP;
 951   TopExp_Explorer aExp;
 952   //
 953   bRet=Standard_False;
 954   aExp.Init(theF1, TopAbs_EDGE);
 955   for (; aExp.More(); aExp.Next()) {
 956     const TopoDS_Edge& aE=(*(TopoDS_Edge*)(&aExp.Current()));
 957     if (!BRep_Tool::Degenerated(aE)) {
 958       aTolF1=BRep_Tool::Tolerance(theF1);
 959       aTolF2=BRep_Tool::Tolerance(theF2);
 960       aTol=aTolF1+aTolF2;
 961       BOPTools_AlgoTools3D::PointNearEdge(aE, theF1, aP2D, aP, theContext);
 962       bRet=theContext->IsValidPointForFace(aP, theF2, aTol);
 963       break;
 964     }
 965   }
 966   return bRet;
 967 }
 968 //=======================================================================
 969 // function: Sense
 970 // purpose:
 971 //=======================================================================
 972 Standard_Integer BOPTools_AlgoTools::Sense (const TopoDS_Face& theF1,
 973                                             const TopoDS_Face& theF2)
 974 {
 975   Standard_Integer iSense=0;
 976   gp_Dir aDNF1, aDNF2;
 977   TopoDS_Edge aE1, aE2;
 978   TopExp_Explorer aExp;
 979   //
 980   aExp.Init(theF1, TopAbs_EDGE);
 981   for (; aExp.More(); aExp.Next()) {
 982     aE1=(*(TopoDS_Edge*)(&aExp.Current()));
 983     if (!BRep_Tool::Degenerated(aE1)) {
 984       if (!BRep_Tool::IsClosed(aE1, theF1)) {
 985         break;
 986       }
 987     }
 988   }
 989   //
 990   aExp.Init(theF2, TopAbs_EDGE);
 991   for (; aExp.More(); aExp.Next()) {
 992     aE2=(*(TopoDS_Edge*)(&aExp.Current()));
 993     if (!BRep_Tool::Degenerated(aE2)) {
 994       if (!BRep_Tool::IsClosed(aE2, theF2)) {
 995         if (aE2.IsSame(aE1)) {
 996           iSense=1;
 997           break;
 998         }
 999       }
1000     }
1001   }
1002   //
1003   if (!iSense) {
1004     return iSense;
1005   }
1006   //
1007   BOPTools_AlgoTools3D::GetNormalToFaceOnEdge(aE1, theF1, aDNF1);
1008   BOPTools_AlgoTools3D::GetNormalToFaceOnEdge(aE2, theF2, aDNF2);
1009   //
1010   iSense=BOPTools_AlgoTools3D::SenseFlag(aDNF1, aDNF2);
1011   //
1012   return iSense;
1013 }
1014 //=======================================================================
1015 // function: IsSplitToReverse
1016 // purpose:
1017 //=======================================================================
1018 Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Shape& theSp,
1019                                                       const TopoDS_Shape& theSr,
1020                                                       Handle(BOPInt_Context)& theContext)
1021 {
1022   Standard_Boolean bRet;
1023   TopAbs_ShapeEnum aType;
1024   //
1025   bRet=Standard_False;
1026   //
1027   aType=theSp.ShapeType();
1028   switch (aType) {
1029     case TopAbs_EDGE: {
1030       const TopoDS_Edge& aESp=(*(TopoDS_Edge*)(&theSp));
1031       const TopoDS_Edge& aESr=(*(TopoDS_Edge*)(&theSr));
1032       bRet=BOPTools_AlgoTools::IsSplitToReverse(aESp, aESr, theContext);
1033     }
1034       break;
1035       //
1036     case TopAbs_FACE: {
1037       const TopoDS_Face& aFSp=(*(TopoDS_Face*)(&theSp));
1038       const TopoDS_Face& aFSr=(*(TopoDS_Face*)(&theSr));
1039       bRet=BOPTools_AlgoTools::IsSplitToReverse(aFSp, aFSr, theContext);
1040     }
1041       break;
1042       //
1043     default:
1044       break;
1045   }
1046   return bRet;
1047 }
1048 //=======================================================================
1049 //function :IsSplitToReverse
1050 //purpose  :
1051 //=======================================================================
1052 Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Face& theFSp,
1053                                                       const TopoDS_Face& theFSr,
1054                                                       Handle(BOPInt_Context)& theContext)
1055 {
1056   Standard_Boolean bRet, bFound, bInFace;
1057   Standard_Real aT1, aT2, aT, aU, aV, aScPr;
1058   gp_Pnt aPFSp, aPFSr;
1059   gp_Dir aDNFSp;
1060   gp_Vec aD1U, aD1V;
1061   Handle(Geom_Surface) aSr, aSp;
1062   TopAbs_Orientation aOrSr, aOrSp;
1063   TopExp_Explorer anExp;
1064   TopoDS_Edge aESp;
1065   //
1066   bRet=Standard_False;
1067   //
1068   aSr=BRep_Tool::Surface(theFSr);
1069   aSp=BRep_Tool::Surface(theFSp);
1070   if (aSr==aSp) {
1071     aOrSr=theFSr.Orientation();
1072     aOrSp=theFSp.Orientation();
1073     bRet=(aOrSr!=aOrSp);
1074     return bRet;
1075   }
1076   //
1077   bFound=Standard_False;
1078   anExp.Init(theFSp, TopAbs_EDGE);
1079   for (; anExp.More(); anExp.Next()) {
1080     aESp=(*(TopoDS_Edge*)(&anExp.Current()));
1081     if (!BRep_Tool::Degenerated(aESp)) {
1082       if (!BRep_Tool::IsClosed(aESp, theFSp)) {
1083         bFound=!bFound;
1084         break;
1085       }
1086     }
1087   }
1088   if (!bFound) {
1089     Standard_Boolean bFlag;
1090     Standard_Integer iErr;
1091     gp_Pnt2d aP2DFSp;
1092     //
1093     iErr=BOPTools_AlgoTools3D::PointInFace(theFSp, aPFSp, aP2DFSp, theContext);
1094     if (iErr) {
1095       return bRet;
1096     }
1097     //
1098     aP2DFSp.Coord(aU, aV);
1099     bFlag=BOPTools_AlgoTools3D::GetNormalToSurface(aSp, aU, aV, aDNFSp);
1100     if (!bFlag) {
1101       return bRet;
1102     }
1103   }
1104   else {
1105     BRep_Tool::Range(aESp, aT1, aT2);
1106     aT=BOPTools_AlgoTools2D::IntermediatePoint(aT1, aT2);
1107     BOPTools_AlgoTools3D::GetApproxNormalToFaceOnEdge(aESp, theFSp, aT, aPFSp, aDNFSp, theContext);
1108   }
1109   //
1110   // Parts of theContext->ComputeVS(..)
1111   GeomAPI_ProjectPointOnSurf& aProjector=theContext->ProjPS(theFSr);
1112   aProjector.Perform(aPFSp);
1113   if (!aProjector.IsDone()) {
1114     return bRet;
1115   }
1116   //
1117   aProjector.LowerDistanceParameters(aU, aV);
1118   gp_Pnt2d aP2D(aU, aV);
1119   bInFace=theContext->IsPointInFace (theFSr, aP2D);
1120   if (!bInFace) {
1121     return bRet;
1122   }
1123   //
1124   aSr->D1(aU, aV, aPFSr, aD1U, aD1V);
1125   gp_Dir aDD1U(aD1U);
1126   gp_Dir aDD1V(aD1V);
1127   gp_Dir aDNFSr=aDD1U^aDD1V;
1128   if (theFSr.Orientation()==TopAbs_REVERSED){
1129     aDNFSr.Reverse();
1130   }
1131   //
1132   aScPr=aDNFSp*aDNFSr;
1133   bRet=(aScPr<0.);
1134   //
1135   return bRet;
1136 }
1137 //=======================================================================
1138 //function :IsSplitToReverse
1139 //purpose  :
1140 //=======================================================================
1141 Standard_Boolean BOPTools_AlgoTools::IsSplitToReverse(const TopoDS_Edge& aEF1,
1142                                                       const TopoDS_Edge& aEF2,
1143                                                       Handle(BOPInt_Context)& theContext)
1144 {
1145   Standard_Boolean bRet, bIsDegenerated;
1146   //
1147   bRet=Standard_False;
1148   bIsDegenerated=(BRep_Tool::Degenerated(aEF1) ||
1149                   BRep_Tool::Degenerated(aEF2));
1150   if (bIsDegenerated) {
1151     return bRet;
1152   }
1153   //
1154   Standard_Real a, b;
1155   TopAbs_Orientation aOrE, aOrSp;
1156   Handle(Geom_Curve)aC1, aC2;
1157   //
1158   aC2=BRep_Tool::Curve(aEF2, a, b);
1159   aC1=BRep_Tool::Curve(aEF1, a, b);
1160   //
1161   if (aC1==aC2) {
1162     aOrE=aEF2.Orientation();
1163     aOrSp=aEF1.Orientation();
1164     bRet=(aOrE!=aOrSp);
1165     return bRet;
1166   }
1167   //
1168   Standard_Real aT1, aT2, aScPr;
1169   gp_Vec aV1, aV2;
1170   gp_Pnt aP;
1171   //
1172   aT1=BOPTools_AlgoTools2D::IntermediatePoint(a, b);
1173   aC1->D0(aT1, aP);
1174   BOPTools_AlgoTools2D::EdgeTangent(aEF1, aT1, aV1);
1175   gp_Dir aDT1(aV1);
1176   //
1177   theContext->ProjectPointOnEdge(aP, aEF2, aT2);
1178   //
1179   BOPTools_AlgoTools2D::EdgeTangent(aEF2, aT2, aV2);
1180   gp_Dir aDT2(aV2);
1181   //
1182   aScPr=aDT1*aDT2;
1183   bRet=(aScPr<0.);
1184   //
1185   return bRet;
1186 }
1187
1188 //=======================================================================
1189 //function : IsHole
1190 //purpose  :
1191 //=======================================================================
1192 Standard_Boolean BOPTools_AlgoTools::IsHole(const TopoDS_Shape& aW,
1193                                             const TopoDS_Shape& aFace)
1194 {
1195   Standard_Boolean bIsHole;
1196   Standard_Integer i, aNbS;
1197   Standard_Real aT1, aT2, aS;
1198   Standard_Real aU1, aU2, aU, dU;
1199   Standard_Real aX1, aY1, aX0, aY0;
1200   TopAbs_Orientation aOr;
1201
1202   gp_Pnt2d aP2D0, aP2D1;
1203   Handle(Geom2d_Curve) aC2D;
1204   TopoDS_Face aF, aFF;
1205   TopoDS_Iterator aItW;
1206   //
1207   bIsHole=Standard_False;
1208   //
1209   aF=(*(TopoDS_Face *)(&aFace));
1210   aFF=aF;
1211   aFF.Orientation(TopAbs_FORWARD);
1212   //
1213   aS=0.;
1214   aItW.Initialize(aW);
1215   for (; aItW.More(); aItW.Next()) {
1216     const TopoDS_Edge& aE=(*(TopoDS_Edge *)(&aItW.Value()));
1217     aOr=aE.Orientation();
1218     if (!(aOr==TopAbs_FORWARD ||
1219           aOr==TopAbs_REVERSED)) {
1220       continue;
1221     }
1222     //
1223     aC2D=BRep_Tool::CurveOnSurface(aE, aFF, aT1, aT2);
1224     if (aC2D.IsNull()) {
1225       break; //xx
1226     }
1227     //
1228     BRepAdaptor_Curve2d aBAC2D(aE, aFF);
1229     aNbS=Geom2dInt_Geom2dCurveTool::NbSamples(aBAC2D);
1230     if (aNbS>2) {
1231       aNbS*=4;
1232     }
1233     //
1234     dU=(aT2-aT1)/(Standard_Real)(aNbS-1);
1235     aU =aT1;
1236     aU1=aT1;
1237     aU2=aT2;
1238     if (aOr==TopAbs_REVERSED) {
1239       aU =aT2;
1240       aU1=aT2;
1241       aU2=aT1;
1242       dU=-dU;
1243     }
1244     //
1245     aC2D->D0(aU, aP2D0);
1246     for(i=2; i<=aNbS; i++) {
1247       aU=aU1+(i-1)*dU;
1248       aC2D->D0(aU, aP2D1);
1249       aP2D0.Coord(aX0, aY0);
1250       aP2D1.Coord(aX1, aY1);
1251       //
1252       aS=aS+(aY0+aY1)*(aX1-aX0);
1253       //
1254       aP2D0=aP2D1;
1255     }
1256   }//for (; aItW.More(); aItW.Next()) {
1257   bIsHole=(aS>0.);
1258   return bIsHole;
1259 }
1260
1261 //=======================================================================
1262 // function: MakeContainer
1263 // purpose:
1264 //=======================================================================
1265 void BOPTools_AlgoTools::MakeContainer(const TopAbs_ShapeEnum theType,
1266                                        TopoDS_Shape& theC)
1267 {
1268   BRep_Builder aBB;
1269   //
1270   switch(theType) {
1271     case TopAbs_COMPOUND:{
1272       TopoDS_Compound aC;
1273       aBB.MakeCompound(aC);
1274       theC=aC;
1275     }
1276       break;
1277       //
1278     case TopAbs_COMPSOLID:{
1279       TopoDS_CompSolid aCS;
1280       aBB.MakeCompSolid(aCS);
1281       theC=aCS;
1282     }
1283       break;
1284       //
1285     case TopAbs_SOLID:{
1286       TopoDS_Solid aSolid;
1287       aBB.MakeSolid(aSolid);
1288       theC=aSolid;
1289     }
1290       break;
1291       //
1292       //
1293     case TopAbs_SHELL:{
1294       TopoDS_Shell aShell;
1295       aBB.MakeShell(aShell);
1296       theC=aShell;
1297     }
1298       break;
1299       //
1300     case TopAbs_WIRE: {
1301       TopoDS_Wire aWire;
1302       aBB.MakeWire(aWire);
1303       theC=aWire;
1304     }
1305       break;
1306       //
1307     default:
1308       break;
1309   }
1310 }
1311 //=======================================================================
1312 // function: MakePCurve
1313 // purpose:
1314 //=======================================================================
1315 void BOPTools_AlgoTools::MakePCurve(const TopoDS_Edge& aE,
1316                                     const TopoDS_Face& aF1,
1317                                     const TopoDS_Face& aF2,
1318                                     const IntTools_Curve& aIC,
1319                                     const Standard_Boolean bPC1,
1320                                     const Standard_Boolean bPC2)
1321
1322 {
1323   Standard_Integer i;
1324   Standard_Real aTolE, aT1, aT2, aOutFirst, aOutLast, aOutTol;
1325   Handle(Geom2d_Curve) aC2D, aC2DA, aC2Dx1;
1326   TopoDS_Face aFFWD;
1327   BRep_Builder aBB;
1328   Standard_Boolean bPC;
1329   //
1330   aTolE=BRep_Tool::Tolerance(aE);
1331   //
1332   const Handle(Geom_Curve)& aC3DE=BRep_Tool::Curve(aE, aT1, aT2);
1333   Handle(Geom_TrimmedCurve)aC3DETrim=new Geom_TrimmedCurve(aC3DE, aT1, aT2);
1334   //
1335   for (i=0; i<2; ++i) {
1336     bPC = !i ? bPC1 : bPC2;
1337     if (!bPC) {
1338       continue;
1339     }
1340     //
1341     if (!i) {
1342       aFFWD=aF1;
1343       aC2Dx1=aIC.FirstCurve2d();
1344     }
1345     else {
1346       aFFWD=aF2;
1347       aC2Dx1=aIC.SecondCurve2d();
1348     }
1349     //
1350     aFFWD.Orientation(TopAbs_FORWARD);
1351     //
1352     aC2D=aC2Dx1;
1353     if (aC2D.IsNull()) {
1354       BOPTools_AlgoTools2D::BuildPCurveForEdgeOnFace(aE, aFFWD);
1355       BOPTools_AlgoTools2D::CurveOnSurface(aE, aFFWD, aC2D,
1356                                        aOutFirst, aOutLast,
1357                                        aOutTol);
1358       }
1359     //
1360     if (aC3DE->IsPeriodic()) {
1361       BOPTools_AlgoTools2D::AdjustPCurveOnFace(aFFWD, aT1, aT2,  aC2D, aC2DA);
1362     }
1363     else {
1364       BOPTools_AlgoTools2D::AdjustPCurveOnFace(aFFWD, aC3DETrim, aC2D, aC2DA);
1365     }
1366     //
1367     aBB.UpdateEdge(aE, aC2DA, aFFWD, aTolE);
1368     //BRepLib::SameParameter(aE);
1369   }
1370   BRepLib::SameParameter(aE);
1371 }
1372 //=======================================================================
1373 // function: MakeEdge
1374 // purpose:
1375 //=======================================================================
1376 void BOPTools_AlgoTools::MakeEdge(const IntTools_Curve& theIC,
1377                                   const TopoDS_Vertex& theV1,
1378                                   const Standard_Real theT1,
1379                                   const TopoDS_Vertex& theV2,
1380                                   const Standard_Real theT2,
1381                                   const Standard_Real theTolR3D,
1382                                   TopoDS_Edge& theE)
1383 {
1384   Standard_Real aTolV;
1385   BRep_Builder aBB;
1386   //
1387   BOPTools_AlgoTools::MakeSectEdge (theIC, theV1, theT1, theV2, theT2, theE);
1388   //
1389   aBB.UpdateEdge(theE, theTolR3D);
1390   //
1391   aTolV=BRep_Tool::Tolerance(theV1);
1392   if (aTolV<theTolR3D) {
1393     aBB.UpdateVertex(theV1, theTolR3D);
1394   }
1395   //
1396   aTolV=BRep_Tool::Tolerance(theV2);
1397   if (aTolV<theTolR3D) {
1398     aBB.UpdateVertex(theV2, theTolR3D);
1399   }
1400 }
1401 //=======================================================================
1402 // function: ComputeVV
1403 // purpose:
1404 //=======================================================================
1405 Standard_Integer BOPTools_AlgoTools::ComputeVV(const TopoDS_Vertex& aV1,
1406                                                const gp_Pnt& aP2,
1407                                                const Standard_Real aTolP2)
1408 {
1409   Standard_Real aTolV1, aTolSum, aTolSum2, aD2;
1410   gp_Pnt aP1;
1411   //
1412   aTolV1=BRep_Tool::Tolerance(aV1);
1413
1414   aTolSum=aTolV1+aTolP2;
1415   aTolSum2=aTolSum*aTolSum;
1416   //
1417   aP1=BRep_Tool::Pnt(aV1);
1418   //
1419   aD2=aP1.SquareDistance(aP2);
1420   if (aD2>aTolSum2) {
1421     return 1;
1422   }
1423   return 0;
1424 }
1425 //=======================================================================
1426 // function: ComputeVV
1427 // purpose:
1428 //=======================================================================
1429 Standard_Integer BOPTools_AlgoTools::ComputeVV(const TopoDS_Vertex& aV1,
1430                                                const TopoDS_Vertex& aV2)
1431 {
1432   Standard_Real aTolV1, aTolV2, aTolSum, aTolSum2, aD2;
1433   gp_Pnt aP1, aP2;
1434   //
1435   aTolV1=BRep_Tool::Tolerance(aV1);
1436   aTolV2=BRep_Tool::Tolerance(aV2);
1437   aTolSum=aTolV1+aTolV2;
1438   aTolSum2=aTolSum*aTolSum;
1439   //
1440   aP1=BRep_Tool::Pnt(aV1);
1441   aP2=BRep_Tool::Pnt(aV2);
1442   //
1443   aD2=aP1.SquareDistance(aP2);
1444   if (aD2>aTolSum2) {
1445     return 1;
1446   }
1447   return 0;
1448 }
1449 //=======================================================================
1450 // function: MakeVertex
1451 // purpose :
1452 //=======================================================================
1453 void BOPTools_AlgoTools::MakeVertex(BOPCol_ListOfShape& aLV,
1454                                     TopoDS_Vertex& aVnew)
1455 {
1456   Standard_Integer aNb;
1457   Standard_Real aTi, aDi, aDmax;
1458   gp_Pnt aPi, aP;
1459   gp_XYZ aXYZ(0.,0.,0.), aXYZi;
1460   BOPCol_ListIteratorOfListOfShape aIt;
1461   //
1462   aNb=aLV.Extent();
1463   if (aNb) {
1464     aIt.Initialize(aLV);
1465     for (; aIt.More(); aIt.Next()) {
1466       TopoDS_Vertex& aVi=*((TopoDS_Vertex*)(&aIt.Value()));
1467       aPi=BRep_Tool::Pnt(aVi);
1468       aXYZi=aPi.XYZ();
1469       aXYZ=aXYZ+aXYZi;
1470     }
1471     //
1472     aXYZ.Divide((Standard_Real)aNb);
1473     aP.SetXYZ(aXYZ);
1474     //
1475     aDmax=-1.;
1476     aIt.Initialize(aLV);
1477     for (; aIt.More(); aIt.Next()) {
1478       TopoDS_Vertex& aVi=*((TopoDS_Vertex*)(&aIt.Value()));
1479       aPi=BRep_Tool::Pnt(aVi);
1480       aTi=BRep_Tool::Tolerance(aVi);
1481       aDi=aP.SquareDistance(aPi);
1482       aDi=fsqrt(aDi);
1483       aDi=aDi+aTi;
1484       if (aDi > aDmax) {
1485         aDmax=aDi;
1486       }
1487     }
1488     //
1489     BRep_Builder aBB;
1490     aBB.MakeVertex (aVnew, aP, aDmax);
1491   }
1492 }
1493 //=======================================================================
1494 //function : GetEdgeOnFace
1495 //purpose  :
1496 //=======================================================================
1497 Standard_Boolean BOPTools_AlgoTools::GetEdgeOnFace(const TopoDS_Edge& theE1,
1498                                                    const TopoDS_Face& theF2,
1499                                                    TopoDS_Edge& theE2)
1500 {
1501   Standard_Boolean bFound;
1502   TopoDS_Iterator aItF, aItW;
1503   //
1504   bFound=Standard_False;
1505   //
1506   aItF.Initialize(theF2);
1507   for (; aItF.More(); aItF.Next()) {
1508     const TopoDS_Shape& aW=aItF.Value();
1509     aItW.Initialize(aW);
1510     for (; aItW.More(); aItW.Next()) {
1511       const TopoDS_Shape& aE=aItW.Value();
1512       if (aE.IsSame(theE1)) {
1513         theE2=(*(TopoDS_Edge*)(&aE));
1514         bFound=!bFound;
1515         return bFound;
1516       }
1517     }
1518   }
1519   return bFound;
1520 }
1521 //=======================================================================
1522 //function : FindFacePairs
1523 //purpose  :
1524 //=======================================================================
1525 Standard_Boolean FindFacePairs (const TopoDS_Edge& theE,
1526                                 const BOPCol_ListOfShape& thLF,
1527                                 BOPTools_ListOfCoupleOfShape& theLCFF,
1528                                 Handle(BOPInt_Context)& theContext)
1529 {
1530   Standard_Boolean bFound;
1531   Standard_Integer i, aNbCEF;
1532   TopAbs_Orientation aOr, aOrC = TopAbs_FORWARD;
1533   BOPCol_MapOfShape aMFP;
1534   TopoDS_Face aF1, aF2;
1535   TopoDS_Edge aEL, aE1;
1536   BOPCol_ListIteratorOfListOfShape aItLF;
1537   BOPTools_CoupleOfShape aCEF, aCFF;
1538   BOPTools_ListOfCoupleOfShape aLCEF, aLCEFx;
1539   BOPTools_ListIteratorOfListOfCoupleOfShape aIt;
1540   //
1541   bFound=Standard_True;
1542   //
1543   // Preface aLCEF
1544   aItLF.Initialize(thLF);
1545   for (; aItLF.More(); aItLF.Next()) {
1546     const TopoDS_Face& aFL=(*(TopoDS_Face*)(&aItLF.Value()));
1547     //
1548     bFound=BOPTools_AlgoTools::GetEdgeOnFace(theE, aFL, aEL);
1549     if (!bFound) {
1550       return bFound; // it can not be so
1551     }
1552     //
1553     aCEF.SetShape1(aEL);
1554     aCEF.SetShape2(aFL);
1555     aLCEF.Append(aCEF);
1556   }
1557   //
1558   aNbCEF=aLCEF.Extent();
1559   while(aNbCEF) {
1560     //
1561     // aLCEFx
1562     aLCEFx.Clear();
1563     aIt.Initialize(aLCEF);
1564     for (i=0; aIt.More(); aIt.Next(), ++i) {
1565       const BOPTools_CoupleOfShape& aCSx=aIt.Value();
1566       const TopoDS_Shape& aEx=aCSx.Shape1();
1567       const TopoDS_Shape& aFx=aCSx.Shape2();
1568       //
1569       aOr=aEx.Orientation();
1570       //
1571       if (!i) {
1572         aOrC=TopAbs::Reverse(aOr);
1573         aE1=(*(TopoDS_Edge*)(&aEx));
1574         aF1=(*(TopoDS_Face*)(&aFx));
1575         aMFP.Add(aFx);
1576         continue;
1577       }
1578       //
1579       if (aOr==aOrC) {
1580         aLCEFx.Append(aCSx);
1581         aMFP.Add(aFx);
1582       }
1583     }
1584     //
1585     // F2
1586     BOPTools_AlgoTools::GetFaceOff(aE1, aF1, aLCEFx, aF2, theContext);
1587     //
1588     aCFF.SetShape1(aF1);
1589     aCFF.SetShape2(aF2);
1590     theLCFF.Append(aCFF);
1591     //
1592     aMFP.Add(aF1);
1593     aMFP.Add(aF2);
1594     //
1595     // refine aLCEF
1596     aLCEFx.Clear();
1597     aLCEFx=aLCEF;
1598     aLCEF.Clear();
1599     aIt.Initialize(aLCEFx);
1600     for (; aIt.More(); aIt.Next()) {
1601       const BOPTools_CoupleOfShape& aCSx=aIt.Value();
1602       const TopoDS_Shape& aFx=aCSx.Shape2();
1603       if (!aMFP.Contains(aFx)) {
1604         aLCEF.Append(aCSx);
1605       }
1606     }
1607     //
1608     aNbCEF=aLCEF.Extent();
1609   }//while(aNbCEF) {
1610   //
1611   return bFound;
1612 }
1613 //=======================================================================
1614 //function : AngleWithRef
1615 //purpose  :
1616 //=======================================================================
1617 Standard_Real AngleWithRef(const gp_Dir& theD1,
1618                            const gp_Dir& theD2,
1619                            const gp_Dir& theDRef)
1620 {
1621   Standard_Real aCosinus, aSinus, aBeta, aHalfPI, aScPr;
1622   gp_XYZ aXYZ;
1623   //
1624   aHalfPI=0.5*M_PI;
1625   //
1626   const gp_XYZ& aXYZ1=theD1.XYZ();
1627   const gp_XYZ& aXYZ2=theD2.XYZ();
1628   aXYZ=aXYZ1.Crossed(aXYZ2);
1629   aSinus=aXYZ.Modulus();
1630   aCosinus=theD1*theD2;
1631   //
1632   aBeta=0.;
1633   if (aSinus>=0.) {
1634     aBeta=aHalfPI*(1.-aCosinus);
1635   }
1636   else {
1637     aBeta=2.*M_PI-aHalfPI*(3.+aCosinus);
1638   }
1639   //
1640   aScPr=aXYZ.Dot(theDRef.XYZ());
1641   if (aScPr<0.) {
1642     aBeta=-aBeta;
1643   }
1644   return aBeta;
1645 }
1646 //=======================================================================
1647 //function : fsqrt
1648 //purpose  :
1649 //=======================================================================
1650 Standard_Real fsqrt(Standard_Real val)
1651 {
1652   union {
1653     int tmp;
1654     float val;
1655   } u;
1656   //
1657   u.val = (float)val;
1658   u.tmp -= 1<<23; /* Remove last bit so 1.0 gives 1.0 */
1659   /* tmp is now an approximation to logbase2(val) */
1660   u.tmp >>= 1; /* divide by 2 */
1661   u.tmp += 1<<29; /* add 64 to exponent: (e+127)/2 =(e/2)+63, */
1662   /* that represents (e/2)-64 but we want e/2 */
1663   return (double)u.val;
1664 }
1665
1666 //=======================================================================
1667 // function: IsBlockInOnFace
1668 // purpose:
1669 //=======================================================================
1670 Standard_Boolean BOPTools_AlgoTools::IsBlockInOnFace (const IntTools_Range& aShrR,
1671                                                       const TopoDS_Face& aF,
1672                                                       const TopoDS_Edge& aE1,
1673                                                       Handle(BOPInt_Context)& aContext)
1674 {
1675   Standard_Boolean bFlag;
1676   Standard_Real f1, l1, ULD, VLD;
1677   gp_Pnt2d aP2D;
1678   gp_Pnt aP11, aP12;
1679   //
1680   aShrR.Range(f1, l1);
1681   Standard_Real dt=0.0075,  k;//dt=0.001,  k;
1682   k=dt*(l1-f1);
1683   f1=f1+k;
1684   l1=l1-k;
1685   //
1686   // Treatment P11
1687   BOPTools_AlgoTools::PointOnEdge(aE1, f1, aP11);
1688   //
1689   GeomAPI_ProjectPointOnSurf& aProjector=aContext->ProjPS(aF);
1690   aProjector.Perform(aP11);
1691   //
1692   bFlag=aProjector.IsDone();
1693   if (!bFlag) {
1694     return bFlag;
1695   }
1696
1697   aProjector.LowerDistanceParameters(ULD, VLD);
1698   aP2D.SetCoord(ULD, VLD);
1699   //
1700   bFlag=aContext->IsPointInOnFace (aF, aP2D);
1701   //
1702   if (!bFlag) {
1703     return bFlag;
1704   }
1705   //
1706   // Treatment P12
1707   BOPTools_AlgoTools::PointOnEdge(aE1, l1, aP12);
1708   //
1709   aProjector.Perform(aP12);
1710   //
1711   bFlag=aProjector.IsDone();
1712   if (!bFlag) {
1713     return bFlag;
1714   }
1715
1716   aProjector.LowerDistanceParameters(ULD, VLD);
1717   aP2D.SetCoord(ULD, VLD);
1718   //
1719   bFlag=aContext->IsPointInOnFace (aF, aP2D);
1720   //
1721   if (!bFlag) {
1722     return bFlag;
1723   }
1724   //
1725   // Treatment intemediate
1726   Standard_Real m1, aTolF, aTolE, aTol, aDist;
1727   m1=IntTools_Tools::IntermediatePoint(f1, l1);
1728   BOPTools_AlgoTools::PointOnEdge(aE1, m1, aP12);
1729   //
1730   aProjector.Perform(aP12);
1731   //
1732   bFlag=aProjector.IsDone();
1733   if (!bFlag) {
1734     return bFlag;
1735   }
1736   //
1737   aTolE=BRep_Tool::Tolerance(aE1);
1738   aTolF=BRep_Tool::Tolerance(aF);
1739   aTol=aTolE+aTolF;
1740   aDist=aProjector.LowerDistance();
1741   if (aDist > aTol){
1742    return Standard_False;
1743   }
1744
1745   aProjector.LowerDistanceParameters(ULD, VLD);
1746   aP2D.SetCoord(ULD, VLD);
1747   //
1748   bFlag=aContext->IsPointInOnFace (aF, aP2D);
1749   //
1750   if (!bFlag) {
1751     return bFlag;
1752   }
1753   return bFlag;
1754 }
1755
1756 //=======================================================================
1757 //function : IsMicroEdge
1758 //purpose  :
1759 //=======================================================================
1760 Standard_Boolean BOPTools_AlgoTools::IsMicroEdge(const TopoDS_Edge& aE,
1761                                                  const Handle(BOPInt_Context)& aCtx)
1762 {
1763   Standard_Boolean bRet;
1764   Standard_Integer iErr;
1765   Standard_Real aT1, aT2, aTmp;
1766   Handle(Geom_Curve) aC3D;
1767   TopoDS_Vertex aV1, aV2;
1768   //
1769   bRet=(BRep_Tool::Degenerated(aE) ||
1770         !BRep_Tool::IsGeometric(aE));
1771   if (bRet) {
1772     return bRet;
1773   }
1774   //
1775   aC3D=BRep_Tool::Curve(aE, aT1, aT2);
1776   TopExp::Vertices(aE, aV1, aV2);
1777   aT1=BRep_Tool::Parameter(aV1, aE);
1778   aT2=BRep_Tool::Parameter(aV2, aE);
1779   if (aT2<aT1) {
1780     aTmp=aT1;
1781     aT1=aT2;
1782     aT2=aTmp;
1783   }
1784   //
1785   BOPInt_ShrunkRange aSR;
1786   aSR.SetData(aE, aT1, aT2, aV1, aV2, aCtx);
1787   aSR.Perform();
1788   iErr=aSR.ErrorStatus();
1789   bRet = !(iErr==0);
1790   //
1791   return bRet;
1792 }
1793
1794 //=======================================================================
1795 //function : GetFaceDir
1796 //purpose  : Get binormal direction for the face in the point aP
1797 //=======================================================================
1798 void GetFaceDir(const TopoDS_Edge& aE,
1799                 const TopoDS_Face& aF,
1800                 const gp_Pnt& aP,
1801                 const Standard_Real aT,
1802                 const gp_Dir& aDTgt,
1803                 gp_Dir& aDN,
1804                 gp_Dir& aDB,
1805                 Handle(BOPInt_Context)& theContext,
1806                 GeomAPI_ProjectPointOnSurf& aProjPL)
1807 {
1808   BOPTools_AlgoTools3D::GetNormalToFaceOnEdge(aE, aF, aT, aDN);
1809   if (aF.Orientation()==TopAbs_REVERSED){
1810     aDN.Reverse();
1811   }
1812   aDB = aDN^aDTgt;
1813   //
1814   gp_Pnt aPx;
1815   if (!FindPointInFace(aE, aF, aP, aDB, aPx, theContext, aProjPL)) {
1816     BOPTools_AlgoTools3D::GetApproxNormalToFaceOnEdge(aE, aF, aT, aPx, aDN, theContext);
1817     aProjPL.Perform(aPx);
1818     aPx = aProjPL.NearestPoint();
1819     gp_Vec aVec(aP, aPx);
1820     aDB.SetXYZ(aVec.XYZ());
1821   }
1822 }
1823
1824 //=======================================================================
1825 //function : FindPointInFace
1826 //purpose  : Find a point in the face in direction of <aDB>
1827 //=======================================================================
1828 Standard_Boolean FindPointInFace(const TopoDS_Edge& aE,
1829                                  const TopoDS_Face& aF,
1830                                  const gp_Pnt& aP,
1831                                  gp_Dir& aDB,
1832                                  gp_Pnt& aPOut,
1833                                  Handle(BOPInt_Context)& theContext,
1834                                  GeomAPI_ProjectPointOnSurf& aProjPL)
1835 {
1836   Standard_Integer aNbItMax;
1837   Standard_Real aDt, aDtMin, aTolE, aTolF, aDist;
1838   Standard_Boolean bRet;
1839   gp_Pnt aP1;
1840   BRepAdaptor_Surface aBAS;
1841   //
1842   bRet = Standard_False;
1843   aTolE = BRep_Tool::Tolerance(aE);
1844   aTolF = BRep_Tool::Tolerance(aF);
1845   aDt = 2*(aTolE+aTolF);
1846   aBAS.Initialize(aF, Standard_False);
1847   //
1848   aDtMin=5.e-4;
1849   if (aDt < aDtMin) {
1850     Standard_Real aR;
1851     GeomAbs_SurfaceType aSType=aBAS.GetType();
1852     switch (aSType) {
1853     case GeomAbs_Cylinder:
1854       aR = aBAS.Cylinder().Radius();
1855       break;
1856     case GeomAbs_Cone:
1857       aR = aBAS.Cone().RefRadius();
1858       break;
1859     case GeomAbs_Sphere:
1860       aR = aBAS.Sphere().Radius();
1861       break;
1862     case GeomAbs_Torus:
1863       aR = aBAS.Torus().MinorRadius();
1864       break;
1865     case GeomAbs_SurfaceOfRevolution:
1866       aR=1.;
1867       break;
1868     default:
1869       aR=0.001;
1870     }
1871     if (aR < 0.01) {
1872       aDtMin=5.e-6;
1873     }
1874     else if (aR > 100.) {
1875       aDtMin*=10;
1876     }
1877     if (aDt < aDtMin) {
1878       aDt=aDtMin;
1879     }
1880   }
1881   //
1882   GeomAPI_ProjectPointOnSurf& aProj=theContext->ProjPS(aF);
1883   aNbItMax = 15;
1884   //
1885   do {
1886     aP1.SetCoord(aP.X()+aDt*aDB.X(),
1887                  aP.Y()+aDt*aDB.Y(),
1888                  aP.Z()+aDt*aDB.Z());
1889     //
1890     aProj.Perform(aP1);
1891     if (!aProj.IsDone()) {
1892       return bRet;
1893     }
1894     aPOut = aProj.NearestPoint();
1895     aDist = aProj.LowerDistance();
1896     //
1897     aProjPL.Perform(aPOut);
1898     aPOut = aProjPL.NearestPoint();
1899     //
1900     gp_Vec aV(aP, aPOut);
1901     aDB.SetXYZ(aV.XYZ());
1902   } while (aDist>Precision::Angular() && --aNbItMax);
1903   //
1904   bRet = aDist < Precision::Angular();
1905   return bRet;
1906 }