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