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