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