0022873: Minor correction in ShapeConstruct_ProjectCurveOnSurface.cxx
[occt.git] / src / ShapeConstruct / ShapeConstruct_ProjectCurveOnSurface.cxx
1 //:k1 abv 16.12.98 K4L PRO10107, PRO10108, PRO10109
2 //:j8 abv 10.12.98 TR10 r0501_db.stp #9423
3 //:S4030 abv, pdn: new methods - interface to standard ProjLib_CompProjectedCurve
4 //%12 pdn 15.02.99 PRO9234 optimizing
5 //%12 pdn 15.02.99 PRO9234 using improved ProjectDegenerated method
6 //    rln 03.03.99 S4135: bm2_sd_t4-A.stp treatment of Geom_SphericalSurface together with V-closed surfaces
7 //:p9 abv 11.03.99 PRO7226 #489490: make IsAnIsoparametric to find nearest case
8 //:q1 abv 15.03.99 (pdn) PRO7226 #525030: limit NextValueOfUV() by tolerance
9 //:q5 abv 19.03.99 code improvement
10 //:q9 abv 23.03.99 PRO7226.stp #489490: cashe for projecting end points
11 //#78 rln 12.03.99 S4135: checking spatial closure with myPreci
12 //    pdn 12.03.99 S4135: creating pcurve with minimal length in the case of densed points
13 //    abv 29.03.99 IsAnIsoparametric with Precision::Confusion
14 //    pdn 09.04.99 IsAnisoparametric uses already computed parameters (S4030, fix PRO14323)
15 //szv#4 S4163
16 //:s5 abv 22.04.99  Adding debug printouts in catch {} blocks
17 //#1  svv 11.01.00  Porting on DEC
18 #include <ShapeConstruct_ProjectCurveOnSurface.ixx>
19
20 #include <Standard_ErrorHandler.hxx>
21 #include <Standard_Failure.hxx>
22
23 #include <Precision.hxx>
24 #include <TColStd_Array1OfInteger.hxx>
25 #include <TColgp_Array1OfPnt.hxx>
26
27 #include <GeomAPI_PointsToBSpline.hxx>
28 #include <Geom2dAPI_Interpolate.hxx>
29 #include <GeomAPI_Interpolate.hxx>
30 #include <Geom2dAdaptor.hxx>
31 #include <Geom2d_Line.hxx>
32 #include <Geom2d_TrimmedCurve.hxx>
33 #include <Geom2d_BSplineCurve.hxx>
34 #include <Geom_BSplineCurve.hxx>
35 #include <Geom_TrimmedCurve.hxx>
36 #include <Geom_RectangularTrimmedSurface.hxx>
37 #include <Geom_SurfaceOfLinearExtrusion.hxx>
38 #include <Geom_SphericalSurface.hxx>
39 #include <Geom_OffsetSurface.hxx>
40 #include <Geom_Plane.hxx>
41 #include <GeomProjLib.hxx>
42 #include <GeomAdaptor_HSurface.hxx>
43 #include <GeomAdaptor_HCurve.hxx>
44
45 #include <ShapeAnalysis_Curve.hxx>
46 #include <ShapeAnalysis_Surface.hxx>
47 #include <ShapeExtend.hxx>
48
49 #include <ProjLib_ProjectedCurve.hxx>
50 #include <ProjLib_CompProjectedCurve.hxx>
51 #include <ProjLib_HCompProjectedCurve.hxx>
52 #include <Approx_CurveOnSurface.hxx>
53
54 #define NCONTROL 23
55   
56 //=======================================================================
57 //function : ShapeConstruct_ProjectCurveOnSurface
58 //purpose  : 
59 //=======================================================================
60
61 ShapeConstruct_ProjectCurveOnSurface::ShapeConstruct_ProjectCurveOnSurface()
62 {
63   myPreci = Precision::Confusion();
64   myBuild = Standard_False;
65   myAdjustOverDegen = 1; //:c0 //szv#4:S4163:12Mar99 was boolean
66   myNbCashe = 0; //:q9
67 }
68
69 //=======================================================================
70 //function : Init
71 //purpose  : 
72 //=======================================================================
73
74  void ShapeConstruct_ProjectCurveOnSurface::Init(const Handle(Geom_Surface)& surf,const Standard_Real preci) 
75 {
76   Init (new ShapeAnalysis_Surface (surf), preci);
77 }
78
79 //=======================================================================
80 //function : Init
81 //purpose  : 
82 //=======================================================================
83
84  void ShapeConstruct_ProjectCurveOnSurface::Init(const Handle(ShapeAnalysis_Surface)& surf,const Standard_Real preci) 
85 {
86   SetSurface (surf);
87   SetPrecision (preci);
88 }
89
90 //=======================================================================
91 //function : SetSurface
92 //purpose  : 
93 //=======================================================================
94
95  void ShapeConstruct_ProjectCurveOnSurface::SetSurface(const Handle(Geom_Surface)& surf) 
96 {
97   SetSurface (new ShapeAnalysis_Surface (surf));
98 }
99
100 //=======================================================================
101 //function : SetSurface
102 //purpose  : 
103 //=======================================================================
104
105  void ShapeConstruct_ProjectCurveOnSurface::SetSurface(const Handle(ShapeAnalysis_Surface)& surf) 
106 {
107   if ( mySurf == surf ) return;
108   mySurf = surf;
109   myNbCashe = 0; //:q9
110 }
111
112 //=======================================================================
113 //function : SetPrecision
114 //purpose  : 
115 //=======================================================================
116
117  void ShapeConstruct_ProjectCurveOnSurface::SetPrecision(const Standard_Real preci) 
118 {
119   myPreci = preci;
120 }
121
122 //=======================================================================
123 //function : BuildCurveMode
124 //purpose  : 
125 //=======================================================================
126
127  Standard_Boolean& ShapeConstruct_ProjectCurveOnSurface::BuildCurveMode()
128 {
129   return myBuild;
130 }
131
132 //=======================================================================
133 //function : AdjustOverDegenMode
134 //purpose  : 
135 //=======================================================================
136 //:c0
137
138 //szv#4:S4163:12Mar99 was Boolean
139  Standard_Integer& ShapeConstruct_ProjectCurveOnSurface::AdjustOverDegenMode()
140  {
141    return myAdjustOverDegen;
142  }
143
144
145 //=======================================================================
146 //function : NbSurfIntervals
147 //purpose  : work-around of bug in standard method 
148 //           GeomAdaptor_Surface->NbIntervals() (PRO16346)
149 //=======================================================================
150
151 static Standard_Integer NbSurfIntervals(const Handle(GeomAdaptor_HSurface)& GAS, const GeomAbs_Shape cont)
152 {
153   Standard_Integer NbU = 0;
154   if (GAS->GetType() == GeomAbs_SurfaceOfExtrusion) {
155     // extract the surface
156     Handle(Geom_SurfaceOfLinearExtrusion) surf = Handle(Geom_SurfaceOfLinearExtrusion)::DownCast(GAS->ChangeSurface().Surface());
157     // build a 3d adaptor curve
158     GeomAdaptor_Curve Adaptor3dCurve(surf->BasisCurve(), GAS->FirstUParameter(), GAS->LastUParameter());
159     if (Adaptor3dCurve.GetType() == GeomAbs_BSplineCurve)
160       NbU = Adaptor3dCurve.NbIntervals(cont);
161   }
162   if (NbU == 0)
163     NbU = GAS->NbUIntervals(cont);
164   return NbU * (GAS->NbVIntervals(cont));
165 }
166
167 //=======================================================================
168 //function : Status
169 //purpose  : 
170 //=======================================================================
171
172  Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::Status (const ShapeExtend_Status Status) const
173 {
174   return ShapeExtend::DecodeStatus (myStatus, Status);
175 }
176
177 //=======================================================================
178 //function : Perform
179 //purpose  : 
180 //=======================================================================
181
182 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::Perform (Handle(Geom_Curve)& c3d,
183                                                                 const Standard_Real First,
184                                                                 const Standard_Real Last,
185                                                                 Handle(Geom2d_Curve)& c2d,
186                                                                 const GeomAbs_Shape,
187                                                                 const Standard_Integer,
188                                                                 const Standard_Integer)
189 {
190   myStatus = ShapeExtend::EncodeStatus (ShapeExtend_OK);
191   //Standard_Boolean OK = Standard_True; //szv#4:S4163:12Mar99 not needed
192   
193   if (mySurf.IsNull()) {
194     c2d.Nullify();
195     myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_FAIL1);
196     return Standard_False;
197   }
198
199 //  Projection Analytique
200   Handle(Geom_Curve) crv3dtrim = c3d;
201   if ( ! c3d->IsKind(STANDARD_TYPE(Geom_BoundedCurve)) )
202     crv3dtrim = new Geom_TrimmedCurve ( c3d, First, Last );
203   c2d = ProjectAnalytic ( crv3dtrim );
204   if (!c2d.IsNull()) {
205     myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE1);
206     return Standard_True;
207   }
208
209 //  Projection par approximation
210
211   // discretize the 3d curve
212
213   Standard_Integer nbrPnt;
214   
215 //   $$$$    :92 abv 28 Jan 98   see PRO10107, big BSplineCurve C0
216   Standard_Integer nbPini = NCONTROL;  // as in BRepCheck_Edge (RLN/Nijni)
217  // 20; // number of points for interpolation, should be "parametric dependent"
218
219   //:92 abv 28 Jan 98: if curve is BSpline with many intervals,
220   // increase number of points to provide at least Degree()+1 points per interval
221   Handle(Geom_BSplineCurve) bspl;
222   if ( c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve)) ) {
223     Handle(Geom_TrimmedCurve) ctrim = Handle(Geom_TrimmedCurve)::DownCast(c3d);
224     bspl = Handle(Geom_BSplineCurve)::DownCast ( ctrim->BasisCurve() );
225   }
226   else bspl = Handle(Geom_BSplineCurve)::DownCast ( c3d );
227   if ( ! bspl.IsNull() ) {
228     Standard_Integer nint = 0;
229     for ( Standard_Integer i=1; i < bspl->NbKnots(); i++ )
230       if ( bspl->Knot(i+1) > First && bspl->Knot(i) < Last ) nint++;
231     Standard_Integer minPnt = nint * ( bspl->Degree() + 1 );
232     while ( nbPini < minPnt ) nbPini += NCONTROL - 1;
233 #ifdef DEBUG
234     if ( nbPini > NCONTROL ) 
235       cout << "Warning: number of points for projecting is " << nbPini << endl;
236 #endif
237   }
238
239 //    $$$$    end :92 (big BSplineCurve C0)
240   
241   // this number should be "parametric dependent"
242   TColgp_Array1OfPnt   points(1, nbPini);
243   TColStd_Array1OfReal params(1, nbPini);
244
245   Standard_Integer iPnt;
246   gp_Pnt p3d;
247   Standard_Real deltaT, t;
248   deltaT = (Last - First) / (nbPini-1);
249   nbrPnt = nbPini;
250   for (iPnt = 1; iPnt <= nbPini; iPnt ++) {
251     if      (iPnt == 1)      t = First;
252     else if (iPnt == nbPini) t = Last;
253     else                     t = First + (iPnt - 1) * deltaT;
254
255     c3d->D0 (t, p3d);
256     points(iPnt) = p3d;
257     params(iPnt) = t;
258   }
259
260 //  CALCUL par approximation
261
262   TColgp_Array1OfPnt2d pnt2d(1, nbrPnt);
263   ApproxPCurve (nbrPnt,points,params,pnt2d,c2d); //szv#4:S4163:12Mar99 OK not needed
264   if (!c2d.IsNull()) {
265     myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE2);
266     return Standard_True;
267   }// cas particulier d iso
268
269 //  INTERPOLATION du resultat
270
271   if ( myBuild ) {
272     Handle(TColgp_HArray1OfPnt) thePnts = new TColgp_HArray1OfPnt  (1, nbPini);
273     Handle(TColStd_HArray1OfReal) theParams = new TColStd_HArray1OfReal(1, nbPini);
274     for (iPnt = 1; iPnt <= nbPini ; iPnt ++) {
275       thePnts->SetValue(iPnt, points(iPnt));
276       theParams->SetValue(iPnt, params(iPnt));
277     }
278
279     Handle(Geom_Curve) newc3d = InterpolateCurve3d (nbPini,thePnts,theParams, c3d);
280     if ( newc3d.IsNull() ) myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_FAIL2);
281     else {
282       myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE3);
283       c3d = newc3d;
284     }
285   }
286
287   Handle(TColgp_HArray1OfPnt2d) thePnts2d = new TColgp_HArray1OfPnt2d(1, nbPini);
288   Handle(TColStd_HArray1OfReal) theParams2d = new TColStd_HArray1OfReal(1, nbPini);
289   for (iPnt = 1; iPnt <= nbPini ; iPnt ++) {
290     theParams2d->SetValue(iPnt, params(iPnt));
291     thePnts2d->SetValue(iPnt, pnt2d(iPnt));
292   }
293
294   c2d = InterpolatePCurve (nbPini, thePnts2d, theParams2d, c3d);
295 //  c2d = ApproximatePCurve (nbPini, thePnts2d, theParams2d, c3d);
296 //  Faut-il aussi reprendre la C3D ?
297
298   myStatus |= ShapeExtend::EncodeStatus (c2d.IsNull() ? ShapeExtend_FAIL1 : ShapeExtend_DONE2);
299   return Status (ShapeExtend_DONE);
300 }
301
302 //=======================================================================
303 //function : PerformByProjLib
304 //purpose  : 
305 //=======================================================================
306
307 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::PerformByProjLib(Handle(Geom_Curve)& c3d,
308                                                                         const Standard_Real First,
309                                                                         const Standard_Real Last,
310                                                                         Handle(Geom2d_Curve)& c2d,
311                                                                         const GeomAbs_Shape continuity,
312                                                                         const Standard_Integer maxdeg,
313                                                                         const Standard_Integer nbinterval)
314 {
315   //Standard_Boolean OK = Standard_True; //szv#4:S4163:12Mar99 unused
316   c2d.Nullify();
317   if (mySurf.IsNull()) {
318     myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL1);
319     return Standard_False;
320   }
321    
322   try {
323     OCC_CATCH_SIGNALS
324     Handle(GeomAdaptor_HSurface) GAS = mySurf->Adaptor3d();
325     Standard_Real URes = GAS->ChangeSurface().UResolution ( myPreci );
326     Standard_Real VRes = GAS->ChangeSurface().VResolution ( myPreci );
327     Handle(GeomAdaptor_HCurve) GAC = new GeomAdaptor_HCurve (c3d,First,Last);
328     ProjLib_CompProjectedCurve Projector ( GAS, GAC, URes, VRes );
329      
330     Standard_Real ubeg, ufin;
331     Standard_Integer nbSol = Projector.NbCurves();
332     if (nbSol==1) {
333       Projector.Bounds(1, ubeg, ufin);
334       if((ubeg<=First)&&(ufin>=Last)) {
335         Standard_Integer nbintervals = ( nbinterval < 1 ? 
336                                         NbSurfIntervals(GAS, GeomAbs_C3)+GAC->NbIntervals(GeomAbs_C3)+2:
337                                         nbinterval);
338         Handle(ProjLib_HCompProjectedCurve) HProjector = new ProjLib_HCompProjectedCurve();
339         HProjector->Set(Projector);
340         Handle(Adaptor2d_HCurve2d) HPCur = HProjector;
341         Approx_CurveOnSurface appr(HPCur, GAS, First, Last, myPreci,
342                                    continuity, maxdeg,
343                                    nbintervals,
344                                    Standard_False, Standard_True);
345         if ( appr.IsDone() )
346           c2d = appr.Curve2d();
347       }
348 #ifdef DEB
349       else
350         cout<<"Warning: ProjLib cutting pcurve "<< First << " -> " << ubeg <<" ; "<< Last << " -> " << ufin << endl;
351 #endif     
352     }
353 #ifdef DEB
354     else cout<<"Warning: ProjLib "<< nbSol << " curves in ProjLib"<<endl;
355 #endif
356     if(c2d.IsNull()) {
357       myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL2);
358       return Standard_False;
359     }
360     else {
361       myStatus = ShapeExtend::EncodeStatus (ShapeExtend_DONE1);
362       return Standard_True;
363     }
364     
365   }
366   catch(Standard_Failure) {
367 #ifdef DEB
368     cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::PerformByProjLib(): Exception: ";
369     Standard_Failure::Caught()->Print(cout); cout << endl;
370 #endif
371     myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL3);
372     c2d.Nullify();
373   }
374   return Standard_False;
375 }
376
377 //=======================================================================
378 //function : PerformAdvanced
379 //purpose  : 
380 //=======================================================================
381
382 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::PerformAdvanced (Handle(Geom_Curve)& c3d,
383                                                                         const Standard_Real First,
384                                                                         const Standard_Real Last,
385                                                                         Handle(Geom2d_Curve)& c2d)
386 {
387   Standard_Boolean hasResult = Standard_False;
388   Standard_Integer nbintervals;
389   
390   Standard_Boolean isStandard = (mySurf->Adaptor3d()->GetType() != GeomAbs_Cylinder);
391 // &&   (mySurf->Adaptor3d()->GetType() != GeomAbs_SurfaceOfRevolution);
392
393   if (isStandard) isStandard = !mySurf->HasSingularities(myPreci);
394   if (isStandard) {
395     Handle(GeomAdaptor_HSurface) GAS = mySurf->Adaptor3d();
396     Handle(GeomAdaptor_HCurve) GAC = new GeomAdaptor_HCurve (c3d,First,Last);
397     nbintervals = NbSurfIntervals(GAS, GeomAbs_C1);//+GAC->NbIntervals(GeomAbs_C3);
398     isStandard = (nbintervals < 2);
399   }
400   if (isStandard) {
401     hasResult = PerformByProjLib(c3d, First, Last, c2d);
402   }
403   if (!hasResult) hasResult = Perform (c3d, First, Last, c2d);
404   return hasResult;
405 }
406
407 //=======================================================================
408 //function : ProjectAnalytic
409 //purpose  : 
410 //=======================================================================
411
412  Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::ProjectAnalytic(const Handle(Geom_Curve)& c3d) const
413 {
414   Handle(Geom2d_Curve) result;
415
416   //:k1 abv 16 Dec 98: limit analytic cases by Plane surfaces only
417   // This is necessary for K4L since it fails on other surfaces
418   // when general method GeomProjLib::Curve2d() is used
419   // Projection is done as in BRep_Tool and BRepCheck_Edge
420   Handle(Geom_Surface) surf = mySurf->Surface();
421   Handle(Geom_Plane) Plane = Handle(Geom_Plane)::DownCast ( surf );
422   if ( Plane.IsNull() ) {
423     Handle(Geom_RectangularTrimmedSurface) RTS = 
424       Handle(Geom_RectangularTrimmedSurface)::DownCast ( surf );
425     if ( ! RTS.IsNull() ) Plane = Handle(Geom_Plane)::DownCast ( RTS->BasisSurface() );
426     else {
427       Handle(Geom_OffsetSurface) OS = 
428         Handle(Geom_OffsetSurface)::DownCast ( surf );
429       if ( ! OS.IsNull() ) 
430         Plane = Handle(Geom_Plane)::DownCast ( OS->BasisSurface() );
431     }
432   }
433   if ( ! Plane.IsNull() ) {
434     Handle(Geom_Curve) ProjOnPlane = 
435       GeomProjLib::ProjectOnPlane (c3d, Plane, 
436                                    Plane->Position().Direction(), Standard_True);
437     Handle(GeomAdaptor_HCurve) HC = new GeomAdaptor_HCurve ( ProjOnPlane );
438     ProjLib_ProjectedCurve Proj ( mySurf->Adaptor3d(), HC );
439
440     result = Geom2dAdaptor::MakeCurve(Proj);
441     if ( result.IsNull() ) return result;
442     if ( result->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve)) ) {
443       Handle(Geom2d_TrimmedCurve) TC = Handle(Geom2d_TrimmedCurve)::DownCast ( result );
444       result = TC->BasisCurve();
445     }
446 #ifdef DEB
447 //    if ( ! result.IsNull() ) cout << "SC_PCONS: analitic projection on plane" << endl;
448 #endif
449     return result;
450   }
451   
452   return result;
453 }
454
455 //=======================================================================
456 //function : ApproxPCurve
457 //purpose  : 
458 //=======================================================================
459
460  Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::ApproxPCurve(const Standard_Integer nbrPnt,
461                                                                      const TColgp_Array1OfPnt& points,
462                                                                      const TColStd_Array1OfReal& params,
463                                                                      TColgp_Array1OfPnt2d& pnt2d,
464                                                                      Handle(Geom2d_Curve)& c2d) 
465 {
466   Standard_Boolean isDone = Standard_True;
467   
468   // test if the curve 3d is a boundary of the surface 
469   // (only for Bezier or BSpline surface)
470   
471   Standard_Boolean isoParam, isoPar2d3d, isoTypeU, p1OnIso, p2OnIso, isoclosed;
472   gp_Pnt2d valueP1, valueP2;
473   Handle(Geom_Curve) cIso;
474   Standard_Real t1, t2;
475   
476   Handle(Standard_Type) sType = mySurf->Surface()->DynamicType();
477   Standard_Boolean isAnalytic = Standard_True;
478   if (sType == STANDARD_TYPE(Geom_BezierSurface) || sType == STANDARD_TYPE(Geom_BSplineSurface)) isAnalytic = Standard_False; 
479   Standard_Real uf, ul, vf, vl;
480   mySurf->Surface()->Bounds(uf, ul, vf, vl);
481   isoclosed = Standard_False;
482   TColStd_Array1OfReal pout(1, nbrPnt);
483   
484   isoParam = IsAnIsoparametric(nbrPnt, points, params, 
485                                isoTypeU, p1OnIso, valueP1, p2OnIso, valueP2,
486                                isoPar2d3d, cIso, t1, t2, pout);
487   
488   // projection of the points on surfaces
489   
490   gp_Pnt p3d;
491   gp_Pnt2d p2d;
492   Standard_Integer i;
493   Standard_Real isoValue=0., isoPar1=0., isoPar2=0., tPar=0., tdeb,tfin;
494   Standard_Real Cf, Cl, parf, parl; //szv#4:S4163:12Mar99 dist not needed
495   
496   //  Le calcul part-il dans le bon sens, c-a-d deb et fin dans le bon ordre ?
497   //  Si uclosed et iso en V, attention isoPar1 ET/OU 2 peut toucher la fermeture
498   if(isoParam){
499     if(isoTypeU){
500       isoValue = valueP1.X();
501       isoPar1 = valueP1.Y();
502       isoPar2 = valueP2.Y();
503       isoclosed = mySurf->IsVClosed(myPreci);//#78 rln 12.03.99 S4135
504       parf = vf;  parl = vl;
505     }
506     else { 
507       isoValue = valueP1.Y();
508       isoPar1 = valueP1.X();
509       isoPar2 = valueP2.X();
510       isoclosed = mySurf->IsUClosed(myPreci);//#78 rln 12.03.99 S4135
511       parf = uf;  parl = ul;
512     }
513     if (!isoPar2d3d && !isAnalytic) {
514       Cf = cIso->FirstParameter();
515       Cl = cIso->LastParameter();
516       if (Precision::IsInfinite(Cf))    Cf = -1000;
517       if (Precision::IsInfinite(Cl))    Cl = +1000;
518       //pdn S4030 optimizing and fix isopar case on PRO41323
519       tdeb = pout(2);
520       //    dist = ShapeAnalysis_Curve().Project (cIso,points(2),myPreci,pt,tdeb,Cf,Cl);
521       //  Chacun des par1 ou par2 est-il sur un bord. Attention first/last : recaler
522       if (isoclosed && (isoPar1 == parf || isoPar1 == parl)) {
523         if (Abs(tdeb-parf) < Abs(tdeb-parl)) isoPar1 = parf;
524         else isoPar1 = parl;
525         if (isoTypeU) valueP1.SetY (isoPar1);
526         else          valueP1.SetX (isoPar1);
527       }
528       if (isoclosed && (isoPar2 == parf || isoPar2 == parl)) {
529         //pdn S4030 optimizing and fix isopar case on PRO41323
530         tfin = pout(nbrPnt-1);
531         //dist =  ShapeAnalysis_Curve().Project (cIso,points(nbrPnt-1),myPreci,pt,tfin,Cf,Cl);
532         if (Abs(tfin-parf) < Abs(tfin-parl)) isoPar2 = parf;
533         else isoPar2 = parl;
534         if (isoTypeU) valueP2.SetY (isoPar2);
535         else          valueP2.SetX (isoPar2);
536       }
537       
538       //  Interversion Par1/Par2 (ne veut que si les 2 sont sur les bords ...)
539       //  Est-ce encore necessaire apres ce qui vient d etre fait ?
540       
541       // PTV 05.02.02 fix for translation face from 12_hp_mouse (PARASOLID) face 24008
542       // if curve is periodic do not change the points
543       // skl change "if" for pout(nbrPnt-1) 19.11.2003
544       if (!isoclosed) {
545         if( (Abs(tdeb-isoPar1)>Abs(tdeb-isoPar2)) &&
546             (Abs(pout(nbrPnt-1)-isoPar2)>Abs(pout(nbrPnt-1)-isoPar1)) ) {
547           gp_Pnt2d valueTmp = valueP1;
548           valueP1 = valueP2;  valueP2 = valueTmp;
549           if (isoTypeU) {
550             isoValue = valueP1.X();
551             isoPar1 = valueP1.Y();
552             isoPar2 = valueP2.Y();
553           }
554           else { 
555             isoValue = valueP1.Y();
556             isoPar1 = valueP1.X();
557             isoPar2 = valueP2.X();
558           }
559           //  Fin calcul sens de courbe iso
560         }
561       } // end of fix check 05.02.02
562     }
563   }
564   
565   //  Si pas isoParam, on a quand meme du p1OnIso/p2OnIso possible ... !!!
566   //  (utile pour detromper bug de projection). Mais detromper aussi circularite
567   //else {
568     //if (p1OnIso) valueP1 =
569     //BestExtremum (valueP1,points(1),points(2));
570     //if (p2OnIso) valueP2 =
571     //BestExtremum (valueP2,points(nbrPnt),points(nbrPnt-1));
572     //}
573   
574   Standard_Real gap = myPreci; //:q1
575   Standard_Boolean ChangeCycle = Standard_False; //skl for OCC3430
576   if( myNbCashe>0 && myCashe3d[0].Distance(points(1))>myCashe3d[0].Distance(points(nbrPnt)) )
577     //if(myCashe3d[0].Distance(points(nbrPnt))<myPreci)
578     if(myCashe3d[0].Distance(points(nbrPnt))<Precision::Confusion())
579       ChangeCycle = Standard_True;
580   //for( i = 1; i <= nbrPnt; i ++) {
581   for(Standard_Integer ii=1; ii<=nbrPnt; ii++) {
582     if(ChangeCycle) //skl for OCC3430
583       i=nbrPnt-ii+1;
584     else
585       i=ii;
586     p3d = points(i);
587     if (isoParam) {
588       
589       if (isoPar2d3d) {
590         if (isoPar2 > isoPar1) tPar = params(i);
591         else                   tPar = t1 + t2 - params(i);
592       } else if (!isAnalytic) {
593         // projection to iso
594         if      (i==1)      tPar = isoPar1;
595         else if (i==nbrPnt) tPar = isoPar2;
596         else {
597           tPar = pout(i);
598           //:S4030  ShapeAnalysis_Curve().Project (cIso,p3d,myPreci,pt,tPar,Cf,Cl); //szv#4:S4163:12Mar99 `dist=` not needed
599         }
600       }
601       
602       if (!isoPar2d3d && isAnalytic) {
603         if      (i == 1)      p2d = valueP1;
604         else if (i == nbrPnt) p2d = valueP2;
605         else {
606           p2d = mySurf->NextValueOfUV(p2d,p3d, myPreci, //%12 pdn 15.02.99 optimizing
607                                       Precision::Confusion()+1000*gap); //:q1
608           gap = mySurf->Gap();
609         }
610       } else {
611         if(isoTypeU)  {  p2d.SetX(isoValue);  p2d.SetY(tPar);     }
612         else          {  p2d.SetX(tPar);      p2d.SetY(isoValue); }
613       }
614     }
615     
616     else {
617       if     ( (i == 1)      && p1OnIso)  p2d = valueP1;
618       else if( (i == nbrPnt) && p2OnIso)  p2d = valueP2;
619       else  {// general case (not an iso)  mais attention aux singularites !
620         if ( ii==1 ) {
621           //:q9 abv 23 Mar 99: use cashe as 1st approach
622           Standard_Integer j; // svv #1
623           for ( j=0; j < myNbCashe; j++ ) 
624             if ( myCashe3d[j].SquareDistance ( p3d ) < myPreci*myPreci ) {
625               p2d = mySurf->NextValueOfUV (myCashe2d[j], p3d, myPreci, 
626                                            Precision::Confusion()+gap);
627               break;
628             }
629           if ( j >= myNbCashe ) p2d = mySurf->ValueOfUV(p3d, myPreci); 
630         }
631         else {
632           p2d = mySurf->NextValueOfUV (p2d, p3d, myPreci, //:S4030: optimizing
633                                        Precision::Confusion()+1000*gap); //:q1
634         }
635         gap = mySurf->Gap();
636       }
637     }
638     pnt2d (i) = p2d;
639     if ( ii > 1 ) {
640       if(ChangeCycle)
641         p2d.SetXY ( 2. * p2d.XY() - pnt2d(i+1).XY() );
642       else
643         p2d.SetXY ( 2. * p2d.XY() - pnt2d(i-1).XY() );
644     }
645   }
646
647   //pdn %12 11.02.99 PRO9234 entity 15402
648   if (!isoPar2d3d) {
649     mySurf->ProjectDegenerated(nbrPnt,points,pnt2d,myPreci,Standard_True);
650     mySurf->ProjectDegenerated(nbrPnt,points,pnt2d,myPreci,Standard_False);
651   }
652   
653   //  attention aux singularites ... (hors cas iso qui les traite deja)
654   //  if (!isoParam) {
655     //    p2d = pnt2d (1);
656     //    if (mySurf->ProjectDegenerated (points(1),myPreci,pnt2d (2),p2d))
657     //      pnt2d (1) = p2d;
658     //    p2d = pnt2d (nbrPnt);
659     //    if (mySurf->ProjectDegenerated (points(nbrPnt),myPreci,pnt2d (nbrPnt-1),p2d))
660     //      pnt2d (nbrPnt) = p2d;
661     //  }
662   
663   // Si la surface est UCLosed et VClosed, on recadre les points
664   // algo un peu complique, on retarde l implementation
665   Standard_Real Up = ul - uf;
666   Standard_Real Vp = vl - vf;
667   Standard_Real dist2d;
668 #ifdef DEBUG
669   if (mySurf->IsUClosed(myPreci) && mySurf->IsVClosed(myPreci)) {//#78 rln 12.03.99 S4135
670     cout << "WARNING : Recadrage incertain sur U & VClosed" << endl;
671   }
672 #endif
673   // Si la surface est UCLosed, on recadre les points
674   if (mySurf->IsUClosed(myPreci)) {//#78 rln 12.03.99 S4135
675     // Premier point dans le domain [uf, ul]
676     Standard_Real prevX, firstX = pnt2d (1).X();
677     while (firstX < uf)  {  firstX += Up;   pnt2d (1).SetX(firstX);  }
678     while (firstX > ul)  {  firstX -= Up;   pnt2d (1).SetX(firstX);  }
679     prevX = firstX;
680     
681     //:97 abv 1 Feb 98: treat case when curve is whole out of surface bounds
682     Standard_Real minX = firstX, maxX = firstX;
683     
684     // On decalle toujours le suivant
685     for (i = 2; i <= nbrPnt; i++) {
686       //      dist2d = pnt2d (i-1).Distance(pnt2d (i));
687       Standard_Real CurX = pnt2d (i).X();
688       dist2d = Abs (CurX - prevX);
689       if (dist2d > ( Up / 2) ) {
690         if        (CurX > prevX + Up/2) {
691           while (CurX > prevX + Up/2) {  CurX -= Up;  pnt2d (i).SetX (CurX);  }
692         } else if (CurX < prevX - Up/2) {
693           while (CurX < prevX - Up/2) {  CurX += Up;  pnt2d (i).SetX (CurX);  }
694         }
695         
696       }
697       prevX = CurX;
698       if ( minX > CurX ) minX = CurX;      //:97
699       else if ( maxX < CurX ) maxX = CurX; //:97
700     }
701     
702     //:97
703     Standard_Real midX = 0.5 * ( minX + maxX );
704     Standard_Real shiftX=0.;
705     if ( midX > ul ) shiftX = -Up;
706     else if ( midX < uf ) shiftX = Up;
707     if ( shiftX != 0. ) 
708       for ( i=1; i <= nbrPnt; i++ ) pnt2d(i).SetX ( pnt2d(i).X() + shiftX );
709   }
710   // Si la surface est VCLosed, on recadre les points
711   // Same code as UClosed : optimisation souhaitable !!
712   // CKY : d abord un code IDENTIQUE A UClosed; PUIS le special Seam ...
713   // Si la surface est UCLosed, on recadre les points
714   //
715   //#69 rln 01.03.99 S4135 bm2_sd_t4-A.stp entity 30
716   //#78 rln 12.03.99 S4135
717   if (mySurf->IsVClosed(myPreci) || mySurf->Surface()->IsKind (STANDARD_TYPE (Geom_SphericalSurface))) {
718     // Premier point dans le domain [vf, vl]
719     Standard_Real prevY, firstY = pnt2d (1).Y();
720     while (firstY < vf)  {  firstY += Vp;  pnt2d (1).SetY(firstY);  }
721     while (firstY > vl)  {  firstY -= Vp;  pnt2d (1).SetY(firstY);  }
722     prevY = firstY;
723     
724     //:97 abv 1 Feb 98: treat case when curve is whole out of surface bounds
725     Standard_Real minY = firstY, maxY = firstY;
726     
727     // On decalle toujours le suivant
728     for (i = 2; i <= nbrPnt; i ++) {
729       //      dist2d = pnt2d (i-1).Distance(pnt2d (i));
730       Standard_Real CurY = pnt2d (i).Y();
731       dist2d = Abs (CurY - prevY);
732       if (dist2d > ( Vp / 2) ) {
733         if        (CurY > prevY + Vp/2) {
734           while (CurY > prevY + Vp/2) {  CurY -= Vp;  pnt2d (i).SetY (CurY);  }
735         } else if (CurY < prevY - Vp/2) {
736           while (CurY < prevY - Vp/2) {  CurY += Vp;  pnt2d (i).SetY (CurY);  }
737         }
738       }
739       prevY = CurY;
740       if ( minY > CurY ) minY = CurY;      //:97
741       else if ( maxY < CurY ) maxY = CurY; //:97
742     }
743     
744     //:97
745     Standard_Real midY = 0.5 * ( minY + maxY );
746     Standard_Real shiftY=0.;
747     if ( midY > vl ) shiftY = -Vp;
748     else if ( midY < vf ) shiftY = Vp;
749     if ( shiftY != 0. ) 
750       for ( i=1; i <= nbrPnt; i++ ) pnt2d(i).SetY ( pnt2d(i).Y() + shiftY );
751   }
752   
753   //#69 rln 01.03.99 S4135 bm2_sd_t4-A.stp entity 30
754   //#78 rln 12.03.99 S4135
755   if (mySurf->IsVClosed(myPreci) || mySurf->Surface()->IsKind (STANDARD_TYPE (Geom_SphericalSurface))) {
756     for (i = 2; i <= nbrPnt; i++) {
757       //#1 rln 11/02/98 ca_exhaust.stp entity #9869 dist2d = pnt2d (i-1).Distance(pnt2d (i));
758       dist2d = Abs (pnt2d(i).Y() - pnt2d(i - 1).Y());
759       if (dist2d > ( Vp / 2) ) {
760         // ATTENTION : il faut regarder ou le decalage se fait.
761         // si plusieurs points sont decalles, il faut plusieurs passes
762         // pour obtenir un resultat correct.
763         // NOT YET IMPLEMENTED
764         
765         // one of those point is incorrectly placed
766         // i.e on the wrong side of the "seam"
767         // on prend le point le plus pres des bords vf ou vl
768         Standard_Boolean prevOnFirst = Standard_False;
769         Standard_Boolean prevOnLast  = Standard_False;
770         Standard_Boolean currOnFirst = Standard_False;
771         Standard_Boolean currOnLast  = Standard_False;
772         
773         //  .X ?  plutot .Y ,  non ?
774         Standard_Real distPrevVF = Abs(pnt2d (i-1).Y() - vf);
775         Standard_Real distPrevVL = Abs(pnt2d (i-1).Y() - vl);
776         Standard_Real distCurrVF = Abs(pnt2d (i).Y() - vf);
777         Standard_Real distCurrVL = Abs(pnt2d (i).Y() - vl);
778         
779         Standard_Real theMin = distPrevVF;
780         prevOnFirst = Standard_True;
781         if (distPrevVL < theMin) {
782           theMin = distPrevVL;
783           prevOnFirst = Standard_False;
784           prevOnLast  = Standard_True;
785         }
786         if (distCurrVF < theMin) {
787           theMin = distCurrVF;
788           prevOnFirst = Standard_False;
789           prevOnLast  = Standard_False;
790           currOnFirst = Standard_True;
791         }
792         if (distCurrVL < theMin) {
793           theMin = distCurrVL;
794           prevOnFirst = Standard_False;
795           prevOnLast  = Standard_False;
796           currOnFirst = Standard_False;
797           currOnLast  = Standard_True;
798         }
799         //  Modifs RLN/Nijni  3-DEC-1997
800         if (prevOnFirst) {
801           // on decalle le point (i-1) en V Last
802           gp_Pnt2d newPrev(pnt2d (i-1).X(), vf); // instead of  vl RLN/Nijni
803           pnt2d (i-1) = newPrev;
804         }
805         else if (prevOnLast) {
806           // on decalle le point (i-1) en V first
807           gp_Pnt2d newPrev(pnt2d (i-1).X(), vl); // instead of  vf RLN/Nijni
808           pnt2d (i-1) = newPrev;
809         }
810         else if (currOnFirst) {
811           // on decalle le point (i) en V Last
812           gp_Pnt2d newCurr(pnt2d (i).X(),vf);  // instead of vl  RLN/Nijni
813           pnt2d (i) = newCurr;
814         }
815         else if (currOnLast) {
816           // on decalle le point (i) en V First
817           gp_Pnt2d newCurr(pnt2d (i).X(), vl); // instead of vf  RLN/Nijni
818           pnt2d (i) = newCurr;
819         }
820         // on verifie
821 #ifdef DEBUG
822         dist2d = pnt2d (i-1).Distance(pnt2d (i));
823         if (dist2d > ( Vp / 2) ) {
824           cout << "Echec dans le recadrage" << endl;
825         }
826 #endif
827       }
828     }
829   }
830   
831   //:c0 abv 20 Feb 98: treat very special case when 3d curve
832   // go over the pole of, e.g., sphere, and partly lies along seam.
833   // 2d representation of such a curve should consist of 3 parts - one on
834   // regular part of surface (interior), one part along degenerated boundary
835   // and one along seam.
836   // Since it cannot be adjusted later by arranging pcurves (curve is single),
837   // to fix it it is nesessary to have a possibility of adjusting seam
838   // part of such curve either to left or right boundary of surface.
839   // Test is performed only if flag AdjustOverDegen is not -1.
840   // If AdjustOverDegen is True, seam part of curve is adjusted to
841   // the left, and if False - to the right parametric boundary 
842   // If treated case is detected, flag DONE4 is set to status
843   // NOTE: currently, precision is Precision::PConfusion() since it 
844   // is enough on encountered example
845   // (ug_turbine-A.stp from ProSTEP Benchmark #3, entities ##2470 & 5680)
846   // (r1001_ac.stp from Test Rally #10, face #35027 and others)
847   if ( myAdjustOverDegen != -1 ) {
848     if ( mySurf->IsUClosed(myPreci) ) {//#78 rln 12.03.99 S4135
849       mySurf->IsDegenerated ( gp_Pnt(0,0,0), myPreci );  // pour calculer les dgnr
850       if ( mySurf->NbSingularities(myPreci) > 0 ) { //rln S4135
851         // 1st, find gap point (degenerated pole)
852         Standard_Real PrevX=0.;
853         Standard_Integer OnBound=0, PrevOnBound=0;
854         Standard_Integer ind; // svv #1
855         Standard_Boolean start = Standard_True;
856         for ( ind=1; ind <= nbrPnt; ind++ ) {
857           Standard_Real CurX = pnt2d(ind).X();
858           // abv 16 Mar 00: trj3_s1-ug.stp #697: ignore points in singularity
859           if ( mySurf->IsDegenerated ( points(ind), Precision::Confusion() ) )
860             continue;
861           OnBound = ( Abs ( Abs ( CurX - 0.5 * ( ul + uf ) ) - Up/2 ) <=
862                      Precision::PConfusion() );
863           if ( ! start &&  Abs ( Abs ( CurX - PrevX ) - Up/2 ) <= 0.01*Up ) 
864             break;
865           start = Standard_False;
866           PrevX = CurX;
867           PrevOnBound = OnBound;
868         }
869         // if found, adjust seam part
870         if ( ind <= nbrPnt ) {
871           PrevX = ( myAdjustOverDegen ? uf : ul );
872           Standard_Real dU = Up/2 + Precision::PConfusion();
873           if ( PrevOnBound ) {
874             pnt2d(ind-1).SetX ( PrevX );
875             for ( Standard_Integer j=ind-2; j >0; j-- ) {
876               Standard_Real CurX = pnt2d(j).X();
877               while ( CurX < PrevX - dU ) pnt2d(j).SetX ( CurX += Up );
878               while ( CurX > PrevX + dU ) pnt2d(j).SetX ( CurX -= Up );
879             }
880           }
881           else if ( OnBound ) {
882             pnt2d(ind).SetX ( PrevX );
883             for ( Standard_Integer j=ind+1; j <= nbrPnt; j++ ) {
884               Standard_Real CurX = pnt2d(j).X();
885               while ( CurX < PrevX - dU ) pnt2d(j).SetX ( CurX += Up );
886               while ( CurX > PrevX + dU ) pnt2d(j).SetX ( CurX -= Up );
887             }
888           }
889           myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE4);
890         }
891       }
892     }
893     else if ( mySurf->IsVClosed(myPreci) ) {//#78 rln 12.03.99 S4135
894       mySurf->IsDegenerated ( gp_Pnt(0,0,0), myPreci );  // pour calculer les dgnr
895       if ( mySurf->NbSingularities(myPreci) > 0 ) { //rln S4135
896         // 1st, find gap point (degenerated pole)
897         Standard_Real PrevY=0.;
898         Standard_Integer OnBound=0, PrevOnBound=0;
899         Standard_Integer ind; // svv #1
900         Standard_Boolean start = Standard_True;
901         for ( ind=1; ind <= nbrPnt; ind++ ) {
902           Standard_Real CurY = pnt2d(ind).Y();
903           // abv 16 Mar 00: trj3_s1-ug.stp #697: ignore points in singularity
904           if ( mySurf->IsDegenerated ( points(ind), Precision::Confusion() ) )
905             continue;
906           OnBound = ( Abs ( Abs ( CurY - 0.5 * ( vl + vf ) ) - Vp/2 ) <=
907                      Precision::PConfusion() );
908           if ( ! start &&  Abs ( Abs ( CurY - PrevY ) - Vp/2 ) <= 0.01*Vp ) 
909             break;
910           start = Standard_False;
911           PrevY = CurY;
912           PrevOnBound = OnBound;
913         }
914         // if found, adjust seam part
915         if ( ind <= nbrPnt ) {
916           PrevY = ( myAdjustOverDegen ? vf : vl );
917           Standard_Real dV = Vp/2 + Precision::PConfusion();
918           if ( PrevOnBound ) {
919             pnt2d(ind-1).SetY ( PrevY );
920             for ( Standard_Integer j=ind-2; j >0; j-- ) {
921               Standard_Real CurY = pnt2d(j).Y();
922               while ( CurY < PrevY - dV ) pnt2d(j).SetY ( CurY += Vp );
923               while ( CurY > PrevY + dV ) pnt2d(j).SetY ( CurY -= Vp );
924             }
925           }
926           else if ( OnBound ) {
927             pnt2d(ind).SetY ( PrevY );
928             for ( Standard_Integer j=ind+1; j <= nbrPnt; j++ ) {
929               Standard_Real CurY = pnt2d(j).Y();
930               while ( CurY < PrevY - dV ) pnt2d(j).SetY ( CurY += Vp );
931               while ( CurY > PrevY + dV ) pnt2d(j).SetY ( CurY -= Vp );
932             }
933           }
934           myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE4);
935         }
936       }
937     }
938   }
939
940   //:q9: fill cashe
941   myNbCashe = 2;
942   if(ChangeCycle) {  // msv 10.08.04: avoid using of uninitialised field
943   //if(myCashe3d[0].Distance(points(1))>Precision::Confusion() &&
944   //   myCashe3d[1].Distance(points(1))>Precision::Confusion()) {
945     myCashe3d[0] = points(1);
946     myCashe3d[1] = points(nbrPnt);
947     myCashe2d[0] = pnt2d(1);
948     myCashe2d[1] = pnt2d(nbrPnt);
949   }
950   else {
951     myCashe3d[1] = points(1);
952     myCashe3d[0] = points(nbrPnt);
953     myCashe2d[1] = pnt2d(1);
954     myCashe2d[0] = pnt2d(nbrPnt);
955   }
956   
957   if (isoParam && isoPar2d3d) {
958     
959     // create directly isoparametrics (PCurve)
960     
961     gp_Vec2d aVec2d(valueP1, valueP2);
962     gp_Dir2d aDir2d(aVec2d);
963     gp_Pnt2d P0;
964     
965     if (isoTypeU) P0.SetCoord(valueP1.X(), valueP1.Y() - params(1)*aDir2d.Y());
966     else          P0.SetCoord(valueP1.X() - params(1)*aDir2d.X(), valueP1.Y());
967     
968     c2d = new Geom2d_Line(P0, aDir2d);
969   }
970
971   if(c2d.IsNull()) {
972     // try create line using pnt2d
973     Standard_Boolean IsLine = Standard_True;
974     Standard_Integer NbPnt2d = pnt2d.Length();
975     if(NbPnt2d >1) {
976       Standard_Real dist = pnt2d(1).SquareDistance(pnt2d(NbPnt2d));
977       Standard_Real dPreci = Precision::PConfusion()*Precision::PConfusion();
978       if(dist >= dPreci) {
979         Standard_Real tol2 = dPreci*dPreci;
980         gp_Vec2d avec (pnt2d(1),pnt2d(NbPnt2d));
981         gp_Dir2d adir (avec);
982         gp_Lin2d alin (pnt2d(1),adir);
983         for(i = 2; i < NbPnt2d; i++) {
984           Standard_Real devia = alin.SquareDistance(pnt2d(i));
985           Standard_Real dist2 = pnt2d(1).SquareDistance(pnt2d(i));
986           Standard_Real step = pnt2d(1).Distance(pnt2d(NbPnt2d))/(NbPnt2d-1);
987           Standard_Real ddist = Abs(pnt2d(1).Distance(pnt2d(i))-step*(i-1));
988           if( devia>tol2 || (dist2-dist)>dPreci || ddist>1.e-3*step ) {
989             IsLine = Standard_False;
990             i = NbPnt2d;
991           }
992         }
993         if(IsLine) {
994           Handle(Geom2d_Line) g2l = new Geom2d_Line(alin);
995           c2d = new Geom2d_TrimmedCurve(g2l,0,pnt2d(1).Distance(pnt2d(NbPnt2d)));
996         }
997       }  
998     }
999   }
1000   
1001   return isDone;
1002 }
1003
1004 //=======================================================================
1005 //function : ApproximatePCurve
1006 //purpose  : 
1007 //=======================================================================
1008
1009 Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::ApproximatePCurve(const Standard_Integer /*nbrPnt*/,
1010                                                                              Handle(TColgp_HArray1OfPnt2d)& points2d, 
1011                                                                              Handle(TColStd_HArray1OfReal)& params,
1012                                                                              const Handle(Geom_Curve)& /*orig*/) const
1013 {
1014 //  Standard_Real resol = Min(mySurf->Adaptor3d()->VResolution(myPreci), mySurf->Adaptor3d()->UResolution(myPreci));
1015   Standard_Real theTolerance2d = myPreci; // (100*nbrPnt);//resol;
1016   Handle(Geom2d_Curve) C2d;
1017   try {
1018     OCC_CATCH_SIGNALS
1019     CheckPoints2d (points2d, params, theTolerance2d);
1020     Standard_Integer numberPnt = points2d->Length();
1021     
1022     TColgp_Array1OfPnt points3d(1,numberPnt);
1023     gp_Pnt2d pnt2d;
1024     gp_Pnt pnt;
1025     Standard_Integer i; // svv #1 
1026     for( i = 1; i <= numberPnt; i++) {
1027       pnt2d = points2d->Value(i);
1028       pnt.SetCoord(pnt2d.X(),pnt2d.Y(),0);
1029       points3d(i) = pnt;
1030     }
1031     
1032     GeomAPI_PointsToBSpline appr(points3d, params->Array1(), 1, 10, GeomAbs_C1, theTolerance2d);
1033     Handle(Geom_BSplineCurve) crv3d = appr.Curve();
1034     Standard_Integer NbPoles = crv3d->NbPoles();
1035     TColgp_Array1OfPnt poles3d (1, NbPoles);
1036     TColgp_Array1OfPnt2d poles2d (1, NbPoles);
1037     crv3d->Poles(poles3d);
1038     for( i = 1; i <= NbPoles; i++) {
1039       pnt2d.SetCoord(poles3d(i).X(),poles3d(i).Y());
1040       poles2d(i) = pnt2d;
1041     }
1042     TColStd_Array1OfReal weights (1,NbPoles);
1043     TColStd_Array1OfInteger multiplicities (1,crv3d->NbKnots());
1044     TColStd_Array1OfReal knots(1,crv3d->NbKnots());
1045     crv3d->Knots(knots);
1046     crv3d->Weights(weights);
1047     crv3d->Multiplicities(multiplicities);
1048     C2d = new Geom2d_BSplineCurve  ( poles2d, weights, knots, multiplicities, crv3d->Degree(), crv3d->IsPeriodic());
1049     return C2d;
1050   }
1051   catch(Standard_Failure) {
1052 #ifdef DEB //:s5
1053     //    debug ...
1054     Standard_Integer nbp = params->Length();
1055     Standard_Integer nb2 = points2d->Length();
1056     cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::ApproximatePCurve(): Exception: ";
1057     Standard_Failure::Caught()->Print(cout); 
1058     cout<<"Pb Geom2dAPI_Approximate, tol2d="<<theTolerance2d<<" NbParams="<<nbp<<" NbPnts="<<nb2<<endl;
1059 //     if (nb2 > nbp) nb2 = nbp;
1060 //     Standard_Real rbp,rb2; rbp = nbp; rb2 = nb2;
1061 //     //    dbl.AddString ("NbP2d/NbParams puis  X Y Param -> mini");
1062 //     dbl.AddReals (rb2,rbp);
1063 //     for (Standard_Integer i = 1; i <= nb2; i ++) {
1064 //       gp_XYZ quoi (points2d->Value(i).X(),points2d->Value(i).Y(),params->Value(i) );
1065 //       dbl.AddXYZ (quoi);
1066 //     }
1067 #endif
1068      C2d.Nullify();
1069   }
1070   return C2d;
1071 }  
1072
1073 //=======================================================================
1074 //function : InterpolatePCurve
1075 //purpose  : 
1076 //=======================================================================
1077
1078  Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::InterpolatePCurve(const Standard_Integer nbrPnt,
1079                                                                               Handle(TColgp_HArray1OfPnt2d)& points2d, 
1080                                                                               Handle(TColStd_HArray1OfReal)& params,
1081                                                                               const Handle(Geom_Curve)& /*orig*/) const
1082 {
1083   Handle(Geom2d_Curve) C2d;    // NULL si echec
1084   Standard_Real theTolerance2d = myPreci / (100 * nbrPnt);
1085   try {
1086     OCC_CATCH_SIGNALS
1087     // on verifie d abord s il n y a pas de points confondus
1088     // si besoin on retouche les valeurs ...
1089     CheckPoints2d (points2d, params, theTolerance2d);
1090     Geom2dAPI_Interpolate myInterPol2d (points2d, params, 
1091                                         Standard_False, theTolerance2d);
1092     myInterPol2d.Perform();
1093     if (myInterPol2d.IsDone()) C2d = myInterPol2d.Curve();
1094   }
1095   catch(Standard_Failure) {
1096 #ifdef DEB //:s5
1097 // //    debug ...
1098     Standard_Integer nbp = params->Length();
1099     Standard_Integer nb2 = points2d->Length();
1100     cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::InterpolatePCurve(): Exception: ";
1101     Standard_Failure::Caught()->Print(cout); 
1102     cout<<"Pb Geom2dAPI_Interpolate, tol2d="<<theTolerance2d<<" NbParams="<<nbp<<" NbPnts="<<nb2<<endl;
1103 //     if (nb2 > nbp) nb2 = nbp;
1104 //     Standard_Real rbp,rb2; rbp = nbp; rb2 = nb2;
1105 // //    dbl.AddString ("NbP2d/NbParams puis  X Y Param -> mini");
1106 //     dbl.AddReals (rb2,rbp);
1107 //     for (Standard_Integer i = 1; i <= nb2; i ++) {
1108 //       gp_XYZ quoi (points2d->Value(i).X(),points2d->Value(i).Y(),params->Value(i) );
1109 //       dbl.AddXYZ (quoi);
1110 //     }
1111 #endif
1112     C2d.Nullify();
1113   }
1114   return C2d;
1115 }
1116
1117 //=======================================================================
1118 //function : InterpolateCurve3d
1119 //purpose  : 
1120 //=======================================================================
1121
1122 Handle(Geom_Curve) ShapeConstruct_ProjectCurveOnSurface::InterpolateCurve3d(const Standard_Integer,
1123                                                                             Handle(TColgp_HArray1OfPnt)& points, 
1124                                                                             Handle(TColStd_HArray1OfReal)& params,
1125                                                                             const Handle(Geom_Curve)& /*orig*/) const
1126 {
1127   Handle(Geom_Curve) C3d;    // NULL si echec
1128   try {
1129     OCC_CATCH_SIGNALS
1130     Standard_Real Tol = myPreci;
1131     CheckPoints(points, params, Tol);
1132     GeomAPI_Interpolate myInterPol(points, params, Standard_False, Tol);
1133     myInterPol.Perform();
1134     if (myInterPol.IsDone()) C3d = myInterPol.Curve();
1135   }
1136   catch(Standard_Failure) {
1137     C3d.Nullify();
1138 #ifdef DEB //:s5
1139     cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::InterpolateCurve3d(): Exception: ";
1140     Standard_Failure::Caught()->Print(cout); cout << endl;
1141 #endif
1142   }
1143   return C3d;
1144 }
1145
1146 //=======================================================================
1147 //function : CheckPoints
1148 //purpose  : 
1149 //=======================================================================
1150
1151  void ShapeConstruct_ProjectCurveOnSurface::CheckPoints(Handle(TColgp_HArray1OfPnt)& points,Handle(TColStd_HArray1OfReal)& params,Standard_Real& preci) const
1152 {
1153   Standard_Integer firstElem = points->Lower();
1154   Standard_Integer lastElem  = points->Upper();
1155   Standard_Integer i;
1156   Standard_Integer nbPntDropped = 0;
1157   Standard_Integer lastValid = firstElem; // indice of last undropped point
1158
1159   // will store 0 when the point is to be removed, 1 otherwise
1160   TColStd_Array1OfInteger tmpParam(firstElem, lastElem);
1161   for (i = firstElem; i<=lastElem ; i++) tmpParam.SetValue(i,1);
1162   Standard_Real DistMin2 = RealLast();
1163   gp_Pnt Prev = points->Value (lastValid);
1164   gp_Pnt Curr;
1165   for (i = firstElem + 1; i <= lastElem ; i ++) {
1166     Curr = points->Value(i);
1167     Standard_Real CurDist2 = Prev.SquareDistance(Curr);
1168     if (CurDist2 < gp::Resolution()) {  // test 0
1169       nbPntDropped ++;
1170       if ( i == lastElem ) tmpParam.SetValue(lastValid, 0); // last point kept
1171       else tmpParam.SetValue(i, 0);    // current dropped, lastValid unchanged
1172     } else {
1173       if (CurDist2 < DistMin2) 
1174         DistMin2 = CurDist2;
1175       // lastValid becomes the current (i.e. i)
1176       lastValid = i;
1177       Prev = Curr;
1178     }
1179   }
1180   if (DistMin2 < RealLast())
1181     preci = 0.9 * Sqrt (DistMin2); // preci est la distance min entre les points on la reduit un peu
1182   if (nbPntDropped == 0)
1183     return;
1184
1185 #ifdef DEBUG
1186   cout << "Warning : removing 3d points for interpolation" << endl;
1187 #endif
1188   // Build new HArrays
1189   Standard_Integer newLast = lastElem - nbPntDropped;
1190   if ((newLast - firstElem + 1)  < 2) {
1191 #ifdef DEBUG
1192     cout << "Too many degenerated points for 3D interpolation" << endl;
1193 #endif
1194     return;
1195   }
1196   Handle(TColgp_HArray1OfPnt) newPnts = 
1197     new TColgp_HArray1OfPnt(firstElem, newLast);
1198   Handle(TColStd_HArray1OfReal) newParams =
1199     new TColStd_HArray1OfReal(firstElem, newLast);
1200   Standard_Integer newCurr = 1;
1201   for (i = firstElem; i<= lastElem ; i++) {
1202     if (tmpParam.Value(i) == 1) { 
1203       newPnts->SetValue(newCurr, points->Value(i));
1204       newParams->SetValue(newCurr, params->Value(i));
1205       newCurr ++;
1206     }
1207   }
1208   points = newPnts;
1209   params = newParams;
1210   // on la reduit un peu
1211 }
1212
1213 //=======================================================================
1214 //function : CheckPoints2d
1215 //purpose  : 
1216 //=======================================================================
1217
1218  void ShapeConstruct_ProjectCurveOnSurface::CheckPoints2d(Handle(TColgp_HArray1OfPnt2d)& points,
1219                                                           Handle(TColStd_HArray1OfReal)& params,
1220                                                           Standard_Real& preci) const
1221 {
1222   Standard_Integer firstElem = points->Lower();
1223   Standard_Integer lastElem  = points->Upper();
1224   Standard_Integer i;
1225   Standard_Integer nbPntDropped = 0;
1226   Standard_Integer lastValid = firstElem; // indice of last undropped point
1227
1228   // will store 0 when the point is to be removed, 1 otherwise
1229   TColStd_Array1OfInteger tmpParam(firstElem, lastElem);
1230   for (i = firstElem; i<=lastElem ; i++) {
1231     tmpParam.SetValue(i,1);
1232   }
1233   Standard_Real DistMin2 = RealLast();
1234   gp_Pnt2d Prev = points->Value(lastValid);
1235   gp_Pnt2d Curr;
1236   for (i = firstElem + 1; i<=lastElem ; i++) {
1237     Curr = points->Value(i);
1238     Standard_Real CurDist2 = Prev.SquareDistance(Curr);
1239     if (CurDist2 < gp::Resolution()) {  // test 0
1240       nbPntDropped ++;
1241       if ( i == lastElem ) tmpParam.SetValue(lastValid, 0); // last point kept
1242       else tmpParam.SetValue(i, 0);    // current dropped, lastValid unchanged
1243     } else {
1244       if (CurDist2 < DistMin2) 
1245         DistMin2 = CurDist2;
1246       // lastValid becomes the current (i.e. i)
1247       lastValid = i;
1248       Prev = Curr;
1249     }
1250   }
1251   if (DistMin2 < RealLast())
1252     preci = 0.9 * Sqrt (DistMin2);
1253   if (nbPntDropped == 0)
1254     return;
1255
1256 #ifdef DEBUG
1257   cout << "Warning : removing 2d points for interpolation" << endl;
1258 #endif
1259   // Build new HArrays
1260   Standard_Integer newLast = lastElem - nbPntDropped;
1261   if ((newLast - firstElem + 1)  < 2) {
1262 #ifdef DEBUG
1263     cout << "Too many degenerated points for 2D interpolation" << endl;
1264 #endif
1265     //pdn 12.02.99 S4135 Creating pcurve with minimal length.
1266     tmpParam.SetValue(firstElem,1);
1267     tmpParam.SetValue(lastElem,1);
1268     gp_XY  lastPnt = points->Value(lastElem).XY();
1269     lastPnt.Add(gp_XY(preci,preci));
1270     points->SetValue(lastElem,lastPnt);
1271     newLast = firstElem+1;
1272     //return;
1273   }
1274   Handle(TColgp_HArray1OfPnt2d) newPnts = 
1275     new TColgp_HArray1OfPnt2d(firstElem, newLast);
1276   Handle(TColStd_HArray1OfReal) newParams =
1277     new TColStd_HArray1OfReal(firstElem, newLast);
1278   Standard_Integer newCurr = 1;
1279   for (i = firstElem; i <= lastElem ; i++) {
1280     if (tmpParam.Value(i) == 1) { 
1281 #ifdef DEBUG
1282       cout << "Point " << i << " : " << points->Value(i).X() << " " << points->Value(i).Y() << " at param " <<  params->Value(i) << endl;
1283 #endif
1284       newPnts->SetValue(newCurr, points->Value(i));
1285       newParams->SetValue(newCurr, params->Value(i));
1286       newCurr ++;
1287     }
1288     else {
1289 #ifdef DEBUG
1290       cout << "Removed " << i << " : " << points->Value(i).X() << " " << points->Value(i).Y() << " at param " <<  params->Value(i) << endl;
1291 #endif
1292     }
1293   }
1294   points = newPnts;
1295   params = newParams;
1296 }
1297
1298 //=======================================================================
1299 //function : IsAnIsoparametric
1300 //purpose  : 
1301 //=======================================================================
1302 //:S4030: modified for optimization
1303 //:p9 abv 11 Mar 99: PRO7226 #489490: find nearest boundary instead of first one
1304
1305  Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::IsAnIsoparametric(const Standard_Integer nbrPnt,
1306                                                                           const TColgp_Array1OfPnt& points,
1307                                                                           const TColStd_Array1OfReal& params,
1308                                                                           Standard_Boolean& isoTypeU,
1309                                                                           Standard_Boolean& p1OnIso,
1310                                                                           gp_Pnt2d& valueP1,
1311                                                                           Standard_Boolean& p2OnIso,
1312                                                                           gp_Pnt2d& valueP2,
1313                                                                           Standard_Boolean& isoPar2d3d,
1314                                                                           Handle(Geom_Curve)& cIso,
1315                                                                           Standard_Real& t1,
1316                                                                           Standard_Real& t2,
1317                                                                           TColStd_Array1OfReal& pout) const
1318 {
1319   try {    // RAJOUT
1320     OCC_CATCH_SIGNALS
1321     
1322   Standard_Real prec = Precision::Confusion();//myPreci;
1323     
1324   Standard_Boolean isoParam = Standard_False;
1325   isoPar2d3d = Standard_False;
1326   
1327   Standard_Real U1, U2, V1, V2;
1328   mySurf->Bounds(U1, U2, V1, V2);
1329   
1330   if ( mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
1331     Handle(Geom_RectangularTrimmedSurface) sTrim =
1332       Handle(Geom_RectangularTrimmedSurface)::DownCast(mySurf->Surface());
1333     sTrim->Bounds(U1, U2, V1, V2);
1334   }
1335   
1336   gp_Pnt pt;
1337   Standard_Integer mpt[2]; mpt[0] = mpt[1] = 0;
1338   Standard_Real t, tpar[2], isoValue=0.;
1339   Standard_Real mindist2;
1340   Standard_Real mind2[2];
1341   mindist2 = mind2[0] = mind2[1] = 4*prec*prec;
1342   
1343   p1OnIso = Standard_False;
1344   p2OnIso = Standard_False;
1345   const Bnd_Box* aBox = 0;
1346   
1347   for (Standard_Integer j=1; (j<=4) /*&& !isoParam*/; j++) {
1348     Standard_Real isoVal=0.;
1349     Standard_Boolean isoU=Standard_False; //szv#4:S4163:12Mar99 `isoU` must be Standard_Boolean
1350     Handle(Geom_Curve) cI;
1351     Standard_Real tt1, tt2;
1352     
1353     if      (j == 1 ) {
1354       if (Precision::IsInfinite(U1)) continue;
1355       cI = mySurf->UIso(U1);
1356       isoU = Standard_True;
1357       isoVal = U1;
1358       aBox = & mySurf->GetBoxUF();
1359     }
1360     else if (j == 2) {
1361       if (Precision::IsInfinite(U2)) continue;
1362       cI = mySurf->UIso(U2);
1363       isoU = Standard_True;
1364       isoVal = U2;
1365       aBox = & mySurf->GetBoxUL();
1366     }
1367     else if (j == 3) {
1368       if (Precision::IsInfinite(V1)) continue;
1369       cI = mySurf->VIso(V1);
1370       isoU = Standard_False;
1371       isoVal = V1;
1372       aBox = & mySurf->GetBoxVF();
1373     }
1374     else if (j == 4) {
1375       if (Precision::IsInfinite(V2)) continue;
1376       cI = mySurf->VIso(V2);
1377       isoU = Standard_False;
1378       isoVal = V2;
1379       aBox = & mySurf->GetBoxVL();
1380     }
1381     if(cI.IsNull())
1382       continue;
1383     
1384     if (isoU)  {  tt1 = V1;  tt2 = V2;  }
1385     else       {  tt1 = U1;  tt2 = U2;  }
1386
1387     gp_Pnt ext1, ext2;
1388     cI->D0(tt1, ext1);
1389     cI->D0(tt2, ext2);
1390
1391 // PATCH CKY 9-JUL-1998 : protection contre singularite
1392     gp_Pnt extmi;
1393     cI->D0( (tt1+tt2)/2,extmi);
1394     if (ext1.IsEqual(ext2,prec) && ext1.IsEqual(extmi,prec)) continue;
1395
1396     Standard_Boolean PtEQext1 = Standard_False;
1397     Standard_Boolean PtEQext2 = Standard_False;
1398
1399     Standard_Real currd2[2], tp[2];
1400     Standard_Integer mp[2];
1401     
1402     for (Standard_Integer i=0; i<2; i++) {
1403       mp[i] = 0;
1404       Standard_Integer k = (i == 0 ? 1 : nbrPnt);
1405
1406       // si ext1 == ext2 => valueP1 == valueP2 => vect null plus tard
1407       currd2[i] = points(k).SquareDistance ( ext1 );
1408       if ( currd2[i] <= prec*prec && !PtEQext1) {
1409         mp[i] = 1;
1410         tp[i] = tt1;
1411         PtEQext1 = Standard_True;
1412         continue;
1413       } 
1414       
1415       currd2[i] = points(k).SquareDistance ( ext2 );
1416       if ( currd2[i] <= prec*prec && !PtEQext2) {
1417         mp[i] = 2;
1418         tp[i] = tt2;
1419         PtEQext2 = Standard_True;
1420         continue;
1421       }  
1422       
1423       // On evite de projecter sur un iso degenere
1424       // on doit egalement le faire pour l apex du cone
1425       if (mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface)) && !isoU) {
1426         continue;
1427       }
1428       
1429       if(aBox->IsOut(points(k))) continue;
1430       
1431       Standard_Real Cf = cI->FirstParameter();
1432       Standard_Real Cl = cI->LastParameter();
1433       if (Precision::IsInfinite(Cf))  Cf = -1000;
1434       if (Precision::IsInfinite(Cl))  Cl = +1000;
1435
1436       ShapeAnalysis_Curve sac;
1437       Standard_Real dist = sac.Project (cI,points(k),prec,pt,t,Cf,Cl);
1438       currd2[i] = dist * dist;
1439       if ((dist <= prec) && (t>= Cf) && (t<=Cl)) {
1440         mp[i]  = 3;
1441         tp[i] = t;
1442       }
1443     }
1444
1445     //:e7 abv 21 Apr 98: ProSTEP TR8, r0501_pe #56679:
1446     // avoid possible null-length curves
1447     if ( mp[0] >0 && mp[1] >0 &&
1448          Abs ( tp[0] - tp[1] ) < Precision::PConfusion() ) continue;
1449     
1450     
1451     if (mp[0] > 0 && 
1452         ( ! p1OnIso || currd2[0] < mind2[0] ) ) {
1453       p1OnIso = Standard_True;
1454       mind2[0] = currd2[0];
1455       if (isoU) valueP1.SetCoord(isoVal, tp[0]);
1456       else      valueP1.SetCoord(tp[0], isoVal);
1457     }
1458
1459     if (mp[1] > 0 &&
1460         ( ! p2OnIso || currd2[1] < mind2[1] ) ) {
1461       p2OnIso = Standard_True;
1462       mind2[1] = currd2[1];
1463       if (isoU) valueP2.SetCoord(isoVal, tp[1]);
1464       else      valueP2.SetCoord(tp[1], isoVal);
1465     }
1466
1467     if ( mp[0] <=0 || mp[1] <=0 ) continue;
1468
1469     Standard_Real md2 = currd2[0] + currd2[1];
1470     if ( mindist2 <= md2 ) continue;
1471     
1472     mindist2 = md2;
1473     mpt[0] = mp[0];
1474     mpt[1] = mp[1];
1475     tpar[0] = tp[0];
1476     tpar[1] = tp[1];
1477     isoTypeU = isoU;
1478     isoValue = isoVal;
1479     cIso = cI;
1480     t1 = tt1;
1481     t2 = tt2;
1482   }
1483     
1484   // probablely it concerns an isoparametrics
1485   if ( mpt[0] >0 && mpt[1] >0 ) {
1486
1487     p1OnIso = p2OnIso = Standard_True;
1488     if (isoTypeU) {
1489       valueP1.SetCoord(isoValue, tpar[0]);
1490       valueP2.SetCoord(isoValue, tpar[1]);
1491     }
1492     else {
1493       valueP1.SetCoord(tpar[0], isoValue);
1494       valueP2.SetCoord(tpar[1], isoValue);
1495     }
1496
1497     if ( mpt[0] != 3 && mpt[1] != 3 ) {
1498       isoPar2d3d = Standard_True;
1499       for (Standard_Integer i=2; i < nbrPnt && isoPar2d3d; i++){
1500         if (tpar[1] > tpar[0])  t = params(i);
1501         else                    t = t1+t2-params(i);
1502         cIso->D0(t, pt);
1503         if (!points(i).IsEqual(pt, prec)) isoPar2d3d = Standard_False;
1504       }
1505     }
1506
1507     if (isoPar2d3d) isoParam = Standard_True;
1508     else {
1509       Standard_Real prevParam = tpar[0];
1510       Standard_Real Cf, Cl;
1511       Standard_Boolean isoByDistance = Standard_True;
1512       Cf = cIso->FirstParameter();
1513       Cl = cIso->LastParameter();
1514       if (Precision::IsInfinite(Cf))  Cf = -1000;
1515       if (Precision::IsInfinite(Cl))  Cl = +1000;
1516         
1517       ShapeAnalysis_Curve sac;
1518       for (Standard_Integer i=2; i < nbrPnt && isoByDistance; i++) {
1519         Standard_Real dist = sac.NextProject (prevParam,cIso,points(i),
1520                                               prec,pt,t,Cf,Cl,
1521                                               Standard_False); //:j8 abv 10.12.98: TR10 r0501_db.stp #9423: avoid adjusting to ends
1522         prevParam = t;
1523         pout(i)=t;
1524         if( (dist > prec) || (t < Cf) || (t > Cl) ) 
1525           isoByDistance = Standard_False;
1526       }
1527       if (isoByDistance) isoParam = Standard_True;
1528     }
1529   }
1530 /*  if (!isoParam) {    CKY 29-mai-1997 : garder tout ce qu on peut ?
1531     p1OnIso = Standard_False;
1532     p2OnIso = Standard_False;
1533   }  */
1534   return isoParam;
1535   }  // RAJOUT
1536   catch(Standard_Failure) {
1537 //  pb : on affiche ce qu on peut
1538 #ifdef DEBUG
1539     for (Standard_Integer numpnt = 1; numpnt <= nbrPnt; numpnt ++) {
1540       cout<<"["<<numpnt<<"]param="<<params(numpnt)<<" point=("<<
1541         points(numpnt).X()<<"  "<<points(numpnt).Y()<<"  "<<points(numpnt).Z()<<")"<<endl;
1542     }
1543 #endif
1544 #ifdef DEB //:s5
1545     cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::IsAnIsoparametric(): Exception: ";
1546     Standard_Failure::Caught()->Print(cout); cout << endl;
1547 #endif
1548     return Standard_False;
1549   }
1550   return Standard_False;  // ramasse-miette
1551 }
1552
1553 /* S4135 : BestExtremum is commented after IsAnIsoparametric works with Precision::Confusion()
1554 //=======================================================================
1555 //function : BestExtremum
1556 //purpose  : auxiliaire prenant le meilleur extremum si ISO car doute possible
1557 //=======================================================================
1558
1559  gp_Pnt2d ShapeConstruct_ProjectCurveOnSurface::BestExtremum(const gp_Pnt2d& P2iso,const gp_Pnt& P3ext,const gp_Pnt& P3next) const
1560 {
1561 //  P2iso a ete calcule depuis P3ext sur une iso externe de la surface
1562 //  En principe bon mais circularite possible ... et IsU/VClosed faillible
1563 //    (si baillement 1e-4 ou 1e-5, on est dedans !). DONC
1564 //  1/ on privilegie l iso mais a tout hasard on verifie si Surf meilleur
1565 //  2/ si iso, attention a la circularite (cas limite)
1566
1567 //  NB : si isoParam, on suppose que P2iso est bon (car il y en a 2). A voir...
1568
1569 //  D abord, calcul p2ext depuis la surface. choix surface/iso
1570   return P2iso;
1571   Standard_Real prec = Precision::Confusion();//myPreci;
1572   gp_Pnt2d P2cal = mySurf->ValueOfUV(P3ext, prec);
1573   gp_Pnt   P3cal = mySurf->Value (P2cal);
1574   Standard_Real dcal = P3ext.Distance (P3cal);
1575   Standard_Real dnxt = P3ext.Distance (P3next);
1576   if (dcal > dnxt) return P2iso;    // en fait protection sur BUG (PRO8468)
1577
1578 //  On choisit entre P2iso et P2cal, le plus proche de P2next ... !!!
1579   gp_Pnt2d P2next = mySurf->ValueOfUV(P3next, prec);
1580   if (P2next.Distance(P2cal) < P2next.Distance(P2iso)) return P2cal;
1581   return P2iso;
1582 }
1583 */