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