0028556: Invalid result of Fuse operation in the test case bugs moddata_2 bug469
[occt.git] / src / BOPAlgo / BOPAlgo_WireSplitter_1.cxx
1 // Created by: Peter KURNEV
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
3 //
4 // This file is part of Open CASCADE Technology software library.
5 //
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
11 //
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
14
15 #include <BOPAlgo_WireEdgeSet.hxx>
16 #include <BOPAlgo_WireSplitter.hxx>
17 #include <BOPCol_IndexedDataMapOfShapeInteger.hxx>
18 #include <BOPCol_IndexedDataMapOfShapeListOfShape.hxx>
19 #include <BOPCol_ListOfShape.hxx>
20 #include <BOPCol_MapOfShape.hxx>
21 #include <BOPCol_SequenceOfPnt2d.hxx>
22 #include <BOPCol_SequenceOfReal.hxx>
23 #include <BOPCol_SequenceOfShape.hxx>
24 #include <BOPTools_AlgoTools2D.hxx>
25 #include <BRep_Builder.hxx>
26 #include <BRep_Tool.hxx>
27 #include <BRepAdaptor_Surface.hxx>
28 #include <Geom2d_Curve.hxx>
29 #include <Geom2d_Line.hxx>
30 #include <Geom2dAdaptor_Curve.hxx>
31 #include <Geom2dInt_GInter.hxx>
32 #include <Geom_Plane.hxx>
33 #include <Geom_RectangularTrimmedSurface.hxx>
34 #include <Geom_Surface.hxx>
35 #include <GeomAdaptor_Surface.hxx>
36 #include <gp_Dir2d.hxx>
37 #include <gp_Pnt2d.hxx>
38 #include <gp_Vec2d.hxx>
39 #include <IntRes2d_Domain.hxx>
40 #include <IntRes2d_IntersectionPoint.hxx>
41 #include <Precision.hxx>
42 #include <TopExp.hxx>
43 #include <TopExp_Explorer.hxx>
44 #include <TopLoc_Location.hxx>
45 #include <TopoDS.hxx>
46 #include <TopoDS_Edge.hxx>
47 #include <TopoDS_Face.hxx>
48 #include <TopoDS_Iterator.hxx>
49 #include <TopoDS_Vertex.hxx>
50 #include <TopoDS_Wire.hxx>
51 #include <TopTools_ShapeMapHasher.hxx>
52 #include <Geom2dLProp_CLProps2d.hxx>
53 typedef NCollection_DataMap \
54   <TopoDS_Shape, Standard_Boolean, TopTools_ShapeMapHasher> \
55   BOPCol_DataMapOfShapeBoolean; 
56 //
57
58 static
59   Standard_Real Angle (const gp_Dir2d& aDir2D);
60
61 static
62   Standard_Real Angle2D (const TopoDS_Vertex& aV,
63                          const TopoDS_Edge& anEdge,
64                          const TopoDS_Face& myFace,
65                          const GeomAdaptor_Surface& aGAS,
66                          const Standard_Boolean aFlag,
67                          const Handle(IntTools_Context)& theContext);
68
69 static
70   void GetNextVertex(const TopoDS_Vertex& aV,
71                      const TopoDS_Edge& aE,
72                      TopoDS_Vertex& aV1);
73
74 static
75   Standard_Real AngleIn(const TopoDS_Edge& aEIn,
76                         const BOPAlgo_ListOfEdgeInfo& aLEInfo);
77
78 static
79   Standard_Integer NbWaysOut(const BOPAlgo_ListOfEdgeInfo& aLEInfo);
80
81 static
82   gp_Pnt2d Coord2dVf (const TopoDS_Edge& aE,
83                       const TopoDS_Face& aF);
84
85 static
86   gp_Pnt2d Coord2d (const TopoDS_Vertex& aV1,
87                     const TopoDS_Edge& aE1,
88                     const TopoDS_Face& aF);
89
90 static 
91   Standard_Real ClockWiseAngle(const Standard_Real aAngleIn,
92                                const Standard_Real aAngleOut);
93
94 static 
95   void Path (const GeomAdaptor_Surface& aGAS,
96              const TopoDS_Face& myFace,
97              const TopoDS_Vertex& aVa,
98              const TopoDS_Edge& aEOuta,
99              BOPAlgo_EdgeInfo& anEdgeInfo,
100              BOPCol_SequenceOfShape& aLS,
101              BOPCol_SequenceOfShape& aVertVa,
102              BOPCol_SequenceOfPnt2d& aCoordVa,
103              BOPTools_ConnexityBlock& aCB,
104              BOPAlgo_IndexedDataMapOfShapeListOfEdgeInfo& mySmartMap,
105              BOPCol_DataMapOfShapeBoolean aVertMap);
106
107 static
108   Standard_Real Angle (const gp_Dir2d& aDir2D);
109
110 static
111   Standard_Real Tolerance2D (const TopoDS_Vertex& aV,
112                              const GeomAdaptor_Surface& aGAS);
113
114
115
116 static
117   Standard_Real UTolerance2D (const TopoDS_Vertex& aV,
118                               const GeomAdaptor_Surface& aGAS);
119 static
120   Standard_Real VTolerance2D (const TopoDS_Vertex& aV,
121                               const GeomAdaptor_Surface& aGAS);
122
123 static
124   void RefineAngles(const TopoDS_Face& myFace,
125                     const BOPCol_ListOfShape&,
126                     BOPAlgo_IndexedDataMapOfShapeListOfEdgeInfo&,
127                     const Handle(IntTools_Context)&);
128
129
130 static
131   void RefineAngles(const TopoDS_Vertex& ,
132                   const TopoDS_Face& ,
133                   const BOPCol_MapOfShape& ,
134                   BOPAlgo_ListOfEdgeInfo&,
135                   const Handle(IntTools_Context)&);
136
137 static
138   Standard_Boolean RefineAngle2D(const TopoDS_Vertex& ,
139                                  const TopoDS_Edge& ,
140                                  const TopoDS_Face& ,
141                                  const Standard_Real ,
142                                  const Standard_Real ,
143                                  const Standard_Real, 
144                                  Standard_Real& ,
145                                  const Handle(IntTools_Context)& );
146
147 //=======================================================================
148 //function : SplitBlock
149 //purpose  : 
150 //=======================================================================
151 void BOPAlgo_WireSplitter::SplitBlock(const TopoDS_Face& myFace,
152                                       BOPTools_ConnexityBlock& aCB,
153                                       const Handle(IntTools_Context)& theContext)
154 {
155   Standard_Boolean bNothingToDo, bIsClosed, bIsIN;
156   Standard_Integer aIx, aNb, i, aCntIn, aCntOut;
157   Standard_Real aAngle;
158   TopAbs_Orientation aOr;
159   TopoDS_Iterator aItS;
160   TopoDS_Vertex aVV;
161   TopoDS_Shape aV1;
162   BOPCol_ListIteratorOfListOfShape aIt;
163   BOPAlgo_ListIteratorOfListOfEdgeInfo aItLEI;
164   //
165   BOPAlgo_IndexedDataMapOfShapeListOfEdgeInfo mySmartMap(100);
166   BOPCol_DataMapOfShapeBoolean aVertMap;
167   //
168   const BOPCol_ListOfShape& myEdges=aCB.Shapes();
169   //
170   // 1.Filling mySmartMap
171   aIt.Initialize(myEdges);
172   for(; aIt.More(); aIt.Next()) {
173     const TopoDS_Edge& aE=(*(TopoDS_Edge *)&aIt.Value());
174     if (!BOPTools_AlgoTools2D::HasCurveOnSurface (aE, myFace)) {
175       continue;
176     }
177     //
178     bIsClosed = BRep_Tool::Degenerated(aE) || 
179                 BRep_Tool::IsClosed(aE, myFace);
180     //
181     aItS.Initialize(aE);
182     for(i = 0; aItS.More(); aItS.Next(), ++i) {
183       const TopoDS_Shape& aV = aItS.Value();
184       aIx = mySmartMap.FindIndex(aV);
185       if (!aIx) {
186         BOPAlgo_ListOfEdgeInfo aLEIx;
187         aIx = mySmartMap.Add(aV, aLEIx);
188       }
189       //
190       BOPAlgo_ListOfEdgeInfo& aLEI = mySmartMap(aIx);
191       BOPAlgo_EdgeInfo aEI;
192       //
193       aEI.SetEdge(aE);
194       aOr = aV.Orientation();
195       bIsIN = (aOr == TopAbs_REVERSED);
196       aEI.SetInFlag(bIsIN);
197       aLEI.Append(aEI);
198       //
199       if (!i) {
200         aV1 = aV;
201       } else {
202         bIsClosed = bIsClosed || aV1.IsSame(aV);
203       }
204       //
205       if (aVertMap.IsBound(aV)) {
206         if (bIsClosed) {
207           aVertMap.ChangeFind(aV) = bIsClosed;
208         }
209       } else {
210         aVertMap.Bind(aV, bIsClosed);
211       }
212     }
213   }
214   //
215   aNb=mySmartMap.Extent();
216   //
217   bNothingToDo=Standard_True;
218   for (i=1; i<=aNb; i++) {
219     aCntIn=0;
220     aCntOut=0;
221     const BOPAlgo_ListOfEdgeInfo& aLEInfo= mySmartMap(i);
222     BOPAlgo_ListIteratorOfListOfEdgeInfo anIt(aLEInfo);
223     for (; anIt.More(); anIt.Next()) {
224       const BOPAlgo_EdgeInfo& aEI=anIt.Value();
225       if (aEI.IsIn()) {
226         aCntIn++;
227       }
228       else {
229         aCntOut++;
230       }
231     }
232     if (aCntIn!=1 || aCntOut!=1) {
233       bNothingToDo=Standard_False;
234       break;
235     }
236   }
237   //
238   // Each vertex has one edge In and one - Out. Good. But it is not enought
239   // to consider that nothing to do with this. We must check edges on TShape
240   // coinsidence. If there are such edges there is something to do with.
241   if (bNothingToDo) {
242     Standard_Integer aNbE, aNbMapEE;
243     Standard_Boolean bFlag;
244     //
245     BOPCol_IndexedDataMapOfShapeListOfShape aMapEE(100);
246     aNbE=myEdges.Extent();
247     //
248     aIt.Initialize(myEdges);
249     for (; aIt.More(); aIt.Next()) {
250       const TopoDS_Shape& aE = aIt.Value();
251       if (!aMapEE.Contains(aE)) {
252         BOPCol_ListOfShape aLEx;
253         aLEx.Append(aE);
254         aMapEE.Add(aE, aLEx);
255       }
256       else {
257         BOPCol_ListOfShape& aLEx=aMapEE.ChangeFromKey(aE);
258         aLEx.Append(aE);
259       }
260     }
261     //
262     bFlag=Standard_True;
263     aNbMapEE=aMapEE.Extent();
264     for (i=1; i<=aNbMapEE; ++i) {
265       const BOPCol_ListOfShape& aLEx=aMapEE(i);
266       aNbE=aLEx.Extent();
267       if (aNbE==1) {// usual case
268         continue;
269       }
270       else if (aNbE==2){
271         const TopoDS_Shape& aE1=aLEx.First();
272         const TopoDS_Shape& aE2=aLEx.Last();
273         if (aE1.IsSame(aE2)) {
274           bFlag=Standard_False;
275           break;
276         }
277       }
278       else {
279         bFlag=Standard_False;
280         break;
281       }
282     }
283     bNothingToDo=bNothingToDo && bFlag;
284   } // if (bNothingToDo) {
285   if (bNothingToDo) {
286     TopoDS_Wire aW;
287     //
288     BOPCol_ListOfShape& aLECB=aCB.ChangeShapes();
289     BOPAlgo_WireSplitter::MakeWire(aLECB, aW);
290     BOPCol_ListOfShape& aLoops=aCB.ChangeLoops();
291     aLoops.Append(aW);
292     //
293     return;
294   }
295   //
296   // 3. Angles in mySmartMap
297   const BRepAdaptor_Surface& aBAS = theContext->SurfaceAdaptor(myFace);
298   const GeomAdaptor_Surface& aGAS=aBAS.Surface();
299   //
300   for (i=1; i<=aNb; i++) {
301     const TopoDS_Vertex& aV = (*(TopoDS_Vertex *)(&mySmartMap.FindKey(i))); 
302     const BOPAlgo_ListOfEdgeInfo& aLEI= mySmartMap(i);
303     aItLEI.Initialize(aLEI);
304     for (; aItLEI.More(); aItLEI.Next()) {
305       BOPAlgo_EdgeInfo& aEI=aItLEI.ChangeValue();
306       const TopoDS_Edge& aE=aEI.Edge();
307       //
308       aVV = aV;
309       bIsIN = aEI.IsIn();
310       aOr = bIsIN ? TopAbs_REVERSED : TopAbs_FORWARD;
311       aVV.Orientation(aOr);
312       aAngle = Angle2D(aVV, aE, myFace, aGAS, bIsIN, theContext);
313       aEI.SetAngle(aAngle);
314     }
315   }// for (i=1; i<=aNb; i++) {
316   //
317   //Theme: The treatment p-curves convergent in node.
318   //The refining the angles of p-curves taking into account 
319   //bounding curves if exist. 
320   RefineAngles(myFace, myEdges, mySmartMap, theContext);
321   //
322   // 4. Do
323   //
324   Standard_Boolean bIsOut, bIsNotPassed;
325   BOPCol_SequenceOfShape aLS, aVertVa;
326   BOPCol_SequenceOfPnt2d aCoordVa;
327   //
328   for (i=1; i<=aNb; ++i) {
329     const TopoDS_Vertex& aVa=(*(TopoDS_Vertex *)(&mySmartMap.FindKey(i))); 
330     const BOPAlgo_ListOfEdgeInfo& aLEI=mySmartMap(i);
331     aItLEI.Initialize(aLEI);
332     for (; aItLEI.More(); aItLEI.Next()) {
333       BOPAlgo_EdgeInfo& aEI=aItLEI.ChangeValue();
334       const TopoDS_Edge& aEOuta=aEI.Edge();
335       //
336       bIsOut=!aEI.IsIn();
337       bIsNotPassed=!aEI.Passed();
338       if (bIsOut && bIsNotPassed) {
339         //
340         aLS.Clear();
341         aVertVa.Clear();
342         aCoordVa.Clear();
343         //
344         Path(aGAS, myFace, aVa, aEOuta, aEI, aLS, 
345              aVertVa, aCoordVa, aCB, mySmartMap, aVertMap);
346       }
347     }
348   }// for (i=1; i<=aNb; ++i) {
349 }
350 //=======================================================================
351 // function: Path
352 // purpose: 
353 //=======================================================================
354 void Path (const GeomAdaptor_Surface& aGAS,
355            const TopoDS_Face& myFace,
356            const TopoDS_Vertex& aVFirst,
357            const TopoDS_Edge& aEFirst,
358            BOPAlgo_EdgeInfo& aEIFirst,
359            BOPCol_SequenceOfShape& aLS,
360            BOPCol_SequenceOfShape& aVertVa,
361            BOPCol_SequenceOfPnt2d& aCoordVa,
362            BOPTools_ConnexityBlock& aCB,
363            BOPAlgo_IndexedDataMapOfShapeListOfEdgeInfo& mySmartMap,
364            BOPCol_DataMapOfShapeBoolean aVertMap)
365      
366 {
367   Standard_Integer i, j, aNb, aNbj;
368   Standard_Real anAngleIn, anAngleOut, anAngle, aMinAngle;
369   Standard_Real aTol2D, aTol2D2, aD2, aTwoPI;
370   Standard_Boolean anIsSameV2d, anIsSameV, anIsFound, anIsOut, anIsNotPassed;
371   Standard_Boolean bIsClosed;
372   TopoDS_Vertex aVa, aVb;
373   TopoDS_Edge aEOuta;
374   BOPAlgo_ListIteratorOfListOfEdgeInfo anIt;
375   Standard_Real eps = Epsilon(1.);
376   //
377   aVa = aVFirst;
378   aEOuta = aEFirst;
379   BOPAlgo_EdgeInfo* anEdgeInfo = &aEIFirst;
380   //
381   aTwoPI = M_PI + M_PI;
382   //
383   // append block
384   //
385   for (;;) {
386     // Do not escape through edge from which you enter 
387     aNb=aLS.Length();
388     if (aNb==1) {
389       const TopoDS_Shape& anEPrev=aLS(aNb);
390       if (anEPrev.IsSame(aEOuta)) {
391         return;
392       }
393     }
394     //
395     anEdgeInfo->SetPassed(Standard_True);
396     aLS.Append(aEOuta);
397     aVertVa.Append(aVa);
398     
399     TopoDS_Vertex pVa=aVa;
400     pVa.Orientation(TopAbs_FORWARD);
401     gp_Pnt2d aPa=Coord2d(pVa, aEOuta, myFace);
402     aCoordVa.Append(aPa);
403     
404     GetNextVertex (pVa, aEOuta, aVb);
405     
406     gp_Pnt2d aPb=Coord2d(aVb, aEOuta, myFace);
407     
408     const BOPAlgo_ListOfEdgeInfo& aLEInfo=mySmartMap.FindFromKey(aVb);
409     //
410     aTol2D = 2.*Tolerance2D(aVb, aGAS);
411     aTol2D2 = aTol2D * aTol2D;
412     //
413     bIsClosed = aVertMap.Find(aVb);
414     {
415       BOPCol_ListOfShape aBuf;
416       Standard_Boolean bHasEdge = Standard_False;
417       //
418       aNb = aLS.Length();
419       for (i = aNb; i>0; --i) {
420         const TopoDS_Shape& aVPrev=aVertVa(i);
421         const gp_Pnt2d& aPaPrev=aCoordVa(i);
422         const TopoDS_Edge& aEPrev = TopoDS::Edge(aLS(i));
423         //
424         aBuf.Append(aEPrev);
425         if (!bHasEdge) {
426           bHasEdge = !BRep_Tool::Degenerated(aEPrev);
427           // do not create wire from degenerated edges only
428           if (!bHasEdge) {
429             continue;
430           }
431         }
432         //
433         anIsSameV = aVPrev.IsSame(aVb);
434         anIsSameV2d = anIsSameV;
435         if (anIsSameV) {
436           if(bIsClosed) {
437             aD2 = aPaPrev.SquareDistance(aPb);
438             anIsSameV2d = aD2 < aTol2D2;
439             if (anIsSameV2d) {
440               Standard_Real udist = fabs(aPaPrev.X() - aPb.X());
441               Standard_Real vdist = fabs(aPaPrev.Y() - aPb.Y());
442               Standard_Real aTolU = 2.*UTolerance2D(aVb, aGAS);
443               Standard_Real aTolV = 2.*VTolerance2D(aVb, aGAS);
444               //
445               if((udist > aTolU) || (vdist > aTolV)) {
446                 anIsSameV2d = Standard_False;
447               }
448             }
449           }
450         }//if (anIsSameV) {
451         //
452         if (anIsSameV && anIsSameV2d) {
453           Standard_Integer iPriz;
454           iPriz=1;
455           if (aBuf.Extent()==2) {
456             if(aBuf.First().IsSame(aBuf.Last())) {
457               iPriz=0;
458             }
459           }
460           if (iPriz) {
461             TopoDS_Wire aW;
462             BOPAlgo_WireSplitter::MakeWire(aBuf, aW);
463             aCB.ChangeLoops().Append(aW);
464           }
465           //
466           aNbj=i-1;
467           if (aNbj<1) {
468             //
469             aLS.Clear();
470             aVertVa.Clear();
471             aCoordVa.Clear();
472             //
473             return;
474           }
475           //
476           BOPCol_SequenceOfShape aLSt, aVertVat;
477           BOPCol_SequenceOfPnt2d aCoordVat;
478           //
479           aVb=(*(TopoDS_Vertex *)(&aVertVa(i))); 
480           //
481           for (j=1; j<=aNbj; ++j) {
482             aLSt.Append(aLS(j));
483             aVertVat.Append(aVertVa(j));
484             aCoordVat.Append(aCoordVa(j));
485           }
486           //
487           aLS.Clear();
488           aVertVa.Clear();
489           aCoordVa.Clear();
490           
491           aLS=aLSt;
492           aVertVa=aVertVat;
493           aCoordVa=aCoordVat;
494           //
495           break;
496         }
497       }
498     }
499     //
500     // aEOutb
501     BOPAlgo_EdgeInfo *pEdgeInfo=NULL;
502     //
503     anAngleIn = AngleIn(aEOuta, aLEInfo);
504     aMinAngle = 100.;
505     anIsFound = Standard_False;
506     Standard_Integer iCnt = NbWaysOut(aLEInfo);
507     Standard_Integer aCurIndexE = 0;
508     anIt.Initialize(aLEInfo);
509     for (; anIt.More(); anIt.Next()) {
510       BOPAlgo_EdgeInfo& anEI=anIt.ChangeValue();
511       const TopoDS_Edge& aE=anEI.Edge();
512       anIsOut=!anEI.IsIn();
513       anIsNotPassed=!anEI.Passed();
514       //
515       if (anIsOut && anIsNotPassed) {
516         aCurIndexE++;
517         //
518         // Is there one way to go out of the vertex 
519         // we have to use it only.
520         //
521         if (!iCnt) {
522           // no way to go . (Error)
523           return;
524         }
525         //
526         if (iCnt==1) {
527           // the one and only way to go out .
528           pEdgeInfo=&anEI;
529           anIsFound=Standard_True;
530           break;
531         }
532         //
533         if (aE.IsSame(aEOuta)) {
534           anAngle = aTwoPI;
535         } else {
536           //check 2d distance
537           if (bIsClosed) {
538             gp_Pnt2d aP2Dx;
539             //
540             aP2Dx = Coord2dVf(aE, myFace);
541             //
542             aD2 = aP2Dx.SquareDistance(aPb);
543             if (aD2 > aTol2D2){
544               continue;
545             }
546           }
547           //
548           // Look for minimal angle and make the choice.
549           anAngleOut=anEI.Angle();
550           anAngle=ClockWiseAngle(anAngleIn, anAngleOut);
551         }
552         if (anAngle < aMinAngle - eps) {
553           aMinAngle=anAngle;
554           pEdgeInfo=&anEI;
555           anIsFound=Standard_True;
556         }
557       }
558     } // for (; anIt.More(); anIt.Next()) 
559     //
560     if (!anIsFound) {
561       // no way to go . (Error)
562       return;
563     }
564     //
565     aVa = aVb;
566     aEOuta = pEdgeInfo->Edge();
567     anEdgeInfo = pEdgeInfo;
568   }
569 }
570 //=======================================================================
571 // function:  ClockWiseAngle
572 // purpose:
573 //=======================================================================
574  Standard_Real ClockWiseAngle(const Standard_Real aAngleIn,
575                               const Standard_Real aAngleOut)
576 {
577   Standard_Real aTwoPi=M_PI+M_PI;
578   Standard_Real dA, A1, A2, AIn, AOut ;
579
580   AIn=aAngleIn;
581   AOut=aAngleOut;
582   if (AIn >= aTwoPi) {
583     AIn=AIn-aTwoPi;
584   }
585   
586   if (AOut >= aTwoPi) {
587     AOut=AOut-aTwoPi;
588   }
589
590   A1=AIn+M_PI;
591   
592   if (A1 >= aTwoPi) {
593     A1=A1-aTwoPi;
594   }
595   
596   A2=AOut;
597   
598   dA=A1-A2;
599   if (dA <= 0.) {
600     dA=aTwoPi+dA;
601   }
602   //xx
603   //else if (dA <= 1.e-15) {
604   else if (dA <= 1.e-14) {
605     dA=aTwoPi;
606   }
607   return dA;
608 }
609 //=======================================================================
610 // function:  Coord2d
611 // purpose:
612 //=======================================================================
613  gp_Pnt2d Coord2d (const TopoDS_Vertex& aV1,
614                    const TopoDS_Edge& aE1,
615                    const TopoDS_Face& aF)
616 {
617   Standard_Real aT, aFirst, aLast;
618   Handle(Geom2d_Curve) aC2D;
619   gp_Pnt2d aP2D1;
620   //
621   aT=BRep_Tool::Parameter (aV1, aE1, aF);
622   aC2D=BRep_Tool::CurveOnSurface(aE1, aF, aFirst, aLast);
623   aC2D->D0 (aT, aP2D1);
624   //
625   return aP2D1;
626 }
627 //=======================================================================
628 // function:  Coord2dVf
629 // purpose:
630 //=======================================================================
631  gp_Pnt2d Coord2dVf (const TopoDS_Edge& aE,
632                      const TopoDS_Face& aF)
633 {
634   Standard_Real aCoord=99.;
635   gp_Pnt2d aP2D1(aCoord, aCoord);
636   TopoDS_Iterator aIt;
637   //
638   aIt.Initialize(aE);
639   for (; aIt.More(); aIt.Next()) {
640     const TopoDS_Shape& aVx=aIt.Value();
641     if (aVx.Orientation()==TopAbs_FORWARD) {
642       
643       const TopoDS_Vertex& aVxx=(*(TopoDS_Vertex *)(&aVx));// TopoDS::Vertex(aVx);
644       aP2D1=Coord2d(aVxx, aE, aF);
645       return aP2D1;
646     }
647   }
648   return aP2D1;
649 }
650
651 //=======================================================================
652 // function: NbWaysOut
653 // purpose: 
654 //=======================================================================
655 Standard_Integer NbWaysOut(const BOPAlgo_ListOfEdgeInfo& aLEInfo)
656 {
657   Standard_Boolean bIsOut, bIsNotPassed;
658   Standard_Integer iCnt=0;
659   BOPAlgo_ListIteratorOfListOfEdgeInfo anIt;
660   //
661   anIt.Initialize(aLEInfo);
662   for (; anIt.More(); anIt.Next()) {
663     const BOPAlgo_EdgeInfo& anEI=anIt.Value();
664     //
665     bIsOut=!anEI.IsIn();
666     bIsNotPassed=!anEI.Passed();
667     if (bIsOut && bIsNotPassed) {
668       iCnt++;
669     }
670   }
671   return iCnt;
672 }
673
674 //=======================================================================
675 // function:  AngleIn
676 // purpose:
677 //=======================================================================
678  Standard_Real AngleIn(const TopoDS_Edge& aEIn,
679                        const BOPAlgo_ListOfEdgeInfo& aLEInfo)
680 {
681   Standard_Real anAngleIn;
682   Standard_Boolean anIsIn;
683   BOPAlgo_ListIteratorOfListOfEdgeInfo anIt;
684
685   anIt.Initialize(aLEInfo);
686   for (; anIt.More(); anIt.Next()) {
687     const BOPAlgo_EdgeInfo& anEdgeInfo=anIt.Value();
688     const TopoDS_Edge& aE=anEdgeInfo.Edge();
689     anIsIn=anEdgeInfo.IsIn();
690     //
691     if (anIsIn && aE==aEIn) {
692       anAngleIn=anEdgeInfo.Angle();
693       return anAngleIn;
694     }
695   }
696   anAngleIn=0.;
697   return anAngleIn;
698 }
699 //=======================================================================
700 // function: GetNextVertex
701 // purpose: 
702 //=======================================================================
703  void GetNextVertex(const TopoDS_Vertex& aV,
704                     const TopoDS_Edge& aE,
705                     TopoDS_Vertex& aV1)
706 {
707   TopoDS_Iterator aIt;
708   //
709   aIt.Initialize(aE);
710   for (; aIt.More(); aIt.Next()) {
711     const TopoDS_Shape& aVx=aIt.Value();
712     if (!aVx.IsEqual(aV)) {
713       aV1=(*(TopoDS_Vertex *)(&aVx)); 
714       return ;
715     }
716   }
717   aV1=aV;
718 }
719 //=======================================================================
720 // function: Angle2D
721 // purpose: 
722 //=======================================================================
723   Standard_Real Angle2D (const TopoDS_Vertex& aV,
724                          const TopoDS_Edge& anEdge,
725                          const TopoDS_Face& myFace,
726                          const GeomAdaptor_Surface& aGAS,
727                          const Standard_Boolean bIsIN,
728                          const Handle(IntTools_Context)& theContext)
729 {
730   Standard_Real aFirst, aLast, aToler, dt, aTV, aTV1, anAngle, aTX;
731   gp_Pnt2d aPV, aPV1;
732   gp_Vec2d aV2D;
733   Handle(Geom2d_Curve) aC2D;
734   //
735   aTV=BRep_Tool::Parameter (aV, anEdge, myFace);
736   if (Precision::IsInfinite(aTV)) {
737     return 0.;
738   }
739   //
740   BOPTools_AlgoTools2D::CurveOnSurface (anEdge, myFace, aC2D, 
741                                         aFirst, aLast, aToler, theContext);
742   Standard_Real tol2d =2.*Tolerance2D(aV, aGAS);
743   //
744   GeomAbs_CurveType aType;
745   Geom2dAdaptor_Curve aGAC2D(aC2D);
746   //
747   dt = Max(aGAC2D.Resolution(tol2d), Precision::PConfusion());
748   //
749   aType=aGAC2D.GetType();
750   if (aType != GeomAbs_Line ) 
751   {
752     Geom2dLProp_CLProps2d LProp(aC2D,  aTV, 2, Precision::PConfusion());
753     if(LProp.IsTangentDefined())
754     {
755       Standard_Real R = LProp.Curvature();
756       if(R > Precision::PConfusion())
757       {
758         R = 1./R;
759         Standard_Real cosphi = R / (R + tol2d);
760         dt = Max(dt, ACos(cosphi)); //to avoid small dt for big R.
761       }
762     }
763   }
764   //for case chl/927/r9
765   aTX=0.05*(aLast - aFirst);//aTX=0.25*(aLast - aFirst);  
766   if (aTX < 5.e-5) {
767     aTX = 5.e-5;
768   }
769   if(dt > aTX) {
770     // to save direction of the curve as much as it possible
771     // in the case of big tolerances
772     dt = aTX; 
773   }
774
775   if (fabs (aTV-aFirst) < fabs(aTV - aLast)) {
776     aTV1=aTV + dt;
777   }
778   else {
779     aTV1=aTV - dt;
780   }
781   //
782   aGAC2D.D0 (aTV1, aPV1);
783   aGAC2D.D0 (aTV, aPV);
784   //
785   aV2D = bIsIN ? gp_Vec2d(aPV1, aPV) : gp_Vec2d(aPV, aPV1);
786   //
787   gp_Dir2d aDir2D(aV2D);
788   anAngle=Angle(aDir2D);
789   //
790   return anAngle;
791 }
792 //=======================================================================
793 // function: Angle
794 // purpose: 
795 //=======================================================================
796 Standard_Real Angle (const gp_Dir2d& aDir2D)
797 {
798   gp_Dir2d aRefDir(1., 0.);
799   Standard_Real anAngle;
800   
801   anAngle = aRefDir.Angle(aDir2D);
802   if (anAngle < 0.)
803     anAngle += M_PI + M_PI;
804   return anAngle;
805 }
806 //=======================================================================
807 // function:  Tolerance2D
808 // purpose:
809 //=======================================================================
810  Standard_Real Tolerance2D (const TopoDS_Vertex& aV,
811                             const GeomAdaptor_Surface& aGAS)         
812 {
813   Standard_Real aTol2D, anUr, aVr, aTolV3D;
814   GeomAbs_SurfaceType aType;
815   //
816   aType=aGAS.GetType();
817   aTolV3D=BRep_Tool::Tolerance(aV);
818
819   anUr=aGAS.UResolution(aTolV3D);
820   aVr =aGAS.VResolution(aTolV3D);
821   aTol2D=(aVr>anUr) ? aVr : anUr;
822   //
823   if (aTol2D < aTolV3D) {
824     aTol2D=aTolV3D;
825   }
826   if (aType==GeomAbs_BSplineSurface) {
827     aTol2D=1.1*aTol2D;
828   }
829   //
830   return aTol2D;
831 }
832
833 //=======================================================================
834 //function : UTolerance2D
835 //purpose  : 
836 //=======================================================================
837 Standard_Real UTolerance2D (const TopoDS_Vertex& aV,
838                             const GeomAdaptor_Surface& aGAS)
839 {
840   const Standard_Real aTolV3D = BRep_Tool::Tolerance(aV);
841   const Standard_Real anUr = aGAS.UResolution(aTolV3D);
842   //
843   return anUr;
844 }
845
846 //=======================================================================
847 //function : VTolerance2D
848 //purpose  : 
849 //=======================================================================
850 Standard_Real VTolerance2D (const TopoDS_Vertex& aV,
851                             const GeomAdaptor_Surface& aGAS)
852 {
853   const Standard_Real aTolV3D = BRep_Tool::Tolerance(aV);
854   const Standard_Real anVr = aGAS.VResolution(aTolV3D);
855   //
856   return anVr;
857 }
858
859 //=======================================================================
860 //function : RefineAngles
861 //purpose  : 
862 //=======================================================================
863 void RefineAngles(const TopoDS_Face& myFace,
864                   const BOPCol_ListOfShape& myEdges,
865                   BOPAlgo_IndexedDataMapOfShapeListOfEdgeInfo& mySmartMap,
866                   const Handle(IntTools_Context)& theContext)
867 {
868   Standard_Integer aNb, i;
869   BOPCol_IndexedDataMapOfShapeInteger aMSI;
870   BOPCol_MapOfShape aMBE;
871   BOPCol_ListIteratorOfListOfShape aIt;
872   //
873   // 1. Boundary Edges
874   aIt.Initialize(myEdges);
875   for(; aIt.More(); aIt.Next()) {
876     const TopoDS_Shape& aE=aIt.Value();
877     if(aMSI.Contains(aE)) {
878       Standard_Integer& iCnt = aMSI.ChangeFromKey(aE);
879       ++iCnt;
880     }
881     else {
882       Standard_Integer iCnt = 1;
883       aMSI.Add(aE, iCnt);
884     }
885   }
886   //
887   aNb = aMSI.Extent();
888   for (i = 1; i <= aNb; ++i) {
889     Standard_Integer iCnt = aMSI(i);
890     if (iCnt == 1) {
891       const TopoDS_Shape& aE = aMSI.FindKey(i);
892       aMBE.Add(aE);
893     }
894   }
895   //
896   aMSI.Clear();
897   //
898   aNb = mySmartMap.Extent();
899   for (i = 1; i <= aNb; ++i) {
900     const TopoDS_Vertex& aV=*((TopoDS_Vertex*)&mySmartMap.FindKey(i)); 
901     BOPAlgo_ListOfEdgeInfo& aLEI=mySmartMap(i);
902     //
903     RefineAngles(aV, myFace, aMBE, aLEI, theContext);
904   }
905 }
906 //=======================================================================
907 typedef NCollection_DataMap \
908   <TopoDS_Shape, Standard_Real, TopTools_ShapeMapHasher> \
909    BOPCol_DataMapOfShapeReal; 
910 typedef BOPCol_DataMapOfShapeReal::Iterator \
911   BOPCol_DataMapIteratorOfDataMapOfShapeReal; 
912 //
913 //=======================================================================
914 //function : RefineAngles
915 //purpose  : 
916 //=======================================================================
917 void RefineAngles(const TopoDS_Vertex& aV,
918                   const TopoDS_Face& myFace,
919                   const BOPCol_MapOfShape& aMBE,
920                   BOPAlgo_ListOfEdgeInfo& aLEI,
921                   const Handle(IntTools_Context)& theContext)
922 {
923   Standard_Boolean bIsIn, bIsBoundary, bRefined; 
924   Standard_Integer iCntBnd, iCntInt;
925   Standard_Real aA, aA1, aA2;
926   BOPCol_DataMapOfShapeReal aDMSR;
927   BOPAlgo_ListIteratorOfListOfEdgeInfo aItLEI;
928   //
929   aA1=0.;  // angle of outgoing edge
930   aA2=0.;  // angle of incoming edge
931   iCntBnd=0;
932   iCntInt=0;
933   aItLEI.Initialize(aLEI);
934   for (; aItLEI.More(); aItLEI.Next()) {
935     BOPAlgo_EdgeInfo& aEI=aItLEI.ChangeValue();
936     const TopoDS_Edge& aE=aEI.Edge();
937     bIsIn=aEI.IsIn();
938     aA=aEI.Angle();
939     //
940     if (aMBE.Contains(aE)) {
941       ++iCntBnd;
942       if (!bIsIn) {
943         aA1=aA;
944       }
945       else {
946         aA2=aA;
947       }
948     }
949     else {
950       ++iCntInt;
951     }
952   }
953   //
954   if (iCntBnd!=2) {
955     return;
956   }
957   //
958   Standard_Real aDelta = ClockWiseAngle(aA2, aA1);
959   aItLEI.Initialize(aLEI);
960   for (; aItLEI.More(); aItLEI.Next()) {
961     BOPAlgo_EdgeInfo& aEI=aItLEI.ChangeValue();
962     const TopoDS_Edge& aE=aEI.Edge();
963     //
964     bIsBoundary=aMBE.Contains(aE);
965     bIsIn=aEI.IsIn();
966     if (bIsBoundary || bIsIn) {
967       continue;
968     }
969     //
970     aA=aEI.Angle();
971     Standard_Real aDA = ClockWiseAngle(aA2, aA);
972     if (aDA < aDelta) {
973       continue;  // already inside
974     }
975     //
976     bRefined=RefineAngle2D(aV, aE, myFace, aA1, aA2, aDelta, aA, theContext);
977     if (bRefined) {
978       aDMSR.Bind(aE, aA);
979     }
980     else if (iCntInt == 2) {
981       aA = (aA <= aA1) ? (aA1 + Precision::Angular()) :
982                          (aA2 - Precision::Angular());
983       aDMSR.Bind(aE, aA);
984     }
985   }
986   
987   if (aDMSR.IsEmpty()) {
988     return;
989   }
990   //
991   // update Angles
992   aItLEI.Initialize(aLEI);
993   for (; aItLEI.More(); aItLEI.Next()) {
994     BOPAlgo_EdgeInfo& aEI=aItLEI.ChangeValue();
995     const TopoDS_Edge& aE=aEI.Edge();
996     //
997     bIsIn=aEI.IsIn();
998     if (!aDMSR.IsBound(aE)) {
999       continue;
1000     }
1001     //
1002     aA=aDMSR.Find(aE);
1003     if (bIsIn) {
1004       aA=aA+M_PI;
1005     }
1006     //
1007     aEI.SetAngle(aA);
1008   }
1009 }
1010 //=======================================================================
1011 //function : RefineAngle2D
1012 //purpose  : 
1013 //=======================================================================
1014 Standard_Boolean RefineAngle2D(const TopoDS_Vertex& aV,
1015                                const TopoDS_Edge& aE,
1016                                const TopoDS_Face& myFace,
1017                                const Standard_Real aA1,
1018                                const Standard_Real aA2,
1019                                const Standard_Real aDelta,
1020                                Standard_Real& aA,
1021                                const Handle(IntTools_Context)& theContext)
1022 {
1023   Standard_Boolean bRet;
1024   Standard_Integer i, j, aNbP;
1025   Standard_Real aTV, aTol, aT1, aT2, dT, aAngle, aT, aTOp;
1026   Standard_Real aTolInt, aAi, aXi, aYi, aT1j, aT2j, aT1max, aT2max, aCf; 
1027   gp_Pnt2d aPV, aP, aP1, aP2;
1028   Handle(Geom2d_Curve) aC2D;
1029   Handle(Geom2d_Line) aLi;
1030   Geom2dAdaptor_Curve aGAC1, aGAC2;
1031   Geom2dInt_GInter aGInter;
1032   IntRes2d_Domain aDomain1, aDomain2;
1033   //
1034   bRet=Standard_True;
1035   aCf=0.01;
1036   aTolInt=1.e-10;  
1037   //
1038   BOPTools_AlgoTools2D::CurveOnSurface(aE, myFace, aC2D, aT1, aT2, aTol, theContext);
1039   aGAC1.Load(aC2D, aT1, aT2);
1040   //
1041   aTV=BRep_Tool::Parameter (aV, aE, myFace);
1042   aGAC1.D0(aTV, aPV);
1043   //
1044   aTOp = (fabs(aTV-aT1) < fabs(aTV-aT2)) ? aT2 : aT1;
1045   //
1046   const Standard_Real MaxDT = 0.3 * (aT2 - aT1);
1047   aGAC1.D0(aT1, aP1);
1048   aGAC1.D0(aT2, aP2);
1049   aDomain1.SetValues(aP1, aT1, aTolInt, aP2, aT2, aTolInt);
1050   //
1051   for (i=0; i<2; ++i) {
1052     aAi=(!i) ? aA1 : (aA2 + M_PI);
1053     aXi=cos(aAi);
1054     aYi=sin(aAi);
1055     gp_Dir2d aDiri(aXi, aYi);
1056     aLi=new Geom2d_Line(aPV, aDiri);
1057     //
1058     aGAC2.Load(aLi);
1059     //
1060     aGInter.Perform(aGAC1, aDomain1, aGAC2, aDomain2, aTolInt, aTolInt);
1061     if (!aGInter.IsDone()) {
1062       continue;
1063     }
1064     //
1065     aNbP = aGInter.NbPoints();
1066     aT1max = aTV;
1067     aT2max = -1.;
1068     for (j = 1; j <= aNbP; ++j) {
1069       const IntRes2d_IntersectionPoint& aIPj = aGInter.Point(j);
1070       aT1j = aIPj.ParamOnFirst();
1071       aT2j = aIPj.ParamOnSecond();
1072       //
1073       if (aT2j > aT2max && Abs(aT1j - aTV) < MaxDT) {
1074         aT2max = aT2j;
1075         aT1max = aT1j;
1076       }
1077     }
1078     //
1079     if (aT2max > 0) {
1080       dT = aTOp - aT1max;
1081       if (Abs(dT) < aTolInt) {
1082         continue;
1083       }
1084       //
1085       aT = aT1max + aCf*dT;
1086       aGAC1.D0(aT, aP);
1087       gp_Vec2d aV2D(aPV, aP);
1088       gp_Dir2d aDir2D(aV2D);
1089       //
1090       aAngle = Angle(aDir2D);
1091       Standard_Real aDA = ClockWiseAngle(aA2, aAngle);
1092       if (aDA < aDelta) {
1093         aA = aAngle;
1094         return bRet;
1095       }
1096     }
1097   }// for (i=0; i<2; ++i) {
1098   return !bRet;
1099 }