Salome HOME
Changes in the presentations of features
[modules/shaper.git] / src / GeomAlgoAPI / GeomAlgoAPI_SketchBuilder.cpp
1 // File:        GeomAlgoAPI_SketchBuilder.cpp
2 // Created:     02 Jun 2014
3 // Author:      Artem ZHIDKOV
4
5 #include <GeomAlgoAPI_SketchBuilder.h>
6
7 #include <set>
8
9 #include <gp_Dir.hxx>
10 #include <gp_Pln.hxx>
11
12 #include <BOPAlgo_Builder.hxx>
13 #include <BOPAlgo_Operation.hxx>
14 #include <BOPAlgo_PaveFiller.hxx>
15 #include <BOPTools.hxx>
16 #include <BRep_Builder.hxx>
17 #include <BRep_Curve3D.hxx>
18 #include <BRep_TEdge.hxx>
19 #include <BRep_TVertex.hxx>
20 #include <BRepAlgoAPI_Cut.hxx>
21 #include <BRepBuilderAPI_MakeFace.hxx>
22 #include <BRepClass_FaceClassifier.hxx>
23 #include <Geom_Curve.hxx>
24 #include <Geom_Plane.hxx>
25 #include <TopExp.hxx>
26 #include <TopExp_Explorer.hxx>
27 #include <TopoDS_Edge.hxx>
28 #include <TopoDS_Face.hxx>
29 #include <TopoDS_Vertex.hxx>
30 #include <TopoDS_Wire.hxx>
31
32 #include <Precision.hxx>
33
34 #ifndef DBL_MAX
35 #define DBL_MAX 1.7976931348623158e+308 
36 #endif
37
38
39 const double tolerance = Precision::Confusion();
40 // This value helps to find direction on the boundaries of curve.
41 // It is not significant for lines, but is used for circles to avoid 
42 // wrong directions of movement (when two edges are tangent on the certain vertex)
43 const double shift     = acos(1.0 - 2.0 * tolerance);
44
45 /// \brief Search first vertex - the vertex with lowest x coordinate, which is used in 2 edges at least
46 static const TopoDS_Shape& findStartVertex(
47           const BOPCol_IndexedDataMapOfShapeListOfShape& theMapVE,
48           const gp_Dir& theDirX, const gp_Dir& theDirY);
49
50 /// \brief Search the vertex on the sketch candidate to be the next one in the loop
51 static void findNextVertex(
52           const TopoDS_Shape&                            theStartVertex, 
53           const BOPCol_IndexedDataMapOfShapeListOfShape& theVertexEdgeMap,
54           const gp_Dir&                                  theStartDir,
55           const gp_Dir&                                  theNormal,
56                 TopoDS_Shape&                            theNextVertex,
57                 TopoDS_Shape&                            theNextEdge,
58                 gp_Dir&                                  theNextDir);
59
60 /// \brief Create planar face using the edges surrounding it
61 static void createFace(const TopoDS_Shape&                      theStartVertex,
62                        const std::list<TopoDS_Shape>::iterator& theStartEdge,
63                        const std::list<TopoDS_Shape>::iterator& theEndOfEdges,
64                        const gp_Pln&                            thePlane,
65                              TopoDS_Face&                       theResFace);
66
67 /// \bief Create planar wire
68 static void createWireList(const TopoDS_Shape&                      theStartVertex,
69                            const std::list<TopoDS_Shape>::iterator& theStartEdge,
70                            const std::list<TopoDS_Shape>::iterator& theEndOfEdges,
71                            const std::set<Handle(TopoDS_TShape)>&   theEdgesInLoops,
72                                  std::list<TopoDS_Wire>&            theResWires);
73
74 /// \brief Calculate outer tengency on the edge in specified vertex
75 static gp_Dir getOuterEdgeDirection(const TopoDS_Shape& theEdge,
76                                     const TopoDS_Shape& theVertex);
77
78 /// \brief Unnecessary edges will be removed from the map.
79 ///        Positions of iterator will be updated
80 static void removeWasteEdges(
81                 std::list<TopoDS_Shape>::iterator&  theStartVertex,
82                 std::list<TopoDS_Shape>::iterator&  theStartEdge,
83           const std::list<TopoDS_Shape>::iterator&  theEndOfVertexes,
84           const std::list<TopoDS_Shape>::iterator&  theEndOfEdges,
85           BOPCol_IndexedDataMapOfShapeListOfShape&  theMapVE);
86
87
88
89
90 void GeomAlgoAPI_SketchBuilder::createFaces(
91           const boost::shared_ptr<GeomAPI_Pnt>&                theOrigin,
92           const boost::shared_ptr<GeomAPI_Dir>&                theDirX,
93           const boost::shared_ptr<GeomAPI_Dir>&                theDirY,
94           const boost::shared_ptr<GeomAPI_Dir>&                theNorm,
95           const std::list< boost::shared_ptr<GeomAPI_Shape> >& theFeatures,
96                 std::list< boost::shared_ptr<GeomAPI_Shape> >& theResultFaces,
97                 std::list< boost::shared_ptr<GeomAPI_Shape> >& theResultWires)
98 {
99   if (theFeatures.empty())
100     return ;
101
102   // Create the list of edges with shared vertexes
103   BOPAlgo_Builder aBuilder;
104   BOPAlgo_PaveFiller aPF;
105   TopoDS_Shape aFeaturesCompound;
106
107   if (theFeatures.size() == 1)
108   { // If there is only one feature, BOPAlgo_Builder will decline to work. Need to process it anyway
109     aFeaturesCompound = theFeatures.front()->impl<TopoDS_Shape>();
110   }
111   else
112   {
113     std::list< boost::shared_ptr<GeomAPI_Shape> >::const_iterator anIt = theFeatures.begin();
114     for (; anIt != theFeatures.end(); anIt++)
115     {
116       boost::shared_ptr<GeomAPI_Shape> aPreview(*anIt);
117       aBuilder.AddArgument(aPreview->impl<TopoDS_Edge>());
118     }
119     aPF.SetArguments(aBuilder.Arguments());
120     aPF.Perform();
121     int aErr = aPF.ErrorStatus();
122     if (aErr) return ;
123     aBuilder.PerformWithFiller(aPF);
124     aErr = aBuilder.ErrorStatus();
125     if (aErr) return ;
126     aFeaturesCompound = aBuilder.Shape();
127   }
128
129   BOPCol_IndexedDataMapOfShapeListOfShape aMapVE; // map between vertexes and edges
130   BOPTools::MapShapesAndAncestors(aFeaturesCompound, TopAbs_VERTEX, TopAbs_EDGE, aMapVE);
131   if (aMapVE.IsEmpty()) // in case of not-initialized circle
132     return;
133
134   gp_Dir aDirX = theDirX->impl<gp_Dir>();
135   gp_Dir aDirY = theDirY->impl<gp_Dir>();
136   gp_Dir aNorm = theNorm->impl<gp_Dir>();
137
138   gp_Pln aPlane(theOrigin->impl<gp_Pnt>(), aNorm);
139
140   // Set of edges used in loops
141   std::set<Handle(TopoDS_TShape)> anEdgesInLoops;
142   // Lists for processed vertexes and edges
143   std::list<TopoDS_Shape> aProcVertexes;
144   std::list<TopoDS_Shape> aProcEdges;
145
146   // Search the start vertex
147   TopoDS_Shape aStartVertex = findStartVertex(aMapVE, aDirX, aDirY);
148   aProcVertexes.push_back(aStartVertex);
149
150   TopoDS_Shape aCurVertex = aStartVertex;
151   gp_Dir aCurDir = aDirY.Reversed();
152   gp_Dir aCurNorm = aNorm.Reversed();
153
154   // Go through the edges and find loops
155   TopoDS_Shape aNextVertex;
156   TopoDS_Shape aBindingEdge;
157   gp_Dir aNextDir;
158   while (aMapVE.Extent() > 0)
159   {
160     findNextVertex(aCurVertex, aMapVE, aCurDir, aCurNorm, aNextVertex, aBindingEdge, aNextDir);
161     aCurNorm = aNorm;
162
163     // Try to find next vertex in the list of already processed
164     std::list<TopoDS_Shape>::iterator          aVertIter = aProcVertexes.begin();
165     std::list<TopoDS_Shape>::iterator          anEdgeIter = aProcEdges.begin();
166     for ( ; aVertIter != aProcVertexes.end(); aVertIter++)
167     {
168       if (aVertIter->TShape() == aNextVertex.TShape())
169         break;
170       if (anEdgeIter != aProcEdges.end())
171         anEdgeIter++;
172     }
173
174     aProcVertexes.push_back(aNextVertex);
175     aProcEdges.push_back(aBindingEdge);
176
177     // The loop was found
178     if (aVertIter != aProcVertexes.end())
179     {
180       // If the binding edge is a full circle, then the list may be empty before addition. Need to update edge iterator
181       if (aProcEdges.size() == 1)
182         anEdgeIter = aProcEdges.begin();
183
184       if (aVertIter != aProcVertexes.begin())
185       {
186         // Check the orientation of the loop
187         gp_Dir aCN = getOuterEdgeDirection(*anEdgeIter, *aVertIter);
188         gp_Dir aCP = getOuterEdgeDirection(aProcEdges.back(), *aVertIter);
189         aCN.Reverse();
190         aCP.Reverse();
191         if (aCN.DotCross(aCP, aNorm) < -tolerance)
192         {
193           // The found loop has wrong orientation and may contain sub-loops.
194           // Need to check it onle again with another initial direction.
195           aCurVertex = *aVertIter;
196           do {
197             aProcVertexes.pop_back();
198             aProcEdges.pop_back();
199           } while (aCurVertex != aProcVertexes.back());
200           aCurDir = getOuterEdgeDirection(aProcEdges.back(), aCurVertex);
201           aCurNorm = aNorm.Reversed();
202           continue;
203         }
204       }
205
206       // When the orientation is correct or the edges looped through
207       // the first element, create new face and remove unnecessary edges.
208       TopoDS_Face aPatch;
209       createFace(*aVertIter, anEdgeIter, aProcEdges.end(), aPlane, aPatch);
210       if (!aPatch.IsNull())
211       {
212         boost::shared_ptr<GeomAPI_Shape> aFace(new GeomAPI_Shape);
213         aFace->setImpl(new TopoDS_Shape(aPatch));
214         theResultFaces.push_back(aFace);
215       }
216       // push the edges used in the loop to the map
217       std::list<TopoDS_Shape>::iterator anIter;
218       for (anIter = anEdgeIter; anIter != aProcEdges.end(); anIter++)
219         anEdgesInLoops.insert(anIter->TShape());
220       // remove unnecessary edges
221       std::list<TopoDS_Shape>::iterator aCopyVLoop = aVertIter;
222       std::list<TopoDS_Shape>::iterator aCopyELoop = anEdgeIter;
223       removeWasteEdges(aVertIter, anEdgeIter, aProcVertexes.end(), aProcEdges.end(), aMapVE);
224
225       // revert the list of remaining edges
226       std::list<TopoDS_Shape> aRemainVertexes;
227       for ( ; aVertIter != aProcVertexes.end(); aVertIter++)
228         aRemainVertexes.push_front(*aVertIter);
229       std::list<TopoDS_Shape> aRemainEdges;
230       for ( ; anEdgeIter != aProcEdges.end(); anEdgeIter++)
231         aRemainEdges.push_front(*anEdgeIter);
232       // remove edges and vertexes used in the loop and add remaining ones
233       aProcVertexes.erase(aCopyVLoop, aProcVertexes.end());
234       aProcVertexes.insert(aProcVertexes.end(), aRemainVertexes.begin(), aRemainVertexes.end());
235       aProcEdges.erase(aCopyELoop, aProcEdges.end());
236       aProcEdges.insert(aProcEdges.end(), aRemainEdges.begin(), aRemainEdges.end());
237
238       // Recalculate current vertex and current direction
239       if (!aProcVertexes.empty())
240       {
241         aNextVertex = aProcVertexes.back();
242         if (!aProcEdges.empty())
243           aNextDir = getOuterEdgeDirection(aProcEdges.back(), aNextVertex);
244         else aNextDir = aDirY;
245       }
246     }
247
248     // if next vertex connected only to alone edge, this is a part of wire (not a closed loop),
249     // we need to go back through the list of already checked edges to find a branching vertex
250     if (!aMapVE.IsEmpty() && aMapVE.Contains(aNextVertex) && 
251         aMapVE.FindFromKey(aNextVertex).Size() == 1)
252     {
253       std::list<TopoDS_Shape>::reverse_iterator aVRIter = aProcVertexes.rbegin();
254       std::list<TopoDS_Shape>::reverse_iterator aERIter = aProcEdges.rbegin();
255       if (aVRIter != aProcVertexes.rend())
256         aVRIter++;
257       if (aERIter != aProcEdges.rend())
258         aERIter++;
259
260       for ( ; aERIter != aProcEdges.rend(); aERIter++, aVRIter++)
261         if (aMapVE.FindFromKey(*aVRIter).Size() > 2)
262           break;
263       if (aERIter != aProcEdges.rend() || 
264          (aVRIter != aProcVertexes.rend() && aMapVE.FindFromKey(*aVRIter).Size() == 1))
265       { // the branching vertex was found or current list of edges is a wire without branches
266         std::list<TopoDS_Shape>::iterator aEIter;
267         TopoDS_Shape aCurEdge;
268         if (aERIter != aProcEdges.rend())
269         {
270           aEIter = aERIter.base();
271           aCurEdge = *aERIter;
272         }
273         else
274           aEIter = aProcEdges.begin();
275         std::list<TopoDS_Wire> aTail;
276         createWireList(*aVRIter, aEIter, aProcEdges.end(), anEdgesInLoops, aTail);
277         std::list<TopoDS_Wire>::const_iterator aTailIter = aTail.begin();
278         for ( ; aTailIter != aTail.end(); aTailIter++)
279           if (!aTailIter->IsNull())
280           {
281             boost::shared_ptr<GeomAPI_Shape> aWire(new GeomAPI_Shape);
282             aWire->setImpl(new TopoDS_Shape(*aTailIter));
283             theResultWires.push_back(aWire);
284           }
285         std::list<TopoDS_Shape>::iterator aVIter = aVRIter.base();
286         std::list<TopoDS_Shape>::iterator aEItCopy = aEIter;
287         removeWasteEdges(--aVIter, aEItCopy, aProcVertexes.end(), aProcEdges.end(), aMapVE);
288
289         aProcEdges.erase(aEIter, aProcEdges.end());
290         aVIter = aVRIter.base();
291         aProcVertexes.erase(aVIter, aProcVertexes.end());
292
293         if (!aProcVertexes.empty())
294         {
295           aNextVertex = aProcVertexes.back();
296           if (!aCurEdge.IsNull())
297             aNextDir = getOuterEdgeDirection(aCurEdge, aNextVertex);
298         }
299       }
300       else
301       { // there is no branching vertex in the list of proceeded,
302         // so we should revert the list and go the opposite way
303         aProcVertexes.reverse();
304         aProcEdges.reverse();
305         aNextVertex = aProcVertexes.back();
306         aNextDir = aProcEdges.empty() ? aDirY : 
307                    getOuterEdgeDirection(aProcEdges.back(), aNextVertex);
308       }
309     }
310
311     // When all edges of current component of connectivity are processed,
312     // we should repeat procedure for finding initial vertex in the remaining
313     if (!aMapVE.IsEmpty() && !aMapVE.Contains(aNextVertex))
314     {
315       aProcVertexes.clear();
316       aProcEdges.clear();
317
318       TopoDS_Shape aStartEdge;
319       aStartVertex = findStartVertex(aMapVE, aDirX, aDirY);
320       aProcVertexes.push_back(aStartVertex);
321       aNextVertex = aStartVertex;
322       aNextDir = aDirY.Reversed();
323       aCurNorm = aNorm.Reversed();
324     }
325
326     aCurVertex = aNextVertex;
327     aCurDir = aNextDir;
328   }
329
330   if (theResultFaces.size() > 1)
331     fixIntersections(theResultFaces);
332 }
333
334 void GeomAlgoAPI_SketchBuilder::fixIntersections(
335           std::list< boost::shared_ptr<GeomAPI_Shape> >& theFaces)
336 {
337   BRepClass_FaceClassifier aClassifier;
338
339   std::list< boost::shared_ptr<GeomAPI_Shape> >::iterator anIter1 = theFaces.begin();
340   std::list< boost::shared_ptr<GeomAPI_Shape> >::iterator anIter2;
341   for ( ; anIter1 != theFaces.end(); anIter1++)
342   {
343     anIter2 = anIter1;
344     for (++anIter2; anIter2 != theFaces.end(); anIter2++)
345     {
346       const TopoDS_Face& aF1 = (*anIter1)->impl<TopoDS_Face>();
347       TopExp_Explorer aVert2((*anIter2)->impl<TopoDS_Shape>(), TopAbs_VERTEX);
348       for ( ; aVert2.More(); aVert2.Next())
349       {
350         Handle(BRep_TVertex) aV = Handle(BRep_TVertex)::DownCast(aVert2.Current().TShape());
351         aClassifier.Perform(aF1, aV->Pnt(), tolerance);
352         if (aClassifier.State() != TopAbs_IN)
353           break;
354       }
355       if (aVert2.More())
356       { // second shape is not inside first, change the shapes order and repeat comparision
357         const TopoDS_Face& aF2 = (*anIter2)->impl<TopoDS_Face>();
358         TopExp_Explorer aVert1((*anIter1)->impl<TopoDS_Shape>(), TopAbs_VERTEX);
359         for ( ; aVert1.More(); aVert1.Next())
360         {
361           Handle(BRep_TVertex) aV = Handle(BRep_TVertex)::DownCast(aVert1.Current().TShape());
362           aClassifier.Perform(aF2, aV->Pnt(), tolerance);
363           if (aClassifier.State() != TopAbs_IN)
364             break;
365         }
366         if (!aVert1.More())
367         { // first shape should be cut from the second
368           BRepAlgoAPI_Cut aCut((*anIter2)->impl<TopoDS_Shape>(),
369                                (*anIter1)->impl<TopoDS_Shape>());
370           (*anIter2)->setImpl(new TopoDS_Shape(aCut.Shape()));
371         }
372       }
373       else
374       { // second shape should be cut from the first
375         BRepAlgoAPI_Cut aCut((*anIter1)->impl<TopoDS_Shape>(),
376                              (*anIter2)->impl<TopoDS_Shape>());
377         (*anIter1)->setImpl(new TopoDS_Shape(aCut.Shape()));
378       }
379     }
380   }
381 }
382
383
384 // =================== Auxiliary functions ====================================
385 const TopoDS_Shape& findStartVertex(
386           const BOPCol_IndexedDataMapOfShapeListOfShape& theMapVE,
387           const gp_Dir& theDirX, const gp_Dir& theDirY)
388 {
389   int aStartVertexInd = 1;
390   double aMaxX = -DBL_MAX;
391   double aMaxY = -DBL_MAX;
392   int aNbVert = theMapVE.Extent();
393   for (int i = 1; i <= aNbVert; i++)
394   {
395     const TopoDS_Vertex& aV = (const TopoDS_Vertex&)theMapVE.FindKey(i);
396     const Handle(BRep_TVertex)& aVert = (const Handle(BRep_TVertex)&)aV.TShape();
397     const gp_Pnt& aVertPnt = aVert->Pnt();
398
399     double aX = aVertPnt.XYZ().Dot(theDirX.XYZ());
400     double aY = aVertPnt.XYZ().Dot(theDirY.XYZ());
401     if ((aX > aMaxX || (fabs(aX - aMaxX) < tolerance && aY > aMaxY)) && 
402         theMapVE.FindFromIndex(i).Extent() > 1)
403     {
404       aMaxX = aX;
405       aMaxY = aY;
406       aStartVertexInd = i;
407     }
408   }
409   return theMapVE.FindKey(aStartVertexInd);
410 }
411
412 void findNextVertex(
413           const TopoDS_Shape&                            theStartVertex,
414           const BOPCol_IndexedDataMapOfShapeListOfShape& theVertexEdgeMap,
415           const gp_Dir&                                  theStartDir,
416           const gp_Dir&                                  theNormal,
417                 TopoDS_Shape&                            theNextVertex,
418                 TopoDS_Shape&                            theNextEdge,
419                 gp_Dir&                                  theNextDir)
420 {
421   const BOPCol_ListOfShape& anEdgesList = theVertexEdgeMap.FindFromKey(theStartVertex);
422   BOPCol_ListOfShape::Iterator aEdIter(anEdgesList);
423   double aBestEdgeProj = DBL_MAX;
424   for ( ; aEdIter.More(); aEdIter.Next())
425   {
426     gp_Dir aTang = getOuterEdgeDirection(aEdIter.Value(), theStartVertex);
427     aTang.Reverse();
428
429     // The projection is normalized in segment (-1, 1),
430     // where (-1, 0] corresponds to the angles (pi/2, 0] between theStartDir and aTang
431     // and [0, 1) corresponds to the angles [0, -pi/2)
432     double aProj = (aTang.Dot(theStartDir) - 1.0) * 0.5;
433     if (fabs(fabs(aProj) - 1) < tolerance)
434       continue;
435     if (theStartDir.DotCross(aTang, theNormal) < tolerance)
436       aProj *= -1.0;
437
438     if (aProj < aBestEdgeProj)
439     {
440       aBestEdgeProj = aProj;
441       theNextEdge = aEdIter.Value();
442       TopExp_Explorer aVertExp(theNextEdge, TopAbs_VERTEX);
443       for ( ; aVertExp.More(); aVertExp.Next())
444         if (aVertExp.Current().TShape() != theStartVertex.TShape())
445         {
446           theNextVertex = aVertExp.Current();
447           theNextDir = getOuterEdgeDirection(aEdIter.Value(), theNextVertex);
448           break;
449         }
450       if (!aVertExp.More())
451       { // This edge is a full circle
452         TopoDS_Vertex aV1, aV2;
453         TopExp::Vertices(*(const TopoDS_Edge*)(&theNextEdge), aV1, aV2);
454         if (aV1.Orientation() == theStartVertex.Orientation())
455           theNextVertex = aV2;
456         else
457           theNextVertex = aV1;
458         theNextDir = getOuterEdgeDirection(aEdIter.Value(), theNextVertex);
459       }
460     }
461   }
462 }
463
464 static void addEdgeToWire(const TopoDS_Edge&  theEdge,
465                           const BRep_Builder& theBuilder,
466                                 TopoDS_Shape& theSpliceVertex,
467                                 TopoDS_Wire&  theWire)
468 {
469   TopoDS_Edge anEdge = theEdge;
470   bool isCurVertChanged = false;
471   TopoDS_Shape aCurVertChanged;
472
473   TopExp_Explorer aVertExp(theEdge, TopAbs_VERTEX);
474   for ( ; aVertExp.More(); aVertExp.Next())
475   {
476     const TopoDS_Shape& aVertex = aVertExp.Current();
477     if (aVertex.TShape() == theSpliceVertex.TShape() && 
478         aVertex.Orientation() != theEdge.Orientation())
479     { // Current vertex is the last for the edge, so its orientation is wrong, need to revert the edge
480       anEdge.Reverse();
481       break;
482     }
483     if (aVertex.TShape() != theSpliceVertex.TShape())
484     {
485       aCurVertChanged = aVertex;
486       isCurVertChanged = true;
487     }
488   }
489   theSpliceVertex = isCurVertChanged ? aCurVertChanged : aVertExp.Current();
490
491   theBuilder.Add(theWire, anEdge);
492 }
493
494 void createFace(const TopoDS_Shape&                      theStartVertex,
495                 const std::list<TopoDS_Shape>::iterator& theStartEdge,
496                 const std::list<TopoDS_Shape>::iterator& theEndOfEdges,
497                 const gp_Pln&                            thePlane,
498                       TopoDS_Face&                       theResFace)
499 {
500   TopoDS_Wire aResWire;
501   BRep_Builder aBuilder;
502   aBuilder.MakeWire(aResWire);
503
504   TopoDS_Shape aCurVertex = theStartVertex;
505   std::list<TopoDS_Shape>::const_iterator anEdgeIter = theStartEdge;
506   for ( ; anEdgeIter != theEndOfEdges; anEdgeIter++)
507   {
508     TopoDS_Edge anEdge = *((TopoDS_Edge*)(&(*anEdgeIter)));
509     if (!anEdge.IsNull())
510       addEdgeToWire(anEdge, aBuilder, aCurVertex, aResWire);
511   }
512
513   BRepBuilderAPI_MakeFace aFaceBuilder(thePlane, aResWire);
514   if (aFaceBuilder.Error() == BRepBuilderAPI_FaceDone)
515     theResFace = aFaceBuilder.Face();
516 }
517
518 void createWireList(const TopoDS_Shape&                      theStartVertex,
519                     const std::list<TopoDS_Shape>::iterator& theStartEdge,
520                     const std::list<TopoDS_Shape>::iterator& theEndOfEdges,
521                     const std::set<Handle(TopoDS_TShape)>&   theEdgesInLoops,
522                           std::list<TopoDS_Wire>&            theResWires)
523 {
524   BRep_Builder aBuilder;
525   bool needNewWire = true;
526   TopoDS_Shape aCurVertex = theStartVertex;
527
528   std::list<TopoDS_Shape>::iterator anIter = theStartEdge;
529   while (anIter != theEndOfEdges)
530   {
531     while (anIter != theEndOfEdges && needNewWire && 
532            theEdgesInLoops.count(anIter->TShape()) != 0)
533     {
534       TopExp_Explorer aVertExp(*anIter, TopAbs_VERTEX);
535       for ( ; aVertExp.More(); aVertExp.Next())
536         if (aVertExp.Current().TShape() != aCurVertex.TShape())
537         {
538           aCurVertex = aVertExp.Current();
539           break;
540         }
541       anIter++;
542     }
543     if (anIter == theEndOfEdges)
544       break;
545
546     if (needNewWire)
547     { // The new wire should be created
548       TopoDS_Wire aWire;
549       aBuilder.MakeWire(aWire);
550       theResWires.push_back(aWire);
551       needNewWire = false;
552     }
553     else if (theEdgesInLoops.count(anIter->TShape()) != 0)
554     { // There was found the edge already used in loop.
555       // Current wire should be released and new one should started
556       needNewWire = true;
557       continue;
558     }
559
560     TopoDS_Edge anEdge = *((TopoDS_Edge*)(&(*anIter)));
561     addEdgeToWire(anEdge, aBuilder, aCurVertex, theResWires.back());
562     anIter++;
563   }
564 }
565
566
567 gp_Dir getOuterEdgeDirection(const TopoDS_Shape& theEdge,
568                              const TopoDS_Shape& theVertex)
569 {
570   const Handle(BRep_TVertex)& aVertex = (const Handle(BRep_TVertex)&)theVertex.TShape();
571   gp_Pnt aVertexPnt = aVertex->Pnt();
572
573   const Handle(BRep_TEdge)& anEdge = (const Handle(BRep_TEdge)&)theEdge.TShape();
574
575   // Convert the edge to the curve to calculate the tangency.
576   // There should be only one curve in the edge.
577   Handle(BRep_Curve3D) aEdCurve = 
578     Handle(BRep_Curve3D)::DownCast(anEdge->Curves().First());
579   double aFirst, aLast;
580   aEdCurve->Range(aFirst, aLast);
581   Handle(Geom_Curve) aCurve = aEdCurve->Curve3D();
582
583   gp_Pnt aPnt;
584   gp_Vec aTang;
585   aCurve->D1(aFirst + shift, aPnt, aTang);
586   aCurve->D0(aFirst, aPnt);
587   if (aVertexPnt.IsEqual(aPnt, tolerance))
588     return gp_Dir(aTang.Reversed());
589
590   aCurve->D1(aLast - shift, aPnt, aTang);
591   return gp_Dir(aTang);
592 }
593
594 void removeWasteEdges(
595                 std::list<TopoDS_Shape>::iterator&  theStartVertex,
596                 std::list<TopoDS_Shape>::iterator&  theStartEdge,
597           const std::list<TopoDS_Shape>::iterator&  theEndOfVertexes,
598           const std::list<TopoDS_Shape>::iterator&  theEndOfEdges,
599           BOPCol_IndexedDataMapOfShapeListOfShape&  theMapVE)
600 {
601   bool isVertStep = true;
602   while (theStartVertex != theEndOfVertexes && theStartEdge != theEndOfEdges)
603   {
604     BOPCol_ListOfShape& aBunch = theMapVE.ChangeFromKey(*theStartVertex);
605     BOPCol_ListOfShape::Iterator anApprEdge(aBunch);
606     for ( ; anApprEdge.More(); anApprEdge.Next())
607       if (anApprEdge.Value() == *theStartEdge)
608         break;
609     if (anApprEdge.More())
610       aBunch.Remove(anApprEdge);
611
612     if (isVertStep)
613       theStartVertex++;
614     else
615     {
616       theStartEdge++;
617       // check current vertex to be a branching point
618       // if so, it will be a new starting point to find a loop
619       if (aBunch.Size() > 1)
620         break;
621     }
622     isVertStep = !isVertStep;
623   }
624
625   // The map of vertex-edges may be changed
626   BOPCol_IndexedDataMapOfShapeListOfShape aMapVECopy;
627   BOPCol_IndexedDataMapOfShapeListOfShape::Iterator aMapIter(theMapVE);
628   for (int ind = 1; aMapIter.More(); aMapIter.Next(), ind++)
629     if (!aMapIter.Value().IsEmpty())
630       aMapVECopy.Add(theMapVE.FindKey(ind), aMapIter.Value());
631   theMapVE.Clear();
632   theMapVE.Exchange(aMapVECopy);
633 }
634