5885154d6180728fc6ef3f6b65502d9dd751989b
[occt.git] / src / ChFi3d / ChFi3d_Builder_0.cxx
1 // Created on: 1993-12-16
2 // Created by: Isabelle GRIGNON
3 // Copyright (c) 1993-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
5 //
6 // This file is part of Open CASCADE Technology software library.
7 //
8 // This library is free software; you can redistribute it and/or modify it under
9 // the terms of the GNU Lesser General Public License version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
13 //
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
16
17 //  modified by ofv - Thu Feb 26 11:18:16 2004 OCC5246
18 //  Modified by skv - Fri Oct 24 14:24:47 2003 OCC4077
19 //  Modified by skv - Mon Jun 16 15:50:44 2003 OCC615
20
21 #include <ChFi3d.hxx>
22 #include <Precision.hxx>
23
24 #include <Standard_NotImplemented.hxx>
25 #include <Standard_ConstructionError.hxx>
26
27 #include <gp.hxx>
28 #include <gp_Circ.hxx>
29 #include <gp_Elips.hxx>
30 #include <gp_Lin.hxx>
31 #include <gp_Pnt.hxx>
32 #include <gp_Pnt2d.hxx>
33 #include <gp_Lin2d.hxx>
34 #include <ElCLib.hxx>
35 #include <ElSLib.hxx>
36 #include <BSplCLib.hxx>
37 #include <GeomLib.hxx>
38
39 #include <TColgp_Array1OfPnt2d.hxx>
40 #include <TColgp_Array1OfPnt.hxx>
41 #include <TColgp_Array1OfXYZ.hxx>
42 #include <TColStd_Array1OfInteger.hxx>
43 #include <TColStd_Array1OfReal.hxx>
44
45 #include <Geom_TrimmedCurve.hxx>
46 #include <Geom_BSplineCurve.hxx>
47 #include <Geom_Surface.hxx>
48 #include <Geom_CylindricalSurface.hxx>
49 #include <Geom_RectangularTrimmedSurface.hxx>
50 #include <Geom_Plane.hxx>
51 #include <Geom_Line.hxx>
52 #include <Geom_Circle.hxx>
53 #include <Geom_Ellipse.hxx>
54 #include <Geom2d_BezierCurve.hxx>
55 #include <Geom2d_BSplineCurve.hxx>
56 #include <Geom2d_Line.hxx>
57 #include <Geom2d_Circle.hxx>
58 #include <Geom2d_Ellipse.hxx>
59 #include <Geom2d_Hyperbola.hxx>
60 #include <Geom2d_Parabola.hxx>
61 #include <Geom2d_TrimmedCurve.hxx>
62 #include <Geom2d_Line.hxx>
63 #include <Geom2d_OffsetCurve.hxx>
64 #include <Geom2dAdaptor_Curve.hxx>
65 #include <Geom2dAdaptor_HCurve.hxx>
66 #include <Adaptor3d_TopolTool.hxx>
67 #include <Adaptor3d_CurveOnSurface.hxx>
68 #include <Adaptor3d_HCurveOnSurface.hxx>
69 #include <GeomAdaptor_HSurface.hxx>
70
71 #include <FairCurve_Batten.hxx>
72 #include <FairCurve_AnalysisCode.hxx>
73 #include <Convert_ParameterisationType.hxx>
74 #include <GeomConvert_CompCurveToBSplineCurve.hxx>
75 #include <GeomConvert.hxx>
76 #include <GeomLib_Interpolate.hxx>
77 #include <GeomAPI_ProjectPointOnSurf.hxx>
78 #include <GeomAPI_ProjectPointOnCurve.hxx>
79 #include <GC_MakeCircle.hxx>
80 #include <BRepAdaptor_Curve.hxx>
81 #include <BRepAdaptor_HCurve.hxx>
82 #include <BRepAdaptor_HCurve2d.hxx>
83 #include <BRepAdaptor_Surface.hxx>
84 #include <BRepTopAdaptor_HVertex.hxx>
85 #include <BRepTopAdaptor_TopolTool.hxx>
86 #include <LocalAnalysis_SurfaceContinuity.hxx>
87
88 #include <BRep_Tool.hxx>
89 #include <BRep_Builder.hxx>
90 #include <BRepTools.hxx>
91 #include <BRepTools_WireExplorer.hxx>
92 #include <BRepLib.hxx>
93 #include <BRepLib_MakeEdge.hxx>
94 #include <BRepLib_MakeWire.hxx>
95 #include <BRepLib_MakeFace.hxx>
96
97 #include <TopAbs.hxx>
98 #include <TopoDS_Shape.hxx>
99 #include <TopoDS_Edge.hxx>
100 #include <TopoDS_Vertex.hxx>
101 #include <TopoDS_Wire.hxx>
102 #include <TopoDS_Face.hxx>
103 #include <TopExp.hxx>
104 #include <TopExp_Explorer.hxx>
105 #include <TopTools_Array1OfShape.hxx>
106
107
108 #include <GeomAbs_Shape.hxx>
109 #include <Bnd_Box2d.hxx>
110
111 //#include <math_FunctionSample.hxx>
112 //#include <math_FunctionAllRoots.hxx>
113 #include <GCPnts_AbscissaPoint.hxx>
114
115 #include <IntCurveSurface_TheQuadCurvFuncOfTheQuadCurvExactHInter.hxx>
116 #include <IntCurveSurface_HInter.hxx>
117 #include <IntCurveSurface_IntersectionPoint.hxx>
118 #include <IntSurf_Quadric.hxx>
119 #include <IntSurf_PntOn2S.hxx>
120 #include <IntSurf_LineOn2S.hxx>
121 #include <IntAna_QuadQuadGeo.hxx>
122 #include <IntAna2d_AnaIntersection.hxx>
123 #include <IntRes2d_IntersectionPoint.hxx>
124 #include <IntWalk_PWalking.hxx>
125 #include <IntPatch_WLine.hxx>
126 #include <Geom2dInt_GInter.hxx>
127 #include <GeomInt_WLApprox.hxx>
128 #include <GeomInt_IntSS.hxx>
129 #include <AppParCurves_MultiBSpCurve.hxx>
130 #include <Approx_SameParameter.hxx>
131
132 #include <TopAbs.hxx>
133 #include <TopoDS_Shape.hxx>
134 #include <TopoDS_Edge.hxx>
135 #include <TopExp.hxx>
136
137 #include <TopOpeBRepDS.hxx>
138 #include <TopOpeBRepDS_Surface.hxx>
139 #include <TopOpeBRepDS_Point.hxx>
140 #include <TopOpeBRepDS_SolidSurfaceInterference.hxx>
141 #include <TopOpeBRepDS_CurvePointInterference.hxx>
142 #include <TopOpeBRepDS_ListOfInterference.hxx>
143 #include <TopOpeBRepDS_InterferenceIterator.hxx>
144 #include <TopOpeBRepTool_TOOL.hxx>
145 #include <ProjLib_ProjectedCurve.hxx>
146
147 #include <BRepBlend_PointOnRst.hxx>
148
149 #include <ChFiDS_HData.hxx>
150 #include <ChFiDS_SurfData.hxx>
151 #include <ChFiDS_FaceInterference.hxx>
152 #include <ChFiDS_Spine.hxx>
153 #include <ChFiDS_FilSpine.hxx>
154 #include <ChFiDS_SequenceOfSurfData.hxx>
155 #include <ChFiDS_Regul.hxx>
156 #include <Law_Function.hxx>
157 #include <Law_Composite.hxx>
158 #include <GeomAPI_PointsToBSpline.hxx>
159 #include <GeomLProp_CLProps.hxx>
160
161 #include <ChFi3d_Builder_0.hxx>
162
163 #ifdef OCCT_DEBUG
164 extern Standard_Boolean ChFi3d_GetcontextFORCEBLEND(); 
165 extern Standard_Boolean ChFi3d_GettraceDRAWINT();
166 extern Standard_Boolean ChFi3d_GettraceDRAWENLARGE();
167 extern Standard_Boolean ChFi3d_GettraceDRAWSPINE();
168 extern Standard_Real  t_sameparam, t_batten;
169 extern void ChFi3d_SettraceDRAWINT(const Standard_Boolean b);
170 extern void ChFi3d_SettraceDRAWSPINE(const Standard_Boolean b);
171 #endif
172
173 #include <stdio.h>
174
175 #include <GeomAdaptor_HCurve.hxx>
176 #include <BRepAdaptor_HSurface.hxx>
177 #include <TopOpeBRepDS_SurfaceCurveInterference.hxx>
178
179 //=======================================================================
180 //function : ChFi3d_InPeriod
181 //purpose  : 
182 //=======================================================================
183 Standard_Real ChFi3d_InPeriod(const Standard_Real U, 
184   const Standard_Real UFirst, 
185   const Standard_Real ULast,
186   const Standard_Real Eps)
187 {
188   const Standard_Real period = ULast - UFirst;
189   Standard_Real u = U;
190   while (Eps < (UFirst-u)) u += period;
191   while (Eps > (ULast -u)) u -= period;
192   if ( u < UFirst) u = UFirst;
193   return u;
194 }
195 //=======================================================================
196 //function : Box 
197 //purpose  : Calculation of min/max uv of the fillet to intersect.
198 //=======================================================================
199 void ChFi3d_Boite(const gp_Pnt2d& p1,const gp_Pnt2d& p2,
200   Standard_Real& mu,Standard_Real& Mu,
201   Standard_Real& mv,Standard_Real& Mv)
202 {
203   mu = Min(p1.X(),p2.X()); Mu = Max(p1.X(),p2.X());
204   mv = Min(p1.Y(),p2.Y()); Mv = Max(p1.Y(),p2.Y());
205 }
206 //=======================================================================
207 //function : Box
208 //purpose  : Calculation of min/max uv of the fillet to intersect.
209 //=======================================================================
210 void ChFi3d_Boite(const gp_Pnt2d& p1,const gp_Pnt2d& p2,
211   const gp_Pnt2d& p3,const gp_Pnt2d& p4,
212   Standard_Real& Du,Standard_Real& Dv,
213   Standard_Real& mu,Standard_Real& Mu,
214   Standard_Real& mv,Standard_Real& Mv)
215 {
216   Standard_Real a,b;
217   a = Min(p1.X(),p2.X());  b = Min(p3.X(),p4.X()); mu = Min(a,b);
218   a = Max(p1.X(),p2.X());  b = Max(p3.X(),p4.X()); Mu = Max(a,b);
219   a = Min(p1.Y(),p2.Y());  b = Min(p3.Y(),p4.Y()); mv = Min(a,b);
220   a = Max(p1.Y(),p2.Y());  b = Max(p3.Y(),p4.Y()); Mv = Max(a,b);
221   Du = Mu - mu;
222   Dv = Mv - mv;
223 }
224 //=======================================================================
225 //function : EnlargeBox and its friends.
226 //purpose  : 
227 //=======================================================================
228 static Handle(Adaptor3d_HSurface) Geometry(TopOpeBRepDS_DataStructure& DStr,
229   const Standard_Integer      ind)
230 {
231   if(ind == 0) return Handle(Adaptor3d_HSurface)();
232   if(ind > 0) {
233     TopoDS_Face F = TopoDS::Face(DStr.Shape(ind));
234     if(F.IsNull()) return Handle(Adaptor3d_HSurface)();
235     Handle(BRepAdaptor_HSurface) HS = new BRepAdaptor_HSurface();
236     HS->ChangeSurface().Initialize(F,0);
237     return HS;
238   }
239   else{
240     Handle(Geom_Surface) S  = DStr.Surface(-ind).Surface();
241     if(S.IsNull()) return Handle(Adaptor3d_HSurface)();
242     return new GeomAdaptor_HSurface(S);
243   }
244 }
245 //=======================================================================
246 //function : ChFi3d_SetPointTolerance
247 //purpose  : 
248 //=======================================================================
249 void ChFi3d_SetPointTolerance(TopOpeBRepDS_DataStructure& DStr,
250   const Bnd_Box&              box,
251   const Standard_Integer      IP)
252 {
253   Standard_Real a,b,c,d,e,f,vtol;
254   box.Get(a,b,c,d,e,f); 
255   d-=a; e-=b; f-=c; 
256   d*=d; e*=e; f*=f;
257   vtol = sqrt(d + e + f) * 1.5;// on prend un petit rab.
258   DStr.ChangePoint(IP).Tolerance(vtol);
259 }
260 //=======================================================================
261 //function : ChFi3d_EnlargeBox
262 //purpose  : 
263 //=======================================================================
264 void ChFi3d_EnlargeBox(const Handle(Geom_Curve)& C,
265   const Standard_Real       wd,
266   const Standard_Real       wf,
267   Bnd_Box&                  box1,
268   Bnd_Box&                  box2)
269 {
270   box1.Add(C->Value(wd));
271   box2.Add(C->Value(wf));
272 }
273 //=======================================================================
274 //function : ChFi3d_EnlargeBox
275 //purpose  : 
276 //=======================================================================
277 void ChFi3d_EnlargeBox(const Handle(Adaptor3d_HSurface)& S,
278   const Handle(Geom2d_Curve)&     PC,
279   const Standard_Real             wd,
280   const Standard_Real             wf,
281   Bnd_Box&                        box1,
282   Bnd_Box&                        box2)
283 {
284   Standard_Real u,v;
285   PC->Value(wd).Coord(u,v);
286   box1.Add(S->Value(u,v));
287   PC->Value(wf).Coord(u,v);
288   box2.Add(S->Value(u,v));
289 }
290 //=======================================================================
291 //function : ChFi3d_EnlargeBox
292 //purpose  : 
293 //=======================================================================
294 void ChFi3d_EnlargeBox(const TopoDS_Edge&           E,
295   const TopTools_ListOfShape&  LF,
296   const Standard_Real          w,
297   Bnd_Box&                     box)
298
299 {
300   BRepAdaptor_Curve BC(E);
301   box.Add(BC.Value(w));
302   TopTools_ListIteratorOfListOfShape It;
303   for(It.Initialize(LF); It.More(); It.Next()) {
304     TopoDS_Face F = TopoDS::Face(It.Value());
305     if(!F.IsNull()) {
306       BC.Initialize(E,F);
307       box.Add(BC.Value(w));
308     }
309   }
310 }
311 //=======================================================================
312 //function : ChFi3d_EnlargeBox
313 //purpose  : 
314 //=======================================================================
315 void ChFi3d_EnlargeBox(TopOpeBRepDS_DataStructure&    DStr,
316   const Handle(ChFiDS_Stripe)&   st, 
317   const Handle(ChFiDS_SurfData)& sd,
318   Bnd_Box&                       b1,
319   Bnd_Box&                       b2,
320   const Standard_Boolean         isfirst)
321 {
322   Standard_Real u,v;
323   const ChFiDS_CommonPoint& cp1 = sd->Vertex(isfirst,1);
324   const ChFiDS_CommonPoint& cp2 = sd->Vertex(isfirst,2);
325   b1.Add(cp1.Point());
326   b2.Add(cp2.Point());
327   const ChFiDS_FaceInterference& fi1 = sd->InterferenceOnS1();
328   const ChFiDS_FaceInterference& fi2 = sd->InterferenceOnS2();
329   const Handle(Geom_Surface)& S = DStr.Surface(sd->Surf()).Surface();
330   const Handle(Geom2d_Curve)& pcs1 = fi1.PCurveOnSurf();
331   const Handle(Geom2d_Curve)& pcs2 = fi2.PCurveOnSurf();
332   const Handle(Geom_Curve)& c3d1 = DStr.Curve(fi1.LineIndex()).Curve();
333   const Handle(Geom_Curve)& c3d2 = DStr.Curve(fi2.LineIndex()).Curve();
334   Handle(Adaptor3d_HSurface) F1 = Geometry(DStr,sd->IndexOfS1());
335   Handle(Adaptor3d_HSurface) F2 = Geometry(DStr,sd->IndexOfS2());
336   Standard_Real p1 = fi1.Parameter(isfirst);
337   if(!c3d1.IsNull()) b1.Add(c3d1->Value(p1));
338   if(!pcs1.IsNull()) {
339     pcs1->Value(p1).Coord(u,v);
340     b1.Add(S->Value(u,v));
341   }
342   if(!F1.IsNull()) {
343     const Handle(Geom2d_Curve)& pcf1 = fi1.PCurveOnFace();
344     if(!pcf1.IsNull()) {
345       pcf1->Value(p1).Coord(u,v);
346       b1.Add(F1->Value(u,v));
347     }
348   }
349   Standard_Real p2 = fi2.Parameter(isfirst);
350   if(!c3d2.IsNull()) b2.Add(c3d2->Value(p2));
351   if(!pcs2.IsNull()) {
352     pcs2->Value(p2).Coord(u,v);
353     b2.Add(S->Value(u,v));
354   }
355   if(!F2.IsNull()) {
356     const Handle(Geom2d_Curve)& pcf2 = fi2.PCurveOnFace();
357     if(!pcf2.IsNull()) {
358       pcf2->Value(p2).Coord(u,v);
359       b2.Add(F2->Value(u,v));
360     }
361   }
362   if(!st.IsNull()) {
363     const Handle(Geom_Curve)& c3d = DStr.Curve(st->Curve(isfirst)).Curve();
364     const Handle(Geom2d_Curve)& c2d = st->PCurve(isfirst);
365     if(st->Orientation(isfirst) == TopAbs_FORWARD) st->Parameters(isfirst,p1,p2);
366     else st->Parameters(isfirst,p2,p1);
367     if(!c3d.IsNull()) {
368       b1.Add(c3d->Value(p1));
369       b2.Add(c3d->Value(p2));
370     }
371     if(!c2d.IsNull()) {
372       c2d->Value(p1).Coord(u,v);
373       b1.Add(S->Value(u,v)); 
374       c2d->Value(p2).Coord(u,v);
375       b2.Add(S->Value(u,v));
376     }
377   }
378 }
379 //=======================================================================
380 //function : conexfaces
381 //purpose  : 
382 //=======================================================================
383 void ChFi3d_conexfaces(const TopoDS_Edge& E,
384   TopoDS_Face&       F1,
385   TopoDS_Face&       F2,
386   const ChFiDS_Map&  EFMap)
387 {
388   TopTools_ListIteratorOfListOfShape It;
389   F1.Nullify();
390   F2.Nullify();
391   for(It.Initialize(EFMap(E));It.More();It.Next()) {  
392     if (F1.IsNull()) {
393       F1 = TopoDS::Face(It.Value());
394     }
395     else {
396       F2 = TopoDS::Face(It.Value());
397       if(!F2.IsSame(F1) || BRep_Tool::IsClosed(E,F1)) {
398         break;
399       }
400       else F2.Nullify();
401     }
402   }  
403 }
404 //=======================================================================
405 //function : EdgeState
406 //purpose  : check concavities for the tops with 3 edges.
407 //=======================================================================
408 ChFiDS_State ChFi3d_EdgeState(TopoDS_Edge* E,
409   const ChFiDS_Map&  EFMap)
410 {
411   ChFiDS_State sst;
412   Standard_Integer i,j;
413   //TopoDS_Face F[3];
414   TopoDS_Face F1,F2,F3,F4,F5,F6;
415   ChFi3d_conexfaces(E[0],F1,F2,EFMap);
416   ChFi3d_conexfaces(E[1],F3,F4,EFMap);
417   ChFi3d_conexfaces(E[2],F5,F6,EFMap);
418
419   /*
420   if(F1.IsSame(F2)) {
421     F[0] = F[1] = F1;
422     if(F1.IsSame(F3)) F[2] = F4;
423     else F[2] = F3;
424   }
425   else if(F3.IsSame(F4)) {
426     F[0] = F[2] = F3;
427     if(F3.IsSame(F1)) F[1] = F2;
428     else F[1] = F1;
429   }
430   else if(F5.IsSame(F6)) {
431     F[1] = F[2] = F5;
432     if(F5.IsSame(F1)) F[0] = F2;
433     else F[0] = F1;
434   }
435   else{
436     if(F1.IsSame(F3) || F1.IsSame(F4)) F[0] = F1;
437     else F[0] = F2;
438     if(F3.IsSame(F[0])) F[2] = F4;
439     else F[2] = F3;
440     if(F5.IsSame(F[2])) F[1] = F6;
441     else F[1] = F5;
442
443   }
444   */
445   
446   //if(F[0].IsNull() || F[1].IsNull() || F[2].IsNull()) sst = ChFiDS_FreeBoundary;
447   if (F2.IsNull() || F4.IsNull() || F6.IsNull())
448     sst = ChFiDS_FreeBoundary;
449   else{
450     TopAbs_Orientation o01,o02,o11,o12,o21,o22;
451     /*
452     i=ChFi3d::ConcaveSide(F[0],F[1],E[0],o01,o02);
453     i=ChFi3d::ConcaveSide(F[0],F[2],E[1],o11,o12);
454     j=ChFi3d::ConcaveSide(F[1],F[2],E[2],o21,o22);
455     */
456     i=ChFi3d::ConcaveSide(F1, F2, E[0], o01, o02);
457     i=ChFi3d::ConcaveSide(F3, F4, E[1], o11, o12);
458     j=ChFi3d::ConcaveSide(F5, F6, E[2], o21, o22);
459     
460     if(o01==o11 && o02==o21 && o12==o22) sst = ChFiDS_AllSame;
461     else if(o12==o22 || i ==10 || j ==10) sst = ChFiDS_OnDiff;
462     else sst = ChFiDS_OnSame;
463   }
464   return sst;
465 }
466 //=======================================================================
467 //function : evalconti
468 //purpose  : Method very fast to code regularities CN. It is necessary to 
469 //           refine the processing.
470 //=======================================================================
471 GeomAbs_Shape ChFi3d_evalconti(const TopoDS_Edge& /*E*/,
472   const TopoDS_Face& F1,
473   const TopoDS_Face& F2)
474 {
475   GeomAbs_Shape cont = GeomAbs_G1;
476   if(!F1.IsSame(F2)) return cont;
477   TopoDS_Face F = F1;
478   F.Orientation(TopAbs_FORWARD);
479   BRepAdaptor_Surface S(F,Standard_False);
480   GeomAbs_SurfaceType typ = S.GetType();
481   if(typ != GeomAbs_Cone && 
482     typ != GeomAbs_Sphere && 
483     typ != GeomAbs_Torus) return cont;
484   return GeomAbs_CN;
485 }
486 //modified by NIZNHY-PKV Wed Dec 15 11:22:35 2010f
487 //=======================================================================
488 //function : KParticular
489 //purpose  : 
490 //=======================================================================
491 Standard_Boolean ChFi3d_KParticular (const Handle(ChFiDS_Spine)& Spine,
492   const Standard_Integer      IE,
493   const BRepAdaptor_Surface&  S1,
494   const BRepAdaptor_Surface&  S2)
495 {
496   Standard_Boolean bRet;
497   //
498   bRet=Standard_True;
499   //
500   Handle(ChFiDS_FilSpine) fs = Handle(ChFiDS_FilSpine)::DownCast(Spine);
501   if(!fs.IsNull() && !fs->IsConstant(IE)) {
502     return !bRet;
503   }
504   //
505   Standard_Boolean bIsPlane1, bIsPlane2;
506   Standard_Real aPA;
507   GeomAbs_CurveType aCT;
508   GeomAbs_SurfaceType aST1, aST2;
509   //
510   aST1=S1.GetType();
511   aST2=S2.GetType();
512   bIsPlane1=(aST1==GeomAbs_Plane);
513   bIsPlane2=(aST2==GeomAbs_Plane);
514   if (!(bIsPlane1 || bIsPlane2)) {
515     return !bRet;
516   }
517   //
518   const BRepAdaptor_Surface& aS1=(bIsPlane1)? S1 : S2;
519   const BRepAdaptor_Surface& aS2=(bIsPlane1)? S2 : S1;
520   aST1=aS1.GetType();
521   aST2=aS2.GetType();
522   //
523   if (!(aST2==GeomAbs_Plane || aST2==GeomAbs_Cylinder || aST2==GeomAbs_Cone)) {
524     return !bRet;
525   } 
526   //
527   const BRepAdaptor_Curve& bc = Spine->CurrentElementarySpine(IE);
528   aCT = bc.GetType();
529   if (!(aCT==GeomAbs_Line || aCT==GeomAbs_Circle)) {
530     return !bRet;
531   }
532   //
533   aPA=Precision::Angular();
534   //
535   if (aST2==GeomAbs_Plane){
536     if (aCT==GeomAbs_Line) { 
537       return bRet;
538     }
539   }
540   else if (aST2==GeomAbs_Cylinder) {
541     const gp_Dir aD1=aS1.Plane().Axis().Direction();
542     const gp_Dir aD2=aS2.Cylinder().Axis().Direction();
543     //
544     if (aCT==GeomAbs_Line && aD1.IsNormal(aD2, aPA)) {
545       return bRet;
546     }
547     else if (aCT==GeomAbs_Circle && aD1.IsParallel(aD2, aPA)) {
548       return bRet;
549     }
550   }
551   else if(aST2==GeomAbs_Cone) {
552     const gp_Dir aD1=aS1.Plane().Axis().Direction();
553     const gp_Dir aD2=aS2.Cone().Axis().Direction();
554     if (aCT == GeomAbs_Circle && aD1.IsParallel(aD2, aPA)) {
555       return bRet;
556     }
557   }  
558   return !bRet;    
559 }
560 //modified by NIZNHY-PKV Wed Dec 15 11:22:43 2010t
561 //=======================================================================
562 //function : BoundFac
563 //purpose  : Resize the limits of surface adjacent to the given box 
564 //           Useful for intersections with known extremities. 
565 //=======================================================================
566 void ChFi3d_BoundFac(BRepAdaptor_Surface& S,
567   const Standard_Real uumin,
568   const Standard_Real uumax,
569   const Standard_Real vvmin,
570   const Standard_Real vvmax,
571   const Standard_Boolean checknaturalbounds)
572 {
573   ChFi3d_BoundSrf(S.ChangeSurface(), uumin,uumax,vvmin,vvmax,checknaturalbounds);
574 }
575 //=======================================================================
576 //function : ChFi3d_BoundSrf
577 //purpose  : Resize the limits of surface adjacent to the given box 
578 //           Useful for intersections with known extremities.
579 //=======================================================================
580 void ChFi3d_BoundSrf(GeomAdaptor_Surface& S,
581   const Standard_Real uumin,
582   const Standard_Real uumax,
583   const Standard_Real vvmin,
584   const Standard_Real vvmax,
585   const Standard_Boolean checknaturalbounds)
586 {
587   Standard_Real umin = uumin, umax = uumax, vmin = vvmin, vmax = vvmax; 
588   Handle(Geom_Surface) surface = S.Surface();
589   Handle(Geom_RectangularTrimmedSurface) 
590     trs = Handle(Geom_RectangularTrimmedSurface)::DownCast(surface);
591   if(!trs.IsNull()) surface = trs->BasisSurface();
592   Standard_Real u1,u2,v1,v2;
593   surface->Bounds(u1,u2,v1,v2);
594   Standard_Real peru=0, perv=0;
595   if(surface->IsUPeriodic()) {
596     peru = surface->UPeriod();
597   }
598   if(surface->IsVPeriodic()) {
599     perv = surface->VPeriod();
600   }
601   Standard_Real Stepu = umax - umin;
602   Standard_Real Stepv = vmax - vmin;
603
604   //It is supposed that box uv is not null in at least 
605   //one direction.
606   Standard_Real scalu = S.UResolution(1.);
607   Standard_Real scalv = S.VResolution(1.);
608
609   Standard_Real step3du = Stepu/scalu; 
610   Standard_Real step3dv = Stepv/scalv;
611
612   if(step3du > step3dv) Stepv = step3du*scalv;
613   if(step3dv > step3du) Stepu = step3dv*scalu;
614
615   if (peru > 0) Stepu = 0.1 * (peru - (umax - umin));
616   if (perv > 0) Stepv = 0.1 * (perv - (vmax - vmin));
617
618   Standard_Real uu1 = umin - Stepu;
619   Standard_Real uu2 = umax + Stepu;
620   Standard_Real vv1 = vmin - Stepv;
621   Standard_Real vv2 = vmax + Stepv;
622   if(checknaturalbounds) {
623     if(!S.IsUPeriodic()) {uu1 = Max(uu1,u1);  uu2 = Min(uu2,u2);}
624     if(!S.IsVPeriodic()) {vv1 = Max(vv1,v1);  vv2 = Min(vv2,v2);}
625   }
626   S.Load(surface,uu1,uu2,vv1,vv2);
627 }
628 //=======================================================================
629 //function : ChFi3d_InterPlaneEdge
630 //purpose  : 
631 //=======================================================================
632 Standard_Boolean  ChFi3d_InterPlaneEdge (const Handle(Adaptor3d_HSurface)& Plan,
633   const Handle(Adaptor3d_HCurve)&  C,
634   Standard_Real& W,
635   const Standard_Boolean Sens,
636   const Standard_Real tolc)
637
638   IntCurveSurface_HInter Intersection;
639   Standard_Integer isol = 0, nbp ,iip;
640   Standard_Real uf = C->FirstParameter(),ul = C->LastParameter();
641   Standard_Real CW;
642
643   Intersection.Perform(C,Plan);
644
645   if(Intersection.IsDone()) {
646     nbp = Intersection.NbPoints();
647     for (iip = 1; iip <= nbp; iip++) {
648       CW = Intersection.Point(iip).W();
649       if(C->IsPeriodic())
650         CW = ElCLib::InPeriod(CW,uf-tolc,uf-tolc+C->Period());
651       if(uf - tolc <= CW && ul + tolc >= CW) {
652         if (isol == 0) {
653           isol = iip; W = CW;
654         }
655         else {
656           if      ( Sens && CW < W) {
657             W = CW; isol = iip;
658           }
659           else if (!Sens && CW > W) {
660             W = CW; isol = iip;
661           }
662         }
663       }
664     }
665   }
666   if(isol == 0) return Standard_False;
667   return Standard_True;
668 }
669 //=======================================================================
670 //function : ExtrSpineCarac
671 //purpose  : 
672 //=======================================================================
673 void ChFi3d_ExtrSpineCarac(const TopOpeBRepDS_DataStructure& DStr,
674   const Handle(ChFiDS_Stripe)& cd,
675   const Standard_Integer i,
676   const Standard_Real p,
677   const Standard_Integer jf,
678   const Standard_Integer sens,
679   gp_Pnt& P,
680   gp_Vec& V,
681   Standard_Real& R) //check if it is necessary to add D1,D2 and DR
682 {
683   // Attention for approximated surfaces it is assumed that e
684   // the parameters of the pcurve are the same as of  
685   // elspine used for its construction.
686   const Handle(Geom_Surface)& fffil = 
687     DStr.Surface(cd->SetOfSurfData()->Value(i)->Surf()).Surface();
688   gp_Pnt2d pp = cd->SetOfSurfData()->Value(i)->Interference(jf).
689     PCurveOnSurf()->Value(p);
690   GeomAdaptor_Surface gs(fffil);
691   P = fffil->Value(pp.X(),pp.Y());
692   gp_Pnt Pbid; gp_Vec Vbid;
693   switch (gs.GetType()) {
694   case GeomAbs_Cylinder :
695     {
696       gp_Cylinder cyl = gs.Cylinder();
697       R = cyl.Radius();
698       ElSLib::D1(pp.X(),pp.Y(),cyl,Pbid,Vbid,V);
699     }
700     break;
701   case GeomAbs_Torus :
702     {
703       gp_Torus tor = gs.Torus();
704       R = tor.MinorRadius();
705       ElSLib::D1(pp.X(),pp.Y(),tor,Pbid,V,Vbid);
706     }
707     break;
708   default:
709     { Standard_Integer nbelspine;
710     const Handle(ChFiDS_Spine)& sp = cd->Spine();
711     Handle(ChFiDS_FilSpine) fsp = Handle(ChFiDS_FilSpine)::DownCast(sp);
712     nbelspine=sp->NbEdges();
713     Handle(ChFiDS_HElSpine) hels;
714     if   (nbelspine==1) hels = sp->ElSpine(1);
715     else hels = sp->ElSpine(p);
716     if(fsp->IsConstant()) { R = fsp->Radius(); }
717     else { R = fsp->Law(hels)->Value(p); }
718     hels->D1(p,Pbid,V);
719     }
720     break;
721   }
722   V.Normalize();
723   if(sens == 1) V.Reverse();
724 }
725 //=======================================================================
726 //function : ChFi3d_CircularSpine
727 //purpose  : Calculate a cicular guideline for the corner created from  
728 //           tangent points and vectors calculated at the extremities 
729 //           of guidelines of start and end fillets.
730 //=======================================================================
731 Handle(Geom_Circle) ChFi3d_CircularSpine(Standard_Real&      WFirst,
732   Standard_Real&      WLast,
733   const gp_Pnt&       Pdeb,
734   const gp_Vec&       Vdeb,
735   const gp_Pnt&       Pfin,
736   const gp_Vec&       Vfin,
737   const Standard_Real rad)
738 {
739   gp_Circ ccc;
740   gp_Pln Pl1(Pdeb,gp_Dir(Vdeb)),Pl2(Pfin,gp_Dir(Vfin));
741   IntAna_QuadQuadGeo LInt (Pl1,Pl2,Precision::Angular(),
742     Precision::Confusion());
743   gp_Lin li;
744   if (LInt.IsDone()) {
745     li = LInt.Line(1);
746     gp_Pnt cendeb = ElCLib::Value(ElCLib::Parameter(li,Pdeb),li);
747     gp_Pnt cenfin = ElCLib::Value(ElCLib::Parameter(li,Pfin),li);
748     gp_Vec vvdeb(cendeb,Pdeb);
749     gp_Vec vvfin(cenfin,Pfin);
750     gp_Dir dddeb(vvdeb);
751     gp_Dir ddfin(vvfin);
752     if(Vdeb.Crossed(vvdeb).Dot(Vfin.Crossed(vvfin)) > 0.) {
753       return Handle(Geom_Circle)();
754     }
755     gp_Ax2 circax2(cendeb,dddeb^ddfin,dddeb);
756     ccc.SetPosition(circax2);
757     ccc.SetRadius(rad);
758     WFirst = 0.;
759     WLast = dddeb.Angle(ddfin);
760     return new Geom_Circle(ccc);
761   }
762
763   return Handle(Geom_Circle)();
764 }
765 //=======================================================================
766 //function : ChFi3d_Spine
767 //purpose  : Calculates the poles of the guideline for the corner from
768 //           tangent points and vectors calculated at the extremities of
769 //           guidelines of start and end fillets.
770 //=======================================================================
771 Handle(Geom_BezierCurve) ChFi3d_Spine(const gp_Pnt&       pd,
772   gp_Vec&             vd,
773   const gp_Pnt&       pf,
774   gp_Vec&             vf,
775   const Standard_Real R)
776 {     
777   TColgp_Array1OfPnt pol(1,4);
778   const Standard_Real fac = 0.5 * tan((M_PI-vd.Angle(vf)) * 0.5);
779   pol(1) = pd;
780   vd.Multiply(fac*R);
781   pol(2).SetCoord(pd.X()+vd.X(),pd.Y()+vd.Y(),pd.Z()+vd.Z());
782   pol(4) = pf;
783   vf.Multiply(fac*R);
784   pol(3).SetCoord(pf.X()+vf.X(),pf.Y()+vf.Y(),pf.Z()+vf.Z());
785   return new Geom_BezierCurve(pol);
786 }
787 //=======================================================================
788 //function : IsInFront
789 //purpose  : Checks if surfdata i1 and i2 are face to face
790 //=======================================================================
791 Standard_Boolean ChFi3d_IsInFront(TopOpeBRepDS_DataStructure& DStr,
792   const Handle(ChFiDS_Stripe)& cd1, 
793   const Handle(ChFiDS_Stripe)& cd2,
794   const Standard_Integer i1,
795   const Standard_Integer i2,
796   const Standard_Integer sens1,
797   const Standard_Integer sens2,
798   Standard_Real& p1,
799   Standard_Real& p2,
800   TopoDS_Face& face,
801   Standard_Boolean& sameside,
802   Standard_Integer& jf1,
803   Standard_Integer& jf2,
804   Standard_Boolean& visavis,
805   const TopoDS_Vertex& Vtx,
806   const Standard_Boolean Check2dDistance,
807   const Standard_Boolean enlarge)
808 {
809   Standard_Boolean isf1 = (sens1 == 1), isf2 = (sens2 == 1);
810   const Handle(ChFiDS_SurfData)& fd1 = cd1->SetOfSurfData()->Value(i1);
811   const Handle(ChFiDS_SurfData)& fd2 = cd2->SetOfSurfData()->Value(i2);
812
813   TopAbs_Orientation Or,OrSave1,OrSave2,OrFace1,OrFace2;
814   visavis = Standard_False;
815   Standard_Real u1 = 0.,u2 = 0.; 
816   Standard_Boolean ss = 0,ok = 0;
817   Standard_Integer j1 = 0,j2 = 0;
818   TopoDS_Face ff;
819   if(fd1->IndexOfS1() == fd2->IndexOfS1()) { 
820     jf1 = 1; jf2 = 1; 
821     face = TopoDS::Face(DStr.Shape(fd1->Index(jf1)));
822     OrSave1 = cd1->Orientation(jf1);
823     Or = OrFace1 = face.Orientation();
824     OrSave2 = cd2->Orientation(jf2);
825     OrFace2 = DStr.Shape(fd2->Index(jf2)).Orientation();
826     visavis = Standard_True;
827     sameside = ChFi3d::SameSide(Or,OrSave1,OrSave2,OrFace1,OrFace2);
828     // The parameters of the other side are not used for orientation. This would raise problems
829     Standard_Integer kf1 = jf1, kf2 = jf2;
830     Standard_Real pref1 = fd1->Interference(kf1).Parameter(isf1);
831     Standard_Real pref2 = fd2->Interference(kf2).Parameter(isf2);
832     gp_Pnt2d P2d;
833     if (Check2dDistance)
834       P2d = BRep_Tool::Parameters( Vtx, face );
835     if(ChFi3d_IntTraces(fd1,pref1,p1,jf1,sens1,fd2,pref2,p2,jf2,sens2,P2d,Check2dDistance,enlarge)) {
836       u1 = p1; u2 = p2; ss = sameside; j1 = jf1; j2 = jf2; ff = face;
837       ok = 1;
838     }
839   }
840   if(fd1->IndexOfS2() == fd2->IndexOfS1()) { 
841     jf1 = 2; jf2 = 1; 
842     face = TopoDS::Face(DStr.Shape(fd1->Index(jf1)));
843     OrSave1 = cd1->Orientation(jf1);
844     Or = OrFace1 = face.Orientation();
845     OrSave2 = cd2->Orientation(jf2);
846     OrFace2 = DStr.Shape(fd2->Index(jf2)).Orientation();
847     visavis = Standard_True;
848     sameside = ChFi3d::SameSide(Or,OrSave1,OrSave2,OrFace1,OrFace2);
849     // The parameters of the other side are not used for orientation. This would raise problems
850     Standard_Integer kf1 = jf1, kf2 = jf2;
851     Standard_Real pref1 = fd1->Interference(kf1).Parameter(isf1);
852     Standard_Real pref2 = fd2->Interference(kf2).Parameter(isf2);
853     gp_Pnt2d P2d;
854     if (Check2dDistance)
855       P2d = BRep_Tool::Parameters( Vtx, face );
856     if(ChFi3d_IntTraces(fd1,pref1,p1,jf1,sens1,fd2,pref2,p2,jf2,sens2,P2d,Check2dDistance,enlarge)) {
857       Standard_Boolean restore = 
858         ok && ((j1 == jf1 && sens1*(p1 - u1) > 0.) || 
859         (j2 == jf2 && sens2*(p2 - u2) > 0.));
860       ok = 1;
861       if(restore) {
862         p1 = u1; p2 = u2; sameside = ss; jf1 = j1; jf2 = j2; face = ff;
863       }
864       else {
865         u1 = p1; u2 = p2; ss = sameside; j1 = jf1; j2 = jf2; ff = face;
866       }
867     }
868     //the re-initialization is added in case p1,... take wrong values
869     else if (ok) {
870       p1 = u1; p2 = u2; sameside = ss; jf1 = j1; jf2 = j2; face = ff;
871     }
872   }
873   if(fd1->IndexOfS1() == fd2->IndexOfS2()) { 
874     jf1 = 1; jf2 = 2; 
875     face = TopoDS::Face(DStr.Shape(fd1->Index(jf1)));
876     OrSave1 = cd1->Orientation(jf1);
877     Or = OrFace1 = face.Orientation();
878     OrSave2 = cd2->Orientation(jf2);
879     OrFace2 = DStr.Shape(fd2->Index(jf2)).Orientation();
880     visavis = Standard_True;
881     sameside = ChFi3d::SameSide(Or,OrSave1,OrSave2,OrFace1,OrFace2);
882     // The parameters of the other side are not used for orientation.
883     Standard_Integer kf1 = jf1, kf2 = jf2;
884     Standard_Real pref1 = fd1->Interference(kf1).Parameter(isf1);
885     Standard_Real pref2 = fd2->Interference(kf2).Parameter(isf2);
886     gp_Pnt2d P2d;
887     if (Check2dDistance)
888       P2d = BRep_Tool::Parameters( Vtx, face );
889     if(ChFi3d_IntTraces(fd1,pref1,p1,jf1,sens1,fd2,pref2,p2,jf2,sens2,P2d,Check2dDistance,enlarge)) {
890       Standard_Boolean restore = 
891         ok && ((j1 == jf1 && sens1*(p1 - u1) > 0.) || 
892         (j2 == jf2 && sens2*(p2 - u2) > 0.));
893       ok = 1;
894       if(restore) {
895         p1 = u1; p2 = u2; sameside = ss; jf1 = j1; jf2 = j2; face = ff;
896       }
897       else {
898         u1 = p1; u2 = p2; ss = sameside; j1 = jf1; j2 = jf2; ff = face;
899       }
900     }
901     //the re-initialization is added in case p1,... take wrong values
902     else if (ok) {
903       p1 = u1; p2 = u2; sameside = ss; jf1 = j1; jf2 = j2; face = ff;
904     }
905   }
906   if(fd1->IndexOfS2() == fd2->IndexOfS2()) { 
907     jf1 = 2; jf2 = 2; 
908     face = TopoDS::Face(DStr.Shape(fd1->Index(jf1)));
909     OrSave1 = cd1->Orientation(jf1);
910     Or = OrFace1 = face.Orientation();
911     OrSave2 = cd2->Orientation(jf2);
912     OrFace2 = DStr.Shape(fd2->Index(jf2)).Orientation();
913     visavis = Standard_True;
914     sameside = ChFi3d::SameSide(Or,OrSave1,OrSave2,OrFace1,OrFace2);
915     // The parameters of the other side are not used for orientation.
916     Standard_Integer kf1 = jf1, kf2 = jf2;
917     Standard_Real pref1 = fd1->Interference(kf1).Parameter(isf1);
918     Standard_Real pref2 = fd2->Interference(kf2).Parameter(isf2);
919     gp_Pnt2d P2d;
920     if (Check2dDistance)
921       P2d = BRep_Tool::Parameters( Vtx, face );
922     if(ChFi3d_IntTraces(fd1,pref1,p1,jf1,sens1,fd2,pref2,p2,jf2,sens2,P2d,Check2dDistance,enlarge)) {
923       Standard_Boolean restore = 
924         ok && ((j1 == jf1 && sens1*(p1 - u1) > 0.) || 
925         (j2 == jf2 && sens2*(p2 - u2) > 0.));
926       ok = 1;
927       if(restore) {
928         p1 = u1; p2 = u2; sameside = ss; jf1 = j1; jf2 = j2; face = ff;
929       }
930       else {
931         u1 = p1; u2 = p2; ss = sameside; j1 = jf1; j2 = jf2; ff = face;
932       }
933     }
934     //the re-initialization is added in case p1,... take wrong values
935     else if (ok) {
936       p1 = u1; p2 = u2; sameside = ss; jf1 = j1; jf2 = j2; face = ff;
937     }
938   }
939   return ok;
940 }
941 //=======================================================================
942 //function : recadre
943 //purpose  : 
944 //=======================================================================
945 static Standard_Real recadre(const Standard_Real p,
946   const Standard_Real ref,
947   const Standard_Integer sens,
948   const Standard_Real first,
949   const Standard_Real last)
950 {
951   const Standard_Real pp = p + (sens > 0 ? (first - last) : (last - first));
952   return ((Abs(pp - ref) < Abs(p - ref))? pp : p);
953 }
954 //=======================================================================
955 //function : ChFi3d_IntTraces
956 //purpose  : 
957 //=======================================================================
958 Standard_Boolean ChFi3d_IntTraces(const Handle(ChFiDS_SurfData)& fd1,
959   const Standard_Real            pref1,
960   Standard_Real&                 p1,
961   const Standard_Integer         jf1,
962   const Standard_Integer         sens1,
963   const Handle(ChFiDS_SurfData)& fd2,
964   const Standard_Real            pref2,
965   Standard_Real&                 p2,
966   const Standard_Integer         jf2,
967   const Standard_Integer         sens2,
968   const gp_Pnt2d&                RefP2d,
969   const Standard_Boolean         Check2dDistance,
970   const Standard_Boolean         enlarge)
971 {
972   Geom2dAdaptor_Curve C1;
973   Geom2dAdaptor_Curve C2;
974   // pcurves are enlarged to be sure that there is intersection
975   // additionally all periodic curves are taken and points on 
976   // them are filtered using a specific criterion.
977
978   Standard_Real first,last,delta = 0.;
979   first = fd1->Interference(jf1).FirstParameter();
980   last = fd1->Interference(jf1).LastParameter();
981   if ((last-first) < Precision::PConfusion())
982     return Standard_False;
983   if(enlarge) delta = Min(0.1,0.05*(last-first));
984   Handle(Geom2d_Curve) pcf1 = fd1->Interference(jf1).PCurveOnFace();
985   if(pcf1.IsNull()) return Standard_False;
986   Standard_Boolean isper1 = pcf1->IsPeriodic();
987   if(isper1) {
988     Handle(Geom2d_TrimmedCurve) tr1 = Handle(Geom2d_TrimmedCurve)::DownCast(pcf1);
989     if(!tr1.IsNull()) pcf1 = tr1->BasisCurve();
990     C1.Load(pcf1);
991   }
992   else C1.Load(pcf1,first-delta,last+delta);
993   Standard_Real first1 = pcf1->FirstParameter(), last1 = pcf1->LastParameter();
994
995   first = fd2->Interference(jf2).FirstParameter();
996   last = fd2->Interference(jf2).LastParameter();
997   if ((last-first) < Precision::PConfusion())
998     return Standard_False;
999   if(enlarge) delta = Min(0.1,0.05*(last-first));
1000   Handle(Geom2d_Curve) pcf2 = fd2->Interference(jf2).PCurveOnFace();
1001   if(pcf2.IsNull()) return Standard_False;
1002   Standard_Boolean isper2 = pcf2->IsPeriodic();
1003   if(isper2) {
1004     Handle(Geom2d_TrimmedCurve) tr2 = Handle(Geom2d_TrimmedCurve)::DownCast(pcf2);
1005     if(!tr2.IsNull()) pcf2 = tr2->BasisCurve();
1006     C2.Load(pcf2);
1007   }
1008   else C2.Load(fd2->Interference(jf2).PCurveOnFace(),first-delta,last+delta);
1009   Standard_Real first2 = pcf2->FirstParameter(), last2 = pcf2->LastParameter();
1010
1011   IntRes2d_IntersectionPoint int2d;
1012   Geom2dInt_GInter Intersection;
1013   Standard_Integer nbpt,nbseg;
1014   gp_Pnt2d p2d;
1015   if(fd1->Interference(jf1).PCurveOnFace() == fd2->Interference(jf2).PCurveOnFace()) {
1016     Intersection.Perform(C1,
1017       Precision::PIntersection(),
1018       Precision::PIntersection());
1019   }
1020   else{
1021     Intersection.Perform(C1,C2,
1022       Precision::PIntersection(),
1023       Precision::PIntersection());
1024   }
1025   if (Intersection.IsDone()) {
1026     if (!Intersection.IsEmpty()) {
1027       nbseg = Intersection.NbSegments();
1028       if ( nbseg > 0 ) { 
1029       }
1030       nbpt = Intersection.NbPoints();
1031       if ( nbpt >= 1 ) {
1032         // The criteria sets to filter the found points in a strict way 
1033         // are missing. Two different criterions chosen somewhat randomly 
1034         // are used :
1035         // - periodic curves : closest to the border.
1036         // - non-periodic curves : the closest to the left of 2 curves
1037         //                             modulo sens1 and sens2
1038         int2d = Intersection.Point(1);
1039         p2d = int2d.Value();
1040         p1 = int2d.ParamOnFirst();
1041         p2 = int2d.ParamOnSecond();
1042         if(isper1) p1 = recadre(p1,pref1,sens1,first1,last1);
1043         if(isper2) p2 = recadre(p2,pref2,sens2,first2,last2);
1044         for(Standard_Integer i = 2; i<=nbpt; i++) {
1045           int2d = Intersection.Point(i);
1046           if(isper1) {
1047             Standard_Real pp1 = int2d.ParamOnFirst();
1048             pp1 = recadre(pp1,pref1,sens1,first1,last1);
1049             if((Abs(pp1 - pref1) < Abs(p1 - pref1))) {
1050               p1 = pp1;
1051               p2 = int2d.ParamOnSecond();
1052               p2d = int2d.Value();
1053             }
1054             //  Modified by skv - Mon Jun 16 15:51:21 2003 OCC615 Begin
1055             else if (Check2dDistance &&
1056               RefP2d.Distance(int2d.Value()) < RefP2d.Distance(p2d)) {
1057                 Standard_Real pp2 = int2d.ParamOnSecond();
1058
1059                 if(isper2)
1060                   pp2 = recadre(pp2,pref2,sens2,first2,last2);
1061
1062                 p1  = pp1;
1063                 p2  = pp2;
1064                 p2d = int2d.Value();
1065             }
1066             //  Modified by skv - Mon Jun 16 15:51:22 2003 OCC615 End
1067           }
1068           else if(isper2) {
1069             Standard_Real pp2 = int2d.ParamOnSecond();
1070             pp2 = recadre(pp2,pref2,sens2,first2,last2);
1071             if((Abs(pp2 - pref2) < Abs(p2 - pref2))) {
1072               p2 = pp2;
1073               p1 = int2d.ParamOnFirst();
1074               p2d = int2d.Value();
1075             }
1076             //  Modified by skv - Mon Jun 16 15:51:21 2003 OCC615 Begin
1077             else if (Check2dDistance &&
1078               RefP2d.Distance(int2d.Value()) < RefP2d.Distance(p2d)) {
1079                 Standard_Real pp1 = int2d.ParamOnFirst();
1080
1081                 if(isper1)
1082                   pp1 = recadre(pp1,pref1,sens1,first1,last1);
1083
1084                 p1  = pp1;
1085                 p2  = pp2;
1086                 p2d = int2d.Value();
1087             }
1088             //  Modified by skv - Mon Jun 16 15:51:22 2003 OCC615 End
1089           }
1090           else if(((int2d.ParamOnFirst() - p1)*sens1 < 0.) &&
1091             ((int2d.ParamOnSecond() - p2)*sens2 < 0.)) {
1092               p1 = int2d.ParamOnFirst();
1093               p2 = int2d.ParamOnSecond();
1094               p2d = int2d.Value();
1095           }
1096           else if((Abs(int2d.ParamOnFirst() - pref1) < Abs(p1 - pref1)) &&
1097             (Abs(int2d.ParamOnSecond() - pref2) < Abs(p2 - pref2))) {
1098               p1 = int2d.ParamOnFirst();
1099               p2 = int2d.ParamOnSecond();
1100               p2d = int2d.Value();
1101           }
1102           else if (Check2dDistance && RefP2d.Distance(int2d.Value()) < RefP2d.Distance(p2d))
1103           {
1104             p1 = int2d.ParamOnFirst();
1105             p2 = int2d.ParamOnSecond();
1106             p2d = int2d.Value();
1107           }
1108         }
1109         return Standard_True; 
1110       }
1111       return Standard_False; 
1112     }
1113     else { return Standard_False; }
1114   }
1115   else { return Standard_False; }
1116
1117 //=======================================================================
1118 //function : Coefficient
1119 //purpose  : 
1120 //=======================================================================
1121 void ChFi3d_Coefficient(const gp_Vec& V3d,
1122   const gp_Vec& D1u,
1123   const gp_Vec& D1v,
1124   Standard_Real& DU,
1125   Standard_Real& DV) 
1126 {
1127   const Standard_Real AA = D1u.SquareMagnitude();
1128   const Standard_Real BB = D1u.Dot(D1v);
1129   const Standard_Real CC = D1v.SquareMagnitude();
1130   const Standard_Real DD = D1u.Dot(V3d);
1131   const Standard_Real EE = D1v.Dot(V3d);
1132   const Standard_Real Delta = AA*CC-BB*BB;
1133   DU = (DD*CC-EE*BB)/Delta;
1134   DV = (AA*EE-BB*DD)/Delta;
1135 }
1136 //=======================================================================
1137 //function : ReparamPcurv
1138 //purpose  : Dans le cas ou la pcurve est une BSpline on verifie 
1139 //           ses parametres et on la reparametre eventuellement.
1140 //=======================================================================
1141 void ChFi3d_ReparamPcurv(const Standard_Real Uf, 
1142   const Standard_Real Ul, 
1143   Handle(Geom2d_Curve)& Pcurv) 
1144 {
1145   if(Pcurv.IsNull()) return;
1146   Standard_Real upcf = Pcurv->FirstParameter();
1147   Standard_Real upcl = Pcurv->LastParameter();
1148   Handle(Geom2d_Curve) basis = Pcurv;
1149   Handle(Geom2d_TrimmedCurve) trpc = Handle(Geom2d_TrimmedCurve)::DownCast(Pcurv);
1150   if(!trpc.IsNull()) basis = trpc->BasisCurve();
1151   Handle(Geom2d_BSplineCurve) pc = Handle(Geom2d_BSplineCurve)::DownCast(basis);
1152   if(pc.IsNull()) return;
1153   if(Abs(upcf - pc->FirstParameter()) > Precision::PConfusion() ||
1154     Abs(upcl - pc->LastParameter()) > Precision::PConfusion()) {
1155       pc->Segment(upcf,upcl);
1156   }
1157   if(Abs(Uf - pc->FirstParameter()) > Precision::PConfusion() ||
1158     Abs(Ul - pc->LastParameter()) > Precision::PConfusion()) {
1159       TColgp_Array1OfPnt2d pol(1,pc->NbPoles());
1160       pc->Poles(pol);
1161       TColStd_Array1OfReal kn(1,pc->NbKnots());
1162       pc->Knots(kn);
1163       TColStd_Array1OfInteger mu(1,pc->NbKnots());
1164       pc->Multiplicities(mu);
1165       Standard_Integer deg = pc->Degree();
1166       BSplCLib::Reparametrize(Uf,Ul,kn);
1167       pc = new Geom2d_BSplineCurve(pol,kn,mu,deg);
1168   }
1169   Pcurv = pc;
1170 }
1171 //=======================================================================
1172 //function : ProjectPCurv
1173 //purpose  : Calculation of the pcurve corresponding to a line of intersection
1174 //           3d. Should be called only in analytic cases.
1175 //=======================================================================
1176 void ChFi3d_ProjectPCurv(const Handle(Adaptor3d_HCurve)&   HCg, 
1177   const Handle(Adaptor3d_HSurface)& HSg, 
1178   Handle(Geom2d_Curve)&           Pcurv,
1179   const Standard_Real             tol,
1180   Standard_Real&                  tolreached) 
1181 {
1182   if (HSg->GetType() != GeomAbs_BezierSurface &&
1183     HSg->GetType() != GeomAbs_BSplineSurface) {
1184
1185       ProjLib_ProjectedCurve Projc (HSg,HCg,tol);
1186       tolreached = Projc.GetTolerance();
1187       switch (Projc.GetType()) {
1188       case GeomAbs_Line : 
1189         {
1190           Pcurv = new Geom2d_Line(Projc.Line());
1191         }
1192         break;
1193       case GeomAbs_Circle :
1194         {
1195           Pcurv = new Geom2d_Circle(Projc.Circle());
1196         }
1197         break;
1198       case GeomAbs_Ellipse :
1199         {
1200           Pcurv = new Geom2d_Ellipse(Projc.Ellipse());
1201         }
1202         break;
1203       case GeomAbs_Hyperbola :
1204         {
1205           Pcurv = new Geom2d_Hyperbola(Projc.Hyperbola());
1206         }
1207         break;
1208       case GeomAbs_Parabola :
1209         {
1210           Pcurv = new Geom2d_Parabola(Projc.Parabola());
1211         }
1212         break;
1213       case GeomAbs_BezierCurve :
1214         {
1215           Pcurv = Projc.Bezier(); 
1216         }
1217         break;
1218       case GeomAbs_BSplineCurve :
1219         {
1220           Pcurv = Projc.BSpline();
1221         }
1222         break;
1223       default:
1224         throw Standard_NotImplemented("echec approximation de la pcurve ");
1225       }
1226   }
1227 }
1228 //=======================================================================
1229 //function : CheckSameParameter
1230 //purpose  : Controls a posteriori that sameparameter worked well
1231 //=======================================================================
1232 Standard_Boolean ChFi3d_CheckSameParameter (const Handle(Adaptor3d_HCurve)&   C3d,
1233   Handle(Geom2d_Curve)&           Pcurv,
1234   const Handle(Adaptor3d_HSurface)& S,
1235   const Standard_Real             tol3d,
1236   Standard_Real&                  tolreached)
1237 {
1238   tolreached = 0.;
1239   Standard_Real f = C3d->FirstParameter();
1240   Standard_Real l = C3d->LastParameter();
1241   Standard_Integer nbp = 45;
1242   Standard_Real step = 1./(nbp -1);
1243   for(Standard_Integer i = 0; i < nbp; i++) {
1244     Standard_Real t,u,v;
1245     t = step * i;
1246     t = (1-t) * f + t * l;
1247     Pcurv->Value(t).Coord(u,v);
1248     gp_Pnt pS = S->Value(u,v);
1249     gp_Pnt pC = C3d->Value(t);
1250     Standard_Real d2 = pS.SquareDistance(pC);
1251     tolreached = Max(tolreached,d2);
1252   }
1253   tolreached = sqrt(tolreached);
1254   if(tolreached > tol3d) {
1255     tolreached *= 2.;
1256     return Standard_False;
1257   }
1258   tolreached *= 2.;
1259   tolreached = Max(tolreached,Precision::Confusion());
1260   return Standard_True;
1261 }
1262 //=======================================================================
1263 //function : SameParameter
1264 //purpose  : Encapsulation of Sameparameter
1265 //=======================================================================
1266 Standard_Boolean ChFi3d_SameParameter(const Handle(Adaptor3d_HCurve)&   C3d,
1267   Handle(Geom2d_Curve)&           Pcurv,
1268   const Handle(Adaptor3d_HSurface)& S,
1269   const Standard_Real             tol3d,
1270   Standard_Real&                  tolreached)
1271 {
1272   if(ChFi3d_CheckSameParameter(C3d,Pcurv,S,tol3d,tolreached)) return Standard_True;
1273   Approx_SameParameter sp(C3d,Pcurv,S,tol3d);
1274   if(sp.IsDone() && !sp.IsSameParameter()) Pcurv = sp.Curve2d();
1275   else if(!sp.IsDone() && !sp.IsSameParameter()) {
1276     return Standard_False;
1277   }
1278   tolreached = sp.TolReached();
1279   return Standard_True;
1280 }
1281 //=======================================================================
1282 //function : SameParameter
1283 //purpose  : Encapsulation de Sameparameter
1284 //=======================================================================
1285 Standard_Boolean ChFi3d_SameParameter(const Handle(Geom_Curve)&   C3d,
1286   Handle(Geom2d_Curve)&       Pcurv,
1287   const Handle(Geom_Surface)& S,
1288   const Standard_Real         Pardeb,
1289   const Standard_Real         Parfin,
1290   const Standard_Real         tol3d,
1291   Standard_Real&              tolreached)
1292 {
1293   /*szv:static*/ Handle(GeomAdaptor_HSurface) hs(new GeomAdaptor_HSurface(S));
1294   /*szv:static*/ Handle(GeomAdaptor_HCurve) hc(new GeomAdaptor_HCurve(C3d,Pardeb,Parfin));
1295   return ChFi3d_SameParameter(hc,Pcurv,hs,tol3d,tolreached);
1296 }
1297 //=======================================================================
1298 //function : ComputePCurv 
1299 //purpose  : Calculates a straight line in form of BSpline 
1300 //           to guarantee the same range and parameters as of the 
1301 //           reference 3D curve.
1302 //=======================================================================
1303 void ChFi3d_ComputePCurv(const Handle(Adaptor3d_HCurve)&   C3d,
1304   const gp_Pnt2d&                 UV1,
1305   const gp_Pnt2d&                 UV2,
1306   Handle(Geom2d_Curve)&           Pcurv,
1307   const Handle(Adaptor3d_HSurface)& S,
1308   const Standard_Real             Pardeb,
1309   const Standard_Real             Parfin,
1310   const Standard_Real             tol3d,
1311   Standard_Real&                  tolreached,
1312   const Standard_Boolean          reverse)
1313 {
1314   ChFi3d_ComputePCurv(UV1,UV2,Pcurv,Pardeb,Parfin,reverse);
1315   ChFi3d_SameParameter(C3d,Pcurv,S,tol3d,tolreached);
1316 }
1317 //=======================================================================
1318 //function : ComputePCurv 
1319 //purpose  : Calculates a straight line in form of BSpline 
1320 //           to guarantee the same range and parameters as of the 
1321 //           reference 3D curve.
1322 //=======================================================================
1323 void ChFi3d_ComputePCurv(const Handle(Geom_Curve)&   C3d,
1324   const gp_Pnt2d&             UV1,
1325   const gp_Pnt2d&             UV2,
1326   Handle(Geom2d_Curve)&       Pcurv,
1327   const Handle(Geom_Surface)& S,
1328   const Standard_Real         Pardeb,
1329   const Standard_Real         Parfin,
1330   const Standard_Real         tol3d,
1331   Standard_Real&              tolreached,
1332   const Standard_Boolean      reverse)
1333 {
1334   Handle(Adaptor3d_HSurface) hs(new GeomAdaptor_HSurface(S));
1335   Handle(Adaptor3d_HCurve) hc(new GeomAdaptor_HCurve(C3d,Pardeb,Parfin));
1336   ChFi3d_ComputePCurv(hc,UV1,UV2,Pcurv,hs,Pardeb,Parfin,tol3d,tolreached,reverse);
1337 }
1338 //=======================================================================
1339 //function : ComputePCurv 
1340 //purpose  : Calculates a straight line in form of BSpline 
1341 //           to guarantee the same range.
1342 //=======================================================================
1343 void ChFi3d_ComputePCurv(const gp_Pnt2d& UV1,
1344   const gp_Pnt2d& UV2,
1345   Handle(Geom2d_Curve)& Pcurv,
1346   const Standard_Real Pardeb,
1347   const Standard_Real Parfin,
1348   const Standard_Boolean reverse)
1349 {
1350   const Standard_Real tol = Precision::PConfusion();
1351   gp_Pnt2d p1,p2;
1352   if (!reverse) {
1353     p1 = UV1;
1354     p2 = UV2;
1355   }
1356   else {
1357     p1 = UV2;
1358     p2 = UV1;
1359   }
1360
1361   if (Abs(p1.X()-p2.X()) <= tol &&
1362     Abs((p2.Y()-p1.Y())-(Parfin-Pardeb)) <= tol) {
1363       gp_Pnt2d ppp(p1.X(),p1.Y()-Pardeb);
1364       Pcurv = new Geom2d_Line(ppp,gp::DY2d());
1365   }
1366   else if (Abs(p1.X()-p2.X()) <= tol &&
1367     Abs((p1.Y()-p2.Y())-(Parfin-Pardeb)) <= tol) {
1368       gp_Pnt2d ppp(p1.X(),p1.Y()+Pardeb);
1369       Pcurv = new Geom2d_Line(ppp,gp::DY2d().Reversed());
1370   }
1371   else if (Abs(p1.Y()-p2.Y()) <= tol &&
1372     Abs((p2.X()-p1.X())-(Parfin-Pardeb)) <= tol) {
1373       gp_Pnt2d ppp(p1.X()-Pardeb,p1.Y());
1374       Pcurv = new Geom2d_Line(ppp,gp::DX2d());
1375   }
1376   else if (Abs(p1.Y()-p2.Y()) <= tol &&
1377     Abs((p1.X()-p2.X())-(Parfin-Pardeb)) <= tol) {
1378       gp_Pnt2d ppp(p1.X()+Pardeb,p1.Y());
1379       Pcurv = new Geom2d_Line(ppp,gp::DX2d().Reversed());
1380   }
1381   else{
1382     TColgp_Array1OfPnt2d p(1,2);
1383     TColStd_Array1OfReal k(1,2);
1384     TColStd_Array1OfInteger m(1,2);
1385     m.Init(2);
1386     k(1) = Pardeb;
1387     k(2) = Parfin;
1388     p(1) = p1;
1389     p(2) = p2;
1390     Pcurv = new Geom2d_BSplineCurve(p,k,m,1);
1391   }
1392   Pcurv = new Geom2d_TrimmedCurve(Pcurv,Pardeb,Parfin);
1393 }
1394 //=======================================================================
1395 //function : ChFi3d_mkbound
1396 //purpose  : 
1397 //=======================================================================
1398 Handle(GeomFill_Boundary) ChFi3d_mkbound(const Handle(Adaptor3d_HSurface)& Fac,
1399   Handle(Geom2d_Curve)& curv, 
1400   const Standard_Integer sens1,
1401   const gp_Pnt2d& pfac1,
1402   const gp_Vec2d& vfac1,
1403   const Standard_Integer sens2,
1404   const gp_Pnt2d& pfac2,
1405   const gp_Vec2d& vfac2,
1406   const Standard_Real t3d,
1407   const Standard_Real ta)
1408 {
1409   gp_Dir2d v1(vfac1);
1410   if(sens1 == 1) v1.Reverse();
1411   gp_Dir2d v2(vfac2);
1412   if(sens2 == 1) v2.Reverse();
1413   curv = ChFi3d_BuildPCurve(Fac,pfac1,v1,pfac2,v2,Standard_False);
1414   return ChFi3d_mkbound(Fac,curv,t3d,ta);
1415 }
1416 //=======================================================================
1417 //function : ChFi3d_mkbound
1418 //purpose  : 
1419 //=======================================================================
1420 Handle(GeomFill_Boundary) ChFi3d_mkbound(const Handle(Adaptor3d_HSurface)& Surf,
1421   Handle(Geom2d_Curve)& curv,
1422   const Standard_Integer sens1,
1423   const gp_Pnt2d& p1,
1424   gp_Vec&   v1,
1425   const Standard_Integer sens2,
1426   const gp_Pnt2d& p2,
1427   gp_Vec& v2,
1428   const Standard_Real t3d,
1429   const Standard_Real ta)
1430 {
1431   if(sens1 == 1) v1.Reverse();
1432   if(sens2 == 1) v2.Reverse();
1433   curv = ChFi3d_BuildPCurve(Surf,p1,v1,p2,v2);
1434   return ChFi3d_mkbound(Surf,curv,t3d,ta);
1435 }
1436 //=======================================================================
1437 //function : ChFi3d_mkbound
1438 //purpose  : 
1439 //=======================================================================
1440 Handle(GeomFill_Boundary) ChFi3d_mkbound(const Handle(Geom_Surface)& s,
1441   const gp_Pnt2d& p1,
1442   const gp_Pnt2d& p2,
1443   const Standard_Real t3d,
1444   const Standard_Real ta,
1445   const Standard_Boolean isfreeboundary)
1446 {
1447   Handle(Adaptor3d_HSurface) HS = new GeomAdaptor_HSurface(s);
1448   return ChFi3d_mkbound(HS,p1,p2,t3d,ta,isfreeboundary);
1449 }
1450 //=======================================================================
1451 //function : ChFi3d_mkbound
1452 //purpose  : 
1453 //=======================================================================
1454 Handle(GeomFill_Boundary) ChFi3d_mkbound(const Handle(Adaptor3d_HSurface)& HS,
1455   const gp_Pnt2d& p1,
1456   const gp_Pnt2d& p2,
1457   const Standard_Real t3d,
1458   const Standard_Real ta,
1459   const Standard_Boolean isfreeboundary)
1460 {
1461   TColgp_Array1OfPnt2d pol(1,2);
1462   pol(1)=p1;
1463   pol(2)=p2;
1464   Handle(Geom2d_Curve) curv = new Geom2d_BezierCurve(pol);
1465   return ChFi3d_mkbound(HS,curv,t3d,ta,isfreeboundary);
1466 }
1467 //=======================================================================
1468 //function : ChFi3d_mkbound
1469 //purpose  : 
1470 //=======================================================================
1471 Handle(GeomFill_Boundary) ChFi3d_mkbound(const Handle(Adaptor3d_HSurface)& HS,
1472   const Handle(Geom2d_Curve)& curv,
1473   const Standard_Real t3d,
1474   const Standard_Real ta,
1475   const Standard_Boolean isfreeboundary)
1476 {
1477   Handle(Geom2dAdaptor_HCurve) HC = new Geom2dAdaptor_HCurve(curv);
1478   Adaptor3d_CurveOnSurface COnS(HC,HS);
1479   if (isfreeboundary) {
1480     Handle(Adaptor3d_HCurveOnSurface) HCOnS = new Adaptor3d_HCurveOnSurface(COnS);
1481     return new GeomFill_SimpleBound(HCOnS,t3d,ta); 
1482   }
1483   return new GeomFill_BoundWithSurf(COnS,t3d,ta);
1484 }
1485 //=======================================================================
1486 //function : ChFi3d_mkbound
1487 //purpose  : 
1488 //=======================================================================
1489 Handle(GeomFill_Boundary) ChFi3d_mkbound(const Handle(Adaptor3d_HSurface)& Fac,
1490   Handle(Geom2d_Curve)& curv, 
1491   const gp_Pnt2d& p1,
1492   const gp_Pnt2d& p2,
1493   const Standard_Real t3d,
1494   const Standard_Real ta,
1495   const Standard_Boolean isfreeboundary)
1496 {
1497   TColgp_Array1OfPnt2d pol(1,2);
1498   pol(1)=p1;
1499   pol(2)=p2;
1500   curv = new Geom2d_BezierCurve(pol);
1501   return ChFi3d_mkbound(Fac,curv,t3d,ta,isfreeboundary);
1502 }
1503 //=======================================================================
1504 //function : ChFi3d_BuildPCurve
1505 //purpose  : 
1506 //=======================================================================
1507 Handle(Geom2d_Curve) ChFi3d_BuildPCurve(const gp_Pnt2d& p1,
1508   gp_Dir2d& d1,
1509   const gp_Pnt2d& p2,
1510   gp_Dir2d& d2,
1511   const Standard_Boolean redresse)
1512 {
1513   gp_Vec2d vref(p1,p2);
1514   gp_Dir2d dref(vref);
1515   Standard_Real mref = vref.Magnitude();
1516   if(redresse) {
1517     if(d1.Dot(dref) < 0.) d1.Reverse();
1518     if(d2.Dot(dref) > 0.) d2.Reverse();
1519   }
1520   //On fait une cubique a la mords moi le noeud
1521   TColgp_Array1OfPnt2d pol(1,4);
1522   pol(1)=p1;
1523   pol(4)=p2;
1524   Standard_Real Lambda1 = Max(Abs(d2.Dot(d1)),Abs(dref.Dot(d1)));
1525   Lambda1 = Max(0.5*mref*Lambda1,1.e-5);
1526   pol(2) = gp_Pnt2d(p1.XY()+Lambda1*d1.XY());
1527   Standard_Real Lambda2 = Max(Abs(d1.Dot(d2)),Abs(dref.Dot(d2)));
1528   Lambda2 = Max(0.5*mref*Lambda2,1.e-5);
1529   pol(3)=gp_Pnt2d(p2.XY()+Lambda2*d2.XY());
1530   return new Geom2d_BezierCurve(pol);
1531 }
1532 //=======================================================================
1533 //function : ChFi3d_BuildPCurve
1534 //purpose  : 
1535 //=======================================================================
1536 Handle(Geom2d_Curve) ChFi3d_BuildPCurve(const Handle(Adaptor3d_HSurface)& Surf,
1537   const gp_Pnt2d&                 p1,
1538   const gp_Vec2d&                 v1,
1539   const gp_Pnt2d&                 p2,
1540   const gp_Vec2d&                 v2,
1541   const Standard_Boolean          redresse)
1542 {
1543   gp_Pnt2d pp1 = p1, pp2 = p2;
1544   gp_Vec2d vv1 = v1, vv2 = v2;
1545   const Standard_Real ures = Surf->UResolution(1.);
1546   const Standard_Real vres = Surf->VResolution(1.);
1547   const Standard_Real invures = 1./ures;
1548   const Standard_Real invvres = 1./vres;
1549   pp1.SetX(invures*pp1.X()); pp1.SetY(invvres*pp1.Y()); 
1550   pp2.SetX(invures*pp2.X()); pp2.SetY(invvres*pp2.Y()); 
1551   vv1.SetX(invures*vv1.X()); vv1.SetY(invvres*vv1.Y()); 
1552   vv2.SetX(invures*vv2.X()); vv2.SetY(invvres*vv2.Y()); 
1553   gp_Dir2d d1(vv1), d2(vv2);
1554   Handle(Geom2d_Curve) g2dc = ChFi3d_BuildPCurve(pp1,d1,pp2,d2,redresse);
1555   Handle(Geom2d_BezierCurve) pc = Handle(Geom2d_BezierCurve)::DownCast(g2dc);
1556   const Standard_Integer nbp = pc->NbPoles();
1557   for(Standard_Integer ip = 1; ip <= nbp; ip++) {
1558     gp_Pnt2d pol = pc->Pole(ip);
1559     pol.SetX(ures*pol.X()); pol.SetY(vres*pol.Y());
1560     pc->SetPole(ip,pol);
1561   }
1562   return pc;
1563 }
1564 //=======================================================================
1565 //function : ChFi3d_BuildPCurve
1566 //purpose  : 
1567 //=======================================================================
1568 Handle(Geom2d_Curve) ChFi3d_BuildPCurve(const Handle(Adaptor3d_HSurface)& Surf,
1569   const gp_Pnt2d&                 p1,
1570   const gp_Vec&                   v1,
1571   const gp_Pnt2d&                 p2,
1572   const gp_Vec&                   v2,
1573   const Standard_Boolean          redresse)
1574 {
1575   gp_Vec D1u,D1v;
1576   gp_Pnt PP1,PP2;
1577   Standard_Real DU,DV;
1578   Surf->D1(p1.X(),p1.Y(),PP1,D1u,D1v);
1579   ChFi3d_Coefficient(v1,D1u,D1v,DU,DV);
1580   gp_Vec2d vv1(DU,DV);
1581   Surf->D1(p2.X(),p2.Y(),PP2,D1u,D1v);
1582   ChFi3d_Coefficient(v2,D1u,D1v,DU,DV);
1583   gp_Vec2d vv2(DU,DV);
1584   gp_Vec Vref(PP1,PP2);
1585   if(redresse) {
1586     if(Vref.Dot(v1) < 0.) vv1.Reverse();
1587     if(Vref.Dot(v2) > 0.) vv2.Reverse();
1588   }
1589   return ChFi3d_BuildPCurve(Surf,p1,vv1,p2,vv2,0);
1590 }
1591 //=======================================================================
1592 //function : ComputeArete
1593 //purpose  : 
1594 // to fill with s.d. a fillet with pcurves constructed as follows
1595 // firstpoint on S1 -------------edge:curve3d/pcurves--->lastpoint on S1
1596 //  |                                                              |
1597 //  |                                                              |
1598 //  |                                                              |
1599 // edge:curve 3d/pcurves           fillet                         edge
1600 //  |         attention it is necessary to test orientation of the fillet before|
1601 //  |         determining the transitions pcurves/fillet           |
1602 //  |                                                              |
1603 //  \/                                                             \/
1604 // firstpoint sur S2 -------------edge:courbe3d/pcurves--->lastpoint sur S2
1605 //
1606 //=======================================================================
1607 void  ChFi3d_ComputeArete(const ChFiDS_CommonPoint&   P1,
1608   const gp_Pnt2d&             UV1,
1609   const ChFiDS_CommonPoint&   P2,
1610   const gp_Pnt2d&             UV2,
1611   const Handle(Geom_Surface)& Surf,
1612   Handle(Geom_Curve)&         C3d,
1613   Handle(Geom2d_Curve)&       Pcurv,
1614   Standard_Real&              Pardeb,
1615   Standard_Real&              Parfin,
1616   const Standard_Real         tol3d,
1617   const Standard_Real         tol2d,
1618   Standard_Real&              tolreached,
1619   const Standard_Integer      IFlag)
1620   // IFlag=0 pcurve et courbe 3d 
1621   // IFlag>0 pcurve (parametrage impose si IFlag=2)
1622 {
1623   /*szv:static*/ Handle(GeomAdaptor_HSurface) hs(new GeomAdaptor_HSurface());
1624   /*szv:static*/ Handle(GeomAdaptor_HCurve) hc(new GeomAdaptor_HCurve());
1625
1626   tolreached = tol3d;
1627
1628   if (Abs(UV1.X()-UV2.X()) <= tol2d) {
1629     if (IFlag == 0) {
1630       Pardeb = UV1.Y();
1631       Parfin = UV2.Y();
1632       C3d = Surf->UIso(UV1.X());
1633       if(Pardeb > Parfin) {
1634         Pardeb = C3d->ReversedParameter(Pardeb);
1635         Parfin = C3d->ReversedParameter(Parfin);
1636         C3d->Reverse();
1637       }
1638       Handle(Geom_TrimmedCurve) tc = Handle(Geom_TrimmedCurve)::DownCast(C3d);
1639       if(!tc.IsNull()) {
1640         C3d = tc->BasisCurve();
1641         if (C3d->IsPeriodic()) {
1642           ElCLib::AdjustPeriodic(C3d->FirstParameter(),C3d->LastParameter(),
1643             tol2d,Pardeb,Parfin);
1644         }
1645       }
1646     }
1647     if(IFlag != 1) {
1648       hs->ChangeSurface().Load(Surf);
1649       hc->ChangeCurve().Load(C3d,Pardeb,Parfin);
1650       const Handle(Adaptor3d_HCurve)& aHCurve = hc; // to avoid ambiguity
1651       ChFi3d_ComputePCurv(aHCurve,UV1,UV2,Pcurv,hs,Pardeb,Parfin,tol3d,tolreached,Standard_False);
1652     }
1653     else{
1654       Pcurv = new Geom2d_Line(UV1,gp_Vec2d(UV1,UV2));
1655     }
1656   }
1657   else if (Abs(UV1.Y()-UV2.Y())<=tol2d) {
1658     //iso v
1659     if (IFlag == 0) {
1660       Pardeb = UV1.X();
1661       Parfin = UV2.X();
1662       C3d = Surf->VIso(UV1.Y());
1663       if(Pardeb > Parfin) {
1664         Pardeb = C3d->ReversedParameter(Pardeb);
1665         Parfin = C3d->ReversedParameter(Parfin);
1666         C3d->Reverse();
1667       }
1668       Handle(Geom_TrimmedCurve) tc = Handle(Geom_TrimmedCurve)::DownCast(C3d);
1669       if(!tc.IsNull()) {
1670         C3d = tc->BasisCurve();
1671         if (C3d->IsPeriodic()) {
1672           ElCLib::AdjustPeriodic(C3d->FirstParameter(),C3d->LastParameter(),
1673             tol2d,Pardeb,Parfin);
1674         }
1675       }
1676     }
1677     if(IFlag != 1) {
1678       hs->ChangeSurface().Load(Surf);
1679       hc->ChangeCurve().Load(C3d,Pardeb,Parfin);
1680       const Handle(Adaptor3d_HCurve)& aHCurve = hc; // to avoid ambiguity
1681       ChFi3d_ComputePCurv(aHCurve,UV1,UV2,Pcurv,hs,Pardeb,Parfin,tol3d,tolreached,Standard_False);
1682     }
1683     else{
1684       Pcurv = new Geom2d_Line(UV1,gp_Vec2d(UV1,UV2));
1685     }
1686   }
1687   else if (IFlag == 0) {
1688
1689     if (P1.IsVertex() || P2.IsVertex() || !P1.IsOnArc() || !P2.IsOnArc()) {
1690       // A straight line is constructed to avoid  
1691       // arc and tangent.
1692       TColgp_Array1OfPnt2d qoles(1,2);
1693       qoles(1)=UV1;
1694       qoles(2)=UV2;
1695       Pcurv = new Geom2d_BezierCurve(qoles);
1696     }
1697     else {
1698       BRepAdaptor_Curve C1(P1.Arc());
1699       gp_Pnt Pp;
1700       gp_Vec Vv1;
1701       C1.D1(P1.ParameterOnArc(),Pp,Vv1);
1702       C1.Initialize(P2.Arc());
1703       gp_Vec Vv2;
1704       C1.D1(P2.ParameterOnArc(),Pp,Vv2);
1705       hs->ChangeSurface().Load(Surf);
1706       Pcurv = ChFi3d_BuildPCurve(hs,UV1,Vv1,UV2,Vv2,Standard_True); 
1707       // There are some cases when PCurve constructed in this way  
1708       // leaves the surface, in particular if it results from an 
1709       // extension. A posteriori checking is required and if
1710       // the curve leaves the surface it is replaced by straight line UV1 UV2
1711       // non regarding the tangency with neighboring arcs!
1712       Bnd_Box2d bs;
1713       Standard_Real umin,umax,vmin,vmax;
1714       Surf->Bounds(umin,umax,vmin,vmax);
1715       bs.Update(umin,vmin,umax,vmax);
1716       bs.SetGap(Precision::PConfusion());
1717       Standard_Boolean aIN = Standard_True;
1718       for(Standard_Integer ii = 1; ii <= 4 && aIN; ii++) {
1719         if(bs.IsOut(Handle(Geom2d_BezierCurve)::DownCast (Pcurv)->Pole(ii))) {
1720           aIN = Standard_False;
1721           TColgp_Array1OfPnt2d qoles(1,2);
1722           qoles(1)=UV1;
1723           qoles(2)=UV2;
1724           Pcurv = new Geom2d_BezierCurve(qoles);
1725         }
1726       }
1727     }
1728     Geom2dAdaptor_Curve AC(Pcurv);
1729     Handle(Geom2dAdaptor_HCurve) AHC = 
1730       new Geom2dAdaptor_HCurve(AC);
1731     GeomAdaptor_Surface AS(Surf);
1732     Handle(GeomAdaptor_HSurface) AHS = 
1733       new GeomAdaptor_HSurface(AS);
1734     Adaptor3d_CurveOnSurface Cs(AHC,AHS);
1735     Pardeb = Cs.FirstParameter();
1736     Parfin = Cs.LastParameter();
1737     Standard_Real avtol;
1738     GeomLib::BuildCurve3d(tol3d,Cs,Pardeb,Parfin,C3d,tolreached,avtol);
1739   }
1740   else {
1741     hs->ChangeSurface().Load(Surf);
1742     hc->ChangeCurve().Load(C3d,Pardeb,Parfin);
1743     ChFi3d_ProjectPCurv(hc,hs,Pcurv,tol3d,tolreached);
1744     gp_Pnt2d p2d = Pcurv->Value(Pardeb);
1745     if(!UV1.IsEqual(p2d,Precision::PConfusion())) {
1746       gp_Vec2d v2d(p2d,UV1);
1747       Pcurv->Translate(v2d);
1748     }
1749   }
1750 }
1751 //=======================================================================
1752 //function : FilCurveInDS
1753 //purpose  : 
1754 //=======================================================================
1755 Handle(TopOpeBRepDS_SurfaceCurveInterference)  ChFi3d_FilCurveInDS
1756   (const Standard_Integer Icurv,
1757   const Standard_Integer Isurf,
1758   const Handle(Geom2d_Curve)& Pcurv,
1759   const TopAbs_Orientation Et)
1760 {
1761   Handle(TopOpeBRepDS_SurfaceCurveInterference) SC1;
1762   SC1 = new TopOpeBRepDS_SurfaceCurveInterference(TopOpeBRepDS_Transition(Et),
1763     TopOpeBRepDS_SURFACE,
1764     Isurf,TopOpeBRepDS_CURVE,Icurv,
1765     Pcurv);
1766   return SC1;
1767 }
1768 //=======================================================================
1769 //function : TrsfTrans
1770 //purpose  : 
1771 //           
1772 //=======================================================================
1773 TopAbs_Orientation ChFi3d_TrsfTrans(const IntSurf_TypeTrans T1) 
1774 {
1775   switch (T1)  {
1776   case IntSurf_In:  return TopAbs_FORWARD;
1777   case IntSurf_Out: return TopAbs_REVERSED;
1778   default:
1779     break;
1780   }
1781   return TopAbs_INTERNAL;
1782 }
1783 //=======================================================================
1784 //function : FilCommonPoint
1785 //purpose  : Loading of the common point
1786 //           management of the case when it happens on already existing vertex.
1787 //=======================================================================
1788 Standard_EXPORT void ChFi3d_FilCommonPoint(const BRepBlend_Extremity& SP,
1789   const IntSurf_TypeTrans TransLine,
1790   const Standard_Boolean Start,
1791   ChFiDS_CommonPoint& CP,
1792   const Standard_Real Tol)
1793 {
1794   //  BRep_Tool Outil;
1795   Standard_Real Dist, maxtol = Max(Tol,CP.Tolerance());
1796
1797   CP.SetPoint(SP.Value()); // One starts with the point and the vector
1798   if (SP.HasTangent()) {
1799     if (Start) {
1800       CP.SetVector(SP.Tangent().Reversed()); // The tangent is oriented to the exit
1801     }
1802     else {
1803       CP.SetVector(SP.Tangent());
1804     }
1805   }
1806
1807   CP.SetParameter(SP.ParameterOnGuide()); // and the parameter of the spine
1808
1809   if (SP.IsVertex()) { // the Vertex is loaded if required
1810     // (inside of a face)
1811     TopoDS_Vertex V =  
1812       Handle(BRepTopAdaptor_HVertex)::DownCast(SP.Vertex())->Vertex();
1813
1814     CP.SetVertex(V);  
1815     Dist = (SP.Value()).Distance(BRep_Tool::Pnt(V));
1816     //// modified by jgv, 18.09.02 for OCC571 ////
1817     //maxtol += Dist;
1818     maxtol = Max( Dist, maxtol );
1819     //////////////////////////////////////////////
1820     CP.SetPoint(BRep_Tool::Pnt(V));
1821
1822     //the sequence of arcs the information is known by thee vertex (ancestor)
1823     //in this case the transitions are not computed, it is done by this program
1824   }
1825
1826   if (SP.NbPointOnRst() != 0) { //  An arc, and/or a vertex is loaded
1827
1828     const BRepBlend_PointOnRst& PR = SP.PointOnRst(1);
1829     Handle(BRepAdaptor_HCurve2d) 
1830       Harc = Handle(BRepAdaptor_HCurve2d)::DownCast(PR.Arc());
1831     if(!Harc.IsNull()) {
1832
1833       Standard_Real DistF, DistL, LeParamAmoi;
1834       Standard_Integer Index_min;
1835       TopoDS_Edge E = Harc->ChangeCurve2d().Edge();
1836
1837       TopoDS_Vertex V[2];
1838       TopExp::Vertices(E, V[0], V[1]);
1839
1840       DistF = (SP.Value()).Distance(BRep_Tool::Pnt(V[0]));
1841       DistL = (SP.Value()).Distance(BRep_Tool::Pnt(V[1]));
1842       if (DistF<DistL) { Index_min = 0;
1843       Dist = DistF; }
1844       else             { Index_min = 1;
1845       Dist = DistL; }
1846
1847       if (Dist <= maxtol + BRep_Tool::Tolerance(V[Index_min]) ) { 
1848         // a prexisting vertex has been met
1849         CP.SetVertex(V[Index_min]); //the old vertex is loaded
1850         CP.SetPoint( BRep_Tool::Pnt(V[Index_min]) );
1851         maxtol = Max(BRep_Tool::Tolerance(V[Index_min]),maxtol);
1852         //// modified by jgv, 18.09.02 for OCC571 ////
1853         //maxtol += Dist;
1854         maxtol = Max( Dist, maxtol );
1855         //////////////////////////////////////////////
1856         LeParamAmoi = BRep_Tool::Parameter(V[Index_min], E);    
1857       }
1858       else {   // Creation of an arc only
1859         maxtol = Max(BRep_Tool::Tolerance(E),maxtol);
1860         maxtol = Max(SP.Tolerance(),maxtol);
1861         LeParamAmoi = PR.ParameterOnArc();
1862       }
1863
1864       // Definition of the arc
1865       TopAbs_Orientation Tr;
1866       TopAbs_Orientation Or = E.Orientation();
1867       if (Start) {
1868         Tr = TopAbs::Reverse(TopAbs::Compose(ChFi3d_TrsfTrans(TransLine),Or));
1869       }
1870       else {
1871         Tr = TopAbs::Compose(ChFi3d_TrsfTrans(TransLine),Or);
1872       }
1873       CP.SetArc(maxtol, E, LeParamAmoi, Tr);
1874     }
1875   }
1876   CP.SetTolerance(maxtol); // Finally, the tolerance.
1877 }
1878
1879 //=======================================================================
1880 //function : SolidIndex
1881 //purpose  : 
1882 //=======================================================================
1883 Standard_Integer ChFi3d_SolidIndex(const Handle(ChFiDS_Spine)&  sp,
1884   TopOpeBRepDS_DataStructure&  DStr,
1885   ChFiDS_Map&                  MapESo,
1886   ChFiDS_Map&                  MapESh)
1887 {
1888   if(sp.IsNull() || sp->NbEdges() == 0) 
1889     throw Standard_Failure("SolidIndex : Spine incomplete");
1890   TopoDS_Shape edref= sp->Edges(1);
1891   TopoDS_Shape shellousolid;
1892   if(!MapESo(edref).IsEmpty()) shellousolid = MapESo(edref).First();
1893   else shellousolid = MapESh(edref).First();
1894   const Standard_Integer solidindex = DStr.AddShape(shellousolid);
1895   return solidindex;
1896 }
1897 //=======================================================================
1898 //function : IndexPointInDS
1899 //purpose  : 
1900 //=======================================================================
1901 Standard_Integer  ChFi3d_IndexPointInDS(const ChFiDS_CommonPoint& P1,
1902   TopOpeBRepDS_DataStructure& DStr) 
1903 {
1904   if (P1.IsVertex()) {
1905     // --------------------------------->  !*!*!* 
1906     // Attention : it is necessary ti implement a mechanism 
1907     // controlling tolerance.
1908     BRep_Builder B;
1909     B.UpdateVertex(P1.Vertex(), P1.Point(), P1.Tolerance());
1910     return DStr.AddShape(P1.Vertex());
1911   }
1912   return DStr.AddPoint(TopOpeBRepDS_Point(P1.Point(),P1.Tolerance()));
1913 }
1914 //=======================================================================
1915 //function : FilPointInDS
1916 //purpose  : 
1917 //=======================================================================
1918 Handle(TopOpeBRepDS_CurvePointInterference) 
1919   ChFi3d_FilPointInDS(const TopAbs_Orientation Et,
1920   const Standard_Integer Ic,
1921   const Standard_Integer Ip,
1922   const Standard_Real Par,
1923   const Standard_Boolean IsVertex)
1924 {
1925   Handle(TopOpeBRepDS_CurvePointInterference) CP1;
1926   if (IsVertex)    
1927     CP1 = new TopOpeBRepDS_CurvePointInterference (TopOpeBRepDS_Transition(Et),
1928     TopOpeBRepDS_CURVE,Ic,
1929     TopOpeBRepDS_VERTEX,Ip,Par); 
1930   else
1931     CP1 = new TopOpeBRepDS_CurvePointInterference (TopOpeBRepDS_Transition(Et),
1932     TopOpeBRepDS_CURVE,Ic,
1933     TopOpeBRepDS_POINT,Ip,Par);
1934   return CP1;
1935 }
1936 //=======================================================================
1937 //function : FilVertexInDS
1938 //purpose  : 
1939 //=======================================================================
1940 Handle(TopOpeBRepDS_CurvePointInterference) 
1941   ChFi3d_FilVertexInDS(const TopAbs_Orientation Et,
1942   const Standard_Integer Ic,
1943   const Standard_Integer Ip,
1944   const Standard_Real Par)
1945 {
1946
1947   Handle(TopOpeBRepDS_CurvePointInterference) CP1 = new
1948     TopOpeBRepDS_CurvePointInterference (TopOpeBRepDS_Transition(Et),
1949     TopOpeBRepDS_CURVE,Ic,
1950     TopOpeBRepDS_VERTEX,Ip,Par);
1951   return CP1;
1952 }
1953 //=======================================================================
1954 //function : Orientation
1955 //purpose  : returns the orientation of the interference (the first found
1956 //           in the list).
1957 //=======================================================================
1958
1959 static Standard_Boolean
1960   ChFi3d_Orientation(const TopOpeBRepDS_ListOfInterference& LI,
1961   const Standard_Integer                 igros,
1962   const Standard_Integer                 ipetit,
1963   TopAbs_Orientation&                    Or,
1964   const Standard_Boolean                 isvertex = Standard_False,
1965   const Standard_Boolean                 aprendre = Standard_False)
1966 {
1967   //In case, when it is necessary to insert a point/vertex, it should be 
1968   //known if this is a point or a vertex, because their index can be the same.
1969   TopOpeBRepDS_Kind typepetit;
1970   if (isvertex)
1971     typepetit =  TopOpeBRepDS_VERTEX;
1972   else
1973     typepetit =  TopOpeBRepDS_POINT;
1974   TopOpeBRepDS_ListIteratorOfListOfInterference itLI(LI);
1975   for (; itLI.More(); itLI.Next() ) {
1976     const Handle(TopOpeBRepDS_Interference)& cur = itLI.Value(); 
1977     TopOpeBRepDS_Kind GK;
1978     TopOpeBRepDS_Kind SK;
1979     Standard_Integer S;
1980     Standard_Integer G;
1981     cur->GKGSKS(GK,G,SK,S);
1982     if (aprendre) {
1983       if ( S == igros && G == ipetit && GK == typepetit) {
1984         Or = cur->Transition().Orientation(TopAbs_IN);
1985         return Standard_True;
1986       }
1987     }
1988     else  {
1989       if ( S == igros && G == ipetit) {
1990         Or = cur->Transition().Orientation(TopAbs_IN);
1991         return Standard_True;
1992       } 
1993     }
1994   }
1995   return Standard_False;
1996 }
1997
1998 //=======================================================================
1999 //function : Contains
2000 //purpose  : Check if the interference does not already exist.
2001 //====================================================================
2002
2003 static Standard_Boolean ChFi3d_Contains
2004   (const TopOpeBRepDS_ListOfInterference& LI,
2005   const Standard_Integer                 igros,
2006   const Standard_Integer                 ipetit,
2007   const Standard_Boolean                 isvertex = Standard_False,
2008   const Standard_Boolean                 aprendre = Standard_False)                             
2009 {
2010   TopAbs_Orientation bidOr;
2011   return ChFi3d_Orientation(LI,igros,ipetit,bidOr,isvertex,aprendre);
2012 }
2013 //=======================================================================
2014 //function : QueryAddVertexInEdge
2015 //purpose  : 
2016 //=======================================================================
2017 static void QueryAddVertexInEdge(TopOpeBRepDS_ListOfInterference& LI,
2018   const Standard_Integer                 IC,
2019   const Standard_Integer                 IV,
2020   const Standard_Real                    par,
2021   const TopAbs_Orientation               Or)
2022 {
2023   TopOpeBRepDS_ListIteratorOfListOfInterference it(LI);
2024   for (; it.More(); it.Next() ) {
2025     const Handle(TopOpeBRepDS_Interference)& cur = it.Value();
2026     Handle(TopOpeBRepDS_CurvePointInterference) cpi (Handle(TopOpeBRepDS_CurvePointInterference)::DownCast(cur));
2027     if(!cpi.IsNull()) {
2028       Standard_Integer newIV = cpi->Geometry();
2029       TopOpeBRepDS_Kind kv = cpi->GeometryType();
2030       TopAbs_Orientation newOr = cpi->Transition().Orientation(TopAbs_IN);
2031       Standard_Real newpar = cpi->Parameter();
2032       if(IV == newIV && kv == TopOpeBRepDS_VERTEX && 
2033         Or == newOr && Abs(par - newpar) < 1.e-10) {
2034           return;
2035       }
2036     }
2037   }
2038   Handle(TopOpeBRepDS_CurvePointInterference) interf = 
2039     ChFi3d_FilVertexInDS(Or,IC,IV,par);
2040   LI.Append(interf);
2041 }
2042
2043 //=======================================================================
2044 //function : CutEdge
2045 //purpose  : 
2046 //=======================================================================
2047 static void CutEdge(const TopoDS_Vertex&           V,
2048   const Handle(ChFiDS_SurfData)& SD,
2049   TopOpeBRepDS_DataStructure&    DStr,
2050   const Standard_Boolean         ,
2051   const Standard_Integer         ons)
2052 {
2053   if(!SD->IsOnCurve(ons)) return;
2054   Standard_Integer IC = SD->IndexOfC(ons);
2055   Standard_Integer IV = DStr.AddShape(V);
2056   TopOpeBRepDS_ListOfInterference& LI = DStr.ChangeShapeInterferences(IC);
2057   TopoDS_Edge E = TopoDS::Edge(DStr.Shape(IC));
2058   E.Orientation(TopAbs_FORWARD);
2059   TopExp_Explorer ex;
2060
2061   // process them checking that it has not been done already.
2062   for(ex.Init(E,TopAbs_VERTEX);ex.More();ex.Next()) {
2063     const TopoDS_Vertex& vv = TopoDS::Vertex(ex.Current());
2064     if(vv.IsSame(V)) {
2065       TopAbs_Orientation Or = TopAbs::Reverse(vv.Orientation());
2066       Standard_Real par = BRep_Tool::Parameter(vv,E);
2067       QueryAddVertexInEdge(LI,IC,IV,par,Or);
2068     }
2069   }
2070 }
2071 //=======================================================================
2072 //function : findIndexPoint
2073 //purpose  : returns in <ipon> index of point bounding a courve interfering
2074 //           with <Fd> and coinciding with last common point on <OnS> face
2075 //=======================================================================
2076 static Standard_Boolean 
2077   findIndexPoint(const TopOpeBRepDS_DataStructure& DStr,
2078   const Handle(ChFiDS_SurfData)&    Fd,
2079   const Standard_Integer            OnS,
2080   Standard_Integer&                 ipoin)
2081 {
2082   ipoin = 0;
2083   gp_Pnt P = Fd->Vertex(Standard_False,OnS).Point();
2084
2085   TopOpeBRepDS_ListIteratorOfListOfInterference SCIIt, CPIIt;
2086
2087   SCIIt.Initialize (DStr.SurfaceInterferences(Fd->Surf()));
2088   for (; SCIIt.More(); SCIIt.Next()) {
2089     Handle(TopOpeBRepDS_SurfaceCurveInterference) SCI =
2090       Handle(TopOpeBRepDS_SurfaceCurveInterference)::DownCast(SCIIt.Value());
2091     if (SCI.IsNull()) continue;
2092     CPIIt.Initialize (DStr.CurveInterferences(SCI->Geometry()));
2093     for (; CPIIt.More(); CPIIt.Next()) {
2094       Handle(TopOpeBRepDS_CurvePointInterference) CPI =
2095         Handle(TopOpeBRepDS_CurvePointInterference)::DownCast(CPIIt.Value());
2096       if (CPI.IsNull()) continue;
2097       Standard_Integer iPoint = CPI->Geometry();
2098       TopOpeBRepDS_Point tp = DStr.Point(iPoint);
2099       if (P.IsEqual(tp.Point(), tp.Tolerance())) {
2100         ipoin = iPoint;
2101         return Standard_True;
2102       }
2103     }
2104   }
2105   return Standard_False;
2106 }
2107 //=======================================================================
2108 //function : FilDS
2109 //purpose  : 
2110 //=======================================================================
2111 void  ChFi3d_FilDS(const Standard_Integer       SolidIndex,
2112   const Handle(ChFiDS_Stripe)& CorDat,
2113   TopOpeBRepDS_DataStructure&  DStr,
2114   ChFiDS_Regularities&         reglist,
2115   const Standard_Real          tol3d,
2116   const Standard_Real          tol2d) 
2117 {
2118   //  BRep_Tool Outil;
2119   TopExp_Explorer ex;
2120   Handle(ChFiDS_Spine) spine = CorDat->Spine();
2121   Standard_Boolean Closed = Standard_False;
2122   Standard_Boolean Degene = 0, isVertex1 = 0, isVertex2 = 0, Singulier_en_Bout = 0;
2123   if(!spine.IsNull()) {
2124     Closed = spine->IsPeriodic();
2125   }
2126   const ChFiDS_SequenceOfSurfData& SeqFil = 
2127     CorDat->SetOfSurfData()->Sequence();
2128   Standard_Integer Ipoin1 = CorDat->IndexFirstPointOnS1();
2129   Standard_Integer Ipoin2 = CorDat->IndexFirstPointOnS2();
2130   Standard_Integer NumEdge = 1;
2131   TopoDS_Vertex BoutdeVtx; 
2132   Standard_Integer Icurv = 0;
2133   Standard_Integer Iarc1 = 0,Iarc2 = 0;
2134   TopAbs_Orientation trafil1 = TopAbs_FORWARD, trafil2 = TopAbs_FORWARD;
2135   Standard_Integer IcFil1,IcFil2,Isurf,Ishape1,Ishape2;
2136   Standard_Real Pardeb = 0.,Parfin = 0.;
2137   TopAbs_Orientation ET1;
2138   Handle(TopOpeBRepDS_CurvePointInterference) Interfp1,Interfp2;
2139   Handle(TopOpeBRepDS_SurfaceCurveInterference) Interfc1,Interfc2;
2140   Handle(TopOpeBRepDS_SurfaceCurveInterference) Interfc3,Interfc4;
2141   Handle(TopOpeBRepDS_CurvePointInterference) Interfp3,Interfp4;
2142   Handle(TopOpeBRepDS_CurvePointInterference) Interfp5,Interfp6;
2143   TopoDS_Face F;
2144   Handle(Geom2d_Curve) PCurv;
2145   TopOpeBRepDS_Curve Crv;
2146
2147   TopOpeBRepDS_ListOfInterference& SolidInterfs = 
2148     DStr.ChangeShapeInterferences(SolidIndex);
2149
2150   ChFiDS_Regul regcout; // for closed and tangent CD
2151   ChFiDS_Regul regfilfil; // for connections Surf/Surf
2152
2153   ChFiDS_CommonPoint V3;
2154   ChFiDS_CommonPoint V4;
2155
2156   // Nullify degenerated ChFi/Faces interferences, eap occ293
2157   Standard_Integer j;
2158   if (SeqFil.Length() > 1) {
2159     for (j=1; j<=SeqFil.Length(); j++) {
2160       Handle(ChFiDS_SurfData) Fd = SeqFil(j);
2161       Standard_Integer onS;
2162       for (onS=1; onS<=2; onS++) {
2163         const ChFiDS_FaceInterference& Fi = Fd->Interference(onS);
2164         IcFil1 = Fi.LineIndex();
2165         if (!IcFil1) continue;
2166         Standard_Real FiLen = Abs(Fi.FirstParameter()-Fi.LastParameter());
2167         if (FiLen > Precision::PConfusion()) continue;
2168         TopOpeBRepDS_Curve& cc = DStr.ChangeCurve(IcFil1);
2169         cc.ChangeCurve().Nullify();
2170
2171         // care of CommonPoint, eap occ354
2172         if (j!=1 && j!=SeqFil.Length()) continue;
2173         Standard_Boolean isfirst = (j==1);
2174         Standard_Integer i = isfirst ? j+1 : j-1;
2175         ChFiDS_CommonPoint& CP1 = SeqFil(i)->ChangeVertex(isfirst,onS);
2176         if (Fd->Vertex(isfirst,onS).IsOnArc() && CP1.IsOnArc()) {
2177           ChFiDS_CommonPoint& CP2 = Fd->ChangeVertex(!isfirst,onS);
2178           CP1.Reset();
2179           CP1.SetPoint(CP2.Point());
2180           CP2.Reset();
2181           CP2.SetPoint(CP1.Point());
2182         }
2183       }
2184     }
2185   }
2186
2187   for (j=1; j<=SeqFil.Length(); j++) {
2188
2189     const Handle(ChFiDS_SurfData)& Fd = SeqFil(j);
2190     Isurf= Fd->Surf();
2191     Ishape1 = Fd->IndexOfS1();
2192     Ishape2 = Fd->IndexOfS2();
2193
2194     // eap, Apr 29 2002, occ 293
2195     // now IsInDS() returns nb of surfaces at end being in DS;
2196     // vars showing which end is in DS
2197     Standard_Boolean isInDS1 = Standard_False, isInDS2 = Standard_False;
2198     if (j <= CorDat->IsInDS(Standard_True)) {
2199       isInDS1 = Standard_True;
2200       isInDS2 = (j+1 <= CorDat->IsInDS(Standard_True));
2201     }
2202     if (SeqFil.Length()-j < CorDat->IsInDS(Standard_False)) {
2203       isInDS2 = Standard_True;
2204       isInDS1 = isInDS1 || SeqFil.Length()-j+1 < CorDat->IsInDS(Standard_False);
2205     }
2206
2207     // creation of SolidSurfaceInterference
2208
2209     Handle(TopOpeBRepDS_SolidSurfaceInterference) 
2210       SSI = new TopOpeBRepDS_SolidSurfaceInterference
2211       (TopOpeBRepDS_Transition(Fd->Orientation()),
2212       TopOpeBRepDS_SOLID,
2213       SolidIndex,
2214       TopOpeBRepDS_SURFACE,
2215       Isurf);
2216
2217     SolidInterfs.Append(SSI);
2218
2219     const ChFiDS_FaceInterference& Fi1 = Fd->InterferenceOnS1();
2220     const ChFiDS_FaceInterference& Fi2 = Fd->InterferenceOnS2();     
2221     const ChFiDS_CommonPoint& V1 = Fd->VertexFirstOnS1();
2222     const ChFiDS_CommonPoint& V2 = Fd->VertexFirstOnS2();
2223
2224     // Processing to manage double interferences
2225     if (j>1) {
2226       if (V1.IsOnArc() && V3.IsOnArc() && V1.Arc().IsSame(V3.Arc())) {
2227         //Iarc1 is initialized
2228         //Iarc1 = DStr.AddShape(V1.Arc());
2229         if (ChFi3d_Contains(DStr.ShapeInterferences(Iarc1),Iarc1,Ipoin1) && 
2230           (V1.TransitionOnArc() != V3.TransitionOnArc()) ) {
2231             Interfp1= ChFi3d_FilPointInDS(V1.TransitionOnArc(),Iarc1,Ipoin1,
2232               V1.ParameterOnArc());
2233             DStr.ChangeShapeInterferences(V1.Arc()).Append(Interfp1);
2234         }
2235       }
2236
2237       if (V2.IsOnArc() && V4.IsOnArc() && V2.Arc().IsSame(V4.Arc())) {
2238         //Iarc2 is initialized
2239         //Iarc2 = DStr.AddShape(V2.Arc());
2240         if ( ChFi3d_Contains(DStr.ShapeInterferences(Iarc2),Iarc2,Ipoin2)  && 
2241           (V2.TransitionOnArc() != V4.TransitionOnArc()) ) {
2242             Interfp2= ChFi3d_FilPointInDS(V2.TransitionOnArc(),Iarc2,Ipoin2,
2243               V2.ParameterOnArc());
2244             DStr.ChangeShapeInterferences(V2.Arc()).Append(Interfp2);
2245         }
2246       }
2247     }
2248
2249     V3 = Fd->VertexLastOnS1();
2250     V4 = Fd->VertexLastOnS2();
2251
2252     if(Ishape1 != 0) {
2253       if(Ishape1 > 0) {
2254         trafil1 = DStr.Shape(Ishape1).Orientation();
2255       }
2256       else{
2257         ChFi3d_Orientation(SolidInterfs,SolidIndex,-Ishape1,trafil1);
2258       }
2259       trafil1 = TopAbs::Compose(trafil1,Fd->Orientation());
2260       trafil1 = TopAbs::Compose(TopAbs::Reverse(Fi1.Transition()),trafil1);
2261       trafil2 = TopAbs::Reverse(trafil1);
2262     }
2263     else{
2264       if(Ishape2 > 0) {
2265         trafil2 = DStr.Shape(Ishape2).Orientation();
2266       }
2267       else{
2268         ChFi3d_Orientation(SolidInterfs,SolidIndex,-Ishape2,trafil2);
2269       }
2270       trafil2 = TopAbs::Compose(trafil2,Fd->Orientation());
2271       trafil2 = TopAbs::Compose(TopAbs::Reverse(Fi2.Transition()),trafil2);
2272       trafil1 = TopAbs::Reverse(trafil2);
2273     }
2274
2275     ET1 = TopAbs::Reverse(trafil1);
2276
2277     // A small paragraph to process contacts of edges, which touch 
2278     // a vertex of the obstacle.
2279     if(V1.IsVertex() && Fd->IsOnCurve1()) {
2280       const TopoDS_Vertex& vv1 = V1.Vertex();
2281       CutEdge(vv1,Fd,DStr,1,1);
2282     }
2283     if(V2.IsVertex() && Fd->IsOnCurve2()) {
2284       const TopoDS_Vertex& vv2 = V2.Vertex();
2285       CutEdge(vv2,Fd,DStr,1,2);
2286     }
2287     if(V3.IsVertex() && Fd->IsOnCurve1()) {
2288       const TopoDS_Vertex& vv3 = V3.Vertex();
2289       CutEdge(vv3,Fd,DStr,0,1);
2290     }
2291     if(V4.IsVertex() && Fd->IsOnCurve2()) {
2292       const TopoDS_Vertex& vv4 = V4.Vertex();
2293       CutEdge(vv4,Fd,DStr,0,2);
2294     }
2295
2296     if (j == 1) {
2297       isVertex1 = V1.IsVertex();
2298       isVertex2 = V2.IsVertex();
2299       Singulier_en_Bout =  (V1.Point().IsEqual(V2.Point(), 0));
2300
2301       if (Singulier_en_Bout) {
2302         // Queue de Billard
2303         if ((!V1.IsVertex()) || (!V2.IsVertex())) {
2304
2305         }
2306         else {
2307           isVertex1 = isVertex2 = Standard_True; //caution...
2308           // The edge is removed from spine starting on this vertex.
2309           TopoDS_Edge Arcspine = spine->Edges(1);
2310           BoutdeVtx = V1.Vertex();
2311           Standard_Integer IArcspine = DStr.AddShape(Arcspine);
2312           Standard_Integer IVtx = CorDat->IndexFirstPointOnS1();
2313
2314           TopAbs_Orientation OVtx = TopAbs_FORWARD;;
2315
2316           for(ex.Init(Arcspine.Oriented(TopAbs_FORWARD),TopAbs_VERTEX); 
2317             ex.More(); ex.Next()) {
2318               if(BoutdeVtx.IsSame(ex.Current())) {
2319                 OVtx = ex.Current().Orientation();
2320                 break;
2321               }
2322           }
2323           OVtx = TopAbs::Reverse(OVtx);
2324           Standard_Real parVtx = BRep_Tool::Parameter(BoutdeVtx,Arcspine);
2325           Handle(TopOpeBRepDS_CurvePointInterference) 
2326             interfv = ChFi3d_FilVertexInDS(OVtx,IArcspine,IVtx,parVtx);
2327           DStr.ChangeShapeInterferences(IArcspine).Append(interfv);
2328         }
2329       }
2330       else {
2331         if (V1.IsOnArc()) {
2332           Iarc1 = DStr.AddShape(V1.Arc()); 
2333           if ( !ChFi3d_Contains(DStr.ShapeInterferences(Iarc1),Iarc1,Ipoin1) ) {
2334             Interfp1= ChFi3d_FilPointInDS(V1.TransitionOnArc(),Iarc1,Ipoin1,
2335               V1.ParameterOnArc(), isVertex1);
2336             DStr.ChangeShapeInterferences(V1.Arc()).Append(Interfp1);
2337           }
2338         }
2339
2340         if (V2.IsOnArc()) {
2341           Iarc2 = DStr.AddShape(V2.Arc());
2342           if ( !ChFi3d_Contains(DStr.ShapeInterferences(Iarc2),Iarc2,Ipoin2) ) {
2343             Interfp2= ChFi3d_FilPointInDS(V2.TransitionOnArc(),Iarc2,Ipoin2,
2344               V2.ParameterOnArc(),isVertex2);
2345             DStr.ChangeShapeInterferences(V2.Arc()).Append(Interfp2);
2346           }
2347         }
2348       }
2349
2350       if (!isInDS1) {
2351         ET1 = TopAbs::Compose(ET1,CorDat->FirstPCurveOrientation());
2352         Icurv = CorDat->FirstCurve();
2353         if(Closed && !Singulier_en_Bout) {
2354           regcout.SetCurve(Icurv);
2355           regcout.SetS1(Isurf,Standard_False);
2356         }
2357         PCurv = CorDat->FirstPCurve();
2358         CorDat->FirstParameters(Pardeb,Parfin);
2359
2360         TopOpeBRepDS_ListOfInterference& Li = DStr.ChangeCurveInterferences(Icurv);
2361         if (Li.IsEmpty()) {
2362           if(CorDat->FirstPCurveOrientation()==TopAbs_REVERSED) {
2363             Interfp1=ChFi3d_FilPointInDS
2364               (TopAbs_REVERSED,Icurv,Ipoin1,Parfin,isVertex1);
2365             Interfp2=ChFi3d_FilPointInDS
2366               (TopAbs_FORWARD,Icurv,Ipoin2,Pardeb,isVertex2);
2367           }
2368           else{
2369             Interfp1=ChFi3d_FilPointInDS
2370               (TopAbs_FORWARD,Icurv,Ipoin1,Pardeb,isVertex1);
2371             Interfp2=ChFi3d_FilPointInDS
2372               (TopAbs_REVERSED,Icurv,Ipoin2,Parfin,isVertex2);
2373           }
2374           Li.Append(Interfp1);
2375           Li.Append(Interfp2);
2376         }
2377         Interfc1= ChFi3d_FilCurveInDS (Icurv,Isurf,PCurv,ET1);
2378         DStr.ChangeSurfaceInterferences(Isurf).Append(Interfc1);
2379         if (Ipoin1 == Ipoin2) {
2380           TopOpeBRepDS_Curve& TCurv = DStr.ChangeCurve(Icurv);
2381           TCurv.ChangeCurve().Nullify();
2382           Handle(TopOpeBRepDS_Interference) bidinterf;
2383           TCurv.SetSCI(Interfc1,bidinterf);         
2384         }
2385       }
2386     } // End of the Initial Processing (j==1)
2387     else {
2388       // ---- Interference between Fillets ------
2389
2390       if (!isInDS1) {// eap, Apr 29 2002, occ 293 
2391
2392         if (Degene && isVertex1) {
2393           // The edge is removed from the spine starting on this vertex.
2394           NumEdge++; // The previous edge of the vertex has already been found.
2395           TopoDS_Edge Arcspine = spine->Edges(NumEdge);
2396           Standard_Integer IArcspine = DStr.AddShape(Arcspine);
2397           Standard_Integer IVtx = DStr.AddShape(BoutdeVtx);
2398           TopAbs_Orientation OVtx = TopAbs_FORWARD;
2399           for(ex.Init(Arcspine.Oriented(TopAbs_FORWARD),TopAbs_VERTEX); 
2400             ex.More(); ex.Next()) {
2401               if(BoutdeVtx.IsSame(ex.Current())) {
2402                 OVtx = ex.Current().Orientation();
2403                 break;
2404               }
2405           }
2406           OVtx = TopAbs::Reverse(OVtx);
2407           Standard_Real parVtx = BRep_Tool::Parameter(BoutdeVtx,Arcspine);
2408           Handle(TopOpeBRepDS_CurvePointInterference) 
2409             interfv = ChFi3d_FilVertexInDS(OVtx,IArcspine,IVtx,parVtx);
2410           DStr.ChangeShapeInterferences(IArcspine).Append(interfv);
2411         } // End of the removal
2412
2413         gp_Pnt2d UV1 = Fd->InterferenceOnS1().PCurveOnSurf()->
2414           Value(Fd->InterferenceOnS1().FirstParameter());
2415         gp_Pnt2d UV2 = Fd->InterferenceOnS2().PCurveOnSurf()-> 
2416           Value(Fd->InterferenceOnS2().FirstParameter());
2417         TopOpeBRepDS_Curve& TCurv = DStr.ChangeCurve(Icurv);
2418         if (Degene) {
2419           // pcurve is associated via SCI to TopOpeBRepDSCurve.
2420           ChFi3d_ComputePCurv(UV1,UV2,PCurv,Pardeb,Parfin);       
2421           Interfc1= ChFi3d_FilCurveInDS (Icurv,Isurf,PCurv,ET1);
2422           DStr.ChangeSurfaceInterferences(Isurf).Append(Interfc1);
2423           TCurv.ChangeCurve().Nullify();
2424           Handle(TopOpeBRepDS_Interference) bidinterf;
2425           TCurv.SetSCI(Interfc1,bidinterf);           
2426         }
2427         else {
2428           regfilfil.SetS2(Isurf,Standard_False);
2429           reglist.Append(regfilfil);
2430           Standard_Real tolreached;
2431           ChFi3d_ComputePCurv(TCurv.ChangeCurve(),UV1,UV2,PCurv,
2432             DStr.Surface(Fd->Surf()).Surface(),
2433             Pardeb,Parfin,tol3d,tolreached);
2434           TCurv.Tolerance(Max(TCurv.Tolerance(),tolreached));
2435           Interfc1= ChFi3d_FilCurveInDS (Icurv,Isurf,PCurv,ET1);
2436           DStr.ChangeSurfaceInterferences(Isurf).Append(Interfc1);
2437         }
2438       }
2439     } // End of Interference between fillets 
2440
2441     // ---- Interference Fillets / Faces
2442     IcFil1 = Fi1.LineIndex();
2443
2444     if (IcFil1!=0 ) {
2445       Interfc3= ChFi3d_FilCurveInDS (IcFil1,Isurf,
2446         Fi1.PCurveOnSurf(),trafil1);
2447       DStr.ChangeSurfaceInterferences(Isurf).Append(Interfc3);
2448       Ishape1 = Fd->IndexOfS1();
2449       // Case of degenerated edge : pcurve is associated via SCI 
2450       // to TopOpeBRepDSCurve.
2451       TopOpeBRepDS_Curve& cc = DStr.ChangeCurve(IcFil1);
2452       if(cc.Curve().IsNull()) {
2453         Handle(TopOpeBRepDS_Interference) bidinterf;
2454         cc.SetSCI(Interfc3,bidinterf);
2455       }
2456       else{
2457         ChFiDS_Regul regon1;
2458         regon1.SetCurve(IcFil1);
2459         regon1.SetS1(Isurf,Standard_False);
2460         if ( Ishape1 < 0 ) {
2461           Ishape1 = -Ishape1;
2462           regon1.SetS2(Ishape1,Standard_False);
2463           Interfc1=ChFi3d_FilCurveInDS(IcFil1,Ishape1,Fi1.PCurveOnFace(),
2464             Fi1.Transition()); 
2465           DStr.ChangeSurfaceInterferences(Ishape1).Append(Interfc1);      
2466         }
2467         else if ( Ishape1 > 0 ) {
2468           regon1.SetS2(Ishape1,Standard_True);
2469           Interfc1=ChFi3d_FilCurveInDS(IcFil1,Ishape1,Fi1.PCurveOnFace(),
2470             Fi1.Transition()); 
2471           DStr.ChangeShapeInterferences(Ishape1).Append(Interfc1);      
2472         }
2473         reglist.Append(regon1);
2474       }
2475       // Indice and type of the point at End
2476       Standard_Integer ipoin;
2477       Standard_Boolean isVertex = Fd->VertexLastOnS1().IsVertex();
2478       if (j == SeqFil.Length()) ipoin = CorDat->IndexLastPointOnS1();
2479       else if ( j == (SeqFil.Length()-1)  &&  /*Closed &&*/
2480         (DStr.Curve(SeqFil.Last()->InterferenceOnS1().
2481         LineIndex()).Curve().IsNull())) {
2482           if (Closed) {
2483             ipoin = CorDat->IndexFirstPointOnS1();
2484             isVertex = SeqFil(1)->VertexFirstOnS1().IsVertex();
2485           } else {
2486             ipoin = CorDat->IndexLastPointOnS1();
2487             isVertex = SeqFil.Last()->VertexLastOnS1().IsVertex();
2488           }
2489       }
2490       else if(DStr.Curve(IcFil1).Curve().IsNull()) {// Rotation !!
2491         ipoin = Ipoin1;
2492         isVertex = isVertex1;
2493       }
2494       else if ( ((j==1) || (j== SeqFil.Length()-1)) && 
2495         ( (Fd->VertexLastOnS1().Point().IsEqual(
2496         SeqFil(1)->VertexFirstOnS1().Point(), 1.e-7)) ||
2497         (Fd->VertexLastOnS1().Point().IsEqual(
2498         SeqFil(SeqFil.Length())->VertexLastOnS1().Point(), 1.e-7))) )
2499         // Case of SurfData cut in "Triangular" way.   
2500         ipoin=CorDat->IndexLastPointOnS1();
2501
2502       // eap, Apr 29 2002, occ 293
2503       else if (isInDS2 && findIndexPoint(DStr, Fd, 1, ipoin)) {
2504
2505       }
2506       else ipoin = ChFi3d_IndexPointInDS(Fd->VertexLastOnS1(),DStr);
2507
2508       TopOpeBRepDS_ListOfInterference& Li = DStr.ChangeCurveInterferences(IcFil1);
2509
2510       if (!ChFi3d_Contains(Li,IcFil1,Ipoin1)) { 
2511
2512         Interfp1 = ChFi3d_FilPointInDS(TopAbs_FORWARD,IcFil1,Ipoin1,
2513           Fi1.FirstParameter(),isVertex1);
2514         DStr.ChangeCurveInterferences(IcFil1).Append(Interfp1);
2515       }
2516       if (ipoin == Ipoin1 || !ChFi3d_Contains(Li,IcFil1,ipoin)) { 
2517         Interfp3 = ChFi3d_FilPointInDS(TopAbs_REVERSED,IcFil1,ipoin,
2518           Fi1.LastParameter(), isVertex);
2519         DStr.ChangeCurveInterferences(IcFil1).Append(Interfp3);
2520       }
2521       Ipoin1 = ipoin;
2522       isVertex1 = isVertex;
2523     }
2524
2525     IcFil2 = Fi2.LineIndex();
2526     if (IcFil2!=0) {
2527       Interfc4=ChFi3d_FilCurveInDS(IcFil2,Isurf,
2528         Fi2.PCurveOnSurf(),trafil2);
2529       DStr.ChangeSurfaceInterferences(Isurf).Append(Interfc4);
2530       Ishape2 = Fd->IndexOfS2();
2531       // Case of degenerated edge : pcurve is associated via SCI 
2532       // to TopOpeBRepDSCurve.
2533       TopOpeBRepDS_Curve& cc = DStr.ChangeCurve(IcFil2);
2534       if(cc.Curve().IsNull()) {
2535         Handle(TopOpeBRepDS_Interference) bidinterf;
2536         cc.SetSCI(Interfc4,bidinterf);
2537       }
2538       else{
2539         ChFiDS_Regul regon2;
2540         regon2.SetCurve(IcFil2);
2541         regon2.SetS1(Isurf,Standard_False);
2542         if ( Ishape2 < 0 ) {
2543           Ishape2 = -Ishape2;
2544           regon2.SetS2(Ishape2,Standard_False);
2545           Interfc2=ChFi3d_FilCurveInDS(IcFil2,Ishape2,Fi2.PCurveOnFace(),
2546             Fi2.Transition());
2547           DStr.ChangeSurfaceInterferences(Ishape2).Append(Interfc2);      
2548         }
2549         else if ( Ishape2 > 0 ) {
2550           regon2.SetS2(Ishape2,Standard_True);
2551           Interfc2=ChFi3d_FilCurveInDS(IcFil2,Ishape2,Fi2.PCurveOnFace(),
2552             Fi2.Transition());
2553           DStr.ChangeShapeInterferences(Ishape2).Append(Interfc2);      
2554         }
2555         reglist.Append(regon2);
2556       }
2557       // Indice and type of the point in End
2558       Standard_Integer ipoin;
2559       Standard_Boolean isVertex = Fd->VertexLastOnS2().IsVertex();
2560       if (j == SeqFil.Length() ) ipoin = CorDat->IndexLastPointOnS2();
2561       else if ( j == (SeqFil.Length()-1)  && /*Closed &&*/
2562         (DStr.Curve(SeqFil.Last()->InterferenceOnS2().
2563         LineIndex()).Curve().IsNull())) {
2564           if (Closed) {
2565             ipoin = CorDat->IndexFirstPointOnS2();
2566             isVertex = SeqFil(1)->VertexFirstOnS2().IsVertex();
2567           } else {
2568             ipoin = CorDat->IndexLastPointOnS2();
2569             isVertex = SeqFil.Last()->VertexLastOnS2().IsVertex();
2570           }
2571       }
2572       else if(DStr.Curve(IcFil2).Curve().IsNull()) { // Rotation !!
2573         ipoin = Ipoin2;
2574         isVertex = isVertex2;
2575       }
2576       else if(Fd->VertexLastOnS2().Point().IsEqual(
2577         Fd->VertexLastOnS1().Point(), 0) ) {  //Pinch !!
2578           ipoin = Ipoin1;
2579           isVertex = isVertex1;
2580       }
2581       else if ( ((j==1) || (j==SeqFil.Length()-1)) && 
2582         ( (Fd->VertexLastOnS2().Point().IsEqual(
2583         SeqFil(1)->VertexFirstOnS2().Point(), 1.e-7)) ||
2584         (Fd->VertexLastOnS2().Point().IsEqual(
2585         SeqFil(SeqFil.Length())->VertexLastOnS2().Point(), 1.e-7))) )
2586         // Case of SurfData cut in "Triangular" way.   
2587         ipoin=CorDat->IndexLastPointOnS2();
2588
2589       // eap, Apr 29 2002, occ 293
2590       else if (isInDS2 && findIndexPoint(DStr, Fd, 2, ipoin)) {
2591
2592       }
2593       else ipoin = ChFi3d_IndexPointInDS(Fd->VertexLastOnS2(),DStr);
2594
2595       TopOpeBRepDS_ListOfInterference& Li = DStr.ChangeCurveInterferences(IcFil2);
2596
2597       if (!ChFi3d_Contains(Li,IcFil2,Ipoin2)) { 
2598         Interfp2 = ChFi3d_FilPointInDS(TopAbs_FORWARD,IcFil2,Ipoin2,
2599           Fi2.FirstParameter(), isVertex2);
2600         DStr.ChangeCurveInterferences(IcFil2).Append(Interfp2);
2601       }
2602       if (ipoin == Ipoin2 || !ChFi3d_Contains(Li,IcFil2,ipoin)) { 
2603         Interfp4= ChFi3d_FilPointInDS(TopAbs_REVERSED,IcFil2,ipoin,
2604           Fi2.LastParameter(), isVertex );
2605         DStr.ChangeCurveInterferences(IcFil2).Append(Interfp4);
2606       }
2607       Ipoin2 = ipoin;
2608       isVertex2 = isVertex;      
2609     }
2610
2611     ET1 = trafil1;
2612     if (j == SeqFil.Length()) {
2613       if (!isInDS2) {
2614         Icurv = CorDat->LastCurve();
2615         if(Closed && !Singulier_en_Bout && (Ipoin1!=Ipoin2)) {
2616           regcout.SetS2(Isurf,Standard_False);
2617           reglist.Append(regcout);
2618         }
2619         PCurv = CorDat->LastPCurve();
2620         ET1 = TopAbs::Compose(ET1,CorDat->LastPCurveOrientation());
2621         CorDat->LastParameters(Pardeb,Parfin);
2622         TopOpeBRepDS_ListOfInterference& Li = DStr.ChangeCurveInterferences(Icurv);
2623         if (Li.IsEmpty()) {
2624           if(CorDat->LastPCurveOrientation()==TopAbs_REVERSED) {
2625             Interfp5=ChFi3d_FilPointInDS
2626               (TopAbs_REVERSED,Icurv,Ipoin1,Parfin, isVertex1);
2627             Interfp6=ChFi3d_FilPointInDS
2628               (TopAbs_FORWARD,Icurv,Ipoin2,Pardeb, isVertex2);
2629           }
2630           else{
2631             Interfp5=ChFi3d_FilPointInDS
2632               (TopAbs_FORWARD,Icurv,Ipoin1,Pardeb, isVertex1);
2633             Interfp6=ChFi3d_FilPointInDS
2634               (TopAbs_REVERSED,Icurv,Ipoin2,Parfin, isVertex2);
2635           }
2636           Li.Append(Interfp5);
2637           Li.Append(Interfp6);
2638         }
2639         Interfc1= ChFi3d_FilCurveInDS(Icurv,Isurf,PCurv,ET1);
2640         DStr.ChangeSurfaceInterferences(Isurf).Append(Interfc1);
2641         if (Ipoin1 == Ipoin2) {
2642           TopOpeBRepDS_Curve& TCurv = DStr.ChangeCurve(Icurv);
2643           TCurv.ChangeCurve().Nullify();
2644           Handle(TopOpeBRepDS_Interference) bidinterf;
2645           TCurv.SetSCI( Interfc1, bidinterf);
2646           //         bidinterf = TCurv.GetSCI1(); 
2647           //      TCurv.SetSCI(bidinterf, Interfc1);        
2648         }
2649       }
2650     }
2651     else {
2652       //      Degene = (Fd->VertexLastOnS1().Point().IsEqual(
2653       //                Fd->VertexLastOnS2().Point(), 0) );
2654
2655       // eap, Apr 29 2002, occ 293 
2656       if (!isInDS2) {
2657
2658         Handle(Geom_Curve) C3d;
2659         Standard_Real tolreached;
2660         ChFi3d_ComputeArete(Fd->VertexLastOnS1(),
2661           Fd->InterferenceOnS1().PCurveOnSurf()->
2662           Value(Fd->InterferenceOnS1().LastParameter()),
2663           Fd->VertexLastOnS2(),
2664           Fd->InterferenceOnS2().PCurveOnSurf()->
2665           Value(Fd->InterferenceOnS2().LastParameter()),
2666           DStr.Surface(Fd->Surf()).Surface(),C3d,PCurv,
2667           Pardeb,Parfin,tol3d,tol2d,tolreached,0);
2668         Crv = TopOpeBRepDS_Curve(C3d,tolreached);
2669         Icurv = DStr.AddCurve(Crv);
2670         regfilfil.SetCurve(Icurv);
2671         regfilfil.SetS1(Isurf,Standard_False);
2672         Interfp5 = ChFi3d_FilPointInDS(TopAbs_FORWARD,Icurv,Ipoin1,Pardeb, isVertex1);
2673         DStr.ChangeCurveInterferences(Icurv).Append(Interfp5);
2674         Interfp6= ChFi3d_FilPointInDS(TopAbs_REVERSED,Icurv,Ipoin2,Parfin, isVertex2);
2675         DStr.ChangeCurveInterferences(Icurv).Append(Interfp6);
2676         Interfc1= ChFi3d_FilCurveInDS(Icurv,Isurf,PCurv,ET1);
2677         DStr.ChangeSurfaceInterferences(Isurf).Append(Interfc1);      
2678       }
2679     }
2680
2681     Degene = V3.Point().IsEqual(V4.Point(), 0);
2682
2683     // Processing of degenerated case     
2684     if (Degene) {
2685       // Queue de Billard
2686       Standard_Boolean Vertex = (V3.IsVertex()) && (V4.IsVertex());
2687       if (!Vertex) {
2688
2689       }
2690       else {
2691         // The edge of the spine starting on this vertex is removed.
2692         Standard_Boolean Trouve = Standard_False;
2693         TopoDS_Edge Arcspine;
2694         TopAbs_Orientation OVtx = TopAbs_FORWARD;
2695         BoutdeVtx = V3.Vertex();
2696
2697         while (NumEdge<= spine->NbEdges() && !Trouve) { 
2698           Arcspine = spine->Edges(NumEdge);
2699           for(ex.Init(Arcspine.Oriented(TopAbs_FORWARD),TopAbs_VERTEX); 
2700             ex.More() && (!Trouve); ex.Next()) {
2701               if(BoutdeVtx.IsSame(ex.Current())) {
2702                 OVtx = ex.Current().Orientation();
2703                 if (Closed &&  (NumEdge == 1)) 
2704                   Trouve = (spine->NbEdges() == 1);
2705                 else Trouve = Standard_True;
2706               }
2707           }
2708           if (!Trouve) NumEdge++; // Go to the next edge
2709         }
2710         Standard_Integer IArcspine = DStr.AddShape(Arcspine);
2711         Standard_Integer IVtx;
2712         if  (j == SeqFil.Length()) {
2713           IVtx = CorDat->IndexLastPointOnS1();
2714         }
2715         else { IVtx = DStr.AddShape(BoutdeVtx); }
2716         OVtx = TopAbs::Reverse(OVtx);
2717         Standard_Real parVtx = BRep_Tool::Parameter(BoutdeVtx,Arcspine);
2718         Handle(TopOpeBRepDS_CurvePointInterference) 
2719           interfv = ChFi3d_FilVertexInDS(OVtx,IArcspine,IVtx,parVtx);
2720         DStr.ChangeShapeInterferences(IArcspine).Append(interfv);
2721       }
2722     } // end of degenerated case
2723     else if (!(Closed && j == SeqFil.Length())) {
2724       // Processing of interference Point / Edges
2725       if (V3.IsOnArc()) {
2726         if(!(V3.IsVertex() && Fd->IsOnCurve1())) {
2727           Iarc1 = DStr.AddShape(V3.Arc());
2728           if ( !ChFi3d_Contains(DStr.ShapeInterferences(Iarc1),Iarc1,Ipoin1,V3.IsVertex(),Standard_True) ) {
2729             Handle(TopOpeBRepDS_CurvePointInterference) Interfpp = 
2730               ChFi3d_FilPointInDS(V3.TransitionOnArc(),
2731               Iarc1,Ipoin1,V3.ParameterOnArc(), V3.IsVertex());
2732             DStr.ChangeShapeInterferences(V3.Arc()).Append(Interfpp);
2733           }
2734         }
2735       }
2736
2737       if (V4.IsOnArc()) {
2738         if(!(V4.IsVertex() && Fd->IsOnCurve2())) {
2739           Iarc2 = DStr.AddShape(V4.Arc());
2740           if ( !ChFi3d_Contains(DStr.ShapeInterferences(Iarc2),Iarc2,Ipoin2,V4.IsVertex(),Standard_True) ) {
2741             Handle(TopOpeBRepDS_CurvePointInterference) Intfpp=
2742               ChFi3d_FilPointInDS(V4.TransitionOnArc(),
2743               Iarc2,Ipoin2,V4.ParameterOnArc(), V4.IsVertex());
2744             DStr.ChangeShapeInterferences(V4.Arc()).Append(Intfpp);
2745           }
2746         }
2747       }
2748     }
2749   }
2750 }
2751 //=======================================================================
2752 //function : StripeEdgeInter
2753 //purpose  : This function examines two stripes for an intersection 
2754 //           between curves of interference with faces. If the intersection
2755 //           exists, it will cause bad result, so it's better to quit.
2756 //remark   : If someone somewhen computes the interference between stripes, 
2757 //           this function will become useless.
2758 //author   : akm, 06/02/02. Against bug OCC119.
2759 //=======================================================================
2760 void ChFi3d_StripeEdgeInter (const Handle(ChFiDS_Stripe)& theStripe1,
2761   const Handle(ChFiDS_Stripe)& theStripe2,
2762   TopOpeBRepDS_DataStructure&  /*DStr*/,
2763   const Standard_Real          tol2d)
2764 {
2765   // Do not check the stripeshaving common corner points
2766   for (Standard_Integer iSur1=1; iSur1<=2; iSur1++)
2767     for (Standard_Integer iSur2=1; iSur2<=2; iSur2++)
2768       if (theStripe1->IndexPoint(0,iSur1)==theStripe2->IndexPoint(0,iSur2) ||
2769         theStripe1->IndexPoint(0,iSur1)==theStripe2->IndexPoint(1,iSur2) ||
2770         theStripe1->IndexPoint(1,iSur1)==theStripe2->IndexPoint(0,iSur2) ||
2771         theStripe1->IndexPoint(1,iSur1)==theStripe2->IndexPoint(1,iSur2))
2772         return;
2773
2774   Handle(ChFiDS_HData) aSurDat1 = theStripe1->SetOfSurfData();
2775   Handle(ChFiDS_HData) aSurDat2 = theStripe2->SetOfSurfData();
2776
2777   Geom2dInt_GInter anIntersector;
2778   Standard_Integer iPart1, iPart2;
2779   Standard_Integer Ishape11, Ishape12, Ishape21, Ishape22;
2780   // Loop on parts of the first stripe
2781   for (iPart1=1; iPart1<=aSurDat1->Length(); iPart1++) 
2782   {
2783     Handle(ChFiDS_SurfData) aDat1 = aSurDat1->Value(iPart1);
2784     Ishape11 = aDat1->IndexOfS1();
2785     Ishape12 = aDat1->IndexOfS2();
2786     // Loop on parts of the second stripe
2787     for (iPart2=1; iPart2<=aSurDat2->Length(); iPart2++) 
2788     {
2789       Handle(ChFiDS_SurfData) aDat2 = aSurDat2->Value(iPart2);
2790       Ishape21 = aDat2->IndexOfS1();
2791       Ishape22 = aDat2->IndexOfS2();
2792
2793       // Find those FaceInterferences able to intersect
2794       ChFiDS_FaceInterference aFI1, aFI2;
2795       if (Ishape11 == Ishape21)
2796       {
2797         aFI1 = aDat1->InterferenceOnS1();
2798         aFI2 = aDat2->InterferenceOnS1();
2799       } 
2800       else if (Ishape11 == Ishape22)
2801       {
2802         aFI1 = aDat1->InterferenceOnS1();
2803         aFI2 = aDat2->InterferenceOnS2();
2804       } 
2805       else if (Ishape12 == Ishape21)
2806       {
2807         aFI1 = aDat1->InterferenceOnS2();
2808         aFI2 = aDat2->InterferenceOnS1();
2809       } 
2810       else if (Ishape12 == Ishape22)
2811       {
2812         aFI1 = aDat1->InterferenceOnS2();
2813         aFI2 = aDat2->InterferenceOnS2();
2814       }
2815       else
2816       {
2817         // No common faces
2818         continue;
2819       }
2820
2821       if (IsEqual (aFI1.FirstParameter(),aFI1.LastParameter()) ||
2822         IsEqual (aFI2.FirstParameter(),aFI2.LastParameter()) ||
2823         aFI1.PCurveOnFace().IsNull() ||
2824         aFI2.PCurveOnFace().IsNull())
2825         // Do not waste time on degenerates
2826         continue;
2827       // Examine for intersections
2828       Geom2dAdaptor_Curve aPCurve1 (aFI1.PCurveOnFace(),
2829         aFI1.FirstParameter(),
2830         aFI1.LastParameter());
2831       Geom2dAdaptor_Curve aPCurve2 (aFI2.PCurveOnFace(),
2832         aFI2.FirstParameter(),
2833         aFI2.LastParameter());
2834       anIntersector.Perform (aPCurve1,
2835         aPCurve2,
2836         tol2d,
2837         Precision::PConfusion());
2838       if (anIntersector.NbSegments() > 0 ||
2839         anIntersector.NbPoints() > 0)
2840         throw StdFail_NotDone("StripeEdgeInter : fillets have too big radiuses");
2841     }
2842   }
2843 }
2844
2845 //=======================================================================
2846 //function : IndexOfSurfData
2847 //purpose  : 
2848 //=======================================================================
2849 Standard_Integer ChFi3d_IndexOfSurfData(const TopoDS_Vertex& V1,
2850   const Handle(ChFiDS_Stripe)& CD,
2851   Standard_Integer& sens)
2852 {
2853   Handle(ChFiDS_Spine) spine = CD->Spine();
2854   Standard_Integer Index = 0;
2855   sens = 1;
2856   TopoDS_Vertex Vref;
2857   const TopoDS_Edge& E = spine->Edges(1);
2858   if (E.Orientation() == TopAbs_REVERSED) Vref = TopExp::LastVertex(E);
2859   else Vref = TopExp::FirstVertex(E); 
2860   if (Vref.IsSame(V1)) Index =1;
2861   else {
2862     const TopoDS_Edge& E1 = spine->Edges(spine->NbEdges());
2863     if (E1.Orientation() == TopAbs_REVERSED) Vref = TopExp::FirstVertex(E1);
2864     else Vref = TopExp::LastVertex(E1); 
2865     sens = -1;
2866     if(CD->SetOfSurfData().IsNull()) return 0;
2867     else if (Vref.IsSame(V1)) Index = CD->SetOfSurfData()->Length();
2868     else throw Standard_ConstructionError ("ChFi3d_IndexOfSurfData() - wrong construction parameters");
2869   }
2870   return Index; 
2871 }  
2872 //=======================================================================
2873 //function : EdgeFromV1
2874 //purpose  : 
2875 //=======================================================================
2876
2877 TopoDS_Edge ChFi3d_EdgeFromV1(const TopoDS_Vertex& V1,
2878   const Handle(ChFiDS_Stripe)& CD,
2879   Standard_Integer& sens)
2880 {
2881   Handle(ChFiDS_Spine) spine = CD->Spine();
2882   sens = 1;
2883   TopoDS_Vertex Vref;
2884   const TopoDS_Edge& E = spine->Edges(1);
2885   if (E.Orientation() == TopAbs_REVERSED) Vref = TopExp::LastVertex(E);
2886   else Vref = TopExp::FirstVertex(E); 
2887   if (Vref.IsSame(V1)) return E;
2888   else
2889   {
2890     const TopoDS_Edge& E1 = spine->Edges(spine->NbEdges());
2891     if (E1.Orientation() == TopAbs_REVERSED) Vref = TopExp::FirstVertex(E1);
2892     else Vref = TopExp::LastVertex(E1); 
2893     sens = -1;
2894     if (Vref.IsSame(V1)) return E1;
2895     else throw Standard_ConstructionError ("ChFi3d_IndexOfSurfData() - wrong construction parameters");
2896   }
2897 }
2898 //=======================================================================
2899 //function : ConvTol2dToTol3d
2900 //purpose  : Comme son nom l indique.
2901 //=======================================================================
2902
2903 Standard_Real ChFi3d_ConvTol2dToTol3d(const Handle(Adaptor3d_HSurface)& S,
2904   const Standard_Real             tol2d)
2905 {
2906   Standard_Real ures = S->UResolution(1.e-7);
2907   Standard_Real vres = S->VResolution(1.e-7);
2908   Standard_Real uresto3d = 1.e-7*tol2d/ures;
2909   Standard_Real vresto3d = 1.e-7*tol2d/vres;
2910   return Max(uresto3d,vresto3d);
2911 }
2912 //=======================================================================
2913 //function : EvalTolReached
2914 //purpose  : The function above is too hard because 
2915 //           parametrization of surfaces is not homogenous.
2916 //=======================================================================
2917
2918 Standard_Real ChFi3d_EvalTolReached(const Handle(Adaptor3d_HSurface)& S1,
2919   const Handle(Geom2d_Curve)&     pc1,
2920   const Handle(Adaptor3d_HSurface)& S2,
2921   const Handle(Geom2d_Curve)&     pc2,
2922   const Handle(Geom_Curve)&       C)
2923 {
2924   Standard_Real distmax = 0.;
2925
2926   Standard_Real f = C->FirstParameter();
2927   Standard_Real l = C->LastParameter();
2928   Standard_Integer nbp = 45;
2929   Standard_Real step = 1./(nbp -1);
2930   for(Standard_Integer i = 0; i < nbp; i++) {
2931     Standard_Real t,u,v;
2932     t = step * i;
2933     t = (1-t) * f + t * l;
2934     pc1->Value(t).Coord(u,v);
2935     gp_Pnt pS1 = S1->Value(u,v);
2936     pc2->Value(t).Coord(u,v);
2937     gp_Pnt pS2 = S2->Value(u,v);
2938     gp_Pnt pC = C->Value(t);
2939     Standard_Real d = pS1.SquareDistance(pC);
2940     if(d>distmax) distmax = d;
2941     d = pS2.SquareDistance(pC);
2942     if(d>distmax) distmax = d;
2943     d = pS1.SquareDistance(pS2);
2944     if(d>distmax) distmax = d;
2945   }
2946   distmax = 1.5*sqrt(distmax);
2947   distmax = Max(distmax, Precision::Confusion());
2948   return distmax;
2949 }
2950
2951 //=======================================================================
2952 //function : trsfsurf
2953 //purpose  : 
2954 //=======================================================================
2955 Handle(Geom_Surface) trsfsurf(const Handle(Adaptor3d_HSurface)& HS,
2956   Handle(Adaptor3d_TopolTool)&      /*dom*/)
2957 {
2958   //Pour l utilisation des domaines voir avec BUBUCH!!
2959   Handle(Geom_Surface) res;
2960   Handle(BRepAdaptor_HSurface) hbs = Handle(BRepAdaptor_HSurface)::DownCast(HS);
2961   Handle(GeomAdaptor_HSurface) hgs = Handle(GeomAdaptor_HSurface)::DownCast(HS);
2962   if(!hbs.IsNull()) {
2963     res = hbs->ChangeSurface().Surface().Surface();
2964     gp_Trsf trsf = hbs->ChangeSurface().Trsf();
2965     res = Handle(Geom_Surface)::DownCast(res->Transformed(trsf));
2966   }
2967   else if(!hgs.IsNull()) {
2968     res = hgs->ChangeSurface().Surface();
2969   }
2970   Handle(Geom_RectangularTrimmedSurface) 
2971     tr = Handle(Geom_RectangularTrimmedSurface)::DownCast(res);
2972   if(!tr.IsNull()) res = tr->BasisSurface();
2973
2974   Standard_Real U1 = HS->FirstUParameter(), U2 = HS->LastUParameter();
2975   Standard_Real V1 = HS->FirstVParameter(), V2 = HS->LastVParameter();
2976   if(!res.IsNull()) {
2977     // Protection against Construction Errors
2978     Standard_Real u1, u2, v1, v2;
2979     res->Bounds( u1, u2, v1, v2);
2980     if (!res->IsUPeriodic()) {
2981       if (U1 < u1) U1 = u1;
2982       if (U2 > u2) U2 = u2;
2983     }
2984     if (!res->IsVPeriodic()) {
2985       if (V1 < v1) V1 = v1;
2986       if (V2 > v2) V2 = v2;
2987     }
2988     res = new Geom_RectangularTrimmedSurface(res,U1,U2,V1,V2);
2989   }
2990   //  Handle(GeomAdaptor_HSurface) temp = new GeomAdaptor_HSurface(res,U1,U2,V1,V2);
2991   //  dom = new Adaptor3d_TopolTool(temp);
2992   return res;
2993 }
2994 //=======================================================================
2995 //function : CurveCleaner
2996 //purpose  : Makes a BSpline as much continued as possible
2997 //           at a given tolerance
2998 //=======================================================================
2999 static void CurveCleaner(Handle(Geom_BSplineCurve)& BS, 
3000   const Standard_Real Tol,
3001   const Standard_Integer MultMin)
3002
3003 {
3004   Standard_Real tol = Tol;
3005   Standard_Integer Mult, ii;
3006   const Standard_Integer NbK=BS->NbKnots();
3007
3008   for (Mult = BS->Degree(); Mult > MultMin; Mult--) {
3009     tol *= 0.5; // Progressive reduction
3010     for (ii=NbK; ii>1; ii--) {
3011       if (BS->Multiplicity(ii) == Mult)
3012         BS->RemoveKnot(ii, Mult-1, tol);
3013     }
3014   }
3015 }
3016 //=======================================================================
3017 //function : ComputeCurves
3018 //purpose  : Calculates intersection between two HSurfaces.
3019 //           It is necessary to know the extremities of intersection and  
3020 //           the surfaces should be processed at input 
3021 //           to fit as good as possible (neither too close nor too far) 
3022 //           the points of beginning and end of the intersection.
3023 //           The analytic intersections are processed separately.
3024 //           <wholeCurv> means that the resulting curve is restricted by
3025 //           boundaries of input surfaces (eap 30 May occ354)
3026 //=======================================================================
3027 Standard_Boolean ChFi3d_ComputeCurves(const Handle(Adaptor3d_HSurface)&   S1,
3028   const Handle(Adaptor3d_HSurface)&   S2,
3029   const TColStd_Array1OfReal& Pardeb,
3030   const TColStd_Array1OfReal& Parfin,
3031   Handle(Geom_Curve)&         C3d,
3032   Handle(Geom2d_Curve)&       Pc1,
3033   Handle(Geom2d_Curve)&       Pc2,
3034   const Standard_Real         tol3d,
3035   const Standard_Real         tol2d,
3036   Standard_Real&              tolreached,
3037   const Standard_Boolean      wholeCurv) 
3038 {
3039   gp_Pnt pdeb1 = S1->Value(Pardeb(1),Pardeb(2));
3040   gp_Pnt pfin1 = S1->Value(Parfin(1),Parfin(2));
3041   gp_Pnt pdeb2 = S2->Value(Pardeb(3),Pardeb(4));
3042   gp_Pnt pfin2 = S2->Value(Parfin(3),Parfin(4));
3043
3044   Standard_Real distrefdeb = pdeb1.Distance(pdeb2);//checks the worthiness 
3045   Standard_Real distreffin = pfin1.Distance(pfin2);//of input data
3046   if(distrefdeb < tol3d) distrefdeb = tol3d;
3047   if(distreffin < tol3d) distreffin = tol3d;
3048
3049   gp_Pnt pdeb,pfin;
3050   pdeb.SetXYZ(0.5*(pdeb1.XYZ()+pdeb2.XYZ()));
3051   pfin.SetXYZ(0.5*(pfin1.XYZ()+pfin2.XYZ()));
3052
3053   Standard_Real distref = 0.005*pdeb.Distance(pfin);
3054   if(distref < distrefdeb) distref = distrefdeb;
3055   if(distref < distreffin) distref = distreffin;
3056
3057   //Some analytic cases are processed separately.
3058   //To reorientate the result of the analythic intersection,
3059   //it is stated that the beginning of the tangent should be
3060   //in the direction of the start/end line.
3061   gp_Vec Vint, Vref(pdeb,pfin);
3062   gp_Pnt Pbid;
3063   Standard_Real Udeb = 0.,Ufin = 0.;
3064   Standard_Real tolr1,tolr2;
3065   tolr1 = tolr2 = tolreached = tol3d;
3066   if((S1->GetType() == GeomAbs_Cylinder && S2->GetType() == GeomAbs_Plane)||
3067     (S1->GetType() == GeomAbs_Plane && S2->GetType() == GeomAbs_Cylinder)) { 
3068       gp_Pln pl;
3069       gp_Cylinder cyl;
3070       if(S1->GetType() == GeomAbs_Plane) {
3071         pl = S1->Plane();
3072         cyl = S2->Cylinder();
3073       }
3074       else{
3075         pl = S2->Plane();
3076         cyl = S1->Cylinder();
3077       }
3078       IntAna_QuadQuadGeo ImpKK(pl,cyl,Precision::Angular(),tol3d);
3079       Standard_Boolean isIntDone = ImpKK.IsDone();
3080
3081       if(ImpKK.TypeInter() == IntAna_Ellipse)
3082       {
3083         const gp_Elips anEl = ImpKK.Ellipse(1);
3084         const Standard_Real aMajorR = anEl.MajorRadius();
3085         const Standard_Real aMinorR = anEl.MinorRadius();
3086         isIntDone = (aMajorR < 100000.0 * aMinorR);
3087       }
3088
3089       if (isIntDone) {
3090         Standard_Boolean c1line = 0;
3091         switch  (ImpKK.TypeInter()) {
3092         case IntAna_Line:
3093           {
3094             c1line = 1;
3095             Standard_Integer nbsol = ImpKK.NbSolutions();
3096             gp_Lin C1;
3097             for(Standard_Integer ilin = 1; ilin <= nbsol; ilin++) {
3098               C1 = ImpKK.Line(ilin);
3099               Udeb = ElCLib::Parameter(C1,pdeb);
3100               gp_Pnt ptest = ElCLib::Value(Udeb,C1);
3101               if(ptest.Distance(pdeb) < tol3d) break;
3102             }
3103             Ufin = ElCLib::Parameter(C1,pfin);
3104             C3d = new Geom_Line(C1);
3105             ElCLib::D1(Udeb,C1,Pbid,Vint);
3106           }
3107           break;
3108         case IntAna_Circle:
3109           {
3110             gp_Circ C1 = ImpKK.Circle(1);
3111             C3d = new Geom_Circle(C1);
3112             Udeb = ElCLib::Parameter(C1,pdeb);
3113             Ufin = ElCLib::Parameter(C1,pfin);
3114             ElCLib::D1(Udeb,C1,Pbid,Vint);
3115           }
3116           break;
3117         case IntAna_Ellipse:
3118           {
3119             gp_Elips C1 = ImpKK.Ellipse(1);
3120             C3d = new Geom_Ellipse(C1);
3121             Udeb = ElCLib::Parameter(C1,pdeb);
3122             Ufin = ElCLib::Parameter(C1,pfin);
3123             ElCLib::D1(Udeb,C1,Pbid,Vint);
3124           }
3125           break;
3126         default:
3127           break;
3128         }
3129         if (Vint.Dot(Vref)<0) {
3130           C3d->Reverse();
3131           if(c1line) {
3132             Udeb = -Udeb;
3133             Ufin = -Ufin;
3134           }
3135           else{
3136             Udeb = 2*M_PI - Udeb;
3137             Ufin = 2*M_PI - Ufin;
3138           }
3139         }
3140         if(!c1line) ElCLib::AdjustPeriodic(0.,2*M_PI,Precision::Angular(),Udeb,Ufin);
3141         Handle(GeomAdaptor_HCurve) HC = new GeomAdaptor_HCurve();
3142         HC->ChangeCurve().Load(C3d,Udeb,Ufin);
3143         ChFi3d_ProjectPCurv(HC,S1,Pc1,tol3d,tolr1);
3144         if(S1->GetType() == GeomAbs_Cylinder) {
3145           Standard_Real x,y;
3146           Pc1->Value(Udeb).Coord(x,y);
3147           x = Pardeb(1) - x;
3148           y = Pardeb(2) - y;
3149           if(Abs(x) >= tol2d || Abs(y) >= tol2d) Pc1->Translate(gp_Vec2d(x,y));
3150         }
3151         ChFi3d_ProjectPCurv(HC,S2,Pc2,tol3d,tolr2);
3152         if(S2->GetType() == GeomAbs_Cylinder) {
3153           Standard_Real x,y;
3154           Pc2->Value(Udeb).Coord(x,y);
3155           x = Pardeb(3) - x;
3156           y = Pardeb(4) - y;
3157           if(Abs(x) >= tol2d || Abs(y) >= tol2d) Pc2->Translate(gp_Vec2d(x,y));
3158         }
3159         C3d = new Geom_TrimmedCurve(C3d,Udeb,Ufin);
3160         tolreached = 1.5*Max(tolr1,tolr2);
3161         tolreached = Min(tolreached,ChFi3d_EvalTolReached(S1,Pc1,S2,Pc2,C3d));
3162         return Standard_True;
3163       }
3164   }    
3165   else if(S1->GetType() == GeomAbs_Plane && S2->GetType() == GeomAbs_Plane) { 
3166     IntAna_QuadQuadGeo LInt(S1->Plane(),S2->Plane(),Precision::Angular(),tol3d);
3167     if (LInt.IsDone()) {
3168       gp_Lin L = LInt.Line(1);
3169       C3d = new Geom_Line(L);
3170       Udeb = ElCLib::Parameter(L,pdeb);
3171       Ufin = ElCLib::Parameter(L,pfin);
3172       ElCLib::D1(Udeb,L,Pbid,Vint);
3173       if (Vint.Dot(Vref)<0) {
3174         C3d->Reverse();
3175         Udeb = - Udeb;
3176         Ufin = - Ufin;
3177       }
3178       Handle(GeomAdaptor_HCurve) HC = new GeomAdaptor_HCurve();
3179       HC->ChangeCurve().Load(C3d,Udeb,Ufin);
3180       ChFi3d_ProjectPCurv(HC,S1,Pc1,tol3d,tolr1);
3181       ChFi3d_ProjectPCurv(HC,S2,Pc2,tol3d,tolr2);
3182       C3d = new Geom_TrimmedCurve(C3d,Udeb,Ufin);
3183       return Standard_True;
3184     }
3185   }
3186   else {
3187     // here GeomInt is approached.
3188     Handle(Adaptor3d_TopolTool) dom1,dom2;
3189     Handle(Geom_Surface) gs1 = trsfsurf(S1,dom1);
3190     Handle(Geom_Surface) gs2 = trsfsurf(S2,dom2);
3191     Standard_Integer nbl ;
3192     if(!gs1.IsNull() && !gs2.IsNull()) {
3193       GeomInt_IntSS inter;
3194       //  Modified by skv - Fri Oct 24 14:24:47 2003 OCC4077 Begin
3195       //       Standard_Real tolap = 1.e-7;//car l approx de la wline est faite dans [0,1]
3196       // Set the lowest tolerance which is used in new boolean operations.
3197       Standard_Real tolap = 2.e-7;
3198       //  Modified by skv - Fri Oct 24 14:24:48 2003 OCC4077 End
3199       inter.Perform(gs1,gs2,tolap,1,1,1);
3200       if(inter.IsDone()) {
3201         nbl = inter.NbLines();   
3202 #if defined(IRIX) || defined(__sgi)
3203         if(nbl==0) {
3204
3205           //  solution of adjustment for SGI 
3206           //  if the intersection of gs1 with gs2 doesnot worke 
3207           //  then the intersection of gs2 with gs1 is attempted.
3208
3209           inter.Perform(gs2,gs1,tolap,1,1,1);
3210           //          inter.Perform(gs2,dom2,gs1,dom1,tolap,1,1,1);
3211           if(!inter.IsDone()) return Standard_False; 
3212           nbl = inter.NbLines(); 
3213
3214           //  if GeomInt does not make the intersection the solution of adjustment 
3215           //  is not attempted 
3216           if (nbl==0) return Standard_False;
3217         }
3218 #endif
3219         GeomAPI_ProjectPointOnCurve proj;
3220         for(Standard_Integer ilin = 1; ilin <= nbl; ilin++) {
3221           if(inter.HasLineOnS1(ilin) && inter.HasLineOnS2(ilin)) {
3222             C3d = inter.Line(ilin);
3223             Pc1 = inter.LineOnS1(ilin);
3224             Pc2 = inter.LineOnS2(ilin);
3225             gp_Pnt ptestdeb, ptestfin;
3226             Standard_Real Uf=0., Ul=0.;
3227             if (wholeCurv) {
3228               Uf = C3d->FirstParameter();
3229               Ul = C3d->LastParameter();
3230               ptestdeb = C3d->Value(Uf);
3231               ptestfin = C3d->Value(Ul);
3232             }
3233             else {
3234               // find end parameters
3235               Standard_Boolean failedF, failedL;
3236               failedF = failedL = Standard_False;
3237               proj.Init( pdeb1, C3d);
3238               if (proj.NbPoints()==0 && distrefdeb > Precision::Confusion())
3239                 proj.Perform( pdeb2 );
3240               if (proj.NbPoints()==0)
3241                 failedF = Standard_True;
3242               else
3243                 Uf = proj.LowerDistanceParameter();
3244               proj.Perform( pfin1 );
3245               if (proj.NbPoints()==0 && distreffin > Precision::Confusion())
3246                 proj.Perform( pfin2 );
3247               if (proj.NbPoints()==0) 
3248                 failedL = Standard_True;
3249               else
3250                 Ul = proj.LowerDistanceParameter();
3251
3252               if (failedF && failedL) {
3253                 Uf = C3d->FirstParameter();
3254                 Ul = C3d->LastParameter();
3255               }
3256               else if (failedF || failedL) {
3257                 // select right end parameter
3258                 Standard_Real Uok = failedF ? Ul : Uf;
3259                 Standard_Real U1 = C3d->FirstParameter(), U2 = C3d->LastParameter();
3260                 Uok = Abs(Uok-U1) > Abs(Uok-U2) ? U1 : U2;
3261                 if (failedF) Uf = Uok;
3262                 else         Ul = Uok;
3263               }
3264               else { // both projected, but where?
3265                 if (Abs(Uf - Ul) < Precision::PConfusion())
3266                   continue;
3267               }
3268               ptestdeb = C3d->Value(Uf);
3269               ptestfin = C3d->Value(Ul);
3270               if (C3d->IsPeriodic() && !(failedF && failedL)) {
3271                 // assure the same order of ends, otherwise TrimmedCurve will take
3272                 // the other part of C3d
3273                 gp_Pnt Ptmp;
3274                 gp_Vec DirOld, DirNew(ptestdeb,ptestfin);
3275                 C3d->D1(Uf, Ptmp, DirOld);
3276                 if (DirOld * DirNew < 0) {
3277                   Standard_Real Utmp = Uf; Uf = Ul; Ul = Utmp;
3278                   Ptmp = ptestdeb; ptestdeb = ptestfin; ptestfin = Ptmp;
3279                 }
3280               }
3281             }
3282             C3d = new Geom_TrimmedCurve(C3d,Uf,Ul);
3283             Pc1 = new Geom2d_TrimmedCurve(Pc1,Uf,Ul);
3284             Pc2 = new Geom2d_TrimmedCurve(Pc2,Uf,Ul);
3285             //is it necesary to invert ?
3286             Standard_Real DistDebToDeb = ptestdeb.Distance(pdeb);
3287             Standard_Real DistDebToFin = ptestdeb.Distance(pfin);
3288             Standard_Real DistFinToFin = ptestfin.Distance(pfin);
3289             Standard_Real DistFinToDeb = ptestfin.Distance(pdeb);
3290             
3291             if (DistDebToDeb > DistDebToFin &&
3292                 DistFinToFin > DistFinToDeb)
3293             {
3294               C3d->Reverse();
3295               Pc1->Reverse();
3296               Pc2->Reverse();
3297               ptestdeb = C3d->Value(C3d->FirstParameter());
3298               ptestfin = C3d->Value(C3d->LastParameter());
3299               DistDebToDeb = ptestdeb.Distance(pdeb);
3300               DistFinToFin = ptestfin.Distance(pfin);
3301             }
3302             if(DistDebToDeb < distref && DistFinToFin < distref)
3303             {
3304               Uf = C3d->FirstParameter();
3305               Ul = C3d->LastParameter();
3306               ChFi3d_ReparamPcurv(Uf,Ul,Pc1);
3307               ChFi3d_ReparamPcurv(Uf,Ul,Pc2);
3308               Standard_Real x,y;
3309               Pc1->Value(Uf).Coord(x,y);
3310               x = Pardeb(1) - x;
3311               y = Pardeb(2) - y;
3312               if(Abs(x) > tol2d || Abs(y) > tol2d) Pc1->Translate(gp_Vec2d(x,y));
3313               Pc2->Value(Uf).Coord(x,y);
3314               x = Pardeb(3) - x;
3315               y = Pardeb(4) - y;
3316               if(Abs(x) > tol2d || Abs(y) > tol2d) Pc2->Translate(gp_Vec2d(x,y));
3317               tolreached = ChFi3d_EvalTolReached(S1,Pc1,S2,Pc2,C3d);
3318               return Standard_True;
3319             }
3320           }
3321         }
3322       }
3323     }
3324   }
3325
3326   // At this stage : 
3327   // classic intersections have failed, the path is approached in vain.
3328
3329   Standard_Real Step = 0.1;
3330   for(;;) {
3331     //Attention the parameters of arrow for the path and
3332     //the tolerance for the approximation can't be taken as those of the  
3333     //Builder, so they are reestimated as much as possible.
3334     Standard_Real fleche = 1.e-3 * pdeb.Distance(pfin);
3335     Standard_Real tolap = 1.e-7;
3336     IntWalk_PWalking
3337       IntKK(S1,S2,tol3d,tol3d,fleche,Step);
3338
3339     //The extremities of the intersection (Pardeb,Parfin) are known,
3340     //one tries to find the start point at the 
3341     //middle to avoid obstacles on the path.
3342     Standard_Boolean depok = Standard_False;
3343     IntSurf_PntOn2S pintdep;
3344     TColStd_Array1OfReal depart(1,4);
3345     for(Standard_Integer ipdep = 2; ipdep <= 7 && !depok; ipdep++) {
3346       Standard_Real alpha = 0.1 * ipdep;
3347       Standard_Real unmoinsalpha = 1. - alpha;
3348       depart(1) = alpha*Pardeb(1) + unmoinsalpha*Parfin(1);
3349       depart(2) = alpha*Pardeb(2) + unmoinsalpha*Parfin(2);
3350       depart(3) = alpha*Pardeb(3) + unmoinsalpha*Parfin(3);
3351       depart(4) = alpha*Pardeb(4) + unmoinsalpha*Parfin(4);
3352       depok = IntKK.PerformFirstPoint(depart,pintdep);
3353     } 
3354     if(!depok) {
3355       return Standard_False;
3356     }
3357     pintdep.Parameters(depart(1),depart(2),depart(3),depart(4));
3358     IntKK.Perform(depart);
3359     if (!IntKK.IsDone()) return Standard_False;
3360     if (IntKK.NbPoints() <= 30) {
3361       Step *= 0.5;
3362       if (Step <= 0.0001) {
3363         return Standard_False;
3364       }
3365     }
3366     else{
3367       // At this stage there is a presentable LineOn2S, it is truncated  
3368       // between the points closest to known  extremites 
3369       // in fact there is a WLine and the approximation is launched.
3370       // Then the result is corrected to get proper start and end points.
3371       const Handle(IntSurf_LineOn2S)& L2S = IntKK.Line();
3372
3373       gp_Pnt codeb1 = S1->Value(Pardeb(1),Pardeb(2));
3374       gp_Pnt codeb2 = S2->Value(Pardeb(3),Pardeb(4));
3375       Standard_Real tol1 = Max(codeb1.Distance(codeb2),tol3d);
3376       Standard_Boolean bondeb = (tol1 == tol3d); 
3377       gp_Pnt pntd(0.5*(codeb1.Coord() + codeb2.Coord()));
3378
3379       gp_Pnt cofin1 = S1->Value(Parfin(1),Parfin(2));
3380       gp_Pnt cofin2 = S2->Value(Parfin(3),Parfin(4));
3381       Standard_Real tol2 = Max(cofin1.Distance(cofin2),tol3d);
3382       Standard_Boolean bonfin = (tol2 == tol3d); 
3383       gp_Pnt pntf(0.5*(cofin1.Coord() + cofin2.Coord()));
3384
3385       Standard_Integer nbp = L2S->NbPoints(), i;
3386       Standard_Real ddeb = Precision::Infinite(); 
3387       Standard_Real dfin = Precision::Infinite();
3388       Standard_Real dd;
3389       Standard_Integer indd = 0, indf = 0;
3390       for(i = 1; i <= nbp; i++) {
3391         dd = L2S->Value(i).Value().Distance(pntd);
3392         if(dd <= ddeb) { ddeb = dd; indd = i;}
3393         dd = L2S->Value(i).Value().Distance(pntf);
3394         if(dd < dfin) { dfin = dd; indf = i;}
3395       }
3396       if(indd > indf) {
3397         L2S->Reverse();
3398         indd = nbp - indd + 1;
3399         indf = nbp - indf + 1;
3400       }
3401       for (i = 1; i < indd; i++) { L2S->RemovePoint(1); nbp--; indf--; }
3402       for (i = indf + 1; i <= nbp; i++) { L2S->RemovePoint(indf + 1); }
3403       nbp = indf;
3404       if(nbp==1) return Standard_False;
3405       //The extremities are inserted in the line if the extremity points on it 
3406       //are too far and if pardeb and parfin are good.
3407       if(ddeb >= tol3d && bondeb) {
3408         IntSurf_PntOn2S p1 = L2S->Value(1);
3409         IntSurf_PntOn2S p2 = L2S->Value(2);
3410
3411         gp_Vec v1(pntd,p1.Value());
3412         gp_Vec v2(p1.Value(),p2.Value());
3413         gp_Vec v3(pntd,p2.Value());
3414         p1.SetValue(pntd,Pardeb(1),Pardeb(2),Pardeb(3),Pardeb(4));
3415         if(v1.Dot(v3) < 0) {
3416           if(v3.Magnitude() < 0.2*v2.Magnitude()) {
3417             L2S->RemovePoint(1);
3418             nbp--;
3419           }
3420           L2S->Value(1,p1);
3421         }
3422         else if(v1.Magnitude() > 0.2*v2.Magnitude()) {
3423           L2S->InsertBefore(1,p1);
3424           nbp++;
3425         }
3426         else{
3427           L2S->Value(1,p1);
3428         }
3429         ddeb = 0.;
3430       }
3431       if(dfin >= tol3d && bonfin) {
3432         IntSurf_PntOn2S p1 = L2S->Value(nbp);
3433         IntSurf_PntOn2S p2 = L2S->Value(nbp - 1);
3434         gp_Vec v1(pntf,p1.Value());
3435         gp_Vec v2(p1.Value(),p2.Value());
3436         gp_Vec v3(pntf,p2.Value());
3437         p1.SetValue(pntf,Parfin(1),Parfin(2),Parfin(3),Parfin(4));
3438         if(v1.Dot(v3) < 0) {
3439           if(v3.Magnitude() < 0.2*v2.Magnitude()) {
3440             L2S->RemovePoint(nbp);
3441             nbp--;
3442           }
3443           L2S->Value(nbp,p1);
3444         }
3445         else if(v1.Magnitude() > 0.2*v2.Magnitude()) {
3446           L2S->Add(p1);
3447           nbp++;
3448         }
3449         else{
3450           L2S->Value(nbp,p1);
3451         }
3452         dfin = 0.;
3453       }      
3454       //
3455       Handle(IntPatch_WLine)    WL = new IntPatch_WLine(L2S,Standard_False);
3456
3457 #ifdef OCCT_DEBUG
3458       //WL->Dump(0);
3459 #endif
3460
3461       GeomInt_WLApprox approx;
3462       approx.SetParameters(tolap, tol2d, 4, 8, 0, 30, Standard_True);
3463       // manage here the approximations that are not useful on planes!
3464       approx.Perform(S1,S2,WL,
3465         Standard_True,Standard_True,Standard_True,
3466         1,nbp);
3467       if(!approx.IsDone()) return Standard_False;
3468       //      tolreached = approx.TolReached3d();
3469       //      Standard_Real tolr2d = approx.TolReached2d();
3470       //      tolreached = Max(tolreached,ChFi3d_ConvTol2dToTol3d(S1,tolr2d));
3471       //      tolreached = Max(tolreached,ChFi3d_ConvTol2dToTol3d(S2,tolr2d));
3472       const AppParCurves_MultiBSpCurve& mbs = approx.Value(1);
3473       Standard_Integer nbpol = mbs.NbPoles();
3474       TColgp_Array1OfPnt pol3d(1,nbpol);
3475       mbs.Curve(1,pol3d);
3476       TColgp_Array1OfPnt2d pol2d1(1,nbpol);
3477       mbs.Curve(2,pol2d1);
3478       TColgp_Array1OfPnt2d pol2d2(1,nbpol);
3479       mbs.Curve(3,pol2d2);
3480       // The extremities of the intersection are reset on known points.
3481       if(ddeb >= tol1) {
3482         pol3d(1) = pntd;
3483         pol2d1(1).SetCoord(Pardeb(1),Pardeb(2));
3484         pol2d2(1).SetCoord(Pardeb(3),Pardeb(4));
3485         //      tolreached = Max(tolreached,ddeb);
3486       }
3487
3488       if(dfin >= tol2) {
3489         pol3d(nbpol) = pntf;
3490         pol2d1(nbpol).SetCoord(Parfin(1),Parfin(2));
3491         pol2d2(nbpol).SetCoord(Parfin(3),Parfin(4));
3492         //      tolreached = Max(tolreached,dfin);
3493       }
3494       const TColStd_Array1OfReal& knots = mbs.Knots();
3495       const TColStd_Array1OfInteger& mults = mbs.Multiplicities();
3496       Standard_Integer deg = mbs.Degree();
3497       C3d = new Geom_BSplineCurve(pol3d,knots,mults,deg);
3498       Pc1 = new Geom2d_BSplineCurve(pol2d1,knots,mults,deg);
3499       Pc2 = new Geom2d_BSplineCurve(pol2d2,knots,mults,deg);
3500       tolreached = ChFi3d_EvalTolReached(S1,Pc1,S2,Pc2,C3d);
3501       tolreached = Max(tolreached,ddeb);
3502       tolreached = Max(tolreached,dfin);
3503       return Standard_True;
3504     }
3505   }
3506 }
3507
3508 //=======================================================================
3509 //function : IntCS
3510 //purpose  : Fast calculation of the intersection curve surface.
3511 //
3512 //=======================================================================
3513
3514 Standard_Boolean ChFi3d_IntCS(const Handle(Adaptor3d_HSurface)& S,
3515                               const Handle(Adaptor3d_HCurve)& C,
3516                               gp_Pnt2d& p2dS,
3517                               Standard_Real& wc)
3518 {
3519   IntCurveSurface_HInter Intersection;
3520
3521   Standard_Real uf = C->FirstParameter(), ul = C->LastParameter();
3522   Standard_Real u1 = S->FirstUParameter(), u2 = S->LastUParameter();
3523   Standard_Real v1 = S->FirstVParameter(), v2 = S->LastVParameter();
3524   IntCurveSurface_IntersectionPoint pint;
3525   Intersection.Perform(C,S);
3526   Standard_Boolean keepfirst = (wc < -1.e100), keeplast = (wc > 1.e100);
3527   Standard_Real temp = 0.;
3528   if(keepfirst) temp = 1.e100;
3529   if(keeplast) temp = -1.e100;
3530   Standard_Real dist = 2.e100;
3531   if(Intersection.IsDone()) {
3532     Standard_Integer nbp = Intersection.NbPoints(),i,isol = 0;
3533     for (i = 1; i <= nbp; i++) {
3534       pint = Intersection.Point(i);
3535       Standard_Real up = pint.U();
3536       Standard_Real vp = pint.V();
3537       if(S->IsUPeriodic()) up = ChFi3d_InPeriod(up,u1,u1+S->UPeriod(),1.e-8);
3538       if(S->IsVPeriodic()) vp = ChFi3d_InPeriod(vp,v1,v1+S->VPeriod(),1.e-8);
3539       if(uf <= pint.W() && ul >= pint.W() &&
3540         u1 <= up && u2 >= up &&
3541         v1 <= vp && v2 >= vp) {
3542           if(keepfirst && pint.W() < temp) {
3543             temp = pint.W();
3544             isol = i;
3545           }
3546           else if(keeplast && pint.W() > temp) {
3547             temp = pint.W();
3548             isol = i;
3549           }
3550           else if(Abs(pint.W() - wc) < dist) {
3551             dist = Abs(pint.W() - wc);
3552             isol = i;
3553           }
3554       }
3555     }
3556     if(isol == 0) return Standard_False;
3557     pint = Intersection.Point(isol);
3558     Standard_Real up = pint.U();
3559     Standard_Real vp = pint.V();
3560     if(S->IsUPeriodic()) up = ChFi3d_InPeriod(up,u1,u1+S->UPeriod(),1.e-8);
3561     if(S->IsVPeriodic()) vp = ChFi3d_InPeriod(vp,v1,v1+S->VPeriod(),1.e-8);
3562     p2dS.SetCoord(up,vp);
3563     wc = pint.W();
3564     return Standard_True;
3565   }
3566   return Standard_False;
3567 }
3568
3569 //=======================================================================
3570 //function : ComputesIntPC
3571 //purpose  : Intersection of two PCurves of type FaceInterference
3572 //           the parameters of the pcurves at the solution point are 
3573 //           UInt1,UInt2
3574 //=======================================================================
3575
3576 void ChFi3d_ComputesIntPC (const ChFiDS_FaceInterference&      Fi1,
3577   const ChFiDS_FaceInterference&      Fi2,
3578   const Handle(GeomAdaptor_HSurface)& HS1,
3579   const Handle(GeomAdaptor_HSurface)& HS2,
3580   Standard_Real&                      UInt1, 
3581   Standard_Real&                      UInt2)
3582 {
3583   gp_Pnt bid;
3584   ChFi3d_ComputesIntPC(Fi1,Fi2,HS1,HS2,UInt1,UInt2,bid);
3585 }
3586
3587 //=======================================================================
3588 //function : ChFi3d_ComputesIntPC
3589 //purpose  : 
3590 //=======================================================================
3591 void ChFi3d_ComputesIntPC (const ChFiDS_FaceInterference&      Fi1,
3592   const ChFiDS_FaceInterference&      Fi2,
3593   const Handle(GeomAdaptor_HSurface)& HS1,
3594   const Handle(GeomAdaptor_HSurface)& HS2,
3595   Standard_Real&                      UInt1, 
3596   Standard_Real&                      UInt2,
3597   gp_Pnt&                             P)
3598 {    
3599   // Only one intersection to be carried out, however, the effort
3600   // is taken to check the extremities by an extrema c3d/c3d
3601   // created on pcurveonsurf of fillets.
3602
3603   Standard_Real x,y,distref2;
3604   Fi1.PCurveOnSurf()->Value(UInt1).Coord(x,y);
3605   gp_Pnt p3d1 = HS1->Value(x,y);
3606   Fi2.PCurveOnSurf()->Value(UInt2).Coord(x,y);
3607   gp_Pnt p3d2 = HS2->Value(x,y);
3608   distref2 = p3d1.SquareDistance(p3d2);
3609   P.SetXYZ(0.5*(p3d1.XYZ() + p3d2.XYZ()));
3610   // recalculation of the extremums
3611   Standard_Real delt1 = 
3612     Min(0.1,0.05*(Fi1.LastParameter() - Fi1.FirstParameter()));
3613   Handle(Geom2dAdaptor_HCurve) hc2d1 = 
3614     new Geom2dAdaptor_HCurve(Fi1.PCurveOnSurf(),UInt1-delt1,UInt1+delt1);
3615   Adaptor3d_CurveOnSurface cons1(hc2d1,HS1);
3616   Standard_Real delt2 = 
3617     Min(0.1,0.05*(Fi2.LastParameter() - Fi2.FirstParameter()));
3618   Handle(Geom2dAdaptor_HCurve) hc2d2 = 
3619     new Geom2dAdaptor_HCurve(Fi2.PCurveOnSurf(),UInt2-delt2,UInt2+delt2);
3620   Adaptor3d_CurveOnSurface cons2(hc2d2,HS2);
3621   Extrema_LocateExtCC ext(cons1,cons2,UInt1,UInt2);
3622   if(ext.IsDone()) {
3623     Standard_Real dist2 = ext.SquareDistance();
3624     if(dist2<distref2) {
3625       Extrema_POnCurv ponc1,ponc2;
3626       ext.Point(ponc1,ponc2);
3627       UInt1 = ponc1.Parameter();
3628       UInt2 = ponc2.Parameter();
3629       gp_Pnt Pnt1 = ponc1.Value();
3630       gp_Pnt Pnt2 = ponc2.Value();
3631       P.SetXYZ(0.5*(Pnt1.XYZ() + Pnt2.XYZ()));
3632     }
3633   }
3634
3635
3636 //=======================================================================
3637 //function : BoundSurf
3638 //purpose  : computes a GeomAdaptor_Surface from the surface of the 
3639 //           SurfData Fd1 and trims it to allow the intersection computation
3640
3641 //=======================================================================
3642 Handle(GeomAdaptor_HSurface) ChFi3d_BoundSurf(TopOpeBRepDS_DataStructure&    DStr,
3643   const Handle(ChFiDS_SurfData)& Fd1,
3644   const Standard_Integer&        IFaCo1,
3645   const Standard_Integer&        IFaArc1)
3646 {
3647   //rmq : as in fact 2 interferences of Fd1 serve only to set limits 
3648   //      indexes IFaCo1 and IFaArc1 are not useful.
3649   //      They are preserver here as an option in case it will be necessary to set 
3650   //      more restrictive limits (with intersection points as additional argument).
3651
3652   Handle(GeomAdaptor_HSurface) HS1 = new GeomAdaptor_HSurface();
3653   GeomAdaptor_Surface& S1 = HS1->ChangeSurface();
3654   S1.Load(DStr.Surface(Fd1->Surf()).Surface());
3655
3656   if ((IFaCo1 == 0)||(IFaArc1 == 0)) 
3657     return HS1;
3658
3659   const ChFiDS_FaceInterference& FiCo1 = Fd1->Interference(IFaCo1);
3660   const ChFiDS_FaceInterference& FiArc1 = Fd1->Interference(IFaArc1);
3661
3662   Standard_Real Du,Dv,mu,Mu,mv,Mv;
3663   gp_Pnt2d UVf1,UVf2,UVl1,UVl2;
3664
3665   UVf1 = FiCo1.PCurveOnSurf()->Value(FiCo1.FirstParameter());
3666   UVl1 = FiCo1.PCurveOnSurf()->Value(FiCo1.LastParameter());
3667   UVf2 = FiArc1.PCurveOnSurf()->Value(FiArc1.FirstParameter());
3668   UVl2 = FiArc1.PCurveOnSurf()->Value(FiArc1.LastParameter());
3669   ChFi3d_Boite(UVf1,UVf2,UVl1,UVl2,Du,Dv,mu,Mu,mv,Mv);
3670   GeomAbs_SurfaceType styp = S1.GetType();
3671   if (styp == GeomAbs_Cylinder) { 
3672     Dv = Max(0.5*Dv,4.*S1.Cylinder().Radius());
3673     Du = 0.;
3674     S1.Load(DStr.Surface(Fd1->Surf()).Surface(),
3675       mu,Mu,mv-Dv,Mv+Dv);
3676   }
3677   //In the case of a torus or cone, it is not necessary that the bounds create a surface with period more than 2PI. 
3678   else if (styp == GeomAbs_Torus ||
3679     styp == GeomAbs_Cone) {
3680       Du = Min(M_PI-0.5*Du,0.1*Du);
3681       Dv = 0.;
3682       S1.Load(DStr.Surface(Fd1->Surf()).Surface(),
3683         mu-Du,Mu+Du,mv,Mv);
3684   }
3685   else if (styp == GeomAbs_Plane) {
3686     Du = Max(0.5*Du,4.*Dv);
3687     Dv = 0.;
3688     S1.Load(DStr.Surface(Fd1->Surf()).Surface(),
3689       mu-Du,Mu+Du,mv,Mv);
3690   }
3691   return HS1;
3692 }
3693 //=======================================================================
3694 //function : SearchPivot
3695 //purpose  : 
3696 //=======================================================================
3697 Standard_Integer ChFi3d_SearchPivot(Standard_Integer* s,
3698   Standard_Real u[3][3],
3699   const Standard_Real t)
3700 {
3701   //           This function finds as pivot a cd the sections which of
3702   //           do not cross on the opposite face. 
3703   //         - probably there will be cases asymmetric to the point that
3704   //           none of tree fillets will match! To be SEEN.
3705   //         - in case when several fillets match the 
3706   //           first one taken is not inevitably the best 
3707   //           it should be refined by comparing the parameters on 
3708   //           guide lines and (/or) radiuses.
3709
3710   Standard_Boolean bondeb,bonfin;
3711   for(Standard_Integer i = 0; i <= 2; i++) {
3712     if(s[(i+1)%3] == 1) {bondeb = (u[(i+1)%3][i]-u[(i+1)%3][(i+2)%3] >= -t);}
3713     else {bondeb = (u[(i+1)%3][i]-u[(i+1)%3][(i+2)%3] <= t);}
3714     if(s[(i+2)%3] == 1) {bonfin = (u[(i+2)%3][i]-u[(i+2)%3][(i+1)%3] >= -t);}
3715     else {bonfin = (u[(i+2)%3][i]-u[(i+2)%3][(i+1)%3] <= t);}
3716     if (bondeb && bonfin) { return i; }
3717   }
3718   return -1;
3719 }
3720
3721
3722
3723 //=======================================================================
3724 //function : SearchFD
3725 //purpose  : 
3726 //=======================================================================
3727 Standard_Boolean ChFi3d_SearchFD(TopOpeBRepDS_DataStructure& DStr,
3728   const Handle(ChFiDS_Stripe)& cd1, 
3729   const Handle(ChFiDS_Stripe)& cd2,
3730   const Standard_Integer sens1,
3731   const Standard_Integer sens2,
3732   Standard_Integer& i1,
3733   Standard_Integer& i2,
3734   Standard_Real& p1,
3735   Standard_Real& p2,
3736   const Standard_Integer ind1,
3737   const Standard_Integer ind2,
3738   TopoDS_Face& face,
3739   Standard_Boolean& sameside,
3740   Standard_Integer& jf1,
3741   Standard_Integer& jf2)
3742 {
3743   Standard_Boolean found = Standard_False;
3744   Standard_Integer id1 = ind1, id2 = ind2;
3745   Standard_Integer if1 = ind1, if2 = ind2;
3746   Standard_Integer l1 = cd1->SetOfSurfData()->Length();
3747   Standard_Integer l2 = cd2->SetOfSurfData()->Length();
3748   Standard_Integer i;
3749   Standard_Boolean fini1 = Standard_False, fini2 = Standard_False;
3750   Standard_Boolean visavis,visavisok = Standard_False;
3751   TopoDS_Vertex Vtx;
3752   while( !found ) {
3753     for(i = id1; (i*sens1) <= (if1*sens1) && !found && !fini2; i = i+sens1 ) { 
3754       if(ChFi3d_IsInFront(DStr,cd1,cd2,i,if2,sens1,sens2,p1,p2,face,sameside,jf1,jf2,visavis,Vtx,Standard_False,0)) {
3755         i1 = i;
3756         i2 = if2;
3757         found = Standard_True;
3758       }
3759       else if (visavis && !visavisok) {
3760         visavisok = Standard_True;
3761         i1 = i;
3762         i2 = if2;
3763       }
3764     }
3765     if(!fini1) {
3766       if1 = if1 + sens1;
3767       if(if1 < 1 || if1 > l1) { if1 = if1 - sens1; fini1 = Standard_True; }
3768     }
3769
3770     for(i = id2; (i*sens2) <= (if2*sens2) && !found && !fini1; i = i+sens2 ) { 
3771       if(ChFi3d_IsInFront(DStr,cd1,cd2,if1,i,sens1,sens2,p1,p2,face,sameside,jf1,jf2,visavis,Vtx,Standard_False,0)) {
3772         i1 = if1;
3773         i2 = i;
3774         found = Standard_True;
3775       }
3776       else if (visavis && !visavisok) {
3777         visavisok = Standard_True;
3778         i1 = if1;
3779         i2 = i;
3780       }
3781     }
3782     if(!fini2) {
3783       if2 = if2 + sens2;
3784       if(if2 < 1 || if2 > l2) { if2 = if2 - sens2; fini2 = Standard_True; }
3785     }
3786     if(fini1 && fini2) break;
3787   }
3788   return found;
3789 }
3790
3791 //=======================================================================
3792 //function : Parameters
3793 //purpose  : compute the parameters <u> and <v> of the 3d point <p3d> 
3794 //           on the surface <S> if it's an analytic surface
3795 //=======================================================================
3796
3797 void ChFi3d_Parameters(const Handle(Geom_Surface)& S,
3798   const gp_Pnt& p3d,
3799   Standard_Real& u,
3800   Standard_Real& v)
3801 {
3802   GeomAdaptor_Surface gas(S);
3803   switch ( gas.GetType() ) {
3804   case GeomAbs_Plane :
3805     ElSLib::Parameters(gas.Plane(),p3d,u,v);
3806     break;
3807   case GeomAbs_Cylinder :
3808     ElSLib::Parameters(gas.Cylinder(),p3d,u,v);
3809     break;
3810   case GeomAbs_Cone :
3811     ElSLib::Parameters(gas.Cone(),p3d,u,v);
3812     break;
3813   case GeomAbs_Sphere :
3814     ElSLib::Parameters(gas.Sphere(),p3d,u,v);
3815     break;
3816   case GeomAbs_Torus :
3817     ElSLib::Parameters(gas.Torus(),p3d,u,v);
3818     break;
3819   case GeomAbs_BezierSurface :
3820   case GeomAbs_BSplineSurface :  
3821   default :
3822     {
3823       GeomAPI_ProjectPointOnSurf tool(p3d,S);
3824       if ( tool.NbPoints() != 1 )
3825         throw StdFail_NotDone ("ChFi3d_Parameters() - no projection results");
3826       else
3827         tool.Parameters(1,u,v);
3828     }
3829   }
3830 }
3831
3832 //=======================================================================
3833 //function : TrimCurve
3834 //purpose  : trims the curve <gc> between the points <FirstP> and
3835 //           <LastP>. The trimmed curve is <gtc>
3836 //=======================================================================
3837
3838 void ChFi3d_TrimCurve(const Handle(Geom_Curve)& gc,
3839   const gp_Pnt& FirstP,
3840   const gp_Pnt& LastP,
3841   Handle(Geom_TrimmedCurve)& gtc)
3842 {
3843   Standard_Real uf = 0.,ul = 0.;
3844   GeomAdaptor_Curve gac(gc);
3845   switch ( gac.GetType() ) {
3846   case GeomAbs_Line :
3847     {
3848       uf = ElCLib::Parameter(gac.Line(),FirstP);
3849       ul = ElCLib::Parameter(gac.Line(),LastP);
3850     }
3851     break;
3852   case GeomAbs_Circle :
3853     {
3854       uf = ElCLib::Parameter(gac.Circle(),FirstP);
3855       ul = ElCLib::Parameter(gac.Circle(),LastP);
3856     }
3857     break;
3858   case GeomAbs_Ellipse :
3859     {
3860       uf = ElCLib::Parameter(gac.Ellipse(),FirstP);
3861       ul = ElCLib::Parameter(gac.Ellipse(),LastP);
3862     }
3863     break;
3864   case GeomAbs_Hyperbola :
3865     {
3866       uf = ElCLib::Parameter(gac.Hyperbola(),FirstP);
3867       ul = ElCLib::Parameter(gac.Hyperbola(),LastP);
3868     }
3869     break;
3870   case GeomAbs_Parabola :
3871     {
3872       uf = ElCLib::Parameter(gac.Parabola(),FirstP);
3873       ul = ElCLib::Parameter(gac.Parabola(),LastP);
3874     }
3875     break;
3876   default :
3877     {
3878       GeomAPI_ProjectPointOnCurve tool(FirstP,gc);
3879       if ( tool.NbPoints() != 1 )
3880         throw StdFail_NotDone ("ChFi3d_TrimCurve() - no projection results for the first point");
3881       else
3882         uf = tool.Parameter(1);
3883       tool.Init(LastP,gc);
3884       if ( tool.NbPoints() != 1 )
3885         throw StdFail_NotDone ("ChFi3d_TrimCurve() - no projection results for the second point");
3886       else
3887         ul = tool.Parameter(1);
3888     }
3889   }
3890   gtc = new Geom_TrimmedCurve(gc,uf,ul);
3891 }
3892
3893
3894
3895 //=======================================================================
3896 //function : GoodExt
3897 //purpose  : 
3898 //=======================================================================
3899 static Standard_Boolean GoodExt(const Handle(Geom_Curve)& C,
3900   const gp_Vec&             V,
3901   const Standard_Real       f,
3902   const Standard_Real       l,
3903   const Standard_Real       a)
3904 {
3905   for(Standard_Integer i = 0; i < 6; i++) {
3906     gp_Pnt d0; gp_Vec d1;
3907     const Standard_Real t = i * 0.2;
3908     C->D1(((1-t)*f+t*l),d0,d1);
3909     const Standard_Real ang = d1.Angle(V);
3910     const Standard_Real angref = a*t + 0.002;
3911     if(ang > angref) return Standard_False;
3912   }
3913   return Standard_True;
3914 }
3915 //=======================================================================
3916 //function : PerformElSpine
3917 //purpose  : 
3918 //=======================================================================
3919 Standard_EXPORT 
3920 void ChFi3d_PerformElSpine(Handle(ChFiDS_HElSpine)& HES,
3921                            Handle(ChFiDS_Spine)&    Spine,
3922                            const GeomAbs_Shape      continuity,
3923                            const Standard_Real      tol,
3924                            const Standard_Boolean   IsOffset)
3925 {
3926
3927   Standard_Boolean periodic, Bof, checkdeb, cepadur,bIsSmooth;
3928   Standard_Integer IEdge,IF,IL,nbed, iToApproxByC2;
3929   Standard_Real WF, WL, Wrefdeb, Wreffin,nwf,nwl,period,pared = 0.,tolpared;
3930   Standard_Real First, Last, epsV, urefdeb, tolrac;
3931   GeomAbs_Shape aContinuity;
3932   gp_Pnt PDeb, PFin, Bout;
3933   gp_Vec VrefDeb, VrefFin;
3934   Handle(Geom_Curve) Cv;
3935   Handle(Geom_BoundedCurve) TC;
3936   Handle(Geom_BSplineCurve) BS, BSpline;
3937   TopoDS_Edge E, Eold;
3938   TopoDS_Vertex V;
3939   //
3940   ChFiDS_ElSpine& ES = HES->ChangeCurve();
3941   WF = ES.FirstParameter();
3942   WL = ES.LastParameter();
3943   Wrefdeb = WF;
3944   Wreffin = WL;
3945   nwf = WF;
3946   nwl = WL;
3947   nbed = Spine->NbEdges();
3948   periodic = Spine->IsPeriodic();
3949   if(periodic) {
3950     period = Spine->Period();
3951     nwf = ElCLib::InPeriod(WF,-tol,period-tol);
3952     IF = Spine->Index(nwf,1);
3953     nwl = ElCLib::InPeriod(WL,tol,period+tol);
3954     IL = Spine->Index(nwl,0);
3955     if(nwl<nwf+tol) IL += nbed;
3956   }
3957   else{
3958     IF = Spine->Index(WF,1);
3959     IL = Spine->Index(WL,0);
3960     Wrefdeb = Max(Spine->FirstParameter(IF),WF);
3961     Wreffin = Min(Spine->LastParameter(IL),WL);
3962   }
3963   //
3964   Spine->D1(WF,PDeb,VrefDeb);
3965   Spine->D1(WL,PFin,VrefFin);
3966   VrefDeb.Normalize();
3967   VrefFin.Normalize();
3968   //
3969   TColgp_Array1OfPnt ExtrapPole(1, 5);
3970   TColgp_Array1OfPnt ExtraCoeffs(1, 5);
3971   TColgp_Array1OfXYZ Cont(1,5);
3972   // Attention on segmente eventuellement la premiere et la
3973   // derniere arete.
3974   // Traitment de la premiere arete
3975   cepadur = 0;
3976   E = (IsOffset)? Spine->OffsetEdges(IF) : Spine->Edges(IF);
3977   Bof = BRepLib::BuildCurve3d(E);
3978   const BRepAdaptor_Curve& edc = Spine->CurrentElementarySpine(IF);
3979   tolpared = edc.Resolution(tol);
3980   Cv = BRep_Tool::Curve(E, First, Last);
3981   //Add vertex with tangent
3982   if (ES.IsPeriodic())
3983   {
3984     Standard_Real ParForElSpine = (E.Orientation() == TopAbs_FORWARD)? First : Last;
3985     gp_Pnt PntForElSpine;
3986     gp_Vec DirForElSpine;
3987     Cv->D1(ParForElSpine, PntForElSpine, DirForElSpine);
3988     ES.AddVertexWithTangent(gp_Ax1(PntForElSpine, DirForElSpine));
3989   }
3990   /////////////////////////
3991   urefdeb = Spine->FirstParameter(IF);
3992   checkdeb = (nwf > urefdeb);
3993   if(checkdeb) {
3994     Spine->Parameter(IF,nwf,pared,0);
3995   }
3996   //
3997   if(E.Orientation() == TopAbs_REVERSED) {
3998     Standard_Real sov = First;
3999     First = Cv->ReversedParameter(Last);
4000     Last = Cv->ReversedParameter(sov);
4001     if(checkdeb) {
4002       pared = Cv->ReversedParameter(pared);
4003     }
4004     else{
4005       pared = First;
4006     }
4007     if(First < pared) {
4008       First = pared; 
4009     }
4010     if(IL == IF) {
4011       Standard_Real ureffin = Spine->LastParameter(IL);
4012       Standard_Boolean checkfin = (nwl < ureffin);
4013       if(checkfin) {
4014         Spine->Parameter(IL,nwl,pared,0);
4015         pared = Cv->ReversedParameter(pared);
4016       }
4017       else {
4018         pared = Last;
4019       }
4020       if(pared < Last) {
4021         Last = pared; 
4022       }
4023     }
4024     Cv = Cv->Reversed();
4025   }//if(E.Orientation() == TopAbs_REVERSED) 
4026   else {//#1
4027     if(!checkdeb) {
4028       pared = First;
4029     }
4030     if(First < pared) {
4031       First = pared; 
4032     }
4033     if(IL == IF) {
4034       Standard_Real ureffin = Spine->LastParameter(IL);
4035       Standard_Boolean checkfin = (nwl < ureffin);
4036       if(checkfin) {
4037         Spine->Parameter(IL,nwl,pared,0);
4038       }
4039       else {
4040         pared = Last;
4041       }
4042       if(pared < Last) {
4043         Last = pared; 
4044       }
4045     }
4046   }// else {//#1
4047   //
4048   if(Abs(Last-First) < tolpared) {
4049     cepadur = 1;
4050   }
4051   //
4052   //Petite veru pour les cas ou un KPart a bouffe l arete
4053   //sans parvenir a terminer. On tire une droite.
4054   if(cepadur) {
4055     Handle(Geom_Line) L;
4056     gp_Pnt ptemp; gp_Vec vtemp;
4057     if(WL < Spine->FirstParameter(1) + tol) {
4058       ES.LastPointAndTgt(ptemp,vtemp);
4059       gp_Dir d(vtemp);
4060       gp_Pnt olin;
4061       olin.ChangeCoord().SetLinearForm(-WL,d.XYZ(),PFin.XYZ());
4062       L = new Geom_Line(olin,d);
4063       ES.SetCurve(L);
4064     }
4065     else if(WF > Spine->LastParameter(nbed) - tol) {
4066       ES.FirstPointAndTgt(ptemp,vtemp);
4067       gp_Dir d(vtemp);
4068       gp_Pnt olin;
4069       olin.ChangeCoord().SetLinearForm(-WF,d.XYZ(),PDeb.XYZ());
4070       L = new Geom_Line(olin,d);
4071       ES.SetCurve(L);
4072     }
4073     return;// => 
4074   }
4075   //
4076   TC = new (Geom_TrimmedCurve)(Cv, First, Last);
4077   BS=GeomConvert::CurveToBSplineCurve(TC);
4078   CurveCleaner(BS, Abs(WL-WF)*1.e-4, 0);
4079   //
4080   //Smoothing of the curve
4081   iToApproxByC2=0;
4082   aContinuity=TC->Continuity();
4083   bIsSmooth=ChFi3d_IsSmooth(TC);
4084   if (aContinuity < GeomAbs_C2 && !bIsSmooth) {
4085     ++iToApproxByC2;
4086     BS = ChFi3d_ApproxByC2(TC);
4087     TC=BS;
4088   }
4089   //
4090   //  Concatenation des aretes suivantes
4091   GeomConvert_CompCurveToBSplineCurve Concat( TC, Convert_QuasiAngular );
4092   //
4093   Eold = E;
4094   for (IEdge=IF+1; IEdge<=IL; ++IEdge) {
4095     Standard_Integer iloc = IEdge;
4096     if(periodic) {
4097       iloc = (IEdge - 1)%nbed + 1;
4098     }
4099     //
4100     E = (IsOffset)? Spine->OffsetEdges(iloc) : Spine->Edges(iloc);
4101     if (BRep_Tool::Degenerated(E)) {
4102       continue;
4103     }
4104     //  
4105     epsV = tol;
4106     Bof = TopExp::CommonVertex(Eold, E, V);
4107     if (Bof) {
4108       epsV = BRep_Tool::Tolerance(V);
4109     }
4110     //
4111     Bof = BRepLib::BuildCurve3d(E);
4112     if (!Bof) {
4113       throw Standard_ConstructionError("PerformElSpine : BuildCurve3d error");
4114     }
4115     //
4116     Cv = BRep_Tool::Curve(E, First, Last);
4117     //Add vertex with tangent
4118     Standard_Real ParForElSpine = (E.Orientation() == TopAbs_FORWARD)? First : Last;
4119     gp_Pnt PntForElSpine;
4120     gp_Vec DirForElSpine;
4121     Cv->D1(ParForElSpine, PntForElSpine, DirForElSpine);
4122     ES.AddVertexWithTangent(gp_Ax1(PntForElSpine, DirForElSpine));
4123     /////////////////////////
4124     if(IEdge == IL) {
4125       Standard_Real ureffin = Spine->LastParameter(iloc);
4126       Standard_Boolean checkfin = (nwl < ureffin);
4127       if(checkfin) {
4128         Spine->Parameter(iloc,nwl,pared,0);
4129       }
4130       else {
4131         pared = Last;
4132       }
4133       if(E.Orientation() == TopAbs_REVERSED) {
4134         Standard_Real sov = First;
4135         First = Cv->ReversedParameter(Last);
4136         Last = Cv->ReversedParameter(sov);
4137         if(checkfin) {
4138           pared = Cv->ReversedParameter(pared);
4139         }
4140         else{
4141           pared = Last;
4142         }
4143         Cv = Cv->Reversed();
4144       }
4145       if(pared < Last) {
4146         Last = pared; 
4147       }
4148     }
4149     //
4150     TC = new (Geom_TrimmedCurve)(Cv, First, Last);
4151     BS = GeomConvert::CurveToBSplineCurve(TC);
4152     CurveCleaner(BS, Abs(WL-WF)*1.e-4, 0);
4153     //
4154     //Smoothing of the curve
4155     aContinuity=TC->Continuity();
4156     bIsSmooth=ChFi3d_IsSmooth(TC);
4157     if (aContinuity < GeomAbs_C2 && !bIsSmooth) {
4158       ++iToApproxByC2;
4159       BS = ChFi3d_ApproxByC2( TC );
4160       TC = BS;
4161     }
4162     //
4163     tolrac = Min(tol, epsV);
4164     Bof = Concat.Add( TC, 2.*tolrac, Standard_True );
4165     // si l'ajout ne s'est pas bien passe on essai d'augmenter la tolerance
4166     if (!Bof) {
4167       Bof = Concat.Add( TC, 2.*epsV, Standard_True );
4168     }
4169     if (!Bof) {
4170       Bof = Concat.Add( TC, 200.*epsV, Standard_True );
4171       if (!Bof) {
4172         throw Standard_ConstructionError("PerformElSpine: spine merged error");
4173       }
4174     }
4175     Eold = E;
4176   }// for (IEdge=IF+1; IEdge<=IL; ++IEdge) {
4177   //
4178   // On a la portion d elspine calculee sans prolongements sur la partie 
4179   // valide des aretes du chemin.
4180   BSpline = Concat.BSplineCurve();
4181   // There is a reparametrisation to maximally connect the abscissas of edges.
4182   TColStd_Array1OfReal BSNoeuds (1, BSpline->NbKnots());
4183   BSpline->Knots(BSNoeuds);
4184   BSplCLib::Reparametrize (Wrefdeb, Wreffin, BSNoeuds);
4185   BSpline->SetKnots(BSNoeuds);
4186   //
4187   // Traitement des Extremites
4188   Standard_Integer caredeb, carefin;
4189   Standard_Real LocalWL, LocalWF, Angle;
4190   GeomAdaptor_Curve gacurve;
4191   Handle(Geom_BSplineCurve) newc;
4192   //
4193   caredeb = 0;
4194   carefin = 0;
4195   Angle = M_PI*0.75;
4196   LocalWL = WL;
4197   LocalWF = WF;
4198   if (!ES.IsPeriodic() && !PDeb.IsEqual(BSpline->Pole(1), tol) ) {
4199     // Prolongement C3 au debut
4200     // afin d'eviter des pts d'inflexions dans la partie utile de la
4201     // spine le prolongement se fait jusqu'a un point eloigne.
4202     if(BSpline->IsRational()) {
4203       caredeb = 1;
4204     }
4205     //
4206     Standard_Real rabdist = Wrefdeb - WF;
4207     Bout = PDeb.Translated(-20*rabdist * VrefDeb);
4208     Standard_Boolean goodext = 0;
4209     for(Standard_Integer icont = 3; icont>=1 && !goodext; icont--) {
4210       Handle(Geom_BoundedCurve) anExtCurve = BSpline;
4211       GeomLib::ExtendCurveToPoint (anExtCurve, Bout, icont, Standard_False);
4212       newc = Handle(Geom_BSplineCurve)::DownCast (anExtCurve);
4213       gacurve.Load(newc);
4214       GCPnts_AbscissaPoint GCP(gacurve,-rabdist,Wrefdeb,WF);
4215       if(GCP.IsDone()) {
4216         WF = GCP.Parameter();
4217         goodext = GoodExt(newc,VrefDeb,Wrefdeb,WF,Angle);
4218       }
4219     }
4220     if(caredeb) {
4221       caredeb = newc->NbKnots() - BSpline->NbKnots();
4222     }
4223     BSpline = newc;
4224     LocalWF = BSpline->FirstParameter();
4225   }
4226   //
4227   if (!ES.IsPeriodic() && !PFin.IsEqual(BSpline->Pole(BSpline->NbPoles()), tol) ) {
4228     // Prolongement C3 en fin
4229     if(BSpline->IsRational()) {
4230       carefin = 1;
4231     }
4232     Standard_Real rabdist = WL - Wreffin;
4233     Bout = PFin.Translated(20*rabdist * VrefFin);
4234     Standard_Boolean goodext = 0;
4235     for(Standard_Integer icont = 3; icont>=1 && !goodext; icont--) {
4236       Handle(Geom_BoundedCurve) anExtCurve = BSpline;
4237       GeomLib::ExtendCurveToPoint (anExtCurve, Bout, icont, Standard_True);
4238       newc = Handle(Geom_BSplineCurve)::DownCast (anExtCurve);
4239       gacurve.Load(newc);
4240       GCPnts_AbscissaPoint GCP(gacurve,rabdist,Wreffin,WL);
4241       if(GCP.IsDone()) {
4242         WL = GCP.Parameter();
4243         goodext = GoodExt(newc, VrefFin, Wreffin,WL,Angle);
4244       }
4245     }
4246     if(carefin) {
4247       carefin = newc->NbKnots() - BSpline->NbKnots();
4248     }
4249     BSpline = newc;
4250     LocalWL = BSpline->LastParameter(); 
4251   }
4252   //
4253   //Reparametrisation et segmentation sur le domaine de la Spine.
4254   if(Abs(BSpline->FirstParameter() - WF)<tol) {
4255     WF = BSpline->FirstParameter();
4256   }
4257   if(Abs(BSpline->LastParameter() - WL)<tol) {
4258     WL = BSpline->LastParameter();
4259   }
4260   //
4261   if ( (LocalWF<WF) ||  (LocalWL>WL)) {   // pour eviter des pb avec segment!
4262     BSpline->Segment(WF, WL);
4263     ES.FirstParameter(WF);
4264     ES.LastParameter(WL);
4265   }
4266   //
4267   if (BSpline->IsRational()) {
4268     Handle(Geom_BSplineCurve) C1;
4269     C1 =  Handle(Geom_BSplineCurve)::DownCast(BSpline->Copy());
4270     GeomConvert::C0BSplineToC1BSplineCurve(C1, tol, 0.1);
4271     // Il faut s'assurer que l'origine n'a pas bouge (cts21158)
4272     if (C1->FirstParameter() == BSpline->FirstParameter()) {
4273       BSpline = C1;
4274     }
4275     else {
4276       //std::cout << "Attention : Echec de C0BSplineToC1 !" << std::endl;
4277     }
4278   }
4279   //
4280   Standard_Integer fk, lk, MultMax, ii;
4281   // Deformation eventuelle pour rendre la spine C2.
4282   // ou C3 pour des approx C2
4283   if((caredeb || carefin) && BSpline->Degree() < 8) {
4284     BSpline->IncreaseDegree(8);
4285   }
4286   //
4287   fk = 2;
4288   lk = BSpline->NbKnots()-1;
4289   if(BSpline->IsPeriodic()) {
4290     fk = 1;
4291   }
4292   if(caredeb) {
4293     fk += caredeb;
4294   }
4295   if(carefin) {
4296     lk -= carefin;
4297   }
4298   //
4299   if (continuity == GeomAbs_C3) {
4300     if (BSpline->Degree() < 7) {
4301       BSpline->IncreaseDegree(7);
4302     }
4303     MultMax = BSpline->Degree() - 3;
4304   }
4305   else {
4306     if (BSpline->Degree() < 5) {
4307       BSpline->IncreaseDegree(5);
4308     }
4309     MultMax = BSpline->Degree() - 2;
4310   }
4311   // correction C2 or C3 (if possible)
4312   CurveCleaner(BSpline, Abs(WL-WF)*1.e-4, 1);
4313   CurveCleaner(BSpline, Abs(WL-WF)*1.e-2, MultMax);
4314   Standard_Integer MultMin = Max(BSpline->Degree() - 4, 1);
4315   for (ii = fk; ii <= lk; ii++) {
4316     if( BSpline->Multiplicity(ii) > MultMax ) {
4317       Bof = BSpline->RemoveKnot(ii, MultMax, Abs(WL-WF)/10);
4318     }
4319     // See C4
4320     if( BSpline->Multiplicity(ii) > MultMin ) {
4321       Bof = BSpline->RemoveKnot(ii, MultMin, Abs(WL-WF)*1.e-4);
4322     }       
4323   }
4324   // elspine periodic => BSpline Periodic
4325   if(ES.IsPeriodic()) {
4326     if(!BSpline->IsPeriodic()) {
4327       BSpline->SetPeriodic();
4328       //modified by NIZNHY-PKV Fri Dec 10 12:20:22 2010ft
4329       if (iToApproxByC2) {
4330         Bof = BSpline->RemoveKnot(1, MultMax, Abs(WL-WF)/10);
4331       }
4332       //Bof = BSpline->RemoveKnot(1, MultMax, Abs(WL-WF)/10);
4333       //modified by NIZNHY-PKV Mon Dec 13 14:12:54 2010t
4334     }
4335   }
4336   else { 
4337     // Otherwise is it necessary to move the poles to adapt 
4338     // them to new tangents ?
4339     Standard_Boolean adjust = Standard_False; 
4340     gp_Pnt P1, P2;
4341     gp_Vec V1, V2;
4342     BSpline->D1(WF, P1, V1);
4343     V1.Normalize();
4344     ES.FirstPointAndTgt(PDeb,VrefDeb);
4345     Standard_Real scaldeb = VrefDeb.Dot(V1);
4346     Standard_Real disdeb = PDeb.Distance(P1);
4347     if((Abs(WF-LocalWF) < 1.e-12) &&
4348       ((scaldeb <= 0.9999999) ||
4349       disdeb >= tol)) {   
4350         // Yes if there was no extension and the tangent is not the good one.
4351         adjust = Standard_True;
4352     } 
4353     BSpline->D1(WL, P2, V2);
4354     V2.Normalize();
4355     ES.LastPointAndTgt(PFin,VrefFin);
4356     Standard_Real scalfin = VrefFin.Dot(V2); 
4357     Standard_Real disfin = PFin.Distance(P2);
4358     if((Abs(WL-LocalWL) < 1.e-12) && 
4359       ((scalfin <= 0.9999999)||
4360       disfin >= tol)) {
4361         // the same at the end
4362         adjust = Standard_True;
4363     }
4364     if(adjust) {
4365       Handle(Geom_BoundedCurve) anExtCurve = BSpline;
4366       GeomLib::AdjustExtremity(anExtCurve, PDeb, PFin, VrefDeb, VrefFin);
4367       BSpline = Handle(Geom_BSplineCurve)::DownCast (anExtCurve);
4368     }
4369   }
4370
4371   // Le Resultat       
4372   ES.SetCurve(BSpline);
4373
4374   //Temporary
4375   //gp_Pnt ptgui;
4376   //gp_Vec d1gui;
4377   //( HES->Curve() ).D1(HES->FirstParameter(),ptgui,d1gui);
4378 }
4379
4380 //=======================================================================
4381 //function : cherche_face1
4382 //purpose  : find face F different from F1 in the map.
4383 // The map contains two faces adjacent to an edge 
4384 //=======================================================================
4385 void ChFi3d_cherche_face1 (const TopTools_ListOfShape & map,
4386   const TopoDS_Face & F1,
4387   TopoDS_Face &  F)   
4388
4389   TopoDS_Face Fcur;
4390   Standard_Boolean trouve=Standard_False;
4391   TopTools_ListIteratorOfListOfShape It;
4392   for (It.Initialize(map);It.More()&&!trouve;It.Next()) { 
4393     Fcur=TopoDS::Face (It.Value());
4394     if (!Fcur.IsSame(F1)) {
4395       F=Fcur;trouve=Standard_True;}
4396   }
4397
4398 //=======================================================================
4399 //function : cherche_element
4400 //purpose  : find edge E of F1 other than E1 and containing vertex V
4401 // Vtx is the other vertex of E 
4402 //=======================================================================
4403 void ChFi3d_cherche_element(const TopoDS_Vertex & V,
4404   const TopoDS_Edge & E1,
4405   const TopoDS_Face & F1,
4406   TopoDS_Edge & E , 
4407   TopoDS_Vertex  & Vtx )
4408
4409   Standard_Integer ie; 
4410   TopoDS_Vertex V1,V2;
4411   Standard_Boolean trouve=Standard_False;
4412   TopoDS_Edge Ecur;
4413   TopTools_IndexedMapOfShape  MapE;
4414   TopExp::MapShapes( F1,TopAbs_EDGE,MapE);
4415   for ( ie=1; ie<= MapE.Extent()&& !trouve; ie++) {  
4416     Ecur = TopoDS::Edge (MapE(ie));
4417     if (!Ecur.IsSame(E1)) { 
4418       TopTools_IndexedMapOfShape  MapV;
4419       TopExp::MapShapes(Ecur, TopAbs_VERTEX, MapV);
4420       if (MapV.Extent()==2) {
4421         V1 = TopoDS::Vertex (MapV(1));
4422         V2 = TopoDS::Vertex (MapV(2));
4423         if (V1.IsSame(V)) { 
4424           Vtx=V2;
4425           E=Ecur;
4426           trouve=Standard_True;
4427         }
4428         else if (V2.IsSame(V)) {
4429           Vtx=V1;
4430           E=Ecur;
4431           trouve=Standard_True;
4432         }
4433       }
4434     }
4435   }
4436
4437 //=======================================================================
4438 //function : cherche_edge
4439 //purpose  : find edge E of F1 other than the list of edges E1 and
4440 //           containing vertex V  Vtx is the other vertex of E. 
4441 //=======================================================================
4442 void ChFi3d_cherche_edge(const TopoDS_Vertex & V,
4443   const  TopTools_Array1OfShape & E1,
4444   const TopoDS_Face & F1,
4445   TopoDS_Edge & E , 
4446   TopoDS_Vertex  & Vtx )
4447
4448   Standard_Integer ie,i; 
4449   TopoDS_Vertex V1,V2;
4450   Standard_Boolean trouve=Standard_False;
4451   TopoDS_Edge Ecur;
4452   Standard_Boolean same;
4453   TopTools_IndexedMapOfShape  MapE;
4454   TopExp::MapShapes( F1,TopAbs_EDGE,MapE);
4455   for ( ie=1; ie<= MapE.Extent()&& !trouve; ie++) {  
4456     Ecur=TopoDS::Edge (MapE(ie));
4457     same=Standard_False;
4458     for (i=E1.Lower();i<=E1.Upper() ;i++) {
4459       if (Ecur.IsSame(E1.Value(i))) same=Standard_True;
4460     }
4461     if (!same) {
4462       TopTools_IndexedMapOfShape  MapV;
4463       TopExp::MapShapes(Ecur, TopAbs_VERTEX, MapV);
4464       if (MapV.Extent()==2) {
4465         V1 = TopoDS::Vertex (MapV(1));
4466         V2 = TopoDS::Vertex (MapV(2));
4467         if (V1.IsSame(V)) { 
4468           Vtx=V2;
4469           E=Ecur;
4470           trouve=Standard_True;
4471         }
4472         else if (V2.IsSame(V)) {
4473           Vtx=V1;
4474           E=Ecur;
4475           trouve=Standard_True;
4476         }
4477       }
4478     }
4479   }
4480 }
4481
4482 //=======================================================================
4483 //function : nbface
4484 //purpose  : calculates the number of faces common to a vertex
4485 //           
4486 //=======================================================================
4487 Standard_Integer  ChFi3d_nbface (const TopTools_ListOfShape & mapVF )
4488 { Standard_Integer nface=0;     
4489 TopTools_ListIteratorOfListOfShape ItF,JtF;
4490 Standard_Integer  fj = 0;
4491 for (ItF.Initialize(mapVF); ItF.More(); ItF.Next()) {
4492   fj++;
4493   Standard_Integer kf = 1;
4494   const TopoDS_Shape& cur = ItF.Value();
4495   for (JtF.Initialize(mapVF); JtF.More( )&&(kf<fj); JtF.Next(), kf++) {
4496     if(cur.IsSame(JtF.Value())) break;
4497   }
4498   if(kf == fj) nface++;
4499 }
4500 return nface;  
4501 }
4502
4503 //=======================================================================
4504 //function : edge_common_faces 
4505 //purpose  :  determines two faces sharing an edge.
4506 //            F1 = F2 if there is an edge to parce 
4507 //=======================================================================
4508 void ChFi3d_edge_common_faces (const TopTools_ListOfShape & mapEF,
4509   TopoDS_Face & F1,
4510   TopoDS_Face &  F2)   
4511 { TopTools_ListIteratorOfListOfShape It;
4512 TopoDS_Face F;
4513 Standard_Boolean trouve;
4514 It.Initialize(mapEF);  
4515 F1=TopoDS::Face(It.Value());
4516 trouve=Standard_False;
4517 for(It.Initialize(mapEF);It.More()&&!trouve;It.Next()) { 
4518   F=TopoDS::Face (It.Value());
4519   if (!F.IsSame(F1)) {
4520     F2=F;trouve=Standard_True;
4521   }
4522 }
4523 if (!trouve) F2=F1;
4524 }
4525
4526 /***********************************************************/
4527 // gives the angle between edges E1 and E2 . Vtx is the 
4528 // vertex common to the edges
4529 /************************************************************/
4530 Standard_Real ChFi3d_AngleEdge (const TopoDS_Vertex & Vtx,
4531   const TopoDS_Edge&  E1,
4532   const TopoDS_Edge &  E2)
4533 {   Standard_Real angle;
4534 BRepAdaptor_Curve BCurv1(E1);
4535 BRepAdaptor_Curve BCurv2(E2);
4536 Standard_Real parE1,parE2;
4537 gp_Vec dir1,dir2 ;
4538 gp_Pnt P1,P2 ;
4539 parE1=BRep_Tool::Parameter(Vtx,E1);
4540 parE2=BRep_Tool::Parameter(Vtx,E2);
4541 BCurv1.D1(parE1,P1,dir1);
4542 BCurv2.D1(parE2,P2,dir2);
4543 if (!Vtx.IsSame(TopExp::FirstVertex(E1))) dir1.Reverse();
4544 if (!Vtx.IsSame(TopExp::FirstVertex(E2)))  dir2.Reverse();
4545 angle=Abs(dir1.Angle(dir2));
4546 return angle;
4547 }
4548
4549 //==================================================================
4550 // ChercheBordsLibres
4551 // determines if vertex V1 has edges on free borders 
4552 // edgelibre1 and edgelibre2 .
4553 // It is supposed that a top can have only 2 edges on free borders
4554 //===================================================================
4555 void ChFi3d_ChercheBordsLibres(const  ChFiDS_Map & myVEMap,
4556   const TopoDS_Vertex & V1,
4557   Standard_Boolean & bordlibre,
4558   TopoDS_Edge & edgelibre1,
4559   TopoDS_Edge & edgelibre2) 
4560
4561   bordlibre=Standard_False;
4562   TopTools_ListIteratorOfListOfShape ItE,ItE1;
4563   Standard_Integer nboccur;
4564   for (ItE.Initialize(myVEMap(V1)); ItE.More()&&!bordlibre; ItE.Next()) {
4565     nboccur=0;
4566     const TopoDS_Edge& cur = TopoDS::Edge(ItE.Value());
4567     if (!BRep_Tool::Degenerated(cur)) {
4568       for (ItE1.Initialize(myVEMap(V1)); ItE1.More(); ItE1.Next()) {
4569         const TopoDS_Edge& cur1 = TopoDS::Edge(ItE1.Value());
4570         if (cur1.IsSame(cur)) nboccur++;
4571       }
4572     }
4573     if (nboccur==1) {
4574       edgelibre1=cur;
4575       bordlibre=Standard_True;
4576     }
4577   }
4578   if (bordlibre) {
4579     bordlibre=Standard_False;
4580     for (ItE.Initialize(myVEMap(V1)); ItE.More()&&!bordlibre; ItE.Next()) {
4581       nboccur=0;
4582       const TopoDS_Edge& cur = TopoDS::Edge(ItE.Value());
4583       if (!BRep_Tool::Degenerated(cur)&&!cur.IsSame(edgelibre1)) {
4584         for (ItE1.Initialize(myVEMap(V1)); ItE1.More(); ItE1.Next()) {
4585           const TopoDS_Edge& cur1 = TopoDS::Edge(ItE1.Value());
4586           if (cur1.IsSame(cur)) nboccur++;
4587         }
4588       }
4589       if (nboccur==1) {
4590         edgelibre2=cur;
4591         bordlibre=Standard_True;
4592       }
4593     }
4594   }
4595 }
4596
4597 Standard_Boolean ChFi3d_isTangentFaces(const TopoDS_Edge &theEdge,
4598                                        const TopoDS_Face &theFace1,
4599                                        const TopoDS_Face &theFace2,
4600                                        const GeomAbs_Shape Order)
4601 {
4602   if (Order == GeomAbs_G1 &&
4603       BRep_Tool::Continuity( theEdge, theFace1, theFace2 ) != GeomAbs_C0)
4604     return Standard_True;
4605
4606   Standard_Real TolC0 = Max(0.001, 1.5*BRep_Tool::Tolerance(theEdge));
4607
4608   Standard_Real aFirst;
4609   Standard_Real aLast;
4610     
4611 // Obtaining of pcurves of edge on two faces.
4612   const Handle(Geom2d_Curve) aC2d1 = BRep_Tool::CurveOnSurface
4613                                                 (theEdge, theFace1, aFirst, aLast);
4614   const Handle(Geom2d_Curve) aC2d2 = BRep_Tool::CurveOnSurface
4615                                                 (theEdge, theFace2, aFirst, aLast);
4616   if (aC2d1.IsNull() || aC2d2.IsNull())
4617     return Standard_False;
4618
4619 // Obtaining of two surfaces from adjacent faces.
4620   Handle(Geom_Surface) aSurf1 = BRep_Tool::Surface(theFace1);
4621   Handle(Geom_Surface) aSurf2 = BRep_Tool::Surface(theFace2);
4622
4623   if (aSurf1.IsNull() || aSurf2.IsNull())
4624     return Standard_False;
4625
4626 // Computation of the number of samples on the edge.
4627   BRepAdaptor_Surface              aBAS1(theFace1);
4628   BRepAdaptor_Surface              aBAS2(theFace2);
4629   Handle(BRepAdaptor_HSurface)     aBAHS1      = new BRepAdaptor_HSurface(aBAS1);
4630   Handle(BRepAdaptor_HSurface)     aBAHS2      = new BRepAdaptor_HSurface(aBAS2);
4631   Handle(BRepTopAdaptor_TopolTool) aTool1      = new BRepTopAdaptor_TopolTool(aBAHS1);
4632   Handle(BRepTopAdaptor_TopolTool) aTool2      = new BRepTopAdaptor_TopolTool(aBAHS2);
4633   Standard_Integer                 aNbSamples1 =     aTool1->NbSamples();
4634   Standard_Integer                 aNbSamples2 =     aTool2->NbSamples();
4635   Standard_Integer                 aNbSamples  =     Max(aNbSamples1, aNbSamples2);
4636
4637
4638 // Computation of the continuity.
4639   Standard_Real    aPar;
4640   Standard_Real    aDelta = (aLast - aFirst)/(aNbSamples - 1);
4641   Standard_Integer i, nbNotDone = 0;
4642
4643   for (i = 1, aPar = aFirst; i <= aNbSamples; i++, aPar += aDelta) {
4644     if (i == aNbSamples) aPar = aLast;
4645
4646     LocalAnalysis_SurfaceContinuity aCont(aC2d1,  aC2d2,  aPar,
4647                                           aSurf1, aSurf2, Order,
4648                                           0.001, TolC0, 0.1, 0.1, 0.1);
4649     if (!aCont.IsDone())
4650       {
4651         nbNotDone++;
4652         continue;
4653       }
4654     
4655     if (Order == GeomAbs_G1)
4656     {
4657       if (!aCont.IsG1())
4658         return Standard_False;
4659     }
4660     else if (!aCont.IsG2())
4661       return Standard_False;
4662   }
4663   
4664   if (nbNotDone == aNbSamples)
4665     return Standard_False;
4666
4667   //Compare normals of tangent faces in the middle point
4668   Standard_Real MidPar = (aFirst + aLast)/2.;
4669   gp_Pnt2d uv1 = aC2d1->Value(MidPar);
4670   gp_Pnt2d uv2 = aC2d2->Value(MidPar);
4671   gp_Dir normal1, normal2;
4672   TopOpeBRepTool_TOOL::Nt( uv1, theFace1, normal1 );
4673   TopOpeBRepTool_TOOL::Nt( uv2, theFace2, normal2 );
4674   Standard_Real dot = normal1.Dot(normal2);
4675   if (dot < 0.)
4676     return Standard_False;
4677   return Standard_True;
4678 }
4679
4680 //=======================================================================
4681 //function : NbNotDegeneratedEdges
4682 //purpose  : calculate the number of non-degenerated edges of Map VEMap(Vtx)
4683 // Attention the edges of junctions are taken into account twice
4684 //=======================================================================
4685 Standard_Integer ChFi3d_NbNotDegeneratedEdges (const TopoDS_Vertex& Vtx,
4686   const ChFiDS_Map& VEMap)
4687 {
4688   TopTools_ListIteratorOfListOfShape ItE;
4689   Standard_Integer nba=VEMap(Vtx).Extent();
4690   for (ItE.Initialize(VEMap(Vtx)); ItE.More(); ItE.Next()) {
4691     const TopoDS_Edge& cur = TopoDS::Edge(ItE.Value());
4692     if (BRep_Tool::Degenerated(cur)) nba--;
4693   }
4694   return nba;
4695 }
4696
4697 //=======================================================================
4698 //function : NbSharpEdges
4699 //purpose  : calculate the number of sharp edges of Map VEMap(Vtx)
4700 // Attention the edges of junctions are taken into account twice
4701 //=======================================================================
4702 Standard_Integer ChFi3d_NbSharpEdges (const TopoDS_Vertex& Vtx,
4703                                       const ChFiDS_Map& VEMap,
4704                                       const ChFiDS_Map& EFMap)
4705 {
4706   TopTools_ListIteratorOfListOfShape ItE;
4707   Standard_Integer nba=VEMap(Vtx).Extent();
4708   for (ItE.Initialize(VEMap(Vtx)); ItE.More(); ItE.Next()) {
4709     const TopoDS_Edge& cur = TopoDS::Edge(ItE.Value());
4710     if (BRep_Tool::Degenerated(cur)) nba--;
4711     else
4712     {
4713       TopoDS_Face F1, F2;
4714       ChFi3d_conexfaces(cur, F1, F2, EFMap);
4715       if (!F2.IsNull() && ChFi3d_isTangentFaces(cur, F1, F2, GeomAbs_G2))
4716         nba--;
4717     }
4718   }
4719   return nba;
4720 }
4721
4722 //=======================================================================
4723 //function : NumberOfEdges
4724 //purpose  : calculate the number of edges arriving to the top Vtx
4725 // degenerated edges are not taken into account. 
4726 //=======================================================================
4727 Standard_Integer ChFi3d_NumberOfEdges(const TopoDS_Vertex& Vtx,
4728   const ChFiDS_Map& VEMap)
4729 {
4730   Standard_Integer nba;
4731   Standard_Boolean bordlibre;
4732   TopoDS_Edge edgelibre1,edgelibre2;
4733   nba=ChFi3d_NbNotDegeneratedEdges(Vtx, VEMap);
4734   ChFi3d_ChercheBordsLibres(VEMap,Vtx,bordlibre,edgelibre1,edgelibre2);
4735   if (bordlibre) nba=(nba-2)/2 +2;
4736   else  nba=nba/2;
4737   return nba;
4738 }
4739
4740 //=======================================================================
4741 //function : NumberOfSharpEdges
4742 //purpose  : calculate the number of edges arriving to the top Vtx
4743 // degenerated edges are not taken into account. 
4744 //=======================================================================
4745 Standard_Integer ChFi3d_NumberOfSharpEdges(const TopoDS_Vertex& Vtx,
4746                                            const ChFiDS_Map& VEMap,
4747                                            const ChFiDS_Map& EFmap)
4748 {
4749   Standard_Integer nba;
4750   Standard_Boolean bordlibre;
4751   TopoDS_Edge edgelibre1,edgelibre2;
4752   nba=ChFi3d_NbSharpEdges(Vtx, VEMap, EFmap);
4753   ChFi3d_ChercheBordsLibres(VEMap,Vtx,bordlibre,edgelibre1,edgelibre2);
4754   if (bordlibre) nba=(nba-2)/2 +2;
4755   else  nba=nba/2;
4756   return nba;
4757 }
4758
4759 //=====================================================
4760 // function cherche_vertex
4761 // finds common vertex between two edges 
4762 //=====================================================
4763
4764 void ChFi3d_cherche_vertex (const TopoDS_Edge & E1,
4765   const TopoDS_Edge & E2,
4766   TopoDS_Vertex & vertex,
4767   Standard_Boolean & trouve)
4768 { Standard_Integer i,j; 
4769 TopoDS_Vertex Vcur1,Vcur2;
4770 trouve=Standard_False;
4771 TopTools_IndexedMapOfShape  MapV1,MapV2;
4772 TopExp::MapShapes( E1,TopAbs_VERTEX,MapV1);
4773 TopExp::MapShapes( E2,TopAbs_VERTEX,MapV2);
4774 for ( i=1; i<= MapV1.Extent()&&!trouve; i++) {
4775   TopoDS_Shape alocalshape = TopoDS_Shape (MapV1(i));
4776   Vcur1=TopoDS::Vertex(alocalshape);
4777   //    Vcur1=TopoDS::Vertex(TopoDS_Shape (MapV1(i)));
4778   for ( j=1; j<= MapV2.Extent()&&!trouve; j++) {
4779     TopoDS_Shape aLocalShape = TopoDS_Shape (MapV2(j));
4780     Vcur2=TopoDS::Vertex(aLocalShape);
4781     //      Vcur2=TopoDS::Vertex(TopoDS_Shape (MapV2(j)));
4782     if (Vcur2.IsSame(Vcur1)) {
4783       vertex=Vcur1;trouve=Standard_True;
4784     }
4785   }
4786 }
4787 }       
4788 //=======================================================================
4789 //function : ChFi3d_Couture
4790 //purpose  : determine si F a une arete de couture
4791 //=======================================================================
4792 void ChFi3d_Couture( const TopoDS_Face & F,
4793   Standard_Boolean & couture,
4794   TopoDS_Edge & edgecouture)
4795 {   TopoDS_Edge Ecur;
4796 couture=Standard_False;
4797 TopTools_IndexedMapOfShape  MapE1;
4798 TopExp::MapShapes( F,TopAbs_EDGE,MapE1);
4799 TopLoc_Location Loc;
4800 Handle(Geom_Surface) Surf =BRep_Tool::Surface(F,Loc);
4801 for ( Standard_Integer i=1; i<= MapE1.Extent()&&!couture; i++) {
4802   TopoDS_Shape aLocalShape = TopoDS_Shape (MapE1(i));
4803   Ecur=TopoDS::Edge(aLocalShape);
4804   //      Ecur=TopoDS::Edge(TopoDS_Shape (MapE1(i)));
4805   if (BRep_Tool::IsClosed(Ecur,Surf,Loc)) {
4806     couture=Standard_True;
4807     edgecouture=Ecur;
4808   } 
4809 }
4810 }
4811
4812 //=======================================================================
4813 //function : ChFi3d_CoutureOnVertex
4814 //purpose  : 
4815 //=======================================================================
4816 void ChFi3d_CoutureOnVertex( const TopoDS_Face & F,
4817   const TopoDS_Vertex & V,
4818   Standard_Boolean & couture,
4819   TopoDS_Edge & edgecouture)
4820 {   TopoDS_Edge Ecur;
4821 couture = Standard_False;
4822 TopTools_IndexedMapOfShape  MapE1;
4823 TopExp::MapShapes( F,TopAbs_EDGE,MapE1);
4824 TopLoc_Location Loc;
4825 Handle(Geom_Surface) Surf = BRep_Tool::Surface(F,Loc);
4826 for ( Standard_Integer i=1; i <= MapE1.Extent(); i++) {
4827   TopoDS_Shape aLocalShape = TopoDS_Shape (MapE1(i));
4828   Ecur=TopoDS::Edge(aLocalShape);
4829   //      Ecur=TopoDS::Edge(TopoDS_Shape (MapE1(i)));
4830   if (BRep_Tool::IsClosed(Ecur,Surf,Loc)) {
4831     TopoDS_Vertex Vf, Vl;
4832     TopExp::Vertices( Ecur, Vf, Vl );
4833     if (Vf.IsSame(V) || Vl.IsSame(V))
4834     {
4835       couture = Standard_True;
4836       edgecouture = Ecur;
4837       break;
4838     }
4839   } 
4840 }
4841 }
4842 //=======================================================================
4843 //function : ChFi3d_IsPseudoSeam
4844 //purpose  : 
4845 //=======================================================================
4846 Standard_Boolean ChFi3d_IsPseudoSeam( const TopoDS_Edge& E,
4847   const TopoDS_Face& F )
4848 {
4849   if (! BRep_Tool::IsClosed( E, F ))
4850     return Standard_False;
4851
4852   Standard_Boolean NeighborSeamFound = Standard_False;
4853   TopoDS_Vertex Vf, Vl, V1, V2;
4854   TopExp::Vertices( E, Vf, Vl );
4855   TopExp_Explorer Explo( F, TopAbs_EDGE );
4856   for (; Explo.More(); Explo.Next())
4857   {
4858     TopoDS_Edge Ecur = TopoDS::Edge( Explo.Current() );
4859     if (! Ecur.IsSame(E))
4860     {
4861       TopExp::Vertices( Ecur, V1, V2 );
4862       if ((V1.IsSame(Vf) || V1.IsSame(Vl) || V2.IsSame(Vf) || V2.IsSame(Vl)) &&
4863         BRepTools::IsReallyClosed( Ecur, F ))
4864       {
4865         NeighborSeamFound = Standard_True;
4866         break;
4867       }
4868     }
4869   }
4870   return NeighborSeamFound;
4871 }
4872
4873 //=======================================================================
4874 //function : ChFi3d_ApproxByC2
4875 //purpose  : 
4876 //=======================================================================
4877 Handle(Geom_BSplineCurve) ChFi3d_ApproxByC2( const Handle(Geom_Curve)& C )
4878 {
4879   Standard_Real First = C->FirstParameter(), Last = C->LastParameter();
4880   Standard_Integer NbPoints = 101;
4881
4882   TColgp_Array1OfPnt Points( 1, NbPoints );
4883   Standard_Real delta = (Last - First) / (NbPoints-1);
4884   for (Standard_Integer i = 1; i <= NbPoints-1; i++)
4885     Points(i) = C->Value(First + (i-1)*delta);
4886   Points(NbPoints) = C->Value(Last);
4887
4888   GeomAPI_PointsToBSpline Approx( Points , Approx_ChordLength, 3, 8, GeomAbs_C2, 1.000001e-3);
4889   Handle(Geom_BSplineCurve) BS = Approx.Curve();
4890   return BS;
4891 }
4892 //=======================================================================
4893 //function : ChFi3d_IsSmooth
4894 //purpose  : 
4895 //=======================================================================
4896 Standard_Boolean ChFi3d_IsSmooth( const Handle(Geom_Curve)& C )
4897 {
4898   GeomAdaptor_Curve GAC( C );
4899
4900   Standard_Integer ii;
4901   Standard_Integer intrv, nbintv = GAC.NbIntervals(GeomAbs_CN);
4902   TColStd_Array1OfReal TI(1,nbintv+1);
4903   GAC.Intervals(TI,GeomAbs_CN);
4904   Standard_Real Resolution = gp::Resolution(), Curvature;
4905   GeomLProp_CLProps LProp(C, 2, Resolution);
4906   gp_Pnt P1, P2;
4907   Standard_Integer Discretisation = 30;
4908
4909   gp_Vec PrevVec;
4910   Standard_Boolean prevVecFound = Standard_False;
4911   Standard_Integer intrvFound = 0;
4912   for( intrv = 1; intrv <= nbintv; intrv++) {
4913     Standard_Real t = TI(intrv);
4914     Standard_Real step = (TI(intrv+1) - t) / Discretisation;
4915     for (ii = 1; ii <= Discretisation; ii++) {
4916       LProp.SetParameter(t);
4917       if (!LProp.IsTangentDefined())
4918         return Standard_False;
4919       Curvature = Abs(LProp.Curvature());
4920       if (Curvature > Resolution) {
4921         C->D0(t, P1);
4922         LProp.CentreOfCurvature(P2);
4923         PrevVec = gp_Vec(P1, P2);
4924         prevVecFound = Standard_True;
4925         break;
4926       }
4927       t += step;
4928     }
4929     if( prevVecFound ) {
4930       intrvFound = intrv;
4931       break;
4932     }
4933   }
4934
4935   if( !prevVecFound )
4936     return Standard_True;
4937
4938   //for (intrv = 1; intrv <= nbintv; intrv++) {
4939   for (intrv = intrvFound; intrv <= nbintv; intrv++) {
4940     Standard_Real t = TI(intrv);
4941     Standard_Real step = (TI(intrv+1) - t) / Discretisation;
4942     for (ii = 1; ii <= Discretisation; ii++)
4943     {    
4944       LProp.SetParameter(t);
4945       if (!LProp.IsTangentDefined())
4946         return Standard_False;
4947       Curvature = Abs(LProp.Curvature());
4948       if (Curvature > Resolution)
4949       {
4950         C->D0(t, P1);
4951         LProp.CentreOfCurvature(P2);
4952         gp_Vec Vec(P1, P2);
4953         Standard_Real Angle = PrevVec.Angle( Vec );
4954         if (Angle > M_PI/3.)
4955           return Standard_False;
4956         Standard_Real Ratio = Vec.Magnitude() / PrevVec.Magnitude();
4957         if (Ratio < 1.)
4958           Ratio = 1. / Ratio;
4959         if (Ratio > 2. && (intrv != nbintv || ii != Discretisation))
4960           return Standard_False;
4961         PrevVec = Vec;
4962       }
4963       t += step;
4964     }
4965   }
4966
4967   return Standard_True;
4968 }