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