6489c4ec53e4950f69a8a8396900b69d3ce9cb79
[occt.git] / src / ShapeUpgrade / ShapeUpgrade_UnifySameDomain.cxx
1 // Created on: 2012-06-09
2 // Created by: jgv@ROLEX
3 // Copyright (c) 2012-2014 OPEN CASCADE SAS
4 //
5 // This file is part of Open CASCADE Technology software library.
6 //
7 // This library is free software; you can redistribute it and/or modify it under
8 // the terms of the GNU Lesser General Public License version 2.1 as published
9 // by the Free Software Foundation, with special exception defined in the file
10 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
11 // distribution for complete text of the license and disclaimer of any warranty.
12 //
13 // Alternatively, this file may be used under the terms of Open CASCADE
14 // commercial license or contractual agreement.
15
16
17 #include <BRep_Builder.hxx>
18 #include <BRep_Tool.hxx>
19 #include <BRepLib.hxx>
20 #include <BRepLib_MakeEdge.hxx>
21 #include <BRepLib_MakeFace.hxx>
22 #include <BRepTopAdaptor_TopolTool.hxx>
23 #include <GC_MakeCircle.hxx>
24 #include <Geom2d_Line.hxx>
25 #include <GCE2d_MakeLine.hxx>
26 #include <Geom2d_TrimmedCurve.hxx>
27 #include <Geom2dConvert.hxx>
28 #include <Geom2dConvert_CompCurveToBSplineCurve.hxx>
29 #include <Geom_BezierCurve.hxx>
30 #include <Geom_BSplineCurve.hxx>
31 #include <Geom_Circle.hxx>
32 #include <Geom_CylindricalSurface.hxx>
33 #include <Geom_ElementarySurface.hxx>
34 #include <Geom_Line.hxx>
35 #include <Geom_OffsetSurface.hxx>
36 #include <Geom_Plane.hxx>
37 #include <Geom_RectangularTrimmedSurface.hxx>
38 #include <Geom_Surface.hxx>
39 #include <Geom_SurfaceOfLinearExtrusion.hxx>
40 #include <Geom_SurfaceOfRevolution.hxx>
41 #include <Geom_SweptSurface.hxx>
42 #include <Geom_TrimmedCurve.hxx>
43 #include <GeomAdaptor_HSurface.hxx>
44 #include <GeomConvert.hxx>
45 #include <GeomConvert_CompCurveToBSplineCurve.hxx>
46 #include <GeomLib_IsPlanarSurface.hxx>
47 #include <gp_Cylinder.hxx>
48 #include <gp_Dir.hxx>
49 #include <gp_Lin.hxx>
50 #include <IntPatch_ImpImpIntersection.hxx>
51 #include <ShapeAnalysis_Edge.hxx>
52 #include <ShapeAnalysis_WireOrder.hxx>
53 #include <ShapeBuild_Edge.hxx>
54 #include <ShapeBuild_ReShape.hxx>
55 #include <ShapeFix_Edge.hxx>
56 #include <ShapeFix_Face.hxx>
57 #include <ShapeFix_Shell.hxx>
58 #include <ShapeFix_Wire.hxx>
59 #include <ShapeUpgrade_UnifySameDomain.hxx>
60 #include <Standard_Type.hxx>
61 #include <TColGeom2d_Array1OfBSplineCurve.hxx>
62 #include <TColGeom2d_HArray1OfBSplineCurve.hxx>
63 #include <TColGeom2d_SequenceOfBoundedCurve.hxx>
64 #include <TColGeom2d_SequenceOfCurve.hxx>
65 #include <TColGeom_Array1OfBSplineCurve.hxx>
66 #include <TColGeom_HArray1OfBSplineCurve.hxx>
67 #include <TColGeom_SequenceOfSurface.hxx>
68 #include <TColStd_Array1OfReal.hxx>
69 #include <TColStd_MapOfInteger.hxx>
70 #include <TopExp.hxx>
71 #include <TopExp_Explorer.hxx>
72 #include <TopoDS.hxx>
73 #include <TopoDS_Edge.hxx>
74 #include <TopoDS_Face.hxx>
75 #include <TopoDS_Shape.hxx>
76 #include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
77 #include <TopTools_IndexedMapOfShape.hxx>
78 #include <TopTools_ListIteratorOfListOfShape.hxx>
79 #include <TopTools_MapOfShape.hxx>
80 #include <TopTools_SequenceOfShape.hxx>
81 #include <gp_Circ.hxx>
82 #include <BRepAdaptor_Curve.hxx>
83 #include <BRepAdaptor_Curve2d.hxx>
84 #include <BRepAdaptor_Surface.hxx>
85 #include <gp_Vec2d.hxx>
86 #include <Extrema_ExtPS.hxx>
87 #include <BRepTools.hxx>
88 #include <BRepTopAdaptor_FClass2d.hxx>
89
90 IMPLEMENT_STANDARD_RTTIEXT(ShapeUpgrade_UnifySameDomain,Standard_Transient)
91
92 struct SubSequenceOfEdges
93 {
94   TopTools_SequenceOfShape SeqsEdges;
95   TopoDS_Edge UnionEdges;
96 };
97
98 static Standard_Real TrueValueOfOffset(const Standard_Real theValue,
99                                        const Standard_Real thePeriod)
100 {
101   if (theValue > 0.)
102     return thePeriod;
103
104   return (-thePeriod);
105 }
106
107 //=======================================================================
108 //function : UpdateBoundaries
109 //purpose  : auxilary
110 //=======================================================================
111 static void UpdateBoundaries(const Handle(Geom2d_Curve)& thePCurve,
112                              const Standard_Real         theFirst,
113                              const Standard_Real         theLast,
114                              Standard_Real& theUfirst,
115                              Standard_Real& theUlast)
116 {
117   const Standard_Integer NbSamples = 4;
118   Standard_Real Delta = (theLast - theFirst)/NbSamples;
119
120   for (Standard_Integer i = 0; i <= NbSamples; i++)
121   {
122     Standard_Real aParam = theFirst + i*Delta;
123     gp_Pnt2d aPoint = thePCurve->Value(aParam);
124     if (aPoint.X() < theUfirst)
125       theUfirst = aPoint.X();
126     if (aPoint.X() > theUlast)
127       theUlast  = aPoint.X();
128   }
129 }
130
131 static void RemoveEdgeFromMap(const TopoDS_Edge& theEdge,
132                               TopTools_IndexedDataMapOfShapeListOfShape& theVEmap)
133 {
134   TopoDS_Vertex VV [2];
135   TopExp::Vertices(theEdge, VV[0], VV[1]);
136   for (Standard_Integer i = 0; i < 2; i++)
137   {
138     TopTools_ListOfShape& Elist = theVEmap.ChangeFromKey(VV[i]);
139     TopTools_ListIteratorOfListOfShape itl(Elist);
140     while (itl.More())
141     {
142       const TopoDS_Shape& anEdge = itl.Value();
143       if (anEdge.IsSame(theEdge))
144         Elist.Remove(itl);
145       else
146         itl.Next();
147     }
148   }
149 }
150
151 static Standard_Real ComputeMinEdgeSize(const TopTools_SequenceOfShape& theEdges,
152                                         const TopoDS_Face&              theRefFace,
153                                         TopTools_MapOfShape&            theEdgesMap)
154 {
155   Standard_Real MinSize = RealLast();
156   
157   for (Standard_Integer ind = 1; ind <= theEdges.Length(); ind++)
158   {
159     const TopoDS_Edge& anEdge = TopoDS::Edge(theEdges(ind));
160     theEdgesMap.Add(anEdge);
161     TopoDS_Vertex V1, V2;
162     TopExp::Vertices(anEdge, V1, V2);
163     BRepAdaptor_Curve2d BAcurve2d(anEdge, theRefFace);
164     gp_Pnt2d FirstP2d = BAcurve2d.Value(BAcurve2d.FirstParameter());
165     gp_Pnt2d LastP2d  = BAcurve2d.Value(BAcurve2d.LastParameter());
166     Standard_Real aSqDist;
167     if (V1.IsSame(V2) &&
168         !BRep_Tool::Degenerated(anEdge))
169     {
170       gp_Pnt2d MidP2d = BAcurve2d.Value((BAcurve2d.FirstParameter()+BAcurve2d.LastParameter())/2);
171       aSqDist = FirstP2d.SquareDistance(MidP2d);
172     }
173     else
174       aSqDist = FirstP2d.SquareDistance(LastP2d);
175     if (aSqDist < MinSize)
176       MinSize = aSqDist;
177   }
178   MinSize = Sqrt(MinSize);
179   return MinSize;
180 }
181
182 static void FindFaceWithMaxUbound(const TopTools_SequenceOfShape& theFaces,
183                                   const TopoDS_Face&              theRefFace,
184                                   const TopTools_MapOfShape&      theEdgesMap,
185                                   Standard_Real&                  theUmax,
186                                   Standard_Integer&               theIndFaceMax)
187 {
188   theUmax = RealFirst();
189   theIndFaceMax = 0;
190   for (Standard_Integer ii = 1; ii <= theFaces.Length(); ii++)
191   {
192     const TopoDS_Face& face_ii = TopoDS::Face(theFaces(ii));
193     Standard_Real Ufirst = RealLast(), Ulast = RealFirst();
194     TopExp_Explorer Explo(face_ii, TopAbs_EDGE);
195     for (; Explo.More(); Explo.Next())
196     {
197       const TopoDS_Edge& anEdge = TopoDS::Edge(Explo.Current());
198       if (!theEdgesMap.Contains(anEdge))
199         continue;
200       Standard_Real fpar, lpar;
201       Handle(Geom2d_Curve) aPCurve = BRep_Tool::CurveOnSurface(anEdge, theRefFace, fpar, lpar);
202       UpdateBoundaries(aPCurve, fpar, lpar, Ufirst, Ulast);
203     }
204     if (Ulast > theUmax)
205     {
206       theUmax = Ulast;
207       theIndFaceMax = ii;
208     }
209   }
210 }
211
212 static void RelocatePCurvesToNewUorigin(const TopTools_SequenceOfShape& theEdges,
213                                         const TopoDS_Shape&             theFirstFace,
214                                         const TopoDS_Face&              theRefFace,
215                                         const Standard_Real             theCoordTol,
216                                         const Standard_Real             theUperiod,
217                                         TopTools_IndexedDataMapOfShapeListOfShape& theVEmap,
218                                         NCollection_DataMap<TopoDS_Shape, Handle(Geom2d_Curve)>& theEdgeNewPCurve,
219                                         TopTools_MapOfShape&      theUsedEdges)
220 {
221   TopTools_MapOfShape EdgesOfFirstFace;
222   TopExp::MapShapes(theFirstFace, EdgesOfFirstFace);
223   
224   for (;;) //walk by contours
225   {
226     //try to find the start edge that:
227     //1. belongs to outer edges of first face;
228     //2. has minimum U-coord of its start point
229     TopoDS_Edge StartEdge;
230     TopAbs_Orientation anOr = TopAbs_FORWARD;
231     Standard_Real anUmin = RealLast();
232     for (Standard_Integer ii = 1; ii <= theEdges.Length(); ii++)
233     {
234       const TopoDS_Edge& anEdge = TopoDS::Edge(theEdges(ii));
235       if (theUsedEdges.Contains(anEdge))
236         continue;
237       if (EdgesOfFirstFace.Contains(anEdge))
238       {
239         if (StartEdge.IsNull())
240         {
241           StartEdge = anEdge;
242           BRepAdaptor_Curve2d StartBAcurve(StartEdge, theRefFace);
243           Standard_Real aFirstParam, aLastParam;
244           if (StartEdge.Orientation() == TopAbs_FORWARD)
245           {
246             aFirstParam = StartBAcurve.FirstParameter();
247             aLastParam  = StartBAcurve.LastParameter();
248           }
249           else
250           {
251             aFirstParam = StartBAcurve.LastParameter();
252             aLastParam  = StartBAcurve.FirstParameter();
253           }
254           gp_Pnt2d aFirstPoint = StartBAcurve.Value(aFirstParam);
255           gp_Pnt2d aLastPoint  = StartBAcurve.Value(aLastParam);
256           if (aFirstPoint.X() < aLastPoint.X())
257           {
258             anUmin = aFirstPoint.X();
259             anOr = TopAbs_FORWARD;
260           }
261           else
262           {
263             anUmin = aLastPoint.X();
264             anOr = TopAbs_REVERSED;
265           }
266         }
267         else
268         {
269           BRepAdaptor_Curve2d aBAcurve(anEdge, theRefFace);
270           Standard_Real aFirstParam, aLastParam;
271           if (anEdge.Orientation() == TopAbs_FORWARD)
272           {
273             aFirstParam = aBAcurve.FirstParameter();
274             aLastParam  = aBAcurve.LastParameter();
275           }
276           else
277           {
278             aFirstParam = aBAcurve.LastParameter();
279             aLastParam  = aBAcurve.FirstParameter();
280           }
281           gp_Pnt2d aFirstPoint = aBAcurve.Value(aFirstParam);
282           gp_Pnt2d aLastPoint  = aBAcurve.Value(aLastParam);
283           if (aFirstPoint.X() < anUmin)
284           {
285             StartEdge = anEdge;
286             anUmin = aFirstPoint.X();
287             anOr = TopAbs_FORWARD;
288           }
289           if (aLastPoint.X() < anUmin)
290           {
291             StartEdge = anEdge;
292             anUmin = aLastPoint.X();
293             anOr = TopAbs_REVERSED;
294           }
295         }
296       } //if (EdgesOfFirstFace.Contains(anEdge))
297     } //for (Standard_Integer ii = 1; ii <= edges.Length(); ii++)
298     
299     if (StartEdge.IsNull()) //all contours are passed
300       break;
301     
302     TopoDS_Vertex StartVertex = (anOr == TopAbs_FORWARD)?
303       TopExp::FirstVertex(StartEdge, Standard_True) : TopExp::LastVertex(StartEdge, Standard_True);
304     TopoDS_Edge CurEdge = StartEdge;
305     Standard_Real fpar, lpar;
306     Handle(Geom2d_Curve) CurPCurve = BRep_Tool::CurveOnSurface(CurEdge, theRefFace, fpar, lpar);
307     CurPCurve = new Geom2d_TrimmedCurve(CurPCurve, fpar, lpar);
308     theEdgeNewPCurve.Bind(CurEdge, CurPCurve);
309     if (CurEdge.Orientation() == TopAbs_REVERSED)
310     { Standard_Real tmp = fpar; fpar = lpar; lpar = tmp; }
311     //Standard_Real StartParam = (anOr == TopAbs_FORWARD)? fpar : lpar;
312     Standard_Real CurParam = (anOr == TopAbs_FORWARD)? lpar : fpar;
313     //gp_Pnt2d StartPoint = CurPCurve->Value(StartParam);
314     gp_Pnt2d CurPoint = CurPCurve->Value(CurParam);
315     for (;;) //collect pcurves of a contour
316     {
317       RemoveEdgeFromMap(CurEdge, theVEmap);
318       theUsedEdges.Add(CurEdge);
319       TopoDS_Vertex CurVertex = (anOr == TopAbs_FORWARD)?
320         TopExp::LastVertex(CurEdge, Standard_True) : TopExp::FirstVertex(CurEdge, Standard_True);
321       
322       const TopTools_ListOfShape& Elist = theVEmap.FindFromKey(CurVertex);
323       if (Elist.IsEmpty())
324         break; //end of contour in 3d
325       
326       TopTools_ListIteratorOfListOfShape itl(Elist);
327       for (; itl.More(); itl.Next())
328       {
329         const TopoDS_Edge& anEdge = TopoDS::Edge(itl.Value());
330         if (anEdge.IsSame(CurEdge))
331           continue;
332         TopoDS_Vertex aFirstVertex = (anOr == TopAbs_FORWARD)?
333           TopExp::FirstVertex(anEdge, Standard_True) : TopExp::LastVertex(anEdge, Standard_True);
334         if (!aFirstVertex.IsSame(CurVertex)) //may be if CurVertex is deg.vertex
335           continue;
336         
337         Handle(Geom2d_Curve) aPCurve = BRep_Tool::CurveOnSurface(anEdge, theRefFace, fpar, lpar);
338         aPCurve = new Geom2d_TrimmedCurve(aPCurve, fpar, lpar);
339         if (anEdge.Orientation() == TopAbs_REVERSED)
340         { Standard_Real tmp = fpar; fpar = lpar; lpar = tmp; }
341         Standard_Real aParam = (anOr == TopAbs_FORWARD)? fpar : lpar;
342         gp_Pnt2d aPoint = aPCurve->Value(aParam);
343         Standard_Real anOffset = CurPoint.X() - aPoint.X();
344         if (!(Abs(anOffset) < theCoordTol ||
345               Abs(Abs(anOffset) - theUperiod) < theCoordTol))
346           continue; //may be if CurVertex is deg.vertex
347         
348         if (Abs(anOffset) > theUperiod/2)
349         {
350           anOffset = TrueValueOfOffset(anOffset, theUperiod);
351           gp_Vec2d aVec(anOffset, 0.);
352           Handle(Geom2d_Curve) aNewPCurve = Handle(Geom2d_Curve)::DownCast(aPCurve->Copy());
353           aNewPCurve->Translate(aVec);
354           aPCurve = aNewPCurve;
355         }
356         theEdgeNewPCurve.Bind(anEdge, aPCurve);
357         CurEdge = anEdge;
358         TopAbs_Orientation CurOr = TopAbs::Compose(anOr, CurEdge.Orientation());
359         CurParam = (CurOr == TopAbs_FORWARD)?
360           aPCurve->LastParameter() : aPCurve->FirstParameter();
361         CurPoint = aPCurve->Value(CurParam);
362         break;
363       }
364     } //for (;;) (collect pcurves of a contour)
365   } //for (;;) (walk by contours)
366 }
367
368 static void InsertWiresIntoFaces(const TopTools_SequenceOfShape& theWires,
369                                  const TopTools_SequenceOfShape& theFaces,
370                                  const TopoDS_Face&              theRefFace)
371 {
372   BRep_Builder BB;
373   for (Standard_Integer ii = 1; ii <= theWires.Length(); ii++)
374   {
375     const TopoDS_Wire& aWire = TopoDS::Wire(theWires(ii));
376     TopoDS_Iterator iter(aWire);
377     const TopoDS_Edge& anEdge = TopoDS::Edge(iter.Value());
378     BRepAdaptor_Curve2d BAcurve2d(anEdge, theRefFace);
379     gp_Pnt2d aPnt2d = BAcurve2d.Value((BAcurve2d.FirstParameter() + BAcurve2d.LastParameter())/2.);
380     TopoDS_Shape RequiredFace;
381     for (Standard_Integer jj = 1; jj <= theFaces.Length(); jj++)
382     {
383       const TopoDS_Face& aFace = TopoDS::Face(theFaces(jj));
384       BRepTopAdaptor_FClass2d Classifier(aFace, Precision::Confusion());
385       TopAbs_State aStatus = Classifier.Perform(aPnt2d);
386       if (aStatus == TopAbs_IN)
387       {
388         RequiredFace = aFace.Oriented (TopAbs_FORWARD);
389         break;
390       }
391     }
392     if (!RequiredFace.IsNull())
393     {
394       BB.Add(RequiredFace, aWire);
395     }
396     else
397     {
398       Standard_ASSERT_INVOKE ("ShapeUpgrade_UnifySameDomain: wire remains unclassified");
399     }
400   }
401 }
402
403 static TopoDS_Face FindCommonFace(const TopoDS_Edge&  theEdge1,
404                                   const TopoDS_Edge&  theEdge2,
405                                   const TopTools_IndexedDataMapOfShapeListOfShape& theVFmap,
406                                   TopAbs_Orientation& theOrOfE1OnFace,
407                                   TopAbs_Orientation& theOrOfE2OnFace)
408 {
409   TopoDS_Vertex aVertex;
410   TopExp::CommonVertex(theEdge1, theEdge2, aVertex);
411   const TopTools_ListOfShape& Flist = theVFmap.FindFromKey(aVertex);
412   TopTools_ListIteratorOfListOfShape itl(Flist);
413   for (; itl.More(); itl.Next())
414   {
415     TopoDS_Face aFace = TopoDS::Face(itl.Value());
416     aFace.Orientation(TopAbs_FORWARD);
417     Standard_Boolean e1found = Standard_False, e2found = Standard_False;
418     TopExp_Explorer Explo(aFace, TopAbs_EDGE);
419     for (; Explo.More(); Explo.Next())
420     {
421       const TopoDS_Shape& anEdge = Explo.Current();
422       if (anEdge.IsSame(theEdge1))
423       {
424         e1found = Standard_True;
425         theOrOfE1OnFace = anEdge.Orientation();
426       }
427       if (anEdge.IsSame(theEdge2))
428       {
429         e2found = Standard_True;
430         theOrOfE2OnFace = anEdge.Orientation();
431       }
432       if (e1found && e2found)
433         return aFace;
434     }
435   }
436
437   TopoDS_Face NullFace;
438   return NullFace;
439 }
440
441 static Standard_Boolean FindClosestPoints(const TopoDS_Edge& theEdge1,
442                                           const TopoDS_Edge& theEdge2,
443                                           const TopTools_IndexedDataMapOfShapeListOfShape& theVFmap,
444                                           TopoDS_Face& theCommonFace,
445                                           Standard_Real& theMinSqDist,
446                                           TopAbs_Orientation& OrOfE1OnFace,
447                                           TopAbs_Orientation& OrOfE2OnFace,
448                                           Standard_Integer& theIndOnE1,
449                                           Standard_Integer& theIndOnE2,
450                                           gp_Pnt2d* PointsOnEdge1,
451                                           gp_Pnt2d* PointsOnEdge2)
452                                 
453 {
454   theCommonFace = FindCommonFace(theEdge1, theEdge2, theVFmap, OrOfE1OnFace, OrOfE2OnFace);
455   if (theCommonFace.IsNull())
456     return Standard_False;
457   
458   Standard_Real fpar1, lpar1, fpar2, lpar2;
459   Handle(Geom2d_Curve) PCurve1 = BRep_Tool::CurveOnSurface(theEdge1, theCommonFace, fpar1, lpar1);
460   Handle(Geom2d_Curve) PCurve2 = BRep_Tool::CurveOnSurface(theEdge2, theCommonFace, fpar2, lpar2);
461   PointsOnEdge1[0] = PCurve1->Value(fpar1);
462   PointsOnEdge1[1] = PCurve1->Value(lpar1);
463   PointsOnEdge2[0] = PCurve2->Value(fpar2);
464   PointsOnEdge2[1] = PCurve2->Value(lpar2);
465   theMinSqDist = RealLast();
466   theIndOnE1 = -1;
467   theIndOnE2 = -1;
468   for (Standard_Integer ind1 = 0; ind1 < 2; ind1++)
469     for (Standard_Integer ind2 = 0; ind2 < 2; ind2++)
470     {
471       Standard_Real aSqDist = PointsOnEdge1[ind1].SquareDistance(PointsOnEdge2[ind2]);
472       if (aSqDist < theMinSqDist)
473       {
474         theMinSqDist = aSqDist;
475         theIndOnE1 = ind1; theIndOnE2 = ind2;
476       }
477     }
478   return Standard_True;
479 }
480
481 //=======================================================================
482 //function : ReconstructMissedSeam
483 //purpose  : auxilary
484 //=======================================================================
485 static void ReconstructMissedSeam(const TopTools_SequenceOfShape& theEdges,
486                                   const TopTools_MapOfShape&      theUsedEdges,
487                                   const TopoDS_Face&              theFrefFace,
488                                   const TopoDS_Vertex&            theCurVertex,
489                                   const gp_Pnt2d&                 theCurPoint,
490                                   const Standard_Real             theUperiod,
491                                   const Standard_Real             theFaceUmin,
492                                   const Standard_Real             theCoordTol,
493                                   TopoDS_Edge&                    theNextEdge,
494                                   TopoDS_Wire&                    theNewWire,
495                                   gp_Pnt2d&                       theNextPoint,
496                                   gp_Pnt2d&                       theStartOfNextEdge,
497                                   TopoDS_Vertex&                  theLastVertexOfSeam,
498                                   TopTools_IndexedDataMapOfShapeListOfShape& theVEmap)
499                                   
500 {
501   Handle(Geom_Surface) RefSurf = BRep_Tool::Surface(theFrefFace);
502   GeomAbs_Shape aContinuity = (RefSurf->IsUPeriodic())? GeomAbs_CN : GeomAbs_C0;
503   
504   Standard_Real Ydir = 1.; //up
505   if (Abs(theCurPoint.X() - theFaceUmin) <= theCoordTol)
506     Ydir = -1.; //down
507
508   //Consider as the candidate to be next edge:
509   //only the edge that has first point with X-coordinate close to X-coordinate of theCurPoint
510   //Choose from candidates the edge that is closest to theCurPoint in the defined direction Ydir
511   Standard_Real MinDeltaY = RealLast();
512   for (Standard_Integer ind = 1; ind <= theEdges.Length(); ind++)
513   {
514     const TopoDS_Edge& aCandidate = TopoDS::Edge(theEdges(ind));
515     if (theUsedEdges.Contains(aCandidate))
516       continue;
517     BRepAdaptor_Curve2d BAcurve2d(aCandidate, theFrefFace);
518     Standard_Real CandParam = (aCandidate.Orientation() == TopAbs_FORWARD)?
519       BAcurve2d.FirstParameter() : BAcurve2d.LastParameter();
520     gp_Pnt2d CandPoint = BAcurve2d.Value(CandParam);
521     Standard_Real DeltaX = Abs(CandPoint.X() - theCurPoint.X());
522     if (DeltaX > theCoordTol)
523       continue;
524     
525     Standard_Real DeltaY = CandPoint.Y() - theCurPoint.Y();
526     DeltaY *= Ydir;
527     if (DeltaY < 0.) //on the other side from CurPoint
528       continue;
529     if (DeltaY < MinDeltaY)
530     {
531       MinDeltaY = DeltaY;
532       theNextEdge = aCandidate;
533       theStartOfNextEdge = CandPoint;
534     }
535   }
536   
537   //Build missed seam edge
538   theLastVertexOfSeam = TopExp::FirstVertex(theNextEdge, Standard_True); //with orientation
539   Standard_Real CurTol  = BRep_Tool::Tolerance(theCurVertex);
540   Standard_Real LastTol = BRep_Tool::Tolerance(theLastVertexOfSeam);
541   Standard_Real anU = (CurTol < LastTol)? theCurPoint.X() : theStartOfNextEdge.X();
542   Handle(Geom_Curve) Uiso = RefSurf->UIso(anU);
543   TopoDS_Vertex V1, V2;
544   Standard_Real Param1, Param2;
545   if (Ydir > 0)
546   {
547     V1 = theCurVertex; V2 = theLastVertexOfSeam;
548     Param1 = theCurPoint.Y(); Param2 = theStartOfNextEdge.Y();
549   }
550   else
551   {
552     V1 = theLastVertexOfSeam; V2 = theCurVertex;
553     Param1 = theStartOfNextEdge.Y(); Param2 = theCurPoint.Y();
554   }
555   TopoDS_Edge MissedSeam = BRepLib_MakeEdge(Uiso, V1, V2, Param1, Param2);
556   Standard_Real Vorigin = 0.;
557   //Correct Param1 and Param2 if needed:
558   //when Uiso-curve is periodic and Param1 and Param2 do not fit into V-range of surface,
559   //BRepLib_MakeEdge may shift Param1 and Param2
560   Standard_Real InitialParam1 = Param1, InitialParam2 = Param2;
561   Handle(Geom_Curve) MissedCurve = BRep_Tool::Curve(MissedSeam, Param1, Param2);
562   if ((Param1 != InitialParam1 || Param2 != InitialParam2) &&
563       MissedCurve->IsPeriodic())
564   {
565     //Vorigin = -(MissedCurve->Period());
566     Vorigin = -(Param1 - InitialParam1);
567   }
568   /////////////////////////////////////
569   Handle(Geom2d_Line) PC1 = new Geom2d_Line(gp_Pnt2d(anU, Vorigin), gp_Dir2d(0., 1.));
570   gp_Vec2d Offset(theUperiod, 0.);
571   if (Ydir > 0)
572     Offset *= -1;
573   Handle(Geom2d_Curve) PC2 = Handle(Geom2d_Curve)::DownCast(PC1->Copy());
574   PC2->Translate(Offset);
575   
576   BRep_Builder BB;
577   if (Ydir > 0)
578     BB.UpdateEdge(MissedSeam, PC1, PC2, theFrefFace, 0.);
579   else
580     BB.UpdateEdge(MissedSeam, PC2, PC1, theFrefFace, 0.);
581   BB.Continuity(MissedSeam, theFrefFace, theFrefFace, aContinuity);
582   if (Ydir < 0)
583     MissedSeam.Reverse();
584
585   BB.Add(theNewWire, MissedSeam);
586   //add newly created edge into VEmap
587   MissedSeam.Reverse();
588   TopExp::MapShapesAndUniqueAncestors(MissedSeam, TopAbs_VERTEX, TopAbs_EDGE, theVEmap);
589                   
590   BRepAdaptor_Curve2d BAcurve2d(theNextEdge, theFrefFace);
591   Standard_Real ParamOfNextPoint = (theNextEdge.Orientation() == TopAbs_FORWARD)?
592     BAcurve2d.LastParameter() : BAcurve2d.FirstParameter();
593   theNextPoint = BAcurve2d.Value(ParamOfNextPoint);
594 }
595
596 //=======================================================================
597 //function : TransformPCurves
598 //purpose  : auxilary
599 //=======================================================================
600 static void TransformPCurves(const TopoDS_Face& theRefFace,
601                              const TopoDS_Face& theFace,
602                              TopTools_MapOfShape& theMapEdgesWithTemporaryPCurves)
603 {
604   BRepAdaptor_Surface BAsurf(theFace, Standard_False);
605
606   Standard_Real Uperiod = 0., Vperiod = 0.;
607   Handle(Geom_Surface) RefSurf = BRep_Tool::Surface(theRefFace);
608   if (RefSurf->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)))
609     RefSurf = (Handle(Geom_RectangularTrimmedSurface)::DownCast(RefSurf))->BasisSurface();
610   if (RefSurf->IsUPeriodic())
611     Uperiod = RefSurf->UPeriod();
612   if (RefSurf->IsVPeriodic())
613     Vperiod = RefSurf->VPeriod();
614
615   GeomAdaptor_Surface RefGAsurf(RefSurf);
616
617   Standard_Real Ufirst = BAsurf.FirstUParameter(),
618     Ulast = BAsurf.LastUParameter(),
619     Vfirst = BAsurf.FirstVParameter(),
620     Vlast = BAsurf.LastVParameter();
621
622   Standard_Real u_mid = (Ufirst + Ulast)/2, v_mid = (Vfirst + Vlast)/2;
623   gp_Pnt MidPoint = BAsurf.Value(u_mid, v_mid);
624   Extrema_ExtPS ProjPS(MidPoint, RefGAsurf,
625                        Precision::PConfusion(), Precision::PConfusion());
626   Standard_Integer indmin = 1;
627   for (Standard_Integer iext = 2; iext <= ProjPS.NbExt(); iext++)
628     if (ProjPS.SquareDistance(iext) < ProjPS.SquareDistance(indmin))
629       indmin = iext;
630   
631   Standard_Real uu, vv;
632   ProjPS.Point(indmin).Parameter(uu,vv);
633   //Check uu and vv
634   if (Abs(u_mid + Uperiod - uu) <= Precision::PConfusion())
635     uu = u_mid;
636   if (Abs(u_mid - uu) <= Precision::PConfusion())
637     uu = u_mid;
638   if (Abs(v_mid + Vperiod - vv) <= Precision::PConfusion())
639     vv = v_mid;
640   if (Abs(v_mid - vv) <= Precision::PConfusion())
641     vv = v_mid;
642   gp_Vec2d Translation(uu - u_mid, vv - v_mid);
643
644   Standard_Boolean X_Reverse = Standard_False, Y_Reverse = Standard_False;
645   Standard_Real u_dx, v_dx, u_dy, v_dy;
646
647   Standard_Real Delta = (Precision::IsInfinite(Ufirst) || Precision::IsInfinite(Ulast))?
648     1. : (Ulast - Ufirst)/4;
649   Standard_Real Offset = (Uperiod == 0.)? Delta : Min(Uperiod/8, Delta);
650   Standard_Real u1 = u_mid + Offset, v1 = v_mid;
651   gp_Pnt DX = BAsurf.Value(u1, v1);
652   ProjPS.Perform(DX);
653   indmin = 1;
654   for (Standard_Integer iext = 2; iext <= ProjPS.NbExt(); iext++)
655     if (ProjPS.SquareDistance(iext) < ProjPS.SquareDistance(indmin))
656       indmin = iext;
657
658   ProjPS.Point(indmin).Parameter(u_dx, v_dx);
659   if (Uperiod != 0. &&
660       Abs(uu - u_dx) > Uperiod/2)
661   {
662     if (uu   < Uperiod/2 &&
663         u_dx > Uperiod/2)
664       X_Reverse = Standard_True;
665   }
666   else if (u_dx < uu)
667     X_Reverse = Standard_True;
668
669   Delta = (Precision::IsInfinite(Vfirst) || Precision::IsInfinite(Vlast))?
670     1. : (Vlast - Vfirst)/4;
671   Offset = (Vperiod == 0.)? Delta : Min(Vperiod/8, Delta);
672   Standard_Real u2 = u_mid, v2 = v_mid + Offset;
673   gp_Pnt DY = BAsurf.Value(u2, v2);
674   ProjPS.Perform(DY);
675   indmin = 1;
676   for (Standard_Integer iext = 2; iext <= ProjPS.NbExt(); iext++)
677     if (ProjPS.SquareDistance(iext) < ProjPS.SquareDistance(indmin))
678       indmin = iext;
679   
680   ProjPS.Point(indmin).Parameter(u_dy, v_dy);
681   if (Vperiod != 0. &&
682       Abs(vv - v_dy) > Vperiod/2)
683   {
684     if (vv   < Vperiod/2 &&
685         v_dy > Vperiod/2)
686       Y_Reverse = Standard_True;
687   }
688   else if (v_dy < vv)
689     Y_Reverse = Standard_True;
690     
691   gp_Trsf2d aTrsf;
692   if (X_Reverse && Y_Reverse)
693     aTrsf.SetMirror(gp::Origin2d());
694   else if (X_Reverse)
695     aTrsf.SetMirror(gp::OY2d());
696   else if (Y_Reverse)
697     aTrsf.SetMirror(gp::OX2d());
698
699   aTrsf.SetTranslationPart(Translation);
700
701   BRep_Builder BB;
702   TopExp_Explorer Explo(theFace, TopAbs_EDGE);
703   for (; Explo.More(); Explo.Next())
704   {
705     const TopoDS_Edge& anEdge = TopoDS::Edge(Explo.Current());
706     if (BRep_Tool::Degenerated(anEdge) &&
707         aTrsf.Form() != gp_Identity)
708       continue;
709     if (BRepTools::IsReallyClosed(anEdge, theFace))
710       continue;
711
712     Standard_Real fpar, lpar;
713     Handle(Geom2d_Curve) PCurveOnRef = BRep_Tool::CurveOnSurface(anEdge, theRefFace, fpar, lpar);
714     if (!PCurveOnRef.IsNull())
715       continue;
716
717     Handle(Geom2d_Curve) aPCurve = BRep_Tool::CurveOnSurface(anEdge, theFace, fpar, lpar);
718     Handle(Geom2d_Curve) aNewPCurve = Handle(Geom2d_Curve)::DownCast(aPCurve->Copy());
719     if (aTrsf.Form() != gp_Identity)
720       aNewPCurve->Transform(aTrsf);
721
722     Standard_Real tmp_first, tmp_last;
723     Handle(Geom2d_Curve) aPCurveOnRefFace = BRep_Tool::CurveOnSurface(anEdge, theRefFace,
724                                                                       tmp_first, tmp_last);
725     if (aPCurveOnRefFace.IsNull())
726       theMapEdgesWithTemporaryPCurves.Add(anEdge);
727     
728     BB.UpdateEdge(anEdge, aNewPCurve, theRefFace, 0.);
729     BB.Range(anEdge, fpar, lpar);
730   }
731 }
732                              
733 //=======================================================================
734 //function : AddPCurves
735 //purpose  : auxilary
736 //=======================================================================
737 static void AddPCurves(const TopTools_SequenceOfShape& theFaces,
738                        const TopoDS_Face&              theRefFace,
739                        TopTools_MapOfShape& theMapEdgesWithTemporaryPCurves)
740 {
741   BRepAdaptor_Surface RefBAsurf(theRefFace, Standard_False);
742
743   GeomAbs_SurfaceType aType = RefBAsurf.GetType();
744   if (aType == GeomAbs_Plane)
745     return;
746
747   for (Standard_Integer i = 1; i <= theFaces.Length(); i++)
748   {
749     const TopoDS_Face& aFace = TopoDS::Face(theFaces(i));
750     if (aFace.IsSame(theRefFace))
751       continue;
752
753     TransformPCurves(theRefFace, aFace, theMapEdgesWithTemporaryPCurves);
754   }
755 }
756
757 //=======================================================================
758 //function : AddOrdinaryEdges
759 //purpose  : auxilary
760 //=======================================================================
761 // adds edges from the shape to the sequence
762 // seams and equal edges are dropped
763 // Returns true if one of original edges dropped
764 static Standard_Boolean AddOrdinaryEdges(TopTools_SequenceOfShape& edges,
765                                          const TopoDS_Shape aShape,
766                                          Standard_Integer& anIndex)
767 {
768   //map of edges
769   TopTools_IndexedMapOfShape aNewEdges;
770   //add edges without seams
771   for(TopExp_Explorer exp(aShape,TopAbs_EDGE); exp.More(); exp.Next()) {
772     TopoDS_Shape edge = exp.Current();
773     if(aNewEdges.Contains(edge))
774       aNewEdges.RemoveKey(edge);
775     else
776       aNewEdges.Add(edge);
777   }
778
779   Standard_Boolean isDropped = Standard_False;
780   //merge edges and drop seams
781   Standard_Integer i;
782   for (i = 1; i <= edges.Length(); i++) {
783     TopoDS_Shape current = edges(i);
784     if(aNewEdges.Contains(current)) {
785
786       aNewEdges.RemoveKey(current);
787       edges.Remove(i);
788       i--;
789
790       if(!isDropped) {
791         isDropped = Standard_True;
792         anIndex = i;
793       }
794     }
795   }
796
797   //add edges to the sequence
798   for (i = 1; i <= aNewEdges.Extent(); i++)
799     edges.Append(aNewEdges(i));
800
801   return isDropped;
802 }
803
804 //=======================================================================
805 //function : getCylinder
806 //purpose  : auxilary
807 //=======================================================================
808 static Standard_Boolean getCylinder(Handle(Geom_Surface)& theInSurface,
809                                     gp_Cylinder& theOutCylinder)
810 {
811   Standard_Boolean isCylinder = Standard_False;
812
813   if (theInSurface->IsKind(STANDARD_TYPE(Geom_CylindricalSurface))) {
814     Handle(Geom_CylindricalSurface) aGC = Handle(Geom_CylindricalSurface)::DownCast(theInSurface);
815
816     theOutCylinder = aGC->Cylinder();
817     isCylinder = Standard_True;
818   }
819   else if (theInSurface->IsKind(STANDARD_TYPE(Geom_SurfaceOfRevolution))) {
820     Handle(Geom_SurfaceOfRevolution) aRS =
821       Handle(Geom_SurfaceOfRevolution)::DownCast(theInSurface);
822     Handle(Geom_Curve) aBasis = aRS->BasisCurve();
823     if (aBasis->IsKind(STANDARD_TYPE(Geom_Line))) {
824       Handle(Geom_Line) aBasisLine = Handle(Geom_Line)::DownCast(aBasis);
825       gp_Dir aDir = aRS->Direction();
826       gp_Dir aBasisDir = aBasisLine->Position().Direction();
827       if (aBasisDir.IsParallel(aDir, Precision::Angular())) {
828         // basis line is parallel to the revolution axis: it is a cylinder
829         gp_Pnt aLoc = aRS->Location();
830         Standard_Real aR = aBasisLine->Lin().Distance(aLoc);
831         gp_Ax3 aCylAx (aLoc, aDir);
832
833         theOutCylinder = gp_Cylinder(aCylAx, aR);
834         isCylinder = Standard_True;
835       }
836     }
837   }
838   else if (theInSurface->IsKind(STANDARD_TYPE(Geom_SurfaceOfLinearExtrusion))) {
839     Handle(Geom_SurfaceOfLinearExtrusion) aLES =
840       Handle(Geom_SurfaceOfLinearExtrusion)::DownCast(theInSurface);
841     Handle(Geom_Curve) aBasis = aLES->BasisCurve();
842     if (aBasis->IsKind(STANDARD_TYPE(Geom_Circle))) {
843       Handle(Geom_Circle) aBasisCircle = Handle(Geom_Circle)::DownCast(aBasis);
844       gp_Dir aDir = aLES->Direction();
845       gp_Dir aBasisDir = aBasisCircle->Position().Direction();
846       if (aBasisDir.IsParallel(aDir, Precision::Angular())) {
847         // basis circle is normal to the extrusion axis: it is a cylinder
848         gp_Pnt aLoc = aBasisCircle->Location();
849         Standard_Real aR = aBasisCircle->Radius();
850         gp_Ax3 aCylAx (aLoc, aDir);
851
852         theOutCylinder = gp_Cylinder(aCylAx, aR);
853         isCylinder = Standard_True;
854       }
855     }
856   }
857   else {
858   }
859
860   return isCylinder;
861 }
862
863 //=======================================================================
864 //function : ClearRts
865 //purpose  : auxilary
866 //=======================================================================
867 static Handle(Geom_Surface) ClearRts(const Handle(Geom_Surface)& aSurface)
868 {
869   if(aSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
870     Handle(Geom_RectangularTrimmedSurface) rts =
871       Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurface);
872     return rts->BasisSurface();
873   }
874   return aSurface;
875 }
876
877 //=======================================================================
878 //function : GetNormalToSurface
879 //purpose  : Gets the normal to surface by the given parameter on edge.
880 //           Returns True if normal was computed.
881 //=======================================================================
882 static Standard_Boolean GetNormalToSurface(const TopoDS_Face& theFace,
883                                            const TopoDS_Edge& theEdge,
884                                            const Standard_Real theP,
885                                            gp_Dir& theNormal)
886 {
887   Standard_Real f, l;
888   // get 2d curve to get point in 2d
889   const Handle(Geom2d_Curve)& aC2d = BRep_Tool::CurveOnSurface(theEdge, theFace, f, l);
890   if (aC2d.IsNull()) {
891     return Standard_False;
892   }
893   //
894   // 2d point
895   gp_Pnt2d aP2d;
896   aC2d->D0(theP, aP2d);
897   //
898   // get D1
899   gp_Vec aDU, aDV;
900   gp_Pnt aP3d;
901   TopLoc_Location aLoc;
902   const Handle(Geom_Surface)& aS = BRep_Tool::Surface(theFace, aLoc);
903   aS->D1(aP2d.X(), aP2d.Y(), aP3d, aDU, aDV);
904   //
905   // compute normal
906   gp_Vec aVNormal = aDU.Crossed(aDV);
907   if (aVNormal.Magnitude() < Precision::Confusion()) {
908     return Standard_False;
909   }
910   //
911   if (theFace.Orientation() == TopAbs_REVERSED) {
912     aVNormal.Reverse();
913   }
914   //
915   aVNormal.Transform(aLoc.Transformation());
916   theNormal = gp_Dir(aVNormal);
917   return Standard_True;
918 }
919
920 //=======================================================================
921 //function : IsSameDomain
922 //purpose  : 
923 //=======================================================================
924 static Standard_Boolean IsSameDomain(const TopoDS_Face& aFace,
925                                      const TopoDS_Face& aCheckedFace,
926                                      const Standard_Real theLinTol,
927                                      const Standard_Real theAngTol)
928 {
929   //checking the same handles
930   TopLoc_Location L1, L2;
931   Handle(Geom_Surface) S1, S2;
932
933   S1 = BRep_Tool::Surface(aFace,L1);
934   S2 = BRep_Tool::Surface(aCheckedFace,L2);
935
936   if (S1 == S2 && L1 == L2)
937     return Standard_True;
938
939   S1 = BRep_Tool::Surface(aFace);
940   S2 = BRep_Tool::Surface(aCheckedFace);
941
942   S1 = ClearRts(S1);
943   S2 = ClearRts(S2);
944
945   //Handle(Geom_OffsetSurface) aGOFS1, aGOFS2;
946   //aGOFS1 = Handle(Geom_OffsetSurface)::DownCast(S1);
947   //aGOFS2 = Handle(Geom_OffsetSurface)::DownCast(S2);
948   //if (!aGOFS1.IsNull()) S1 = aGOFS1->BasisSurface();
949   //if (!aGOFS2.IsNull()) S2 = aGOFS2->BasisSurface();
950
951   // case of two planar surfaces:
952   // all kinds of surfaces checked, including b-spline and bezier
953   GeomLib_IsPlanarSurface aPlanarityChecker1(S1, theLinTol);
954   if (aPlanarityChecker1.IsPlanar()) {
955     GeomLib_IsPlanarSurface aPlanarityChecker2(S2, theLinTol);
956     if (aPlanarityChecker2.IsPlanar()) {
957       gp_Pln aPln1 = aPlanarityChecker1.Plan();
958       gp_Pln aPln2 = aPlanarityChecker2.Plan();
959
960       if (aPln1.Position().Direction().IsParallel(aPln2.Position().Direction(), theAngTol) &&
961         aPln1.Distance(aPln2) < theLinTol) {
962         return Standard_True;
963       }
964     }
965   }
966
967   // case of two elementary surfaces: use OCCT tool
968   // elementary surfaces: ConicalSurface, CylindricalSurface,
969   //                      Plane, SphericalSurface and ToroidalSurface
970   if (S1->IsKind(STANDARD_TYPE(Geom_ElementarySurface)) &&
971       S2->IsKind(STANDARD_TYPE(Geom_ElementarySurface)))
972   {
973     Handle(GeomAdaptor_HSurface) aGA1 = new GeomAdaptor_HSurface(S1);
974     Handle(GeomAdaptor_HSurface) aGA2 = new GeomAdaptor_HSurface(S2);
975
976     Handle(BRepTopAdaptor_TopolTool) aTT1 = new BRepTopAdaptor_TopolTool();
977     Handle(BRepTopAdaptor_TopolTool) aTT2 = new BRepTopAdaptor_TopolTool();
978
979     try {
980       IntPatch_ImpImpIntersection anIIInt(aGA1, aTT1, aGA2, aTT2, theLinTol, theLinTol);
981       if (!anIIInt.IsDone() || anIIInt.IsEmpty())
982         return Standard_False;
983
984       return anIIInt.TangentFaces();
985     }
986     catch (Standard_Failure const&) {
987       return Standard_False;
988     }
989   }
990
991   // case of two cylindrical surfaces, at least one of which is a swept surface
992   // swept surfaces: SurfaceOfLinearExtrusion, SurfaceOfRevolution
993   if ((S1->IsKind(STANDARD_TYPE(Geom_CylindricalSurface)) ||
994        S1->IsKind(STANDARD_TYPE(Geom_SweptSurface))) &&
995       (S2->IsKind(STANDARD_TYPE(Geom_CylindricalSurface)) ||
996        S2->IsKind(STANDARD_TYPE(Geom_SweptSurface))))
997   {
998     gp_Cylinder aCyl1, aCyl2;
999     if (getCylinder(S1, aCyl1) && getCylinder(S2, aCyl2)) {
1000       if (fabs(aCyl1.Radius() - aCyl2.Radius()) < theLinTol) {
1001         gp_Dir aDir1 = aCyl1.Position().Direction();
1002         gp_Dir aDir2 = aCyl2.Position().Direction();
1003         if (aDir1.IsParallel(aDir2, Precision::Angular())) {
1004           gp_Pnt aLoc1 = aCyl1.Location();
1005           gp_Pnt aLoc2 = aCyl2.Location();
1006           gp_Vec aVec12 (aLoc1, aLoc2);
1007           if (aVec12.SquareMagnitude() < theLinTol*theLinTol ||
1008               aVec12.IsParallel(aDir1, Precision::Angular())) {
1009             return Standard_True;
1010           }
1011         }
1012       }
1013     }
1014   }
1015
1016   return Standard_False;
1017 }
1018
1019 //=======================================================================
1020 //function : UpdateMapOfShapes
1021 //purpose  :
1022 //=======================================================================
1023 static void UpdateMapOfShapes(TopTools_MapOfShape& theMapOfShapes,
1024                               Handle(ShapeBuild_ReShape)& theContext)
1025 {
1026   for (TopTools_MapIteratorOfMapOfShape it(theMapOfShapes); it.More(); it.Next()) {
1027     const TopoDS_Shape& aShape = it.Value();
1028     TopoDS_Shape aContextShape = theContext->Apply(aShape);
1029     if (!aContextShape.IsSame(aShape))
1030       theMapOfShapes.Add(aContextShape);
1031   }
1032 }
1033
1034 //=======================================================================
1035 //function : GlueEdgesWithPCurves
1036 //purpose  : Glues the pcurves of the sequence of edges
1037 //           and glues their 3d curves
1038 //=======================================================================
1039 static TopoDS_Edge GlueEdgesWithPCurves(const TopTools_SequenceOfShape& aChain,
1040                                         const TopoDS_Vertex& FirstVertex,
1041                                         const TopoDS_Vertex& LastVertex)
1042 {
1043   Standard_Integer i, j;
1044
1045   TopoDS_Edge FirstEdge = TopoDS::Edge(aChain(1));
1046   TColGeom_SequenceOfSurface SurfSeq;
1047   NCollection_Sequence<TopLoc_Location> LocSeq;
1048   
1049   for (int aCurveIndex = 0;; aCurveIndex++)
1050   {
1051     Handle(Geom2d_Curve) aCurve;
1052     Handle(Geom_Surface) aSurface;
1053     TopLoc_Location aLocation;
1054     Standard_Real aFirst, aLast;
1055     BRep_Tool::CurveOnSurface (FirstEdge, aCurve, aSurface, aLocation, aFirst, aLast, aCurveIndex);
1056     if (aCurve.IsNull())
1057       break;
1058
1059     SurfSeq.Append(aSurface);
1060     LocSeq.Append(aLocation);
1061   }
1062
1063   Standard_Real fpar, lpar;
1064   BRep_Tool::Range(FirstEdge, fpar, lpar);
1065   TopoDS_Edge PrevEdge = FirstEdge;
1066   TopoDS_Vertex CV;
1067   Standard_Real MaxTol = 0.;
1068   
1069   TopoDS_Edge ResEdge;
1070   BRep_Builder BB;
1071
1072   Standard_Integer nb_curve = aChain.Length();   //number of curves
1073   TColGeom_Array1OfBSplineCurve tab_c3d(0,nb_curve-1);                    //array of the curves
1074   TColStd_Array1OfReal tabtolvertex(0,nb_curve-1); //(0,nb_curve-2);  //array of the tolerances
1075     
1076   TopoDS_Vertex PrevVertex = FirstVertex;
1077   for (i = 1; i <= nb_curve; i++)
1078   {
1079     TopoDS_Edge anEdge = TopoDS::Edge(aChain(i));
1080     TopoDS_Vertex VF, VL;
1081     TopExp::Vertices(anEdge, VF, VL);
1082     Standard_Boolean ToReverse = (!VF.IsSame(PrevVertex));
1083     
1084     Standard_Real Tol1 = BRep_Tool::Tolerance(VF);
1085     Standard_Real Tol2 = BRep_Tool::Tolerance(VL);
1086     if (Tol1 > MaxTol)
1087       MaxTol = Tol1;
1088     if (Tol2 > MaxTol)
1089       MaxTol = Tol2;
1090     
1091     if (i > 1)
1092     {
1093       TopExp::CommonVertex(PrevEdge, anEdge, CV);
1094       Standard_Real Tol = BRep_Tool::Tolerance(CV);
1095       tabtolvertex(i-2) = Tol;
1096     }
1097     
1098     Handle(Geom_Curve) aCurve = BRep_Tool::Curve(anEdge, fpar, lpar);
1099     Handle(Geom_TrimmedCurve) aTrCurve = new Geom_TrimmedCurve(aCurve, fpar, lpar);
1100     tab_c3d(i-1) = GeomConvert::CurveToBSplineCurve(aTrCurve);
1101     GeomConvert::C0BSplineToC1BSplineCurve(tab_c3d(i-1), Precision::Confusion());
1102     if (ToReverse)
1103       tab_c3d(i-1)->Reverse();
1104     PrevVertex = (ToReverse)? VF : VL;
1105     PrevEdge = anEdge;
1106   }
1107   Handle(TColGeom_HArray1OfBSplineCurve)  concatcurve;     //array of the concatenated curves
1108   Handle(TColStd_HArray1OfInteger)        ArrayOfIndices;  //array of the remining Vertex
1109   Standard_Boolean closed_flag = Standard_False;
1110   GeomConvert::ConcatC1(tab_c3d,
1111                         tabtolvertex,
1112                         ArrayOfIndices,
1113                         concatcurve,
1114                         closed_flag,
1115                         Precision::Confusion());   //C1 concatenation
1116   
1117   if (concatcurve->Length() > 1)
1118   {
1119     GeomConvert_CompCurveToBSplineCurve Concat(concatcurve->Value(concatcurve->Lower()));
1120     
1121     for (i = concatcurve->Lower()+1; i <= concatcurve->Upper(); i++)
1122       Concat.Add( concatcurve->Value(i), MaxTol, Standard_True );
1123     
1124     concatcurve->SetValue(concatcurve->Lower(), Concat.BSplineCurve());
1125   }
1126   Handle(Geom_BSplineCurve) ResCurve = concatcurve->Value(concatcurve->Lower());
1127   
1128   TColGeom2d_SequenceOfBoundedCurve ResPCurves;
1129   for (j = 1; j <= SurfSeq.Length(); j++)
1130   {
1131     TColGeom2d_Array1OfBSplineCurve tab_c2d(0,nb_curve-1); //array of the pcurves
1132     
1133     PrevVertex = FirstVertex;
1134     PrevEdge = FirstEdge;
1135     for (i = 1; i <= nb_curve; i++)
1136     {
1137       TopoDS_Edge anEdge = TopoDS::Edge(aChain(i));
1138       TopoDS_Vertex VF, VL;
1139       TopExp::Vertices(anEdge, VF, VL);
1140       Standard_Boolean ToReverse = (!VF.IsSame(PrevVertex));
1141
1142       Handle(Geom2d_Curve) aPCurve =
1143         BRep_Tool::CurveOnSurface(anEdge, SurfSeq(j), LocSeq(j), fpar, lpar);
1144       if (aPCurve.IsNull())
1145         continue;
1146       Handle(Geom2d_TrimmedCurve) aTrPCurve = new Geom2d_TrimmedCurve(aPCurve, fpar, lpar);
1147       tab_c2d(i-1) = Geom2dConvert::CurveToBSplineCurve(aTrPCurve);
1148       Geom2dConvert::C0BSplineToC1BSplineCurve(tab_c2d(i-1), Precision::Confusion());
1149       if (ToReverse)
1150         tab_c2d(i-1)->Reverse();
1151       PrevVertex = (ToReverse)? VF : VL;
1152       PrevEdge = anEdge;
1153     }
1154     Handle(TColGeom2d_HArray1OfBSplineCurve)  concatc2d;     //array of the concatenated curves
1155     Handle(TColStd_HArray1OfInteger)        ArrayOfInd2d;  //array of the remining Vertex
1156     closed_flag = Standard_False;
1157     Geom2dConvert::ConcatC1(tab_c2d,
1158                             tabtolvertex,
1159                             ArrayOfInd2d,
1160                             concatc2d,
1161                             closed_flag,
1162                             Precision::Confusion());   //C1 concatenation
1163     
1164     if (concatc2d->Length() > 1)
1165     {
1166       Geom2dConvert_CompCurveToBSplineCurve Concat2d(concatc2d->Value(concatc2d->Lower()));
1167       
1168       for (i = concatc2d->Lower()+1; i <= concatc2d->Upper(); i++)
1169         Concat2d.Add( concatc2d->Value(i), MaxTol, Standard_True );
1170       
1171       concatc2d->SetValue(concatc2d->Lower(), Concat2d.BSplineCurve());
1172     }
1173     Handle(Geom2d_BSplineCurve) aResPCurve = concatc2d->Value(concatc2d->Lower());
1174     ResPCurves.Append(aResPCurve);
1175   }
1176   
1177   ResEdge = BRepLib_MakeEdge(ResCurve,
1178                              FirstVertex, LastVertex,
1179                              ResCurve->FirstParameter(), ResCurve->LastParameter());
1180   BB.SameRange(ResEdge, Standard_False);
1181   BB.SameParameter(ResEdge, Standard_False);
1182   for (j = 1; j <= ResPCurves.Length(); j++)
1183   {
1184     BB.UpdateEdge(ResEdge, ResPCurves(j), SurfSeq(j), LocSeq(j), MaxTol);
1185     BB.Range(ResEdge, SurfSeq(j), LocSeq(j), ResPCurves(j)->FirstParameter(), ResPCurves(j)->LastParameter());
1186   }
1187
1188   BRepLib::SameParameter(ResEdge, MaxTol, Standard_True);
1189   
1190   return ResEdge;
1191 }
1192
1193 //=======================================================================
1194 //function : MergeSubSeq
1195 //purpose  : Merges a sequence of edges into one edge if possible
1196 //=======================================================================
1197
1198 static Standard_Boolean MergeSubSeq(const TopTools_SequenceOfShape& theChain,
1199                                     const TopTools_IndexedDataMapOfShapeListOfShape& theVFmap,
1200                                     TopoDS_Edge& OutEdge, 
1201                                     double theAngTol, 
1202                                     Standard_Boolean ConcatBSplines,
1203                                     Standard_Boolean isSafeInputMode,
1204                                     Handle(ShapeBuild_ReShape)& theContext)
1205 {
1206   ShapeAnalysis_Edge sae;
1207   BRep_Builder B;
1208   // union edges in chain
1209   int j;
1210   Standard_Real fp1,lp1,fp2,lp2;
1211   Standard_Boolean IsUnionOfLinesPossible = Standard_True;
1212   Standard_Boolean IsUnionOfCirclesPossible = Standard_True;
1213   Handle(Geom_Curve) c3d1, c3d2;
1214   for(j = 1; j < theChain.Length(); j++) 
1215   {
1216     TopoDS_Edge edge1 = TopoDS::Edge(theChain.Value(j));
1217     TopoDS_Edge edge2 = TopoDS::Edge(theChain.Value(j+1));
1218
1219     if (BRep_Tool::Degenerated(edge1) &&
1220         BRep_Tool::Degenerated(edge2))
1221     {
1222       //Find the closest points in 2d
1223       TopoDS_Edge edgeFirst = TopoDS::Edge(theChain.First());
1224       TopoDS_Edge edgeLast  = TopoDS::Edge(theChain.Last());
1225       TopoDS_Face CommonFace;
1226       Standard_Real MinSqDist;
1227       TopAbs_Orientation OrOfE1OnFace, OrOfE2OnFace;
1228       Standard_Integer IndOnE1, IndOnE2;
1229       gp_Pnt2d PointsOnEdge1 [2], PointsOnEdge2 [2];
1230       if (!FindClosestPoints(edgeFirst, edgeLast, theVFmap, CommonFace,
1231                              MinSqDist, OrOfE1OnFace, OrOfE2OnFace,
1232                              IndOnE1, IndOnE2, PointsOnEdge1, PointsOnEdge2))
1233         return Standard_False;
1234       
1235       //Define indices corresponding to extremities of future edge
1236       IndOnE1 = 1 - IndOnE1;
1237       IndOnE2 = 1 - IndOnE2;
1238
1239       //Construct new degenerated edge
1240       gp_Pnt2d StartPoint = PointsOnEdge1[IndOnE1];
1241       gp_Pnt2d EndPoint   = PointsOnEdge2[IndOnE2];
1242       if ((OrOfE1OnFace == TopAbs_FORWARD  && IndOnE1 == 1) ||
1243           (OrOfE1OnFace == TopAbs_REVERSED && IndOnE1 == 0))
1244       { gp_Pnt2d Tmp = StartPoint; StartPoint = EndPoint; EndPoint = Tmp; }
1245       
1246       Handle(Geom2d_Line) aLine = GCE2d_MakeLine(StartPoint, EndPoint);
1247
1248       TopoDS_Vertex aVertex = TopExp::FirstVertex(edgeFirst);
1249       TopoDS_Vertex StartVertex = aVertex, EndVertex = aVertex;
1250       StartVertex.Orientation(TopAbs_FORWARD);
1251       EndVertex.Orientation(TopAbs_REVERSED);
1252       
1253       TopoDS_Edge NewEdge;
1254       B.MakeEdge(NewEdge);
1255       B.UpdateEdge(NewEdge, aLine, CommonFace, Precision::Confusion());
1256       B.Range(NewEdge, 0., StartPoint.Distance(EndPoint));
1257       B.Add (NewEdge, StartVertex);
1258       B.Add (NewEdge, EndVertex);
1259       B.Degenerated(NewEdge, Standard_True);
1260       OutEdge = NewEdge;
1261       return Standard_True;
1262     }
1263     
1264     c3d1 = BRep_Tool::Curve(edge1,fp1,lp1);
1265     c3d2 = BRep_Tool::Curve(edge2,fp2,lp2);
1266
1267     if(c3d1.IsNull() || c3d2.IsNull()) 
1268       return Standard_False;
1269
1270     while(c3d1->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
1271       Handle(Geom_TrimmedCurve) tc =
1272         Handle(Geom_TrimmedCurve)::DownCast(c3d1);
1273       c3d1 = tc->BasisCurve();
1274     }
1275     while(c3d2->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
1276       Handle(Geom_TrimmedCurve) tc =
1277         Handle(Geom_TrimmedCurve)::DownCast(c3d2);
1278       c3d2 = tc->BasisCurve();
1279     }
1280     if( c3d1->IsKind(STANDARD_TYPE(Geom_Line)) && c3d2->IsKind(STANDARD_TYPE(Geom_Line)) ) {
1281       Handle(Geom_Line) L1 = Handle(Geom_Line)::DownCast(c3d1);
1282       Handle(Geom_Line) L2 = Handle(Geom_Line)::DownCast(c3d2);
1283       gp_Dir Dir1 = L1->Position().Direction();
1284       gp_Dir Dir2 = L2->Position().Direction();
1285       if(!Dir1.IsParallel(Dir2,theAngTol))  
1286         IsUnionOfLinesPossible = Standard_False;
1287     }
1288     else
1289       IsUnionOfLinesPossible = Standard_False;
1290     if( c3d1->IsKind(STANDARD_TYPE(Geom_Circle)) && c3d2->IsKind(STANDARD_TYPE(Geom_Circle)) ) {
1291       Handle(Geom_Circle) C1 = Handle(Geom_Circle)::DownCast(c3d1);
1292       Handle(Geom_Circle) C2 = Handle(Geom_Circle)::DownCast(c3d2);
1293       gp_Pnt P01 = C1->Location();
1294       gp_Pnt P02 = C2->Location();
1295       if (P01.Distance(P02) > Precision::Confusion())
1296         IsUnionOfCirclesPossible = Standard_False;
1297     }
1298     else
1299       IsUnionOfCirclesPossible = Standard_False;
1300   }
1301   if (IsUnionOfLinesPossible && IsUnionOfCirclesPossible)
1302     return Standard_False;
1303
1304   //union of lines is possible
1305   if (IsUnionOfLinesPossible)
1306   {
1307     TopoDS_Vertex V[2];
1308     V[0] = sae.FirstVertex(TopoDS::Edge(theChain.First()));
1309     gp_Pnt PV1 = BRep_Tool::Pnt(V[0]);
1310     V[1] = sae.LastVertex(TopoDS::Edge(theChain.Last()));
1311     gp_Pnt PV2 = BRep_Tool::Pnt(V[1]);
1312     gp_Vec Vec(PV1, PV2);
1313     if (isSafeInputMode) {
1314       for (int k = 0; k < 2; k++) {
1315         if (!theContext->IsRecorded(V[k])) {
1316           TopoDS_Vertex Vcopy = TopoDS::Vertex(V[k].EmptyCopied());
1317           theContext->Replace(V[k], Vcopy);
1318           V[k] = Vcopy;
1319         }
1320         else
1321           V[k] = TopoDS::Vertex(theContext->Apply(V[k]));
1322       }
1323     }
1324     Handle(Geom_Line) L = new Geom_Line(gp_Ax1(PV1,Vec));
1325     Standard_Real dist = PV1.Distance(PV2);
1326     Handle(Geom_TrimmedCurve) tc = new Geom_TrimmedCurve(L,0.0,dist);
1327     TopoDS_Edge E;
1328     B.MakeEdge (E, tc ,Precision::Confusion());
1329     B.Add (E,V[0]);  B.Add (E,V[1]);
1330     B.UpdateVertex(V[0], 0., E, 0.);
1331     B.UpdateVertex(V[1], dist, E, 0.);
1332     OutEdge = E;
1333     return Standard_True;
1334   }
1335
1336   if (IsUnionOfCirclesPossible)
1337   {
1338     double f,l;
1339     TopoDS_Edge FE = TopoDS::Edge(theChain.First());
1340     Handle(Geom_Curve) c3d = BRep_Tool::Curve(FE,f,l);
1341
1342     while(c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
1343       Handle(Geom_TrimmedCurve) tc =
1344         Handle(Geom_TrimmedCurve)::DownCast(c3d);
1345       c3d = tc->BasisCurve();
1346     }
1347     Handle(Geom_Circle) Cir = Handle(Geom_Circle)::DownCast(c3d);
1348
1349     TopoDS_Vertex V[2];
1350     V[0] = sae.FirstVertex(FE);
1351     V[1] = sae.LastVertex(TopoDS::Edge(theChain.Last()));
1352     TopoDS_Edge E;
1353     if (V[0].IsSame(V[1])) {
1354       // closed chain
1355       BRepAdaptor_Curve adef(FE);
1356       Handle(Geom_Circle) Cir1;
1357       double FP, LP;
1358       if ( FE.Orientation() == TopAbs_FORWARD)
1359       {
1360         FP = adef.FirstParameter();
1361         LP = adef.LastParameter();
1362       }
1363       else
1364       {
1365         FP = adef.LastParameter();
1366         LP = adef.FirstParameter();
1367       }
1368       if (Abs(FP) < Precision::PConfusion())
1369       {
1370         B.MakeEdge (E,Cir, Precision::Confusion());
1371         B.Add(E,V[0]);
1372         B.Add(E,V[1]);
1373         E.Orientation(FE.Orientation());
1374       }
1375       else
1376       {
1377         GC_MakeCircle MC1 (adef.Value(FP), adef.Value((FP + LP) * 0.5), adef.Value(LP));
1378         if (MC1.IsDone())
1379           Cir1 = MC1.Value();
1380         else
1381           return Standard_False;
1382         B.MakeEdge (E, Cir1, Precision::Confusion());
1383         B.Add(E,V[0]);
1384         B.Add(E,V[1]);
1385       }
1386     }
1387     else {
1388       if (isSafeInputMode) {
1389         for (int k = 0; k < 2; k++) {
1390           if (!theContext->IsRecorded(V[k])) {
1391             TopoDS_Vertex Vcopy = TopoDS::Vertex(V[k].EmptyCopied());
1392             theContext->Replace(V[k], Vcopy);
1393             V[k] = Vcopy;
1394           }
1395           else
1396             V[k] = TopoDS::Vertex(theContext->Apply(V[k]));
1397         }
1398       }
1399       gp_Pnt PV1 = BRep_Tool::Pnt(V[0]);
1400       gp_Pnt PV2 = BRep_Tool::Pnt(V[1]);
1401       TopoDS_Vertex VM = sae.LastVertex(FE);
1402       gp_Pnt PVM = BRep_Tool::Pnt(VM);
1403       GC_MakeCircle MC (PV1,PVM,PV2);
1404       Handle(Geom_Circle) C = MC.Value();
1405       gp_Pnt P0 = C->Location();
1406       gp_Dir D1(gp_Vec(P0,PV1));
1407       gp_Dir D2(gp_Vec(P0,PV2));
1408       Standard_Real fpar = C->XAxis().Direction().Angle(D1);
1409       if(fabs(fpar)>Precision::Confusion()) {
1410         // check orientation
1411         gp_Dir ND =  C->XAxis().Direction().Crossed(D1);
1412         if(ND.IsOpposite(C->Axis().Direction(),Precision::Confusion())) {
1413           fpar = -fpar;
1414         }
1415       }
1416       Standard_Real lpar = C->XAxis().Direction().Angle(D2);
1417       if(fabs(lpar)>Precision::Confusion()) {
1418         // check orientation
1419         gp_Dir ND =  C->XAxis().Direction().Crossed(D2);
1420         if(ND.IsOpposite(C->Axis().Direction(),Precision::Confusion())) {
1421           lpar = -lpar;
1422         }
1423       }
1424       if (lpar < fpar) lpar += 2*M_PI;
1425       Handle(Geom_TrimmedCurve) tc = new Geom_TrimmedCurve(C,fpar,lpar);
1426       B.MakeEdge (E,tc,Precision::Confusion());
1427       B.Add(E,V[0]);
1428       B.Add(E,V[1]);
1429       B.UpdateVertex(V[0], fpar, E, 0.);
1430       B.UpdateVertex(V[1], lpar, E, 0.);
1431     }
1432     OutEdge = E;
1433     return Standard_True;
1434   }
1435   if (theChain.Length() > 1 && ConcatBSplines) {
1436     // second step: union edges with various curves
1437     // skl for bug 0020052 from Mantis: perform such unions
1438     // only if curves are bspline or bezier
1439
1440     TopoDS_Vertex VF = sae.FirstVertex(TopoDS::Edge(theChain.First()));
1441     TopoDS_Vertex VL = sae.LastVertex(TopoDS::Edge(theChain.Last()));
1442     Standard_Boolean NeedUnion = Standard_True;
1443     for(j = 1; j <= theChain.Length(); j++) {
1444       TopoDS_Edge edge = TopoDS::Edge(theChain.Value(j));
1445       TopLoc_Location Loc;
1446       Handle(Geom_Curve) c3d = BRep_Tool::Curve(edge,Loc,fp1,lp1);
1447       if(c3d.IsNull()) continue;
1448       while(c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
1449         Handle(Geom_TrimmedCurve) tc =
1450           Handle(Geom_TrimmedCurve)::DownCast(c3d);
1451         c3d = tc->BasisCurve();
1452       }
1453       if( ( c3d->IsKind(STANDARD_TYPE(Geom_BSplineCurve)) ||
1454             c3d->IsKind(STANDARD_TYPE(Geom_BezierCurve)) ) ) continue;
1455       NeedUnion = Standard_False;
1456       break;
1457     }
1458     if(NeedUnion) {
1459 #ifdef OCCT_DEBUG
1460       std::cout<<"can not make analitical union => make approximation"<<std::endl;
1461 #endif
1462       TopoDS_Edge E = GlueEdgesWithPCurves(theChain, VF, VL);
1463       OutEdge = E;
1464       return Standard_True;
1465     }
1466     else {
1467 #ifdef OCCT_DEBUG
1468       std::cout<<"can not make approximation for such types of curves"<<std::endl;
1469 #endif
1470       return Standard_False;
1471     }
1472   }
1473   return Standard_False;
1474 }
1475
1476 //=======================================================================
1477 //function : IsMergingPossible
1478 //purpose  : Checks if merging of two edges is possible
1479 //=======================================================================
1480
1481 static Standard_Boolean IsMergingPossible(const TopoDS_Edge& edge1, const TopoDS_Edge& edge2, 
1482                                           double theAngTol, double theLinTol, 
1483                                           const TopTools_MapOfShape& AvoidEdgeVrt, const bool theLineDirectionOk,
1484                                           const gp_Pnt& theFirstPoint, const gp_Vec& theDirectionVec,
1485                                           const TopTools_IndexedDataMapOfShapeListOfShape& theVFmap)
1486 {
1487   Standard_Boolean IsDegE1 = BRep_Tool::Degenerated(edge1);
1488   Standard_Boolean IsDegE2 = BRep_Tool::Degenerated(edge2);
1489   
1490   if (IsDegE1 && IsDegE2)
1491   {
1492     //Find connstion point in 2d
1493     TopoDS_Face CommonFace;
1494     Standard_Real MinSqDist;
1495     TopAbs_Orientation OrOfE1OnFace, OrOfE2OnFace;
1496     Standard_Integer IndOnE1, IndOnE2;
1497     gp_Pnt2d PointsOnEdge1 [2], PointsOnEdge2 [2];
1498     if (!FindClosestPoints(edge1, edge2, theVFmap, CommonFace,
1499                            MinSqDist, OrOfE1OnFace, OrOfE2OnFace,
1500                            IndOnE1, IndOnE2, PointsOnEdge1, PointsOnEdge2))
1501       return Standard_False;
1502     
1503     if (MinSqDist <= Precision::SquareConfusion())
1504       return Standard_True;
1505     
1506     return Standard_False;
1507   }
1508   else if (IsDegE1 || IsDegE2)
1509     return Standard_False;
1510   
1511   TopoDS_Vertex CV = TopExp::LastVertex(edge1, Standard_True);
1512   if (CV.IsNull() || AvoidEdgeVrt.Contains(CV))
1513     return Standard_False;
1514
1515   BRepAdaptor_Curve ade1(edge1);
1516   BRepAdaptor_Curve ade2(edge2);
1517
1518   GeomAbs_CurveType t1 = ade1.GetType();
1519   GeomAbs_CurveType t2 = ade2.GetType();
1520
1521   if( t1 == GeomAbs_Circle && t2 == GeomAbs_Circle)
1522   {
1523     if (ade1.Circle().Location().Distance(ade2.Circle().Location()) > Precision::Confusion())
1524       return Standard_False;
1525   }
1526
1527   if( ( (t1 != GeomAbs_BezierCurve && t1 != GeomAbs_BSplineCurve) ||
1528       (t2 != GeomAbs_BezierCurve && t2 != GeomAbs_BSplineCurve)) && t1 != t2)
1529     return Standard_False;
1530
1531   gp_Vec Diff1, Diff2;
1532   gp_Pnt P1, P2;
1533   if (edge1.Orientation() == TopAbs_FORWARD)
1534     ade1.D1(ade1.LastParameter(), P1, Diff1);
1535   else
1536   {
1537     ade1.D1(ade1.FirstParameter(), P1, Diff1);
1538     Diff1 = -Diff1;
1539   }
1540
1541   if (edge2.Orientation() == TopAbs_FORWARD)
1542     ade2.D1(ade2.FirstParameter(), P2, Diff2);
1543   else
1544   {
1545     ade2.D1(ade2.LastParameter(), P2, Diff2);
1546     Diff2 = -Diff2;
1547   }
1548
1549   if (Diff1.Angle(Diff2) > theAngTol)
1550     return Standard_False;
1551
1552   if (theLineDirectionOk && t2 == GeomAbs_Line)
1553   {
1554     // Check that the accumulated deflection does not exceed the linear tolerance
1555     Standard_Real aLast = (edge2.Orientation() == TopAbs_FORWARD) ?
1556       ade2.LastParameter() : ade2.FirstParameter();
1557     gp_Vec aCurV(theFirstPoint, ade2.Value(aLast));
1558     Standard_Real aDD = theDirectionVec.CrossSquareMagnitude(aCurV);
1559     if (aDD > theLinTol*theLinTol)
1560       return Standard_False;
1561
1562     // Check that the accumulated angle does not exceed the angular tolerance.
1563     // For symmetry, check the angle between vectors of:
1564     // - first edge and resulting curve, and
1565     // - the last edge and resulting curve.
1566     if (theDirectionVec.Angle(aCurV) > theAngTol || Diff2.Angle(aCurV) > theAngTol)
1567       return Standard_False;
1568   }
1569
1570   return Standard_True;
1571 }
1572
1573 //=======================================================================
1574 //function : GetLineEdgePoints
1575 //purpose  : 
1576 //=======================================================================
1577 static Standard_Boolean GetLineEdgePoints(const TopoDS_Edge& theInpEdge, gp_Pnt& theFirstPoint, gp_Vec& theDirectionVec)
1578 {
1579   double f, l;
1580   Handle(Geom_Curve) aCur = BRep_Tool::Curve(theInpEdge, f, l);
1581   if(aCur.IsNull()) 
1582     return Standard_False;
1583
1584   Handle(Geom_TrimmedCurve) aTC = Handle(Geom_TrimmedCurve)::DownCast(aCur);
1585   if (!aTC.IsNull())
1586     aCur = aTC->BasisCurve();
1587
1588   if (aCur->DynamicType() != STANDARD_TYPE(Geom_Line))
1589     return Standard_False;
1590
1591   if (theInpEdge.Orientation() == TopAbs_REVERSED) {
1592     Standard_Real tmp = f;
1593     f = l;
1594     l = tmp;
1595   }
1596   theFirstPoint = aCur->Value(f);
1597   gp_Pnt aLP = aCur->Value(l);
1598   theDirectionVec = aLP.XYZ().Subtracted(theFirstPoint.XYZ());
1599   theDirectionVec.Normalize();
1600   return Standard_True;
1601 }
1602
1603 //=======================================================================
1604 //function : GenerateSubSeq
1605 //purpose  : Generates sub-sequences of edges from sequence of edges
1606 //Edges from each subsequences can be merged into the one edge  
1607 //=======================================================================
1608
1609 static void GenerateSubSeq (const TopTools_SequenceOfShape& anInpEdgeSeq,
1610                             NCollection_Sequence<SubSequenceOfEdges>& SeqOfSubSeqOfEdges,
1611                             Standard_Boolean IsClosed, double theAngTol, double theLinTol, 
1612                             const TopTools_MapOfShape& AvoidEdgeVrt,
1613                             const TopTools_IndexedDataMapOfShapeListOfShape& theVFmap)
1614 {
1615   Standard_Boolean isOk = Standard_False;
1616   TopoDS_Edge edge1, edge2;
1617
1618   SubSequenceOfEdges SubSeq;
1619   TopoDS_Edge RefEdge = TopoDS::Edge(anInpEdgeSeq(1));
1620   SubSeq.SeqsEdges.Append(RefEdge);
1621   SeqOfSubSeqOfEdges.Append(SubSeq);
1622
1623   gp_Pnt aFirstPoint;
1624   gp_Vec aDirectionVec;
1625   Standard_Boolean isLineDirectionOk = GetLineEdgePoints(RefEdge, aFirstPoint, aDirectionVec);  
1626   
1627   for (int i = 1; i < anInpEdgeSeq.Length(); i++)
1628   {
1629     edge1 = TopoDS::Edge(anInpEdgeSeq(i));
1630     edge2 = TopoDS::Edge(anInpEdgeSeq(i+1));
1631     isOk = IsMergingPossible(edge1, edge2, theAngTol, theLinTol,
1632                              AvoidEdgeVrt, isLineDirectionOk, aFirstPoint, aDirectionVec, theVFmap);
1633     if (!isOk)
1634     {
1635       SubSequenceOfEdges aSubSeq;
1636       aSubSeq.SeqsEdges.Append(edge2);
1637       SeqOfSubSeqOfEdges.Append(aSubSeq);
1638       isLineDirectionOk = GetLineEdgePoints(edge2, aFirstPoint, aDirectionVec);
1639     }
1640     else
1641       SeqOfSubSeqOfEdges.ChangeLast().SeqsEdges.Append(edge2);
1642   }
1643   /// check first and last chain segments
1644   if (IsClosed && SeqOfSubSeqOfEdges.Length() > 1)
1645   {
1646     edge1 = TopoDS::Edge(anInpEdgeSeq.Last());
1647     edge2 = TopoDS::Edge(anInpEdgeSeq.First());
1648     if (IsMergingPossible(edge1, edge2, theAngTol, theLinTol,
1649                           AvoidEdgeVrt, Standard_False, aFirstPoint, aDirectionVec, theVFmap))
1650     {
1651       SeqOfSubSeqOfEdges.ChangeLast().SeqsEdges.Append(SeqOfSubSeqOfEdges.ChangeFirst().SeqsEdges);
1652       SeqOfSubSeqOfEdges.Remove(1);
1653     }
1654   }
1655 }
1656
1657 //=======================================================================
1658 //function : MergeEdges
1659 //purpose  : auxilary
1660 //=======================================================================
1661 static Standard_Boolean MergeEdges(TopTools_SequenceOfShape& SeqEdges,
1662                                    const TopTools_IndexedDataMapOfShapeListOfShape& theVFmap,
1663                                    const Standard_Real theAngTol,
1664                                    const Standard_Real theLinTol,
1665                                    const Standard_Boolean ConcatBSplines,
1666                                    const Standard_Boolean isSafeInputMode,
1667                                    Handle(ShapeBuild_ReShape)& theContext,
1668                                    NCollection_Sequence<SubSequenceOfEdges>& SeqOfSubSeqOfEdges,
1669                                    const TopTools_MapOfShape& NonMergVrt)
1670 {
1671   TopTools_IndexedDataMapOfShapeListOfShape aMapVE;
1672   Standard_Integer j;
1673   TopTools_MapOfShape VerticesToAvoid;
1674   for (j = 1; j <= SeqEdges.Length(); j++)
1675   {
1676     TopoDS_Edge anEdge = TopoDS::Edge(SeqEdges(j));
1677     // fill in the map V-E
1678     for (TopoDS_Iterator it(anEdge.Oriented(TopAbs_FORWARD)); it.More(); it.Next())
1679     {
1680       TopoDS_Shape aV = it.Value();
1681       if (aV.Orientation() == TopAbs_FORWARD || aV.Orientation() == TopAbs_REVERSED)
1682       {
1683         if (!aMapVE.Contains(aV))
1684           aMapVE.Add(aV, TopTools_ListOfShape());
1685         aMapVE.ChangeFromKey(aV).Append(anEdge);
1686       }
1687     }
1688   }
1689   VerticesToAvoid.Unite(NonMergVrt);
1690
1691   // do loop while there are unused edges
1692   TopTools_MapOfShape aUsedEdges;
1693   for (;;)
1694   {
1695     TopoDS_Edge edge;
1696     for(j=1; j <= SeqEdges.Length(); j++)
1697     {
1698       edge = TopoDS::Edge(SeqEdges.Value(j));
1699       if (!aUsedEdges.Contains(edge))
1700         break;
1701     }
1702     if (j > SeqEdges.Length())
1703       break; // all edges have been used
1704
1705     // make chain for unite
1706     TopTools_SequenceOfShape aChain;
1707     aChain.Append(edge);
1708     aUsedEdges.Add(edge);
1709     TopoDS_Vertex V[2];
1710     TopExp::Vertices(edge, V[0], V[1], Standard_True);
1711
1712     // connect more edges to the chain in both directions
1713     for (j = 0; j < 2; j++)
1714     {
1715       Standard_Boolean isAdded = Standard_True;
1716       while (isAdded)
1717       {
1718         isAdded = Standard_False;
1719         if (V[j].IsNull())
1720           break;
1721         const TopTools_ListOfShape& aLE = aMapVE.FindFromKey(V[j]);
1722         for (TopTools_ListIteratorOfListOfShape itL(aLE); itL.More(); itL.Next())
1723         {
1724           edge = TopoDS::Edge(itL.Value());
1725           if (!aUsedEdges.Contains(edge))
1726           {
1727             TopoDS_Vertex V2[2];
1728             TopExp::Vertices(edge, V2[0], V2[1], Standard_True);
1729             // the neighboring edge must have V[j] reversed and located on the opposite end
1730             if (V2[1 - j].IsEqual(V[j].Reversed()))
1731             {
1732               if (j == 0)
1733                 aChain.Prepend(edge);
1734               else
1735                 aChain.Append(edge);
1736               aUsedEdges.Add(edge);
1737               V[j] = V2[j];
1738               isAdded = Standard_True;
1739               break;
1740             }
1741           }
1742         }
1743       }
1744     }
1745
1746     if (aChain.Length() < 2)
1747       continue;
1748
1749     Standard_Boolean IsClosed = Standard_False;
1750     if (V[0].IsSame ( V[1] ))
1751       IsClosed = Standard_True;
1752
1753     // split chain by vertices at which merging is not possible
1754     NCollection_Sequence<SubSequenceOfEdges> aOneSeq;
1755     GenerateSubSeq(aChain, aOneSeq, IsClosed, theAngTol, theLinTol, VerticesToAvoid, theVFmap);
1756
1757     // put sub-chains in the result
1758     SeqOfSubSeqOfEdges.Append(aOneSeq);
1759   }
1760
1761   for (int i = 1; i <= SeqOfSubSeqOfEdges.Length(); i++)
1762   {
1763     TopoDS_Edge UE;
1764     if (SeqOfSubSeqOfEdges(i).SeqsEdges.Length() < 2)
1765       continue;
1766     if (MergeSubSeq(SeqOfSubSeqOfEdges(i).SeqsEdges, theVFmap,
1767                     UE, theAngTol, 
1768                     ConcatBSplines, isSafeInputMode, theContext))
1769       SeqOfSubSeqOfEdges(i).UnionEdges = UE;
1770   }
1771   return Standard_True;
1772 }
1773
1774 //=======================================================================
1775 //function : MergeSeq
1776 //purpose  : Tries to unify the sequence of edges with the set of
1777 //           another edges which lies on the same geometry
1778 //=======================================================================
1779 static Standard_Boolean MergeSeq (TopTools_SequenceOfShape& SeqEdges,
1780                                   const TopTools_IndexedDataMapOfShapeListOfShape& theVFmap,
1781                                   const Standard_Real theAngTol,
1782                                   const Standard_Real theLinTol,
1783                                   const Standard_Boolean ConcatBSplines,
1784                                   const Standard_Boolean isSafeInputMode,
1785                                   Handle(ShapeBuild_ReShape)& theContext,
1786                                   const TopTools_MapOfShape& nonMergVert)
1787 {
1788   NCollection_Sequence<SubSequenceOfEdges> SeqOfSubsSeqOfEdges;
1789   if (MergeEdges(SeqEdges, theVFmap, theAngTol, theLinTol, ConcatBSplines, isSafeInputMode,
1790                  theContext, SeqOfSubsSeqOfEdges, nonMergVert))
1791   {
1792     for (Standard_Integer i = 1; i <= SeqOfSubsSeqOfEdges.Length(); i++ )
1793     {
1794       if (SeqOfSubsSeqOfEdges(i).UnionEdges.IsNull())
1795         continue;
1796
1797       theContext->Merge(SeqOfSubsSeqOfEdges(i).SeqsEdges,
1798         SeqOfSubsSeqOfEdges(i).UnionEdges);
1799     }
1800     return Standard_True;
1801   }
1802   return Standard_False;
1803 }
1804
1805 //=======================================================================
1806 //function : CheckSharedVertices
1807 //purpose  : Checks the sequence of edges on the presence of shared vertex 
1808 //=======================================================================
1809
1810 static void CheckSharedVertices(const TopTools_SequenceOfShape& theSeqEdges, 
1811                                 const TopTools_IndexedDataMapOfShapeListOfShape& theMapEdgesVertex,
1812                                 const TopTools_MapOfShape& theMapKeepShape,
1813                                 TopTools_MapOfShape& theShareVertMap)
1814 {
1815   ShapeAnalysis_Edge sae;
1816   TopTools_SequenceOfShape SeqVertexes;
1817   TopTools_MapOfShape MapVertexes;
1818   for (Standard_Integer k = 1; k <= theSeqEdges.Length(); k++ )
1819   {
1820     TopoDS_Vertex aV1 = sae.FirstVertex(TopoDS::Edge(theSeqEdges(k)));
1821     TopoDS_Vertex aV2 = sae.LastVertex(TopoDS::Edge(theSeqEdges(k)));
1822     if (!MapVertexes.Add(aV1))
1823       SeqVertexes.Append(aV1);
1824     if (!MapVertexes.Add(aV2))
1825       SeqVertexes.Append(aV2);
1826   }
1827
1828   for (Standard_Integer k = 1; k <= SeqVertexes.Length()/* && !IsSharedVertexPresent*/; k++ )
1829   {
1830     const TopTools_ListOfShape& ListEdgesV1 = theMapEdgesVertex.FindFromKey(SeqVertexes(k));
1831     if (ListEdgesV1.Extent() > 2 || theMapKeepShape.Contains(SeqVertexes(k)))
1832       theShareVertMap.Add(SeqVertexes(k));
1833   }
1834   //return theShareVertMap.IsEmpty() ? false : true;
1835 }
1836
1837 //=======================================================================
1838 //function : ShapeUpgrade_UnifySameDomain
1839 //purpose  : Constructor
1840 //=======================================================================
1841
1842 ShapeUpgrade_UnifySameDomain::ShapeUpgrade_UnifySameDomain()
1843   : myLinTol(Precision::Confusion()),
1844     myAngTol(Precision::Angular()),
1845     myUnifyFaces(Standard_True),
1846     myUnifyEdges (Standard_True),
1847     myConcatBSplines (Standard_False),
1848     myAllowInternal (Standard_False),
1849     mySafeInputMode(Standard_True),
1850     myHistory(new BRepTools_History)
1851 {
1852   myContext = new ShapeBuild_ReShape;
1853 }
1854
1855 //=======================================================================
1856 //function : ShapeUpgrade_UnifySameDomain
1857 //purpose  : Constructor
1858 //=======================================================================
1859
1860 ShapeUpgrade_UnifySameDomain::ShapeUpgrade_UnifySameDomain(const TopoDS_Shape& aShape,
1861                                                            const Standard_Boolean UnifyEdges,
1862                                                            const Standard_Boolean UnifyFaces,
1863                                                            const Standard_Boolean ConcatBSplines)
1864   : myInitShape (aShape),
1865     myLinTol(Precision::Confusion()),
1866     myAngTol(Precision::Angular()),
1867     myUnifyFaces(UnifyFaces),
1868     myUnifyEdges (UnifyEdges),
1869     myConcatBSplines (ConcatBSplines),
1870     myAllowInternal (Standard_False),
1871     mySafeInputMode (Standard_True),
1872     myShape (aShape),
1873     myHistory(new BRepTools_History)
1874 {
1875   myContext = new ShapeBuild_ReShape;
1876 }
1877
1878 //=======================================================================
1879 //function : Initialize
1880 //purpose  : 
1881 //=======================================================================
1882
1883 void ShapeUpgrade_UnifySameDomain::Initialize(const TopoDS_Shape& aShape,
1884                                               const Standard_Boolean UnifyEdges,
1885                                               const Standard_Boolean UnifyFaces,
1886                                               const Standard_Boolean ConcatBSplines)
1887 {
1888   myInitShape = aShape;
1889   myShape = aShape;
1890   myUnifyEdges = UnifyEdges;
1891   myUnifyFaces = UnifyFaces;
1892   myConcatBSplines = ConcatBSplines;
1893
1894   myContext->Clear();
1895   myKeepShapes.Clear();
1896   myHistory->Clear();
1897 }
1898
1899 //=======================================================================
1900 //function : AllowInternalEdges
1901 //purpose  : 
1902 //=======================================================================
1903
1904 void ShapeUpgrade_UnifySameDomain::AllowInternalEdges (const Standard_Boolean theValue)
1905 {
1906   myAllowInternal = theValue;
1907 }
1908
1909 //=======================================================================
1910 //function : SetSafeInputMode
1911 //purpose  : 
1912 //=======================================================================
1913
1914 void ShapeUpgrade_UnifySameDomain::SetSafeInputMode(Standard_Boolean theValue)
1915 {
1916   mySafeInputMode = theValue;
1917 }
1918
1919 //=======================================================================
1920 //function : KeepShape
1921 //purpose  : 
1922 //=======================================================================
1923
1924 void ShapeUpgrade_UnifySameDomain::KeepShape(const TopoDS_Shape& theShape)
1925 {
1926   if (theShape.ShapeType() == TopAbs_EDGE || theShape.ShapeType() == TopAbs_VERTEX)
1927     myKeepShapes.Add(theShape);
1928 }
1929
1930 //=======================================================================
1931 //function : KeepShapes
1932 //purpose  : 
1933 //=======================================================================
1934
1935 void ShapeUpgrade_UnifySameDomain::KeepShapes(const TopTools_MapOfShape& theShapes)
1936 {
1937   for (TopTools_MapIteratorOfMapOfShape it(theShapes); it.More(); it.Next()) {
1938     if (it.Value().ShapeType() == TopAbs_EDGE || it.Value().ShapeType() == TopAbs_VERTEX)
1939       myKeepShapes.Add(it.Value());
1940   }
1941 }
1942
1943 //=======================================================================
1944 //function : UnifyFaces
1945 //purpose  : 
1946 //=======================================================================
1947
1948 void ShapeUpgrade_UnifySameDomain::UnifyFaces()
1949 {
1950   // creating map of edge faces for the whole shape
1951   TopTools_IndexedDataMapOfShapeListOfShape aGMapEdgeFaces;
1952   TopExp::MapShapesAndAncestors(myShape, TopAbs_EDGE, TopAbs_FACE, aGMapEdgeFaces);
1953   
1954   // unify faces in each shell separately
1955   TopExp_Explorer exps;
1956   for (exps.Init(myShape, TopAbs_SHELL); exps.More(); exps.Next())
1957     IntUnifyFaces(exps.Current(), aGMapEdgeFaces);
1958
1959   // gather all faces out of shells in one compound and unify them at once
1960   BRep_Builder aBB;
1961   TopoDS_Compound aCmp;
1962   aBB.MakeCompound(aCmp);
1963   Standard_Integer nbf = 0;
1964   for (exps.Init(myShape, TopAbs_FACE, TopAbs_SHELL); exps.More(); exps.Next(), nbf++)
1965     aBB.Add(aCmp, exps.Current());
1966
1967   if (nbf > 0)
1968     IntUnifyFaces(aCmp, aGMapEdgeFaces);
1969   
1970   myShape = myContext->Apply(myShape);
1971 }
1972
1973 //=======================================================================
1974 //function : SetFixWireModes
1975 //purpose  : 
1976 //=======================================================================
1977
1978 static void SetFixWireModes(ShapeFix_Face& theSff)
1979 {
1980   Handle(ShapeFix_Wire) aFixWire = theSff.FixWireTool();
1981   aFixWire->FixSelfIntersectionMode() = 0;
1982   aFixWire->FixNonAdjacentIntersectingEdgesMode() = 0;
1983   aFixWire->FixLackingMode() = 0;
1984   aFixWire->FixNotchedEdgesMode() = 0;
1985   aFixWire->ModifyTopologyMode() = Standard_False;
1986   aFixWire->ModifyRemoveLoopMode() = 0;
1987   aFixWire->FixGapsByRangesMode() = Standard_False;
1988   aFixWire->FixSmallMode() = 0;
1989 }
1990
1991 //=======================================================================
1992 //function : IntUnifyFaces
1993 //purpose  : 
1994 //=======================================================================
1995
1996 void ShapeUpgrade_UnifySameDomain::IntUnifyFaces(const TopoDS_Shape& theInpShape,
1997    TopTools_IndexedDataMapOfShapeListOfShape& theGMapEdgeFaces)
1998 {
1999   // creating map of edge faces for the shape
2000   TopTools_IndexedDataMapOfShapeListOfShape aMapEdgeFaces;
2001   TopExp::MapShapesAndAncestors(theInpShape, TopAbs_EDGE, TopAbs_FACE, aMapEdgeFaces);
2002
2003   // map of processed shapes
2004   TopTools_MapOfShape aProcessed;
2005
2006   // processing each face
2007   TopExp_Explorer exp;
2008   for (exp.Init(theInpShape, TopAbs_FACE); exp.More(); exp.Next()) {
2009     
2010     TopoDS_Face aFace = TopoDS::Face(exp.Current());
2011
2012     if (aProcessed.Contains(aFace))
2013       continue;
2014
2015     // Boundary edges for the new face
2016     TopTools_SequenceOfShape edges;
2017
2018     Standard_Integer dummy;
2019     AddOrdinaryEdges(edges, aFace, dummy);
2020
2021     // Faces to get unified with the current faces
2022     TopTools_SequenceOfShape faces;
2023
2024     // Add the current face for unification
2025     faces.Append(aFace);
2026
2027     // surface and location to construct result
2028     TopLoc_Location aBaseLocation;
2029     Handle(Geom_Surface) aBaseSurface = BRep_Tool::Surface(aFace,aBaseLocation);
2030     aBaseSurface = ClearRts(aBaseSurface);
2031     TopAbs_Orientation RefFaceOrientation = aFace.Orientation();
2032
2033     //Take original surface
2034     TopoDS_Face RefFace;
2035     BRep_Builder BB;
2036     BB.MakeFace(RefFace, aBaseSurface, aBaseLocation, 0.);
2037     RefFace.Orientation(RefFaceOrientation);
2038     TopTools_MapOfShape MapEdgesWithTemporaryPCurves; //map of edges not lying on RefFace
2039     //these edges may be updated by temporary pcurves
2040     
2041     Standard_Real Uperiod = (aBaseSurface->IsUPeriodic())? aBaseSurface->UPeriod() : 0.;
2042
2043     // find adjacent faces to union
2044     Standard_Integer i;
2045     for (i = 1; i <= edges.Length(); i++) {
2046       TopoDS_Edge edge = TopoDS::Edge(edges(i));
2047       if (BRep_Tool::Degenerated(edge))
2048         continue;
2049
2050       // get connectivity of the edge in the global shape
2051       const TopTools_ListOfShape& aGList = theGMapEdgeFaces.FindFromKey(edge);
2052       if (!myAllowInternal && (aGList.Extent() != 2 || myKeepShapes.Contains(edge))) {
2053         // non manifold case is not processed unless myAllowInternal
2054         continue;
2055       }
2056       //
2057       // Get the faces connected through the edge in the current shape
2058       const TopTools_ListOfShape& aList = aMapEdgeFaces.FindFromKey(edge);
2059       if (aList.Extent() < 2) {
2060         continue;
2061       }
2062
2063       // for a planar face create and store pcurve of edge on face
2064       // to speed up all operations
2065       if (!mySafeInputMode && aBaseSurface->IsKind(STANDARD_TYPE(Geom_Plane)))
2066         BRepLib::BuildPCurveForEdgeOnPlane(edge, aFace);
2067
2068       // get normal of the face to compare it with normals of other faces
2069       gp_Dir aDN1;
2070       //
2071       // take intermediate point on edge to compute the normal
2072       Standard_Real f, l;
2073       BRep_Tool::Range(edge, f, l);
2074       Standard_Real aTMid = (f + l) * .5;
2075       //
2076       Standard_Boolean bCheckNormals = GetNormalToSurface(aFace, edge, aTMid, aDN1);
2077       //
2078       // Process the faces
2079       TopTools_ListIteratorOfListOfShape anIter(aList);
2080       for (; anIter.More(); anIter.Next()) {
2081         
2082         TopoDS_Face aCheckedFace = TopoDS::Face(anIter.Value());
2083         if (aCheckedFace.IsSame(aFace))
2084           continue;
2085
2086         if (aProcessed.Contains(aCheckedFace))
2087           continue;
2088
2089         if (bCheckNormals) {
2090           // get normal of checked face using the same parameter on edge
2091           gp_Dir aDN2;
2092           if (GetNormalToSurface(aCheckedFace, edge, aTMid, aDN2)) {
2093             // and check if the adjacent faces are having approximately same normals
2094             Standard_Real anAngle = aDN1.Angle(aDN2);
2095             if (anAngle > myAngTol) {
2096               continue;
2097             }
2098           }
2099         }
2100         //
2101         if (IsSameDomain(aFace,aCheckedFace, myLinTol, myAngTol)) {
2102
2103           if (AddOrdinaryEdges(edges,aCheckedFace,dummy)) {
2104             // sequence edges is modified
2105             i = dummy;
2106           }
2107
2108           faces.Append(aCheckedFace);
2109           aProcessed.Add(aCheckedFace);
2110           break;
2111         }
2112       }
2113     }
2114
2115     if (faces.Length() > 1) {
2116       //Add correct pcurves for the reference surface to the edges of other faces
2117       AddPCurves(faces, RefFace, MapEdgesWithTemporaryPCurves);
2118       
2119       // fill in the connectivity map for selected faces
2120       TopTools_IndexedDataMapOfShapeListOfShape aMapEF;
2121       for (i = 1; i <= faces.Length(); i++) {
2122         TopExp::MapShapesAndAncestors(faces(i), TopAbs_EDGE, TopAbs_FACE, aMapEF);
2123       }
2124       // Collect keep edges and multi-connected edges, i.e. edges that are internal to
2125       // the set of selected faces and have connections to other faces.
2126       TopTools_ListOfShape aKeepEdges;
2127       for (i = 1; i <= aMapEF.Extent(); i++) {
2128         const TopTools_ListOfShape& aLF = aMapEF(i);
2129         if (aLF.Extent() == 2) {
2130           const TopoDS_Shape& aE = aMapEF.FindKey(i);
2131           const TopTools_ListOfShape& aGLF = theGMapEdgeFaces.FindFromKey(aE);
2132           if (aGLF.Extent() > 2 || myKeepShapes.Contains(aE)) {
2133             aKeepEdges.Append(aE);
2134           }
2135         }
2136       } 
2137       if (!aKeepEdges.IsEmpty()) {
2138         if  (!myAllowInternal) {
2139           // Remove from the selection the faces which have no other connect edges 
2140           // and contain multi-connected edges and/or keep edges.
2141           TopTools_MapOfShape anAvoidFaces;
2142           TopTools_ListIteratorOfListOfShape it(aKeepEdges);
2143           for (; it.More(); it.Next()) {
2144             const TopoDS_Shape& aE = it.Value();
2145             const TopTools_ListOfShape& aLF = aMapEF.FindFromKey(aE);
2146             anAvoidFaces.Add(aLF.First());
2147             anAvoidFaces.Add(aLF.Last());
2148           }
2149           for (i = 1; i <= faces.Length(); i++) {
2150             if (anAvoidFaces.Contains(faces(i))) {
2151               // update the boundaries of merged area, for that
2152               // remove from 'edges' the edges of this face and add to 'edges' 
2153               // the edges of this face that were not present in 'edges' before
2154               Standard_Boolean hasConnectAnotherFaces = Standard_False;
2155               TopExp_Explorer ex(faces(i), TopAbs_EDGE);
2156               for (; ex.More() && !hasConnectAnotherFaces; ex.Next()) {
2157                 TopoDS_Shape aE = ex.Current();
2158                 const TopTools_ListOfShape& aLF = aMapEF.FindFromKey(aE);
2159                 if (aLF.Extent() > 1) {
2160                   for (it.Init(aLF); it.More() && !hasConnectAnotherFaces; it.Next()) {
2161                     if (!anAvoidFaces.Contains(it.Value()))
2162                       hasConnectAnotherFaces = Standard_True;
2163                   }
2164                 }
2165               }
2166               if (!hasConnectAnotherFaces) {
2167                 AddOrdinaryEdges(edges, faces(i), dummy);
2168                 faces.Remove(i);
2169                 i--;
2170               }
2171             }
2172           }
2173           // check if the faces with keep edges contained in 
2174           // already updated the boundaries of merged area
2175           if (!faces.IsEmpty()) {
2176             TopTools_MapOfShape aMapFaces;
2177             for (i = 1; i <= faces.Length(); i++) {
2178               aMapFaces.Add(faces(i));
2179             }
2180             for (it.Init(aKeepEdges); it.More(); it.Next()) {
2181               const TopoDS_Shape& aE = it.Value();
2182               const TopTools_ListOfShape& aLF = aMapEF.FindFromKey(aE);
2183               if (aLF.Extent() < 2)
2184                 continue;
2185               if (aMapFaces.Contains(aLF.First()) && 
2186                   aMapFaces.Contains(aLF.Last())) {
2187                 for (i = 1; i <= faces.Length(); i++) {
2188                   if (faces(i).IsEqual(aLF.First()) ||
2189                       faces(i).IsEqual(aLF.Last())) {
2190                     AddOrdinaryEdges(edges, faces(i), dummy);
2191                     faces.Remove(i);
2192                     i--;
2193                   }
2194                 }
2195               }
2196             }
2197           }
2198         } //if (!myAllowInternal)
2199         else { //internal edges are allowed
2200           // add multi-connected and keep edges as internal in new face
2201           TopTools_ListIteratorOfListOfShape it(aKeepEdges);
2202           for (; it.More(); it.Next()) {
2203             const TopoDS_Shape& aE = it.Value();
2204             edges.Append(aE.Oriented(TopAbs_INTERNAL));
2205           }
2206         }
2207       } //if (!aKeepEdges.IsEmpty())
2208     } //if (faces.Length() > 1)
2209
2210     TopTools_IndexedDataMapOfShapeListOfShape aMapEF;
2211     for (i = 1; i <= faces.Length(); i++)
2212       TopExp::MapShapesAndUniqueAncestors(faces(i), TopAbs_EDGE, TopAbs_FACE, aMapEF);
2213     
2214     //Correct orientation of edges
2215     for (Standard_Integer ii = 1; ii <= edges.Length(); ii++)
2216     {
2217       const TopoDS_Shape& anEdge = edges(ii);
2218       const TopTools_ListOfShape& aLF = aMapEF.FindFromKey(anEdge);
2219       if (myAllowInternal &&
2220           myKeepShapes.Contains(anEdge) &&
2221           aLF.Extent() == 2)
2222         edges(ii).Orientation(TopAbs_INTERNAL);
2223       
2224       if (anEdge.Orientation() != TopAbs_INTERNAL)
2225         for (Standard_Integer jj = 1; jj <= faces.Length(); jj++)
2226         {
2227           const TopoDS_Shape aCurFace = faces(jj);
2228           Standard_Boolean found = Standard_False;
2229           TopExp_Explorer Explo(aCurFace, TopAbs_EDGE);
2230           for (; Explo.More(); Explo.Next())
2231           {
2232             const TopoDS_Shape& aCurEdge = Explo.Current();
2233             if (anEdge.IsSame(aCurEdge))
2234             {
2235               edges(ii) = aCurEdge;
2236               found = Standard_True;
2237               break;
2238             }
2239           }
2240           if (found)
2241             break;
2242         }
2243     }
2244
2245     //Exclude internal edges
2246     TopTools_IndexedMapOfShape InternalEdges;
2247     Standard_Integer ind_e = 1;
2248     while (ind_e <= edges.Length())
2249     {
2250       const TopoDS_Shape& anEdge = edges(ind_e);
2251       if (anEdge.Orientation() == TopAbs_INTERNAL)
2252       {
2253         InternalEdges.Add(anEdge);
2254         edges.Remove(ind_e);
2255       }
2256       else
2257         ind_e++;
2258     }    
2259
2260     if (RefFaceOrientation == TopAbs_REVERSED)
2261       for (Standard_Integer ii = 1; ii <= edges.Length(); ii++)
2262         edges(ii).Reverse();
2263     TopoDS_Face F_RefFace = RefFace;
2264     F_RefFace.Orientation(TopAbs_FORWARD);
2265     
2266     // all faces collected in the sequence. Perform union of faces
2267     if (faces.Length() > 1)
2268     {
2269       Standard_Real CoordTol = Precision::Confusion();
2270
2271       TopTools_IndexedDataMapOfShapeListOfShape VEmap;
2272       for (Standard_Integer ind = 1; ind <= edges.Length(); ind++)
2273         TopExp::MapShapesAndUniqueAncestors(edges(ind), TopAbs_VERTEX, TopAbs_EDGE, VEmap);
2274       
2275       //Perform relocating to new U-origin
2276       //Define boundaries in 2d space of RefFace
2277       if (Uperiod != 0.)
2278       {
2279         TopTools_MapOfShape edgesMap;
2280         CoordTol = ComputeMinEdgeSize(edges, F_RefFace, edgesMap);
2281         CoordTol /= 10.;
2282         CoordTol = Max(CoordTol, Precision::Confusion());
2283
2284         //try to find a real seam edge - if it exists, do nothing
2285         Standard_Boolean SeamFound = Standard_False;
2286         for (Standard_Integer ii = 1; ii <= faces.Length(); ii++)
2287         {
2288           const TopoDS_Face& face_ii = TopoDS::Face(faces(ii));
2289           TopoDS_Wire anOuterWire = BRepTools::OuterWire(face_ii);
2290           TopoDS_Iterator itw(anOuterWire);
2291           for (; itw.More(); itw.Next())
2292           {
2293             const TopoDS_Edge& anEdge = TopoDS::Edge(itw.Value());
2294             if (BRepTools::IsReallyClosed(anEdge, face_ii))
2295             {
2296               SeamFound = Standard_True;
2297               break;
2298             }
2299           }
2300         }
2301         
2302         if (!SeamFound)
2303         {
2304           //try to find the origin of U in 2d space
2305           //so that all the faces are in [origin, origin + Uperiod]
2306           Standard_Real Umax;
2307           Standard_Integer i_face_max;
2308           FindFaceWithMaxUbound(faces, F_RefFace, edgesMap, Umax, i_face_max);
2309           
2310           TopTools_MapOfShape UsedEdges;
2311           NCollection_DataMap<TopoDS_Shape, Handle(Geom2d_Curve)> EdgeNewPCurve;
2312
2313           //Relocate pcurves to new U-origin
2314           RelocatePCurvesToNewUorigin(edges, faces(i_face_max), F_RefFace, CoordTol, Uperiod,
2315                                       VEmap, EdgeNewPCurve, UsedEdges);
2316
2317           //PCurves from unused edges (may be degenerated edges)
2318           for (Standard_Integer ind = 1; ind <= edges.Length(); ind++)
2319           {
2320             const TopoDS_Edge& anEdge = TopoDS::Edge(edges(ind));
2321             if (!UsedEdges.Contains(anEdge))
2322             {
2323               Standard_Real fpar, lpar;
2324               Handle(Geom2d_Curve) aPCurve = BRep_Tool::CurveOnSurface(anEdge, F_RefFace, fpar, lpar);
2325               aPCurve = new Geom2d_TrimmedCurve(aPCurve, fpar, lpar);
2326               EdgeNewPCurve.Bind(anEdge, aPCurve);
2327             }
2328           }
2329           
2330           //Restore VEmap
2331           VEmap.Clear();
2332           for (Standard_Integer ind = 1; ind <= edges.Length(); ind++)
2333             TopExp::MapShapesAndUniqueAncestors(edges(ind), TopAbs_VERTEX, TopAbs_EDGE, VEmap);
2334
2335           //Find NewUmin and NewUmax
2336           Standard_Real NewUmin = RealLast(), NewUmax = RealFirst();
2337           for (Standard_Integer ii = 1; ii <= edges.Length(); ii++)
2338           {
2339             const Handle(Geom2d_Curve)& aPCurve = EdgeNewPCurve(edges(ii));
2340             UpdateBoundaries(aPCurve, aPCurve->FirstParameter(), aPCurve->LastParameter(),
2341                              NewUmin, NewUmax);
2342           }
2343           
2344           if (NewUmax - NewUmin < Uperiod - CoordTol &&
2345               !(-Precision::Confusion() < NewUmin && NewUmin < Uperiod+Precision::Confusion() &&
2346                 -Precision::Confusion() < NewUmax && NewUmax < Uperiod+Precision::Confusion()))
2347           {
2348             //we can build a face without seam edge:
2349             //update the edges with earlier computed relocated pcurves
2350             //fitting into (NewUorigin, NewUorigin + Uperiod)
2351             Standard_Real umin, umax, vmin, vmax;
2352             aBaseSurface->Bounds(umin, umax, vmin, vmax);
2353             Standard_Real RestSpaceInU = Uperiod - (NewUmax - NewUmin);
2354             Standard_Real NewUorigin = NewUmin - RestSpaceInU/2;
2355             if (NewUorigin < umin)
2356               NewUorigin = umin;
2357             Handle(Geom_Surface) NewSurf;
2358             if (NewUorigin == umin)
2359               NewSurf = aBaseSurface;
2360             else
2361               NewSurf = new Geom_RectangularTrimmedSurface(aBaseSurface,
2362                                                            NewUorigin, NewUorigin + Uperiod,
2363                                                            Standard_True); //trim U
2364             TopoDS_Face OldRefFace = RefFace;
2365             Handle(Geom2d_Curve) NullPCurve;
2366             RefFace.Nullify();
2367             BB.MakeFace(RefFace, NewSurf, aBaseLocation, 0.);
2368             for (Standard_Integer ii = 1; ii <= edges.Length(); ii++)
2369             {
2370               TopoDS_Edge anEdge = TopoDS::Edge(edges(ii));
2371               if (MapEdgesWithTemporaryPCurves.Contains(anEdge))
2372                 BB.UpdateEdge(anEdge, NullPCurve, OldRefFace, 0.);
2373               const Handle(Geom2d_Curve)& aPCurve = EdgeNewPCurve(anEdge);
2374               BB.UpdateEdge(anEdge, aPCurve, RefFace, 0.);
2375             }
2376           }
2377         } //if (!SeamFound)
2378       } //if (Uperiod != 0.)
2379       ////////////////////////////////////
2380       F_RefFace = RefFace;
2381       F_RefFace.Orientation(TopAbs_FORWARD);
2382       
2383       TopTools_SequenceOfShape NewFaces, NewWires;
2384       
2385       if (Uperiod == 0)
2386       {
2387         //Set the "period" for closed non-periodic surface
2388         TopLoc_Location aLoc;
2389         Handle(Geom_Surface) aSurf = BRep_Tool::Surface(RefFace, aLoc);
2390         if (aSurf->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)))
2391           aSurf = (Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurf))->BasisSurface();
2392         Standard_Real Ufirst, Ulast, Vfirst, Vlast;
2393         aSurf->Bounds(Ufirst, Ulast, Vfirst, Vlast);
2394         if (aSurf->IsUClosed())
2395           Uperiod = Ulast - Ufirst;
2396       }
2397
2398       TopTools_MapOfShape UsedEdges;
2399
2400       Standard_Real FaceUmin = RealLast();
2401       for (Standard_Integer ii = 1; ii <= edges.Length(); ii++)
2402       {
2403         const TopoDS_Edge& anEdge = TopoDS::Edge(edges(ii));
2404         BRepAdaptor_Curve2d aBAcurve(anEdge, F_RefFace);
2405         gp_Pnt2d aFirstPoint = aBAcurve.Value(aBAcurve.FirstParameter());
2406         gp_Pnt2d aLastPoint  = aBAcurve.Value(aBAcurve.LastParameter());
2407         if (aFirstPoint.X() < FaceUmin)
2408           FaceUmin = aFirstPoint.X();
2409         if (aLastPoint.X() < FaceUmin)
2410           FaceUmin = aLastPoint.X();
2411       }
2412
2413       //Building new wires from <edges>
2414       //and build faces
2415       while (!edges.IsEmpty())
2416       {
2417         //try to find non-degenerated edge
2418         TopoDS_Edge StartEdge = TopoDS::Edge(edges(1));
2419         Standard_Integer istart = 1;
2420         while (BRep_Tool::Degenerated(StartEdge) &&
2421                istart < edges.Length())
2422         {
2423           istart++;
2424           StartEdge = TopoDS::Edge(edges(istart));
2425         }
2426
2427         TopoDS_Wire aNewWire;
2428         BB.MakeWire(aNewWire);
2429         BB.Add(aNewWire, StartEdge);
2430         RemoveEdgeFromMap(StartEdge, VEmap);
2431         
2432         Standard_Real fpar, lpar;
2433         Handle(Geom2d_Curve) StartPCurve = BRep_Tool::CurveOnSurface(StartEdge, F_RefFace, fpar, lpar);
2434         TopoDS_Vertex StartVertex, CurVertex;
2435         TopExp::Vertices(StartEdge, StartVertex, CurVertex, Standard_True); //with orientation
2436         Standard_Real StartParam, CurParam;
2437         if (StartEdge.Orientation() == TopAbs_FORWARD)
2438         {
2439           StartParam = fpar; CurParam = lpar;
2440         }
2441         else
2442         {
2443           StartParam = lpar; CurParam = fpar;
2444         }
2445         gp_Pnt2d StartPoint = StartPCurve->Value(StartParam);
2446         gp_Pnt2d CurPoint   = StartPCurve->Value(CurParam);
2447         
2448         TopoDS_Edge CurEdge = StartEdge;
2449         for (;;) //loop till the end of current new wire
2450         {
2451           TopoDS_Edge NextEdge;
2452           gp_Pnt2d NextPoint;
2453           
2454           const TopTools_ListOfShape& Elist = VEmap.FindFromKey(CurVertex);
2455           TopTools_ListIteratorOfListOfShape itl(Elist);
2456           if (Elist.IsEmpty())
2457           {
2458             if (CurVertex.IsSame(StartVertex))
2459             {
2460               //Points of two vertices coincide in 3d but may be not in 2d
2461               if (Uperiod != 0. &&
2462                   Abs(StartPoint.X() - CurPoint.X()) > Uperiod/2) //end of parametric space
2463               {
2464                 //<edges> do not contain seams => we must reconstruct the seam up to <NextEdge>
2465                 gp_Pnt2d StartOfNextEdge;
2466                 TopoDS_Vertex LastVertexOfSeam;
2467                 ReconstructMissedSeam(edges, UsedEdges, F_RefFace, CurVertex,
2468                                       CurPoint, Uperiod, FaceUmin, CoordTol,
2469                                       NextEdge, aNewWire, NextPoint,
2470                                       StartOfNextEdge, LastVertexOfSeam, VEmap);
2471               }
2472               else
2473               {
2474                 break; //end of wire
2475               }
2476             }
2477           }
2478           
2479           if (NextEdge.IsNull())
2480           {
2481             Standard_Boolean EndOfWire = Standard_False;
2482             
2483             TopTools_ListOfShape TmpElist, TrueElist;
2484             //<TrueElist> will be the list of candidates to become <NextEdge>
2485             for (itl.Initialize(Elist); itl.More(); itl.Next())
2486             {
2487               const TopoDS_Edge& anEdge = TopoDS::Edge(itl.Value());
2488               if (UsedEdges.Contains(anEdge))
2489                 continue;
2490               TopoDS_Vertex aFirstVertex = TopExp::FirstVertex(anEdge, Standard_True);
2491               if (!aFirstVertex.IsSame(CurVertex))
2492                 continue;
2493               TmpElist.Append(anEdge);
2494             }
2495             if (TmpElist.Extent() <= 1 ||
2496                 Uperiod != 0.)
2497               TrueElist.Assign(TmpElist);
2498             else
2499             {
2500               //we must choose the closest direction - the biggest angle
2501               Standard_Real MaxAngle = RealFirst();
2502               TopoDS_Edge TrueEdge;
2503               Handle(Geom2d_Curve) CurPCurve = BRep_Tool::CurveOnSurface(CurEdge, F_RefFace, fpar, lpar);
2504               CurParam = (CurEdge.Orientation() == TopAbs_FORWARD)? lpar : fpar;
2505               gp_Vec2d CurDir;
2506               CurPCurve->D1(CurParam, CurPoint, CurDir);
2507               CurDir.Normalize();
2508               if (CurEdge.Orientation() == TopAbs_REVERSED)
2509                 CurDir.Reverse();
2510               for (itl.Initialize(TmpElist); itl.More(); itl.Next())
2511               {
2512                 const TopoDS_Edge& anEdge = TopoDS::Edge(itl.Value());
2513                 Handle(Geom2d_Curve) aPCurve = BRep_Tool::CurveOnSurface(anEdge, F_RefFace, fpar, lpar);
2514                 Standard_Real aParam = (anEdge.Orientation() == TopAbs_FORWARD)? fpar : lpar;
2515                 gp_Pnt2d aPoint;
2516                 gp_Vec2d aDir;
2517                 aPCurve->D1(aParam, aPoint, aDir);
2518                 aDir.Normalize();
2519                 if (anEdge.Orientation() == TopAbs_REVERSED)
2520                   aDir.Reverse();
2521                 Standard_Real anAngle = CurDir.Angle(aDir);
2522                 if (anAngle > MaxAngle)
2523                 {
2524                   MaxAngle = anAngle;
2525                   TrueEdge = anEdge;
2526                 }
2527               }
2528               TrueElist.Append(TrueEdge);
2529             }
2530
2531             //Find next edge in TrueElist
2532             for (itl.Initialize(TrueElist); itl.More(); itl.Next())
2533             {
2534               const TopoDS_Edge& anEdge = TopoDS::Edge(itl.Value());
2535               
2536               Handle(Geom2d_Curve) aPCurve = BRep_Tool::CurveOnSurface(anEdge, F_RefFace, fpar, lpar);
2537               Standard_Real aParam = (anEdge.Orientation() == TopAbs_FORWARD)? fpar : lpar;
2538               gp_Pnt2d aPoint = aPCurve->Value(aParam);
2539               Standard_Real Diff = Abs(aPoint.X() - CurPoint.X());
2540               if (Uperiod != 0. &&
2541                   Diff > CoordTol &&
2542                   Abs(Diff - Uperiod) > CoordTol) //may be is is a deg.vertex
2543                 continue;
2544               
2545               //Check: may be <CurPoint> and <aPoint> are on Period from each other
2546               if (Uperiod != 0. && Diff > Uperiod/2) //end of parametric space
2547               {
2548                 //<edges> do not contain seams => we must reconstruct the seam up to <NextEdge>
2549                 gp_Pnt2d StartOfNextEdge;
2550                 TopoDS_Vertex LastVertexOfSeam;
2551                 ReconstructMissedSeam(edges, UsedEdges, F_RefFace, CurVertex,
2552                                       CurPoint, Uperiod, FaceUmin, CoordTol,
2553                                       NextEdge, aNewWire, NextPoint,
2554                                       StartOfNextEdge, LastVertexOfSeam, VEmap);
2555                 
2556                 //Check: may be it is the end
2557                 if (LastVertexOfSeam.IsSame(StartVertex) &&
2558                     Abs(StartPoint.X() - StartOfNextEdge.X()) < Uperiod/2)
2559                   EndOfWire = Standard_True;
2560                 
2561                 break;
2562               } //if (Uperiod != 0. && Abs(aPoint.X() - CurPoint.X()) > Uperiod/2)
2563               else
2564               {
2565                 NextEdge = anEdge;
2566                 Standard_Real LastParam = (NextEdge.Orientation() == TopAbs_FORWARD)? lpar : fpar;
2567                 NextPoint = aPCurve->Value(LastParam);
2568                 break;
2569               }
2570             } //for (itl.Initialize(TrueElist); itl.More(); itl.Next())
2571             
2572             if (EndOfWire)
2573               break;
2574           }
2575           
2576           if (NextEdge.IsNull())
2577           {
2578             //throw Standard_ConstructionError("Construction of unified wire failed: no proper edge found");
2579             return;
2580           }
2581           else
2582           {
2583             CurPoint = NextPoint;
2584             CurEdge = NextEdge;
2585             CurVertex = TopExp::LastVertex(CurEdge, Standard_True); //with orientation
2586             BB.Add(aNewWire, CurEdge);
2587             UsedEdges.Add(CurEdge);
2588             RemoveEdgeFromMap(CurEdge, VEmap);
2589           }
2590         } //for (;;)
2591         UsedEdges.Add(StartEdge);
2592         
2593         //Remove used edges from sequence
2594         Standard_Integer ind = 1;
2595         while (ind <= edges.Length())
2596         {
2597           if (UsedEdges.Contains(edges(ind)))
2598             edges.Remove(ind);
2599           else
2600             ind++;
2601         }
2602         
2603         //add just built wire to current face or save it in the sequence of wires
2604         Standard_Boolean EdgeOnBoundOfSurfFound = Standard_False;
2605         TopoDS_Iterator itw(aNewWire);
2606         for (; itw.More(); itw.Next())
2607         {
2608           const TopoDS_Edge& anEdge = TopoDS::Edge(itw.Value());
2609           if (BRep_Tool::IsClosed(anEdge, RefFace))
2610           {
2611             EdgeOnBoundOfSurfFound = Standard_True;
2612             break;
2613           }
2614         }
2615         if (EdgeOnBoundOfSurfFound) //this wire can not be a hole
2616         {
2617           TopLoc_Location aLoc;
2618           Handle(Geom_Surface) aSurf = BRep_Tool::Surface(RefFace, aLoc);
2619           TopoDS_Face aResult;
2620           BB.MakeFace(aResult,aSurf,aLoc,0);
2621           BB.Add(aResult, aNewWire);
2622           aResult.Orientation(RefFaceOrientation);
2623           NewFaces.Append(aResult);
2624         }
2625         else //may be this wire is a hole
2626         {
2627           NewWires.Append(aNewWire);
2628         }
2629       } //while (!edges.IsEmpty())
2630
2631       //Build wires from internal edges
2632       TopTools_IndexedDataMapOfShapeListOfShape IntVEmap;
2633       for (Standard_Integer ii = 1; ii <= InternalEdges.Extent(); ii++)
2634         TopExp::MapShapesAndAncestors(InternalEdges(ii), TopAbs_VERTEX, TopAbs_EDGE, IntVEmap);
2635       TopTools_SequenceOfShape InternalWires;
2636       while (!InternalEdges.IsEmpty())
2637       {
2638         TopoDS_Edge aFirstEdge = TopoDS::Edge(InternalEdges(1));
2639         InternalEdges.RemoveFromIndex(1);
2640         TopoDS_Wire anInternalWire;
2641         BB.MakeWire(anInternalWire);
2642         BB.Add(anInternalWire, aFirstEdge);
2643         TopoDS_Edge EndEdges [2];
2644         EndEdges[0] = EndEdges[1] = aFirstEdge;
2645         TopoDS_Vertex VV [2];
2646         TopExp::Vertices(aFirstEdge, VV[0], VV[1]);
2647         for (;;)
2648         {
2649           if (VV[0].IsSame(VV[1])) //closed wire
2650             break;
2651           Standard_Boolean found = Standard_False;
2652           for (Standard_Integer ii = 0; ii < 2; ii++)
2653           {
2654             const TopTools_ListOfShape& Elist = IntVEmap.FindFromKey(VV[ii]);
2655             TopTools_ListIteratorOfListOfShape itl(Elist);
2656             for (; itl.More(); itl.Next())
2657             {
2658               TopoDS_Edge anEdge = TopoDS::Edge(itl.Value());
2659               if (anEdge.IsSame(EndEdges[ii]))
2660                 continue;
2661               found = Standard_True;
2662               InternalEdges.RemoveKey(anEdge);
2663               BB.Add(anInternalWire, anEdge);
2664               TopoDS_Vertex V1, V2;
2665               TopExp::Vertices(anEdge, V1, V2);
2666               VV[ii] = (V1.IsSame(VV[ii]))? V2 : V1;
2667               EndEdges[ii] = anEdge;
2668               break;
2669             }
2670           }
2671           if (!found) //end of open wire
2672             break;
2673         }
2674         InternalWires.Append(anInternalWire);
2675       }
2676
2677       //Insert new faces instead of old ones
2678       if (NewFaces.IsEmpty())
2679       {
2680         //one face without seam
2681         TopLoc_Location aLoc;
2682         Handle(Geom_Surface) aSurf = BRep_Tool::Surface(RefFace, aLoc);
2683         TopoDS_Face aResult;
2684         BB.MakeFace(aResult,aSurf,aLoc,0.);
2685         for (Standard_Integer ii = 1; ii <= NewWires.Length(); ii++)
2686           BB.Add(aResult, NewWires(ii));
2687         for (Standard_Integer ii = 1; ii <= InternalWires.Length(); ii++)
2688           BB.Add(aResult, InternalWires(ii));
2689         aResult.Orientation(RefFaceOrientation);
2690         myContext->Merge(faces, aResult);
2691       }
2692       else if (NewFaces.Length() == 1)
2693       {
2694         for (Standard_Integer ii = 1; ii <= NewWires.Length(); ii++)
2695           BB.Add(NewFaces(1), NewWires(ii));
2696         for (Standard_Integer ii = 1; ii <= InternalWires.Length(); ii++)
2697           BB.Add(NewFaces(1), InternalWires(ii));
2698         myContext->Merge(faces, NewFaces(1));
2699       }
2700       else
2701       {
2702         //Insert new wires and internal wires into correspondent faces
2703         InsertWiresIntoFaces(NewWires, NewFaces, RefFace);
2704         InsertWiresIntoFaces(InternalWires, NewFaces, RefFace);
2705         
2706         NCollection_Sequence<TopTools_MapOfShape> Emaps;
2707         for (Standard_Integer ii = 1; ii <= faces.Length(); ii++)
2708         {
2709           TopTools_MapOfShape aEmap;
2710           TopExp::MapShapes(faces(ii), aEmap);
2711           Emaps.Append(aEmap);
2712         }
2713         for (Standard_Integer ii = 1; ii <= NewFaces.Length(); ii++)
2714         {
2715           TopTools_SequenceOfShape facesForThisFace;
2716           TopTools_MapOfShape UsedFaces;
2717           TopExp_Explorer Explo(NewFaces(ii), TopAbs_EDGE);
2718           for (; Explo.More(); Explo.Next())
2719           {
2720             const TopoDS_Edge& anEdge = TopoDS::Edge(Explo.Current());
2721             if (BRep_Tool::Degenerated(anEdge) ||
2722                 BRep_Tool::IsClosed(anEdge, RefFace))
2723               continue;
2724             Standard_Integer jj;
2725             for (jj = 1; jj <= Emaps.Length(); jj++)
2726               if (Emaps(jj).Contains(anEdge))
2727                 break;
2728             if (UsedFaces.Add(faces(jj)))
2729               facesForThisFace.Append(faces(jj));
2730           }
2731           myContext->Merge(facesForThisFace, NewFaces(ii));
2732         }
2733       }
2734     } //if (faces.Length() > 1)
2735   } // end processing each face
2736 }
2737
2738 //=======================================================================
2739 //function : UnifyEdges
2740 //purpose  : 
2741 //=======================================================================
2742 void ShapeUpgrade_UnifySameDomain::UnifyEdges()
2743 {
2744   TopoDS_Shape aRes = myContext->Apply(myShape);
2745   // creating map of edge faces
2746   TopTools_IndexedDataMapOfShapeListOfShape aMapEdgeFaces;
2747   TopExp::MapShapesAndAncestors(aRes, TopAbs_EDGE, TopAbs_FACE, aMapEdgeFaces);
2748   // creating map of vertex edges
2749   TopTools_IndexedDataMapOfShapeListOfShape aMapEdgesVertex;
2750   TopExp::MapShapesAndUniqueAncestors(aRes, TopAbs_VERTEX, TopAbs_EDGE, aMapEdgesVertex);
2751   // creating map of vertex faces
2752   TopTools_IndexedDataMapOfShapeListOfShape aVFmap;
2753   TopExp::MapShapesAndUniqueAncestors(aRes, TopAbs_VERTEX, TopAbs_FACE, aVFmap);
2754
2755   if (mySafeInputMode)
2756     UpdateMapOfShapes(myKeepShapes, myContext);
2757
2758   // Sequence of the edges of the shape
2759   TopTools_SequenceOfShape aSeqEdges;
2760   const Standard_Integer aNbE = aMapEdgeFaces.Extent();
2761   for (Standard_Integer i = 1; i <= aNbE; ++i)
2762     aSeqEdges.Append(aMapEdgeFaces.FindKey(i));
2763
2764   // Prepare map of shared vertices (with the number of connected edges greater then 2)
2765   TopTools_MapOfShape aSharedVert;
2766   CheckSharedVertices(aSeqEdges, aMapEdgesVertex, myKeepShapes, aSharedVert);
2767   // Merge the edges avoiding removal of the shared vertices
2768   Standard_Boolean isMerged = MergeSeq(aSeqEdges, aVFmap, myAngTol, myLinTol, myConcatBSplines,
2769                                        mySafeInputMode, myContext, aSharedVert);
2770   // Collect faces to rebuild
2771   TopTools_IndexedMapOfShape aChangedFaces;
2772   if (isMerged)
2773   {
2774     for (Standard_Integer i = 1; i <= aNbE; ++i)
2775     {
2776       const TopoDS_Shape& aE = aMapEdgeFaces.FindKey(i);
2777       if (myContext->IsRecorded(aE))
2778       {
2779         TopTools_ListIteratorOfListOfShape it(aMapEdgeFaces(i));
2780         for (; it.More(); it.Next())
2781           aChangedFaces.Add(it.Value());
2782       }
2783     }
2784   }
2785
2786   // fix changed faces and replace them in the local context
2787   Standard_Real aPrec = Precision::Confusion();
2788   for (Standard_Integer i = 1; i <= aChangedFaces.Extent(); i++) {
2789     TopoDS_Face aFace = TopoDS::Face(myContext->Apply(aChangedFaces.FindKey(i)));
2790     if (aFace.IsNull())
2791       continue;
2792
2793     // for a planar face create and store pcurve of edge on face
2794     // to speed up all operations; but this is allowed only when non-safe mode in force
2795     if (!mySafeInputMode)
2796     {
2797       TopLoc_Location aLoc;
2798       Handle(Geom_Surface) aSurface = BRep_Tool::Surface(aFace, aLoc);
2799       aSurface = ClearRts(aSurface);
2800       if (aSurface->IsKind(STANDARD_TYPE(Geom_Plane)))
2801       {
2802         TopTools_ListOfShape aLE;
2803         for (TopExp_Explorer anEx(aFace, TopAbs_EDGE); anEx.More(); anEx.Next())
2804           aLE.Append(anEx.Current());
2805         BRepLib::BuildPCurveForEdgesOnPlane(aLE, aFace);
2806       }
2807     }
2808
2809     ShapeFix_Face sff(aFace);
2810     if (mySafeInputMode)
2811       sff.SetContext(myContext);
2812     sff.SetPrecision(aPrec);
2813     sff.SetMinTolerance(aPrec);
2814     sff.SetMaxTolerance(Max(1., aPrec*1000.));
2815     sff.FixOrientationMode() = 0;
2816     sff.FixAddNaturalBoundMode() = 0;
2817     sff.FixIntersectingWiresMode() = 0;
2818     sff.FixLoopWiresMode() = 0;
2819     sff.FixSplitFaceMode() = 0;
2820     sff.FixPeriodicDegeneratedMode() = 0;
2821     SetFixWireModes(sff);
2822     sff.Perform();
2823     TopoDS_Shape aNewFace = sff.Face();
2824     myContext->Replace(aFace,aNewFace);
2825   }
2826
2827   if (aChangedFaces.Extent() > 0) {
2828     // fix changed shell and replace it in the local context
2829     TopoDS_Shape aRes1 = myContext->Apply(aRes);
2830     Standard_Boolean isChanged = Standard_False;
2831     TopExp_Explorer expsh;
2832     for (expsh.Init(aRes1, TopAbs_SHELL); expsh.More(); expsh.Next()) {
2833       TopoDS_Shell aShell = TopoDS::Shell(expsh.Current());
2834       Handle(ShapeFix_Shell) sfsh = new ShapeFix_Shell;
2835       sfsh->FixFaceOrientation(aShell);
2836       TopoDS_Shape aNewShell = sfsh->Shell();
2837       if (!aNewShell.IsSame(aShell)) {
2838         myContext->Replace(aShell, aNewShell);
2839         isChanged = Standard_True;
2840       }
2841     }
2842     if (isChanged)
2843       aRes1 = myContext->Apply(aRes1);
2844     myContext->Replace(myShape, aRes1);
2845   }
2846
2847   myShape = myContext->Apply(myShape);
2848 }
2849
2850 //=======================================================================
2851 //function : Build
2852 //purpose  : builds the resulting shape
2853 //=======================================================================
2854 void ShapeUpgrade_UnifySameDomain::Build() 
2855 {
2856   if (myUnifyFaces)
2857     UnifyFaces();
2858   if (myUnifyEdges)
2859     UnifyEdges();
2860
2861   // Fill the history of modifications during the operation
2862   FillHistory();
2863 }
2864
2865 //=======================================================================
2866 //function : FillHistory
2867 //purpose  : Fill the history of modifications during the operation
2868 //=======================================================================
2869 void ShapeUpgrade_UnifySameDomain::FillHistory()
2870 {
2871   if (myHistory.IsNull())
2872     // History is not requested
2873     return;
2874
2875   // Only Vertices, Edges and Faces can be modified during unification.
2876   // Thus, only these kind of shapes should be checked.
2877
2878   // Get history from the context.
2879   // It contains all modifications of the operation. Some of these
2880   // modifications become not relevant and should be filtered.
2881   Handle(BRepTools_History) aCtxHistory = myContext->History();
2882
2883   // Explore the history of the context and fill
2884   // the history of UnifySameDomain algorithm
2885   Handle(BRepTools_History) aUSDHistory = new BRepTools_History();
2886
2887   // Map all Vertices, Edges, Faces and Solids in the input shape
2888   TopTools_IndexedMapOfShape aMapInputShape;
2889   TopExp::MapShapes(myInitShape, TopAbs_VERTEX, aMapInputShape);
2890   TopExp::MapShapes(myInitShape, TopAbs_EDGE  , aMapInputShape);
2891   TopExp::MapShapes(myInitShape, TopAbs_FACE  , aMapInputShape);
2892   TopExp::MapShapes(myInitShape, TopAbs_SOLID , aMapInputShape);
2893
2894   // Map all Vertices, Edges, Faces and Solids in the result shape
2895   TopTools_IndexedMapOfShape aMapResultShapes;
2896   TopExp::MapShapes(myShape, TopAbs_VERTEX, aMapResultShapes);
2897   TopExp::MapShapes(myShape, TopAbs_EDGE  , aMapResultShapes);
2898   TopExp::MapShapes(myShape, TopAbs_FACE  , aMapResultShapes);
2899   TopExp::MapShapes(myShape, TopAbs_SOLID , aMapResultShapes);
2900
2901   // Iterate on all input shapes and get their modifications
2902   Standard_Integer i, aNb = aMapInputShape.Extent();
2903   for (i = 1; i <= aNb; ++i)
2904   {
2905     const TopoDS_Shape& aS = aMapInputShape(i);
2906
2907     // Check the shape itself to be present in the result
2908     if (aMapResultShapes.Contains(aS))
2909     {
2910       // The shape is present in the result as is, thus has not been modified
2911       continue;
2912     }
2913
2914     // Check if the shape has been modified during the operation
2915     const TopTools_ListOfShape& aLSImages = aCtxHistory->Modified(aS);
2916     if (aLSImages.IsEmpty())
2917     {
2918       // The shape has not been modified and not present in the result,
2919       // thus it has been removed
2920       aUSDHistory->Remove(aS);
2921       continue;
2922     }
2923
2924     // Check the images of the shape to be present in the result
2925     Standard_Boolean bRemoved = Standard_True;
2926     TopTools_ListIteratorOfListOfShape aItLSIm(aLSImages);
2927     for (; aItLSIm.More(); aItLSIm.Next())
2928     {
2929       const TopoDS_Shape& aSIm = aItLSIm.Value();
2930       if (aMapResultShapes.Contains(aSIm))
2931       {
2932         if (!aSIm.IsSame(aS))
2933           // Image is found in the result, thus the shape has been modified
2934           aUSDHistory->AddModified(aS, aSIm);
2935         bRemoved = Standard_False;
2936       }
2937     }
2938
2939     if (bRemoved)
2940     {
2941       // No images are found in the result, thus the shape has been removed
2942       aUSDHistory->Remove(aS);
2943     }
2944   }
2945
2946   // Merge the history of the operation into global history
2947   myHistory->Merge(aUSDHistory);
2948 }