projects / occt.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
0028966: Coding Rules - remove Adaptor2d_HCurve2d, Adaptor3d_HCurve and Adaptor3d_HSu...
[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
7fd59977 31
42cf5bc1 32#include <Approx_CurveOnSurface.hxx>
33#include <Geom2d_BSplineCurve.hxx>
9a6ea9c4 34#include <Geom2d_Circle.hxx>
42cf5bc1 35#include <Geom2d_Curve.hxx>
9a6ea9c4 36#include <Geom2d_Ellipse.hxx>
37#include <Geom2d_Hyperbola.hxx>
42cf5bc1 38#include <Geom2d_Line.hxx>
9a6ea9c4 39#include <Geom2d_Parabola.hxx>
7fd59977 40#include <Geom2d_TrimmedCurve.hxx>
42cf5bc1 41#include <Geom2dAdaptor.hxx>
42#include <Geom2dAPI_Interpolate.hxx>
43#include <Geom_BezierSurface.hxx>
44#include <Geom_BoundedCurve.hxx>
7fd59977 45#include <Geom_BSplineCurve.hxx>
42cf5bc1 46#include <Geom_BSplineSurface.hxx>
47#include <Geom_Curve.hxx>
7fd59977 48#include <Geom_OffsetSurface.hxx>
49#include <Geom_Plane.hxx>
42cf5bc1 50#include <Geom_RectangularTrimmedSurface.hxx>
51#include <Geom_SphericalSurface.hxx>
52#include <Geom_Surface.hxx>
53#include <Geom_SurfaceOfLinearExtrusion.hxx>
54#include <Geom_TrimmedCurve.hxx>
c22b52d6 55#include <GeomAdaptor_Curve.hxx>
56#include <GeomAdaptor_Surface.hxx>
42cf5bc1 57#include <GeomAPI_Interpolate.hxx>
58#include <GeomAPI_PointsToBSpline.hxx>
59#include <GeomProjLib.hxx>
60#include <gp_Pnt2d.hxx>
2c26a53d 61#include <ElCLib.hxx>
42cf5bc1 62#include <NCollection_Sequence.hxx>
63#include <Precision.hxx>
64#include <ProjLib_CompProjectedCurve.hxx>
65#include <ProjLib_HCompProjectedCurve.hxx>
66#include <ProjLib_ProjectedCurve.hxx>
02fd709b 67#include <ShapeAnalysis.hxx>
7fd59977 68#include <ShapeAnalysis_Curve.hxx>
69#include <ShapeAnalysis_Surface.hxx>
42cf5bc1 70#include <ShapeConstruct_ProjectCurveOnSurface.hxx>
7fd59977 71#include <ShapeExtend.hxx>
42cf5bc1 72#include <Standard_ErrorHandler.hxx>
73#include <Standard_Failure.hxx>
74#include <Standard_Type.hxx>
75#include <TColgp_Array1OfPnt.hxx>
76#include <TColStd_Array1OfInteger.hxx>
2c26a53d 77#include <IntRes2d_Domain.hxx>
78#include <IntCurve_IntConicConic.hxx>
7fd59977 79
42cf5bc1 80#include <algorithm>
25e59720 81IMPLEMENT_STANDARD_RTTIEXT(ShapeConstruct_ProjectCurveOnSurface,Standard_Transient)
92efcf78 82
7fd59977 83#define NCONTROL 23
b311480e 84
2c26a53d 85
86static void AdjustSecondPointToFirstPoint(const gp_Pnt2d& theFirstPoint,
87 gp_Pnt2d& theSecondPoint,
88 const Handle(Geom_Surface)& theSurf)
89{
90 if (theSurf->IsUPeriodic())
91 {
92 Standard_Real UPeriod = theSurf->UPeriod();
93 Standard_Real NewU = ElCLib::InPeriod(theSecondPoint.X(),
94 theFirstPoint.X() - UPeriod/2,
95 theFirstPoint.X() + UPeriod/2);
96 theSecondPoint.SetX(NewU);
97 }
98 if (theSurf->IsVPeriodic())
99 {
100 Standard_Real VPeriod = theSurf->VPeriod();
101 Standard_Real NewV = ElCLib::InPeriod(theSecondPoint.Y(),
102 theFirstPoint.Y() - VPeriod/2,
103 theFirstPoint.Y() + VPeriod/2);
104 theSecondPoint.SetY(NewV);
105 }
106}
107
108
7fd59977 109//=======================================================================
110//function : ShapeConstruct_ProjectCurveOnSurface
111//purpose :
112//=======================================================================
113
114ShapeConstruct_ProjectCurveOnSurface::ShapeConstruct_ProjectCurveOnSurface()
115{
116 myPreci = Precision::Confusion();
117 myBuild = Standard_False;
118 myAdjustOverDegen = 1; //:c0 //szv#4:S4163:12Mar99 was boolean
119 myNbCashe = 0; //:q9
120}
121
122//=======================================================================
123//function : Init
124//purpose :
125//=======================================================================
126
127 void ShapeConstruct_ProjectCurveOnSurface::Init(const Handle(Geom_Surface)& surf,const Standard_Real preci)
128{
129 Init (new ShapeAnalysis_Surface (surf), preci);
130}
131
132//=======================================================================
133//function : Init
134//purpose :
135//=======================================================================
136
137 void ShapeConstruct_ProjectCurveOnSurface::Init(const Handle(ShapeAnalysis_Surface)& surf,const Standard_Real preci)
138{
139 SetSurface (surf);
140 SetPrecision (preci);
141}
142
143//=======================================================================
144//function : SetSurface
145//purpose :
146//=======================================================================
147
148 void ShapeConstruct_ProjectCurveOnSurface::SetSurface(const Handle(Geom_Surface)& surf)
149{
150 SetSurface (new ShapeAnalysis_Surface (surf));
151}
152
153//=======================================================================
154//function : SetSurface
155//purpose :
156//=======================================================================
157
158 void ShapeConstruct_ProjectCurveOnSurface::SetSurface(const Handle(ShapeAnalysis_Surface)& surf)
159{
160 if ( mySurf == surf ) return;
161 mySurf = surf;
162 myNbCashe = 0; //:q9
163}
164
165//=======================================================================
166//function : SetPrecision
167//purpose :
168//=======================================================================
169
170 void ShapeConstruct_ProjectCurveOnSurface::SetPrecision(const Standard_Real preci)
171{
172 myPreci = preci;
173}
174
175//=======================================================================
176//function : BuildCurveMode
177//purpose :
178//=======================================================================
179
180 Standard_Boolean& ShapeConstruct_ProjectCurveOnSurface::BuildCurveMode()
181{
182 return myBuild;
183}
184
185//=======================================================================
186//function : AdjustOverDegenMode
187//purpose :
188//=======================================================================
189//:c0
190
191//szv#4:S4163:12Mar99 was Boolean
192 Standard_Integer& ShapeConstruct_ProjectCurveOnSurface::AdjustOverDegenMode()
193 {
194 return myAdjustOverDegen;
195 }
196
7fd59977 197//=======================================================================
198//function : Status
199//purpose :
200//=======================================================================
201
9fd2d2c3 202 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::Status (const ShapeExtend_Status theStatus) const
7fd59977 203{
9fd2d2c3 204 return ShapeExtend::DecodeStatus (myStatus, theStatus);
7fd59977 205}
206
207//=======================================================================
208//function : Perform
209//purpose :
210//=======================================================================
7fd59977 211Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::Perform (Handle(Geom_Curve)& c3d,
212 const Standard_Real First,
213 const Standard_Real Last,
214 Handle(Geom2d_Curve)& c2d,
2c26a53d 215 const Standard_Real TolFirst,
216 const Standard_Real TolLast)
7fd59977 217{
218 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_OK);
219 //Standard_Boolean OK = Standard_True; //szv#4:S4163:12Mar99 not needed
220
221 if (mySurf.IsNull()) {
222 c2d.Nullify();
223 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_FAIL1);
224 return Standard_False;
225 }
7fd59977 226// Projection Analytique
227 Handle(Geom_Curve) crv3dtrim = c3d;
228 if ( ! c3d->IsKind(STANDARD_TYPE(Geom_BoundedCurve)) )
229 crv3dtrim = new Geom_TrimmedCurve ( c3d, First, Last );
230 c2d = ProjectAnalytic ( crv3dtrim );
231 if (!c2d.IsNull()) {
232 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE1);
233 return Standard_True;
234 }
235
236// Projection par approximation
237
238 // discretize the 3d curve
239
240 Standard_Integer nbrPnt;
241
242// $$$$ :92 abv 28 Jan 98 see PRO10107, big BSplineCurve C0
243 Standard_Integer nbPini = NCONTROL; // as in BRepCheck_Edge (RLN/Nijni)
244 // 20; // number of points for interpolation, should be "parametric dependent"
245
246 //:92 abv 28 Jan 98: if curve is BSpline with many intervals,
247 // increase number of points to provide at least Degree()+1 points per interval
248 Handle(Geom_BSplineCurve) bspl;
249 if ( c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve)) ) {
250 Handle(Geom_TrimmedCurve) ctrim = Handle(Geom_TrimmedCurve)::DownCast(c3d);
251 bspl = Handle(Geom_BSplineCurve)::DownCast ( ctrim->BasisCurve() );
252 }
253 else bspl = Handle(Geom_BSplineCurve)::DownCast ( c3d );
254 if ( ! bspl.IsNull() ) {
255 Standard_Integer nint = 0;
256 for ( Standard_Integer i=1; i < bspl->NbKnots(); i++ )
2c26a53d 257 {
258 if ( bspl->Knot(i+1) > First && bspl->Knot(i) < Last )
259 nint++;
260 }
7fd59977 261 Standard_Integer minPnt = nint * ( bspl->Degree() + 1 );
262 while ( nbPini < minPnt ) nbPini += NCONTROL - 1;
0797d9d3 263#ifdef OCCT_DEBUG
7fd59977 264 if ( nbPini > NCONTROL )
04232180 265 std::cout << "Warning: number of points for projecting is " << nbPini << std::endl;
7fd59977 266#endif
267 }
268
269// $$$$ end :92 (big BSplineCurve C0)
270
271 // this number should be "parametric dependent"
2c26a53d 272 TColgp_SequenceOfPnt points;
273 TColStd_SequenceOfReal params;
9a6ea9c4 274 NCollection_Sequence<Standard_Real> aKnotCoeffs;
7fd59977 275 gp_Pnt p3d;
9a6ea9c4 276 Standard_Integer iPnt;
277
278 // In case of bspline compute parametrization speed on each
279 // knot interval inside [aFirstParam, aLastParam].
280 // If quotient = (MaxSpeed / MinSpeed) >= aMaxQuotientCoeff then
281 // use PerformByProjLib algorithm.
282 if(!bspl.IsNull())
283 {
9a6ea9c4 284 Standard_Real aFirstParam = First; // First parameter of current interval.
285 Standard_Real aLastParam = Last; // Last parameter of current interval.
286
287 // First index computation.
288 Standard_Integer anIdx = 1;
289 for(; anIdx <= bspl->NbKnots() && aFirstParam < Last; anIdx++)
290 {
291 if(bspl->Knot(anIdx) > First)
292 {
293 break;
294 }
295 }
296
6a2ee094 297 GeomAdaptor_Curve aC3DAdaptor(c3d);
67285996 298 Standard_Real aMinParSpeed = Precision::Infinite(); // Minimal parameterization speed.
9a6ea9c4 299 for(; anIdx <= bspl->NbKnots() && aFirstParam < Last; anIdx++)
300 {
301 // Fill current knot interval.
302 aLastParam = Min(Last, bspl->Knot(anIdx));
6a2ee094 303 Standard_Integer aNbIntPnts = NCONTROL;
304 // Number of inner points is adapted according to the length of the interval
305 // to avoid a lot of calculations on small range of parameters.
306 if (anIdx > 1)
307 {
308 const Standard_Real aLenThres = 1.e-2;
309 const Standard_Real aLenRatio =
310 (aLastParam - aFirstParam) / (bspl->Knot(anIdx) - bspl->Knot(anIdx - 1));
311 if (aLenRatio < aLenThres)
312 {
313 aNbIntPnts = Standard_Integer(aLenRatio / aLenThres * aNbIntPnts);
314 if (aNbIntPnts < 2)
315 aNbIntPnts = 2;
316 }
317 }
9a6ea9c4 318 Standard_Real aStep = (aLastParam - aFirstParam) / (aNbIntPnts - 1);
319 Standard_Integer anIntIdx;
320 gp_Pnt p3d1, p3d2;
6a2ee094 321 // Start filling from first point.
322 aC3DAdaptor.D0(aFirstParam, p3d1);
323
9a6ea9c4 324 Standard_Real aLength3d = 0.0;
6a2ee094 325 for(anIntIdx = 1; anIntIdx < aNbIntPnts; anIntIdx++)
9a6ea9c4 326 {
6a2ee094 327 Standard_Real aParam = aFirstParam + aStep * anIntIdx;
328 aC3DAdaptor.D0 (aParam, p3d2);
67285996 329 const Standard_Real aDist = p3d2.Distance(p3d1);
330
331 aLength3d += aDist;
6a2ee094 332 p3d1 = p3d2;
67285996 333
334 aMinParSpeed = Min(aMinParSpeed, aDist / aStep);
9a6ea9c4 335 }
6a2ee094 336 const Standard_Real aCoeff = aLength3d / (aLastParam - aFirstParam);
337 if (Abs(aCoeff) > gp::Resolution())
338 aKnotCoeffs.Append(aCoeff);
9a6ea9c4 339 aFirstParam = aLastParam;
340 }
341
342 Standard_Real anEvenlyCoeff = 0;
343 if (aKnotCoeffs.Size() > 0)
344 {
345 anEvenlyCoeff = *std::max_element(aKnotCoeffs.begin(), aKnotCoeffs.end()) /
346 *std::min_element(aKnotCoeffs.begin(), aKnotCoeffs.end());
347 }
348
349 const Standard_Real aMaxQuotientCoeff = 1500.0;
67285996 350 if (anEvenlyCoeff > aMaxQuotientCoeff &&
351 aMinParSpeed > Precision::Confusion() )
9a6ea9c4 352 {
353 PerformByProjLib(c3d, First, Last, c2d);
354 // PerformByProjLib fail detection:
355 if (!c2d.IsNull())
356 {
357 return Status (ShapeExtend_DONE);
358 }
359 }
360 }
361
7fd59977 362 Standard_Real deltaT, t;
363 deltaT = (Last - First) / (nbPini-1);
364 nbrPnt = nbPini;
9a6ea9c4 365 for (iPnt = 1; iPnt <= nbPini; iPnt ++)
366 {
7fd59977 367 if (iPnt == 1) t = First;
368 else if (iPnt == nbPini) t = Last;
369 else t = First + (iPnt - 1) * deltaT;
370
371 c3d->D0 (t, p3d);
2c26a53d 372 points.Append(p3d);
373 params.Append(t);
7fd59977 374 }
375
9a6ea9c4 376 // CALCUL par approximation
2c26a53d 377 TColgp_SequenceOfPnt2d pnt2d;
378 ApproxPCurve (nbrPnt,c3d,TolFirst,TolLast,
379 points,params,pnt2d,c2d); //szv#4:S4163:12Mar99 OK not needed
380 nbPini = points.Length();
7fd59977 381 if (!c2d.IsNull()) {
382 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE2);
383 return Standard_True;
384 }// cas particulier d iso
385
386// INTERPOLATION du resultat
387
388 if ( myBuild ) {
389 Handle(TColgp_HArray1OfPnt) thePnts = new TColgp_HArray1OfPnt (1, nbPini);
390 Handle(TColStd_HArray1OfReal) theParams = new TColStd_HArray1OfReal(1, nbPini);
391 for (iPnt = 1; iPnt <= nbPini ; iPnt ++) {
392 thePnts->SetValue(iPnt, points(iPnt));
393 theParams->SetValue(iPnt, params(iPnt));
394 }
395
396 Handle(Geom_Curve) newc3d = InterpolateCurve3d (nbPini,thePnts,theParams, c3d);
397 if ( newc3d.IsNull() ) myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_FAIL2);
398 else {
399 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE3);
400 c3d = newc3d;
401 }
402 }
403
404 Handle(TColgp_HArray1OfPnt2d) thePnts2d = new TColgp_HArray1OfPnt2d(1, nbPini);
405 Handle(TColStd_HArray1OfReal) theParams2d = new TColStd_HArray1OfReal(1, nbPini);
406 for (iPnt = 1; iPnt <= nbPini ; iPnt ++) {
407 theParams2d->SetValue(iPnt, params(iPnt));
408 thePnts2d->SetValue(iPnt, pnt2d(iPnt));
409 }
7fd59977 410 c2d = InterpolatePCurve (nbPini, thePnts2d, theParams2d, c3d);
b7353545 411
412 // Faut-il aussi reprendre la C3D ?
7fd59977 413 myStatus |= ShapeExtend::EncodeStatus (c2d.IsNull() ? ShapeExtend_FAIL1 : ShapeExtend_DONE2);
414 return Status (ShapeExtend_DONE);
415}
416
417//=======================================================================
418//function : PerformByProjLib
419//purpose :
420//=======================================================================
421
422Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::PerformByProjLib(Handle(Geom_Curve)& c3d,
423 const Standard_Real First,
424 const Standard_Real Last,
425 Handle(Geom2d_Curve)& c2d,
9a6ea9c4 426 const GeomAbs_Shape /*continuity*/,
427 const Standard_Integer /*maxdeg */,
428 const Standard_Integer /*nbinterval */)
7fd59977 429{
430 //Standard_Boolean OK = Standard_True; //szv#4:S4163:12Mar99 unused
431 c2d.Nullify();
432 if (mySurf.IsNull()) {
433 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL1);
434 return Standard_False;
435 }
436
9a6ea9c4 437 try
438 {
7fd59977 439 OCC_CATCH_SIGNALS
c22b52d6 440 Handle(GeomAdaptor_Surface) GAS = mySurf->Adaptor3d();
441 Handle(GeomAdaptor_Curve) GAC = new GeomAdaptor_Curve (c3d,First,Last);
9a6ea9c4 442 ProjLib_ProjectedCurve Projector(GAS, GAC);
443
444 switch (Projector.GetType())
445 {
446 case GeomAbs_Line :
447 c2d = new Geom2d_Line(Projector.Line());
448 break;
449 case GeomAbs_Circle :
450 c2d = new Geom2d_Circle(Projector.Circle());
451 break;
452 case GeomAbs_Ellipse :
453 c2d = new Geom2d_Ellipse(Projector.Ellipse());
454 break;
455 case GeomAbs_Parabola :
456 c2d = new Geom2d_Parabola(Projector.Parabola());
457 break;
458 case GeomAbs_Hyperbola :
459 c2d = new Geom2d_Hyperbola(Projector.Hyperbola());
460 break;
461 case GeomAbs_BSplineCurve :
462 c2d = Projector.BSpline();
463 break;
1aec3320 464 default:
9a6ea9c4 465 // Not possible, handling added to avoid gcc warning.
466 break;
7fd59977 467 }
9a6ea9c4 468
469 if(c2d.IsNull())
470 {
7fd59977 471 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL2);
472 return Standard_False;
473 }
9a6ea9c4 474 else
475 {
7fd59977 476 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_DONE1);
477 return Standard_True;
478 }
479
480 }
9775fa61 481 catch(Standard_Failure const& anException) {
0797d9d3 482#ifdef OCCT_DEBUG
04232180 483 std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::PerformByProjLib(): Exception: ";
484 anException.Print(std::cout); std::cout << std::endl;
7fd59977 485#endif
9775fa61 486 (void)anException;
7fd59977 487 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL3);
488 c2d.Nullify();
489 }
490 return Standard_False;
491}
492
7fd59977 493//=======================================================================
494//function : ProjectAnalytic
495//purpose :
496//=======================================================================
497
498 Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::ProjectAnalytic(const Handle(Geom_Curve)& c3d) const
499{
500 Handle(Geom2d_Curve) result;
501
502 //:k1 abv 16 Dec 98: limit analytic cases by Plane surfaces only
503 // This is necessary for K4L since it fails on other surfaces
504 // when general method GeomProjLib::Curve2d() is used
505 // Projection is done as in BRep_Tool and BRepCheck_Edge
506 Handle(Geom_Surface) surf = mySurf->Surface();
507 Handle(Geom_Plane) Plane = Handle(Geom_Plane)::DownCast ( surf );
508 if ( Plane.IsNull() ) {
509 Handle(Geom_RectangularTrimmedSurface) RTS =
510 Handle(Geom_RectangularTrimmedSurface)::DownCast ( surf );
511 if ( ! RTS.IsNull() ) Plane = Handle(Geom_Plane)::DownCast ( RTS->BasisSurface() );
512 else {
513 Handle(Geom_OffsetSurface) OS =
514 Handle(Geom_OffsetSurface)::DownCast ( surf );
515 if ( ! OS.IsNull() )
516 Plane = Handle(Geom_Plane)::DownCast ( OS->BasisSurface() );
517 }
518 }
519 if ( ! Plane.IsNull() ) {
520 Handle(Geom_Curve) ProjOnPlane =
521 GeomProjLib::ProjectOnPlane (c3d, Plane,
522 Plane->Position().Direction(), Standard_True);
c22b52d6 523 Handle(GeomAdaptor_Curve) HC = new GeomAdaptor_Curve ( ProjOnPlane );
7fd59977 524 ProjLib_ProjectedCurve Proj ( mySurf->Adaptor3d(), HC );
525
526 result = Geom2dAdaptor::MakeCurve(Proj);
527 if ( result.IsNull() ) return result;
528 if ( result->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve)) ) {
529 Handle(Geom2d_TrimmedCurve) TC = Handle(Geom2d_TrimmedCurve)::DownCast ( result );
530 result = TC->BasisCurve();
531 }
63c629aa 532
7fd59977 533 return result;
534 }
535
536 return result;
537}
538
15b54261 539 //! Fix possible period jump and handle walking period parameter.
540 static Standard_Boolean fixPeriodictyTroubles(gp_Pnt2d *thePnt, // pointer to gp_Pnt2d[4] beginning
541 Standard_Integer theIdx, // Index of objective coord: 1 ~ X, 2 ~ Y
02fd709b 542 Standard_Real thePeriod, // Period on objective coord
543 Standard_Integer theSavedPoint, // Point number to choose period
544 Standard_Real theSavedParam) // Param from cashe to choose period
545{
546 Standard_Real aSavedParam;
547 Standard_Integer aSavedPoint;
548 Standard_Real aMinParam = 0.0, aMaxParam = thePeriod;
549 if (theSavedPoint < 0) {
550 // normalize to first period by default
551 aSavedParam = 0.5 * thePeriod;
552 aSavedPoint = 0;
553 }
554 else {
555 aSavedParam = theSavedParam;
556 aSavedPoint = theSavedPoint;
557 while (aMinParam > aSavedParam) {
558 aMinParam -= thePeriod;
559 aMaxParam -= thePeriod;
560 }
561 while (aMaxParam < aSavedParam) {
562 aMinParam += thePeriod;
563 aMaxParam += thePeriod;
564 }
565 }
15b54261 566
02fd709b 567 Standard_Real aFixIsoParam = aMinParam;
568 Standard_Boolean isIsoLine = Standard_False;
569 if (aMaxParam - aSavedParam < Precision::PConfusion() ||
570 aSavedParam - aMinParam < Precision::PConfusion()) {
571 aFixIsoParam = aSavedParam;
572 isIsoLine = Standard_True;
573 }
574 // normalize all coordinates to [aMinParam, aMaxParam)
575 for (Standard_Integer i = 0; i < 4; i++) {
576 Standard_Real aParam = thePnt[i].Coord(theIdx);
577 Standard_Real aShift = ShapeAnalysis::AdjustToPeriod(aParam, aMinParam, aMaxParam);
578 aParam += aShift;
579 // Walk over period coord -> not walking on another isoline in parameter space.
580 if (isIsoLine) {
581 if (aMaxParam - aParam < Precision::PConfusion() || aParam - aMinParam < Precision::PConfusion())
582 aParam = aFixIsoParam;
583 }
584 else {
585 if (aMaxParam - aParam < Precision::PConfusion())
586 aParam = aMaxParam;
587 if (aParam - aMinParam < Precision::PConfusion())
588 aParam = aMinParam;
589 }
15b54261 590
02fd709b 591 thePnt[i].SetCoord(theIdx, aParam);
592 }
15b54261 593
02fd709b 594 // find possible period jump and increasing or decreasing coordinates vector
595 Standard_Boolean isJump = Standard_False;
596 Standard_Real aPrevDiff = 0.0;
597 Standard_Real aSumDiff = 1.0;
598 for (Standard_Integer i = 0; i < 3; i++) {
599 Standard_Real aDiff = thePnt[i + 1].Coord(theIdx) - thePnt[i].Coord(theIdx);
600 if (aDiff < -Precision::PConfusion()) {
601 aSumDiff *= -1.0;
602 }
603 //if first derivative changes its sign then period jump may exists in this place
604 if (aDiff * aPrevDiff < -Precision::PConfusion()) {
605 isJump = Standard_True;
606 }
607 aPrevDiff = aDiff;
608 }
15b54261 609
02fd709b 610 if (!isJump)
611 return Standard_False;
612
613 if (aSumDiff > 0) { // decreasing sequence (parameters decrease twice(--) and one period jump(+))
614 for (Standard_Integer i = aSavedPoint; i > 0; i--)
615 if (thePnt[i].Coord(theIdx) > thePnt[i - 1].Coord(theIdx)) {
616 thePnt[i - 1].SetCoord(theIdx, thePnt[i - 1].Coord(theIdx) + thePeriod);
617 }
618 for (Standard_Integer i = aSavedPoint; i < 3; i++)
619 if (thePnt[i].Coord(theIdx) < thePnt[i + 1].Coord(theIdx)) {
620 thePnt[i + 1].SetCoord(theIdx, thePnt[i + 1].Coord(theIdx) - thePeriod);
621 }
622 }
623 else {// increasing sequence (parameters increase twice(++) and one period jump(-))
624 for (Standard_Integer i = aSavedPoint; i > 0; i--)
625 if (thePnt[i].Coord(theIdx) < thePnt[i - 1].Coord(theIdx)) {
626 thePnt[i - 1].SetCoord(theIdx, thePnt[i - 1].Coord(theIdx) - thePeriod);
627 }
628 for (Standard_Integer i = aSavedPoint; i < 3; i++)
629 if (thePnt[i].Coord(theIdx) > thePnt[i + 1].Coord(theIdx)) {
630 thePnt[i + 1].SetCoord(theIdx, thePnt[i + 1].Coord(theIdx) + thePeriod);
631 }
632 }
15b54261 633
02fd709b 634 // Do not return false, because for nonlinear 2d curves vector of parameters
635 // may change its first derivative and shifted parameters will be broken for this case.
636 return Standard_True;
15b54261 637 }
638
639//=======================================================================
640//function : getLine
641//purpose :
642//=======================================================================
643
644 Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::getLine(
2c26a53d 645 const TColgp_SequenceOfPnt& thepoints,
646 const TColStd_SequenceOfReal& theparams,
647 TColgp_SequenceOfPnt2d& thePnt2ds,
15b54261 648 Standard_Real theTol,
02fd709b 649 Standard_Boolean &isRecompute,
650 Standard_Boolean &isFromCashe) const
15b54261 651 {
652 Standard_Integer nb = thepoints.Length();
653 gp_Pnt aP[4];
654 aP[0] = thepoints(1);
655 aP[1] = thepoints(2);
656 aP[2] = thepoints(nb - 1);
657 aP[3] = thepoints(nb);
658 gp_Pnt2d aP2d[4];
659 Standard_Integer i = 0;
660
661 Standard_Real aTol2 = theTol * theTol;
662 Standard_Boolean isPeriodicU = mySurf->Surface()->IsUPeriodic();
663 Standard_Boolean isPeriodicV = mySurf->Surface()->IsVPeriodic();
664
665 // Workaround:
666 // Protection against bad "tolerance" shapes.
667 if (aTol2 > 1.0)
668 {
669 theTol = Precision::Confusion();
670 aTol2 = theTol * theTol;
671 }
02fd709b 672 if (aTol2 < Precision::SquareConfusion())
673 aTol2 = Precision::SquareConfusion();
15b54261 674 Standard_Real anOldTol2 = aTol2;
02fd709b 675 // auxiliary variables to choose period for connection with previous 2dcurve (if exist)
676 Standard_Integer aSavedPointNum = -1;
677 gp_Pnt2d aSavedPoint;
15b54261 678
679 // project first and last points
680 for( ; i < 4; i +=3)
681 {
682 Standard_Integer j;
2c26a53d 683 for (j = 0; j < myNbCashe; j++)
c48e2889 684 {
15b54261 685 if ( myCashe3d[j].SquareDistance (aP[i] ) < aTol2)
686 {
c48e2889 687 aP2d[i] = mySurf->NextValueOfUV (myCashe2d[j], aP[i], theTol, theTol);
02fd709b 688 aSavedPointNum = i;
689 aSavedPoint = myCashe2d[j];
690 if (i == 0)
691 isFromCashe = Standard_True;
15b54261 692 break;
693 }
c48e2889 694 }
2c26a53d 695
c48e2889 696 if (j >= myNbCashe)
c48e2889 697 aP2d[i] = mySurf->ValueOfUV(aP[i], theTol);
2c26a53d 698
c48e2889 699 Standard_Real aDist = mySurf->Gap();
700 Standard_Real aCurDist = aDist * aDist;
2c26a53d 701 if( aTol2 < aDist * aDist)
c48e2889 702 aTol2 = aCurDist;
15b54261 703 }
2c26a53d 704
15b54261 705 if ( isPeriodicU || isPeriodicV )
706 {
707 // Compute second and last but one c2d points.
708 for(i = 1; i < 3; i++)
709 {
710 Standard_Integer j;
2c26a53d 711 for (j = 0; j < myNbCashe; j++)
c48e2889 712 {
15b54261 713 if ( myCashe3d[j].SquareDistance (aP[i] ) < aTol2)
714 {
715 aP2d[i] = mySurf->NextValueOfUV (myCashe2d[j], aP[i], theTol, theTol);
02fd709b 716 aSavedPointNum = i;
717 aSavedPoint = myCashe2d[j];
15b54261 718 break;
719 }
c48e2889 720 }
2c26a53d 721
c48e2889 722 if (j >= myNbCashe)
c48e2889 723 aP2d[i] = mySurf->ValueOfUV(aP[i], theTol);
2c26a53d 724
c48e2889 725 Standard_Real aDist = mySurf->Gap();
726 Standard_Real aCurDist = aDist * aDist;
2c26a53d 727 if( aTol2 < aDist * aDist)
c48e2889 728 aTol2 = aCurDist;
15b54261 729 }
730
731 if (isPeriodicU)
732 {
02fd709b 733 isRecompute = fixPeriodictyTroubles(&aP2d[0], 1 /* X Coord */, mySurf->Surface()->UPeriod(), aSavedPointNum, aSavedPoint.X());
15b54261 734 }
735
736 if (isPeriodicV)
737 {
02fd709b 738 isRecompute = fixPeriodictyTroubles(&aP2d[0], 2 /* Y Coord */, mySurf->Surface()->VPeriod(), aSavedPointNum, aSavedPoint.Y());
15b54261 739 }
740 }
741
02fd709b 742 if (isRecompute && mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface))) {
743 // Do not try to make line, because in this case may be very special case when 3d curve
744 // go over the pole of, e.g., sphere, and partly lies along seam
745 // (see ApproxPCurve() for more information).
746 return 0;
747 }
748
15b54261 749 thePnt2ds.SetValue(1, aP2d[0]);
750 thePnt2ds.SetValue(nb, aP2d[3]);
751
752 // Restore old tolerance in 2d space to avoid big gap cases.
753 aTol2 = anOldTol2;
754 // Check that straight line in 2d with parameterisation as in 3d will fit
755 // fit 3d curve at all points.
756 Standard_Real dPar = theparams(nb) - theparams(1);
757 if ( Abs(dPar) < Precision::PConfusion() )
758 return 0;
759 gp_Vec2d aVec0 (aP2d[0], aP2d[3]);
760 gp_Vec2d aVec = aVec0 / dPar;
02fd709b 761 Handle(Geom_Surface) aSurf = mySurf->Surface();
762 Standard_Boolean isNormalCheck = aSurf->IsCNu(1) && aSurf->IsCNv(1);
763 if (isNormalCheck) {
764 for(i = 1; i <= nb; i++)
765 {
766 gp_XY aCurPoint = aP2d[0].XY() + aVec.XY() * (theparams(i) - theparams(1));
767 gp_Pnt aCurP;
768 gp_Vec aNormalVec, aDu, aDv;
769 aSurf->D1(aCurPoint.X(), aCurPoint.Y(), aCurP, aDu, aDv);
770 aNormalVec = aDu ^ aDv;
771 if (aNormalVec.SquareMagnitude() < Precision::SquareConfusion()) {
772 isNormalCheck = Standard_False;
773 break;
774 }
775 gp_Lin aNormalLine(aCurP, gp_Dir(aNormalVec));
776 Standard_Real aDist = aNormalLine.Distance(thepoints(i));
777 if (aDist > theTol)
778 return 0;
779 }
780 }
781 if (!isNormalCheck) {
782 Standard_Real aFirstPointDist = mySurf->Surface()->Value(aP2d[0].X(), aP2d[0].Y()).
783 SquareDistance(thepoints(1));
784 aTol2 = Max(aTol2, aTol2 * 2 * aFirstPointDist);
785 for(i = 2; i < nb; i++)
786 {
787 gp_XY aCurPoint = aP2d[0].XY() + aVec.XY() * (theparams(i) - theparams(1));
788 gp_Pnt aCurP;
789 aSurf->D0(aCurPoint.X(), aCurPoint.Y(), aCurP);
790 Standard_Real aDist1 = aCurP.SquareDistance(thepoints(i));
791
792 if(Abs (aFirstPointDist - aDist1) > aTol2)
793 return 0;
794 }
15b54261 795 }
796
797 // check if pcurve can be represented by Geom2d_Line (parameterised by length)
798 Standard_Real aLLength = aVec0.Magnitude();
799 if ( Abs (aLLength - dPar) <= Precision::PConfusion() )
800 {
801 gp_XY aDirL = aVec0.XY() / aLLength;
802 gp_Pnt2d aPL (aP2d[0].XY() - theparams(1) * aDirL);
803 return new Geom2d_Line (aPL, gp_Dir2d(aDirL));
804 }
805
806 // create straight bspline
807 TColgp_Array1OfPnt2d aPoles(1, 2);
808 aPoles(1) = aP2d[0];
809 aPoles(2) = aP2d[3];
810
811 TColStd_Array1OfReal aKnots(1,2);
812 aKnots(1) = theparams(1);
813 aKnots(2) = theparams(theparams.Length());
814
815 TColStd_Array1OfInteger aMults(1,2);
816 aMults(1) = 2;
817 aMults(2) = 2;
818 Standard_Integer aDegree = 1;
819 Handle(Geom2d_BSplineCurve) abspl2d =
820 new Geom2d_BSplineCurve (aPoles, aKnots, aMults, aDegree);
821 return abspl2d;
822 }
823
7fd59977 824//=======================================================================
825//function : ApproxPCurve
826//purpose :
827//=======================================================================
828
15b54261 829 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::ApproxPCurve(const Standard_Integer nbrPnt,
2c26a53d 830 const Handle(Geom_Curve)& c3d,
831 const Standard_Real TolFirst,
832 const Standard_Real TolLast,
833 TColgp_SequenceOfPnt& points,
834 TColStd_SequenceOfReal& params,
835 TColgp_SequenceOfPnt2d& pnt2d,
836 Handle(Geom2d_Curve)& c2d)
7fd59977 837{
15b54261 838 // for performance, first try to handle typical case when pcurve is straight
839 Standard_Boolean isRecompute = Standard_False;
02fd709b 840 Standard_Boolean isFromCasheLine = Standard_False;
2c26a53d 841 for (Standard_Integer iseq = 1; iseq <= nbrPnt; iseq++)
842 {
843 gp_Pnt2d aP2d(0.,0.);
844 pnt2d.Append(aP2d);
845 }
02fd709b 846 c2d = getLine(points, params, pnt2d, myPreci, isRecompute, isFromCasheLine);
15b54261 847 if(!c2d.IsNull())
848 {
02fd709b 849 // fill cashe
850 Standard_Boolean ChangeCycle = Standard_False;
851 if(myNbCashe>0 && myCashe3d[0].Distance(points(1)) > myCashe3d[0].Distance(points(nbrPnt)) &&
852 myCashe3d[0].Distance(points(nbrPnt))<Precision::Confusion())
853 ChangeCycle = Standard_True;
854 myNbCashe = 2;
855 if(ChangeCycle) {
856 myCashe3d[0] = points(1);
857 myCashe3d[1] = points(nbrPnt);
858 myCashe2d[0] = pnt2d(1);
859 myCashe2d[1] = pnt2d(nbrPnt);
860 }
861 else {
862 myCashe3d[1] = points(1);
863 myCashe3d[0] = points(nbrPnt);
864 myCashe2d[1] = pnt2d(1);
865 myCashe2d[0] = pnt2d(nbrPnt);
866 }
15b54261 867 return Standard_True;
868 }
869 Standard_Boolean isDone = Standard_True;
7fd59977 870 // test if the curve 3d is a boundary of the surface
871 // (only for Bezier or BSpline surface)
872
873 Standard_Boolean isoParam, isoPar2d3d, isoTypeU, p1OnIso, p2OnIso, isoclosed;
874 gp_Pnt2d valueP1, valueP2;
875 Handle(Geom_Curve) cIso;
876 Standard_Real t1, t2;
877
878 Handle(Standard_Type) sType = mySurf->Surface()->DynamicType();
879 Standard_Boolean isAnalytic = Standard_True;
880 if (sType == STANDARD_TYPE(Geom_BezierSurface) || sType == STANDARD_TYPE(Geom_BSplineSurface)) isAnalytic = Standard_False;
881 Standard_Real uf, ul, vf, vl;
882 mySurf->Surface()->Bounds(uf, ul, vf, vl);
883 isoclosed = Standard_False;
884 TColStd_Array1OfReal pout(1, nbrPnt);
885
886 isoParam = IsAnIsoparametric(nbrPnt, points, params,
887 isoTypeU, p1OnIso, valueP1, p2OnIso, valueP2,
888 isoPar2d3d, cIso, t1, t2, pout);
889
890 // projection of the points on surfaces
891
892 gp_Pnt p3d;
893 gp_Pnt2d p2d;
7fd59977 894 Standard_Real isoValue=0., isoPar1=0., isoPar2=0., tPar=0., tdeb,tfin;
895 Standard_Real Cf, Cl, parf, parl; //szv#4:S4163:12Mar99 dist not needed
896
897 // Le calcul part-il dans le bon sens, c-a-d deb et fin dans le bon ordre ?
898 // Si uclosed et iso en V, attention isoPar1 ET/OU 2 peut toucher la fermeture
899 if(isoParam){
900 if(isoTypeU){
901 isoValue = valueP1.X();
902 isoPar1 = valueP1.Y();
903 isoPar2 = valueP2.Y();
904 isoclosed = mySurf->IsVClosed(myPreci);//#78 rln 12.03.99 S4135
905 parf = vf; parl = vl;
906 }
907 else {
908 isoValue = valueP1.Y();
909 isoPar1 = valueP1.X();
910 isoPar2 = valueP2.X();
911 isoclosed = mySurf->IsUClosed(myPreci);//#78 rln 12.03.99 S4135
912 parf = uf; parl = ul;
913 }
914 if (!isoPar2d3d && !isAnalytic) {
915 Cf = cIso->FirstParameter();
916 Cl = cIso->LastParameter();
917 if (Precision::IsInfinite(Cf)) Cf = -1000;
918 if (Precision::IsInfinite(Cl)) Cl = +1000;
919 //pdn S4030 optimizing and fix isopar case on PRO41323
920 tdeb = pout(2);
921 // dist = ShapeAnalysis_Curve().Project (cIso,points(2),myPreci,pt,tdeb,Cf,Cl);
922 // Chacun des par1 ou par2 est-il sur un bord. Attention first/last : recaler
923 if (isoclosed && (isoPar1 == parf || isoPar1 == parl)) {
924 if (Abs(tdeb-parf) < Abs(tdeb-parl)) isoPar1 = parf;
925 else isoPar1 = parl;
926 if (isoTypeU) valueP1.SetY (isoPar1);
927 else valueP1.SetX (isoPar1);
928 }
929 if (isoclosed && (isoPar2 == parf || isoPar2 == parl)) {
930 //pdn S4030 optimizing and fix isopar case on PRO41323
931 tfin = pout(nbrPnt-1);
932 //dist = ShapeAnalysis_Curve().Project (cIso,points(nbrPnt-1),myPreci,pt,tfin,Cf,Cl);
933 if (Abs(tfin-parf) < Abs(tfin-parl)) isoPar2 = parf;
934 else isoPar2 = parl;
935 if (isoTypeU) valueP2.SetY (isoPar2);
936 else valueP2.SetX (isoPar2);
937 }
938
939 // Interversion Par1/Par2 (ne veut que si les 2 sont sur les bords ...)
940 // Est-ce encore necessaire apres ce qui vient d etre fait ?
941
942 // PTV 05.02.02 fix for translation face from 12_hp_mouse (PARASOLID) face 24008
943 // if curve is periodic do not change the points
944 // skl change "if" for pout(nbrPnt-1) 19.11.2003
945 if (!isoclosed) {
946 if( (Abs(tdeb-isoPar1)>Abs(tdeb-isoPar2)) &&
947 (Abs(pout(nbrPnt-1)-isoPar2)>Abs(pout(nbrPnt-1)-isoPar1)) ) {
948 gp_Pnt2d valueTmp = valueP1;
949 valueP1 = valueP2; valueP2 = valueTmp;
950 if (isoTypeU) {
951 isoValue = valueP1.X();
952 isoPar1 = valueP1.Y();
953 isoPar2 = valueP2.Y();
954 }
955 else {
956 isoValue = valueP1.Y();
957 isoPar1 = valueP1.X();
958 isoPar2 = valueP2.X();
959 }
960 // Fin calcul sens de courbe iso
961 }
962 } // end of fix check 05.02.02
963 }
964 }
965
966 // Si pas isoParam, on a quand meme du p1OnIso/p2OnIso possible ... !!!
967 // (utile pour detromper bug de projection). Mais detromper aussi circularite
968 //else {
969 //if (p1OnIso) valueP1 =
970 //BestExtremum (valueP1,points(1),points(2));
971 //if (p2OnIso) valueP2 =
972 //BestExtremum (valueP2,points(nbrPnt),points(nbrPnt-1));
973 //}
974
b7353545 975 Standard_Real Up = ul - uf;
976 Standard_Real Vp = vl - vf;
7fd59977 977 Standard_Real gap = myPreci; //:q1
978 Standard_Boolean ChangeCycle = Standard_False; //skl for OCC3430
02fd709b 979 // auxiliaruy variables to shift 2dcurve, according to previous
980 Standard_Boolean isFromCashe = Standard_False;
981 gp_Pnt2d aSavedPoint;
7fd59977 982 if( myNbCashe>0 && myCashe3d[0].Distance(points(1))>myCashe3d[0].Distance(points(nbrPnt)) )
983 //if(myCashe3d[0].Distance(points(nbrPnt))<myPreci)
984 if(myCashe3d[0].Distance(points(nbrPnt))<Precision::Confusion())
985 ChangeCycle = Standard_True;
b7353545 986 Standard_Boolean needResolveUJump = Standard_False;
987 Standard_Boolean needResolveVJump = Standard_False;
988 gp_Pnt prevP3d;
989 gp_Pnt2d prevP2d;
7fd59977 990 //for( i = 1; i <= nbrPnt; i ++) {
991 for(Standard_Integer ii=1; ii<=nbrPnt; ii++) {
2c26a53d 992 const Standard_Integer aPntIndex = ChangeCycle ? (nbrPnt - ii + 1) : ii;
993 p3d = points (aPntIndex);
7fd59977 994 if (isoParam) {
b7353545 995
7fd59977 996 if (isoPar2d3d) {
b7353545 997 if (isoPar2 > isoPar1) tPar = params(aPntIndex);
998 else tPar = t1 + t2 - params(aPntIndex);
7fd59977 999 }
b7353545 1000 else if (!isAnalytic) {
1001 // projection to iso
1002 if (aPntIndex == 1) tPar = isoPar1;
1003 else if (aPntIndex == nbrPnt) tPar = isoPar2;
1004 else {
1005 tPar = pout(aPntIndex);
1006 //:S4030 ShapeAnalysis_Curve().Project (cIso,p3d,myPreci,pt,tPar,Cf,Cl); //szv#4:S4163:12Mar99 `dist=` not needed
1007 }
1008 }
1009
7fd59977 1010 if (!isoPar2d3d && isAnalytic) {
b7353545 1011 if (aPntIndex == 1) p2d = valueP1;
1012 else if (aPntIndex == nbrPnt) p2d = valueP2;
7fd59977 1013 else {
b7353545 1014 p2d = mySurf->NextValueOfUV(p2d, p3d, myPreci, //%12 pdn 15.02.99 optimizing
1015 Precision::Confusion() + 1000 * gap); //:q1
1016 gap = mySurf->Gap();
1017 }
1018 }
1019 else {
1020 if (isoTypeU) { p2d.SetX(isoValue); p2d.SetY(tPar); }
1021 else { p2d.SetX(tPar); p2d.SetY(isoValue); }
7fd59977 1022 }
1023 }
b7353545 1024
7fd59977 1025 else {
b7353545 1026 if ((aPntIndex == 1) && p1OnIso) p2d = valueP1;
1027 else if ((aPntIndex == nbrPnt) && p2OnIso) p2d = valueP2;
1028 else {// general case (not an iso) mais attention aux singularites !
15b54261 1029 // first and last points are already computed by getLine()
2c26a53d 1030 if (aPntIndex == 1 || aPntIndex == nbrPnt)
15b54261 1031 {
1032 if (!isRecompute)
1033 {
b7353545 1034 p2d = pnt2d(aPntIndex);
15b54261 1035 gap = mySurf->Gap();
2c26a53d 1036 if (aPntIndex == 1) {
02fd709b 1037 isFromCashe = isFromCasheLine;
1038 aSavedPoint = p2d;
1039 }
15b54261 1040 continue;
1041 }
1042 else
1043 {
1044 //:q9 abv 23 Mar 99: use cashe as 1st approach
1045 Standard_Integer j; // svv #1
c48e2889 1046 for (j = 0; j < myNbCashe; ++j)
1047 {
b7353545 1048 if (myCashe3d[j].SquareDistance(p3d) < myPreci*myPreci)
15b54261 1049 {
b7353545 1050 p2d = mySurf->NextValueOfUV(myCashe2d[j], p3d, myPreci, Precision::Confusion() + gap);
2c26a53d 1051 if (aPntIndex == 1)
c48e2889 1052 {
02fd709b 1053 isFromCashe = Standard_True;
1054 aSavedPoint = myCashe2d[j];
1055 }
15b54261 1056 break;
1057 }
c48e2889 1058 }
1059 if (j >= myNbCashe)
1060 {
1061 p2d = mySurf->ValueOfUV(p3d, myPreci);
1062 }
15b54261 1063 }
1064 }
7fd59977 1065 else {
b7353545 1066 p2d = mySurf->NextValueOfUV(p2d, p3d, myPreci, //:S4030: optimizing
1067 Precision::Confusion() + 1000 * gap); //:q1
7fd59977 1068 }
b7353545 1069 gap = mySurf->Gap();
7fd59977 1070 }
1071 }
b7353545 1072 pnt2d(aPntIndex) = p2d;
1073 if (nbrPnt > 23 && ii > 2 && ii < nbrPnt)
1074 {
1075 // additional check for possible invalid jump by U or V parameter
1076 if (fabs(p2d.X() - prevP2d.X()) > 0.95*Up && prevP3d.Distance(p3d) < myPreci && !mySurf->IsUClosed(myPreci)
1077 && mySurf->NbSingularities(myPreci) > 0 && mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_BSplineSurface)))
1078 {
1079 needResolveUJump = Standard_True;
1080 }
1081 if (fabs(p2d.Y() - prevP2d.Y()) > 0.95*Vp && prevP3d.Distance(p3d) < myPreci && !mySurf->IsVClosed(myPreci)
1082 && mySurf->NbSingularities(myPreci) > 0 && mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_BSplineSurface)))
1083 {
1084 needResolveVJump = Standard_True;
1085 }
1086 }
1087 prevP3d = p3d;
1088 prevP2d = p2d;
1089 if (ii > 1) {
1090 if (ChangeCycle)
1091 p2d.SetXY(2. * p2d.XY() - pnt2d(aPntIndex + 1).XY());
7fd59977 1092 else
b7353545 1093 p2d.SetXY(2. * p2d.XY() - pnt2d(aPntIndex - 1).XY());
7fd59977 1094 }
1095 }
1096
1097 //pdn %12 11.02.99 PRO9234 entity 15402
1098 if (!isoPar2d3d) {
1099 mySurf->ProjectDegenerated(nbrPnt,points,pnt2d,myPreci,Standard_True);
1100 mySurf->ProjectDegenerated(nbrPnt,points,pnt2d,myPreci,Standard_False);
1101 }
2c26a53d 1102
1103 //Check the extremities of 3d curve for coinciding with singularities of surf
1104 //Standard_Integer NbSing = mySurf->NbSingularities(Precision::Confusion());
1105 gp_Pnt PointFirst = points.First(), PointLast = points.Last();
1106 Standard_Real aTolFirst = (TolFirst == -1)? Precision::Confusion() : TolFirst;
1107 Standard_Real aTolLast = (TolLast == -1)? Precision::Confusion() : TolLast;
1108 for (Standard_Integer i = 1; ; i++)
1109 {
1110 Standard_Real aPreci, aFirstPar, aLastPar;
1111 gp_Pnt aP3d;
1112 gp_Pnt2d aFirstP2d, aLastP2d;
1113 Standard_Boolean IsUiso;
1114 if (!mySurf->Singularity(i, aPreci, aP3d, aFirstP2d, aLastP2d, aFirstPar, aLastPar, IsUiso))
1115 break;
1116 if (aPreci <= Precision::Confusion() &&
1117 PointFirst.Distance(aP3d) <= aTolFirst)
1118 {
1119 CorrectExtremity(c3d, params, pnt2d,
1120 Standard_True, //first point
1121 aFirstP2d,
1122 IsUiso);
1123 }
1124 if (aPreci <= Precision::Confusion() &&
1125 PointLast.Distance(aP3d) <= aTolLast)
1126 {
1127 CorrectExtremity(c3d, params, pnt2d,
1128 Standard_False, //last point
1129 aFirstP2d,
1130 IsUiso);
1131 }
1132 }
7fd59977 1133
1134 // attention aux singularites ... (hors cas iso qui les traite deja)
1135 // if (!isoParam) {
1136 // p2d = pnt2d (1);
1137 // if (mySurf->ProjectDegenerated (points(1),myPreci,pnt2d (2),p2d))
1138 // pnt2d (1) = p2d;
1139 // p2d = pnt2d (nbrPnt);
1140 // if (mySurf->ProjectDegenerated (points(nbrPnt),myPreci,pnt2d (nbrPnt-1),p2d))
1141 // pnt2d (nbrPnt) = p2d;
1142 // }
1143
1144 // Si la surface est UCLosed et VClosed, on recadre les points
1145 // algo un peu complique, on retarde l implementation
b7353545 1146 //Standard_Real Up = ul - uf;
1147 //Standard_Real Vp = vl - vf;
7fd59977 1148 Standard_Real dist2d;
2c26a53d 1149 const Standard_Real TolOnUPeriod = Precision::Confusion() * Up;
1150 const Standard_Real TolOnVPeriod = Precision::Confusion() * Vp;
0797d9d3 1151#ifdef OCCT_DEBUG
7fd59977 1152 if (mySurf->IsUClosed(myPreci) && mySurf->IsVClosed(myPreci)) {//#78 rln 12.03.99 S4135
04232180 1153 std::cout << "WARNING : Recadrage incertain sur U & VClosed" << std::endl;
7fd59977 1154 }
1155#endif
1156 // Si la surface est UCLosed, on recadre les points
b7353545 1157 if (mySurf->IsUClosed(myPreci) || needResolveUJump) {//#78 rln 12.03.99 S4135
7fd59977 1158 // Premier point dans le domain [uf, ul]
1159 Standard_Real prevX, firstX = pnt2d (1).X();
02fd709b 1160 if (!isFromCashe) {
1161 // do not shift 2dcurve, if it connects to previous
1162 while (firstX < uf) { firstX += Up; pnt2d (1).SetX(firstX); }
1163 while (firstX > ul) { firstX -= Up; pnt2d (1).SetX(firstX); }
1164 }
1165 // shift first point, according to cashe
1166 if (mySurf->Surface()->IsUPeriodic() && isFromCashe) {
1167 Standard_Real aMinParam = uf, aMaxParam = ul;
1168 while (aMinParam > aSavedPoint.X()) {
1169 aMinParam -= Up;
1170 aMaxParam -= Up;
1171 }
1172 while (aMaxParam < aSavedPoint.X()) {
1173 aMinParam += Up;
1174 aMaxParam += Up;
1175 }
1176 Standard_Real aShift = ShapeAnalysis::AdjustToPeriod(firstX, aMinParam, aMaxParam);
1177 firstX += aShift;
1178 pnt2d(1).SetX(firstX);
1179 }
7fd59977 1180 prevX = firstX;
1181
1182 //:97 abv 1 Feb 98: treat case when curve is whole out of surface bounds
1183 Standard_Real minX = firstX, maxX = firstX;
2c26a53d 1184 Standard_Boolean ToAdjust = Standard_False;
7fd59977 1185
1186 // On decalle toujours le suivant
2c26a53d 1187 for (Standard_Integer aPntIter = 2; aPntIter <= pnt2d.Length(); ++aPntIter)
1188 {
1189 // dist2d = pnt2d (aPntIter - 1).Distance(pnt2d (aPntIter));
1190 Standard_Real CurX = pnt2d (aPntIter).X();
7fd59977 1191 dist2d = Abs (CurX - prevX);
2c26a53d 1192 if (dist2d > ( Up / 2) )
1193 {
1194 InsertAdditionalPointOrAdjust(ToAdjust, 1, Up, TolOnUPeriod,
1195 CurX, prevX,
1196 c3d,
1197 aPntIter,
1198 points, params, pnt2d);
7fd59977 1199 }
1200 prevX = CurX;
1201 if ( minX > CurX ) minX = CurX; //:97
1202 else if ( maxX < CurX ) maxX = CurX; //:97
1203 }
1204
1205 //:97
02fd709b 1206 if (!isFromCashe) {
1207 // do not shift 2dcurve, if it connects to previous
1208 Standard_Real midX = 0.5 * ( minX + maxX );
1209 Standard_Real shiftX=0.;
1210 if ( midX > ul ) shiftX = -Up;
1211 else if ( midX < uf ) shiftX = Up;
1212 if ( shiftX != 0. )
2c26a53d 1213 for (Standard_Integer aPntIter = 1; aPntIter <= pnt2d.Length(); ++aPntIter)
1214 pnt2d (aPntIter).SetX ( pnt2d (aPntIter).X() + shiftX );
02fd709b 1215 }
7fd59977 1216 }
1217 // Si la surface est VCLosed, on recadre les points
1218 // Same code as UClosed : optimisation souhaitable !!
1219 // CKY : d abord un code IDENTIQUE A UClosed; PUIS le special Seam ...
1220 // Si la surface est UCLosed, on recadre les points
1221 //
1222 //#69 rln 01.03.99 S4135 bm2_sd_t4-A.stp entity 30
1223 //#78 rln 12.03.99 S4135
b7353545 1224 if (mySurf->IsVClosed(myPreci) || needResolveVJump ||
1225 mySurf->Surface()->IsKind (STANDARD_TYPE (Geom_SphericalSurface))) {
7fd59977 1226 // Premier point dans le domain [vf, vl]
1227 Standard_Real prevY, firstY = pnt2d (1).Y();
02fd709b 1228 if (!isFromCashe) {
1229 // do not shift 2dcurve, if it connects to previous
1230 while (firstY < vf) { firstY += Vp; pnt2d (1).SetY(firstY); }
1231 while (firstY > vl) { firstY -= Vp; pnt2d (1).SetY(firstY); }
1232 }
1233 // shift first point, according to cashe
1234 if (mySurf->Surface()->IsVPeriodic() && isFromCashe) {
1235 Standard_Real aMinParam = vf, aMaxParam = vl;
1236 while (aMinParam > aSavedPoint.Y()) {
1237 aMinParam -= Vp;
1238 aMaxParam -= Vp;
1239 }
1240 while (aMaxParam < aSavedPoint.Y()) {
1241 aMinParam += Vp;
1242 aMaxParam += Vp;
1243 }
1244 Standard_Real aShift = ShapeAnalysis::AdjustToPeriod(firstY, aMinParam, aMaxParam);
1245 firstY += aShift;
1246 pnt2d(1).SetY(firstY);
1247 }
7fd59977 1248 prevY = firstY;
1249
1250 //:97 abv 1 Feb 98: treat case when curve is whole out of surface bounds
1251 Standard_Real minY = firstY, maxY = firstY;
2c26a53d 1252 Standard_Boolean ToAdjust = Standard_False;
7fd59977 1253
1254 // On decalle toujours le suivant
2c26a53d 1255 for (Standard_Integer aPntIter = 2; aPntIter <= pnt2d.Length(); ++aPntIter)
1256 {
7fd59977 1257 // dist2d = pnt2d (i-1).Distance(pnt2d (i));
2c26a53d 1258 Standard_Real CurY = pnt2d (aPntIter).Y();
7fd59977 1259 dist2d = Abs (CurY - prevY);
2c26a53d 1260 if (dist2d > ( Vp / 2) )
1261 {
1262 InsertAdditionalPointOrAdjust(ToAdjust, 2, Vp, TolOnVPeriod,
1263 CurY, prevY,
1264 c3d,
1265 aPntIter,
1266 points, params, pnt2d);
7fd59977 1267 }
1268 prevY = CurY;
1269 if ( minY > CurY ) minY = CurY; //:97
1270 else if ( maxY < CurY ) maxY = CurY; //:97
1271 }
1272
1273 //:97
02fd709b 1274 if (!isFromCashe) {
1275 // do not shift 2dcurve, if it connects to previous
1276 Standard_Real midY = 0.5 * ( minY + maxY );
1277 Standard_Real shiftY=0.;
1278 if ( midY > vl ) shiftY = -Vp;
1279 else if ( midY < vf ) shiftY = Vp;
1280 if ( shiftY != 0. )
2c26a53d 1281 for (Standard_Integer aPntIter = 1; aPntIter <= pnt2d.Length(); ++aPntIter)
1282 pnt2d(aPntIter).SetY ( pnt2d(aPntIter).Y() + shiftY );
02fd709b 1283 }
7fd59977 1284 }
1285
1286 //#69 rln 01.03.99 S4135 bm2_sd_t4-A.stp entity 30
1287 //#78 rln 12.03.99 S4135
1288 if (mySurf->IsVClosed(myPreci) || mySurf->Surface()->IsKind (STANDARD_TYPE (Geom_SphericalSurface))) {
2c26a53d 1289 for (Standard_Integer aPntIter = 2; aPntIter <= pnt2d.Length(); ++aPntIter) {
7fd59977 1290 //#1 rln 11/02/98 ca_exhaust.stp entity #9869 dist2d = pnt2d (i-1).Distance(pnt2d (i));
2c26a53d 1291 dist2d = Abs (pnt2d (aPntIter).Y() - pnt2d (aPntIter - 1).Y());
7fd59977 1292 if (dist2d > ( Vp / 2) ) {
1293 // ATTENTION : il faut regarder ou le decalage se fait.
1294 // si plusieurs points sont decalles, il faut plusieurs passes
1295 // pour obtenir un resultat correct.
1296 // NOT YET IMPLEMENTED
1297
1298 // one of those point is incorrectly placed
1299 // i.e on the wrong side of the "seam"
1300 // on prend le point le plus pres des bords vf ou vl
1301 Standard_Boolean prevOnFirst = Standard_False;
1302 Standard_Boolean prevOnLast = Standard_False;
1303 Standard_Boolean currOnFirst = Standard_False;
1304 Standard_Boolean currOnLast = Standard_False;
1305
1306 // .X ? plutot .Y , non ?
2c26a53d 1307 Standard_Real distPrevVF = Abs(pnt2d (aPntIter - 1).Y() - vf);
1308 Standard_Real distPrevVL = Abs(pnt2d (aPntIter - 1).Y() - vl);
1309 Standard_Real distCurrVF = Abs(pnt2d (aPntIter).Y() - vf);
1310 Standard_Real distCurrVL = Abs(pnt2d (aPntIter).Y() - vl);
7fd59977 1311
1312 Standard_Real theMin = distPrevVF;
1313 prevOnFirst = Standard_True;
1314 if (distPrevVL < theMin) {
1315 theMin = distPrevVL;
1316 prevOnFirst = Standard_False;
1317 prevOnLast = Standard_True;
1318 }
1319 if (distCurrVF < theMin) {
1320 theMin = distCurrVF;
1321 prevOnFirst = Standard_False;
1322 prevOnLast = Standard_False;
1323 currOnFirst = Standard_True;
1324 }
1325 if (distCurrVL < theMin) {
1326 theMin = distCurrVL;
1327 prevOnFirst = Standard_False;
1328 prevOnLast = Standard_False;
1329 currOnFirst = Standard_False;
1330 currOnLast = Standard_True;
1331 }
1332 // Modifs RLN/Nijni 3-DEC-1997
1333 if (prevOnFirst) {
2c26a53d 1334 // on decalle le point (aPntIter - 1) en V Last
1335 gp_Pnt2d newPrev(pnt2d (aPntIter - 1).X(), vf); // instead of vl RLN/Nijni
1336 pnt2d (aPntIter - 1) = newPrev;
7fd59977 1337 }
1338 else if (prevOnLast) {
2c26a53d 1339 // on decalle le point (aPntIter - 1) en V first
1340 gp_Pnt2d newPrev(pnt2d (aPntIter - 1).X(), vl); // instead of vf RLN/Nijni
1341 pnt2d (aPntIter - 1) = newPrev;
7fd59977 1342 }
1343 else if (currOnFirst) {
2c26a53d 1344 // on decalle le point (aPntIter) en V Last
1345 gp_Pnt2d newCurr(pnt2d (aPntIter).X(),vf); // instead of vl RLN/Nijni
1346 pnt2d (aPntIter) = newCurr;
7fd59977 1347 }
1348 else if (currOnLast) {
2c26a53d 1349 // on decalle le point (aPntIter) en V First
1350 gp_Pnt2d newCurr(pnt2d (aPntIter).X(), vl); // instead of vf RLN/Nijni
1351 pnt2d (aPntIter) = newCurr;
7fd59977 1352 }
1353 // on verifie
0797d9d3 1354#ifdef OCCT_DEBUG
2c26a53d 1355 dist2d = pnt2d (aPntIter - 1).Distance(pnt2d (aPntIter));
7fd59977 1356 if (dist2d > ( Vp / 2) ) {
04232180 1357 std::cout << "Echec dans le recadrage" << std::endl;
7fd59977 1358 }
1359#endif
1360 }
1361 }
1362 }
1363
1364 //:c0 abv 20 Feb 98: treat very special case when 3d curve
1365 // go over the pole of, e.g., sphere, and partly lies along seam.
1366 // 2d representation of such a curve should consist of 3 parts - one on
1367 // regular part of surface (interior), one part along degenerated boundary
1368 // and one along seam.
1369 // Since it cannot be adjusted later by arranging pcurves (curve is single),
1370 // to fix it it is nesessary to have a possibility of adjusting seam
1371 // part of such curve either to left or right boundary of surface.
1372 // Test is performed only if flag AdjustOverDegen is not -1.
1373 // If AdjustOverDegen is True, seam part of curve is adjusted to
1374 // the left, and if False - to the right parametric boundary
1375 // If treated case is detected, flag DONE4 is set to status
1376 // NOTE: currently, precision is Precision::PConfusion() since it
1377 // is enough on encountered example
1378 // (ug_turbine-A.stp from ProSTEP Benchmark #3, entities ##2470 & 5680)
1379 // (r1001_ac.stp from Test Rally #10, face #35027 and others)
1380 if ( myAdjustOverDegen != -1 ) {
1381 if ( mySurf->IsUClosed(myPreci) ) {//#78 rln 12.03.99 S4135
1382 mySurf->IsDegenerated ( gp_Pnt(0,0,0), myPreci ); // pour calculer les dgnr
1383 if ( mySurf->NbSingularities(myPreci) > 0 ) { //rln S4135
1384 // 1st, find gap point (degenerated pole)
1385 Standard_Real PrevX=0.;
1386 Standard_Integer OnBound=0, PrevOnBound=0;
1387 Standard_Integer ind; // svv #1
1388 Standard_Boolean start = Standard_True;
2c26a53d 1389 for ( ind=1; ind <= pnt2d.Length(); ind++ ) {
7fd59977 1390 Standard_Real CurX = pnt2d(ind).X();
1391 // abv 16 Mar 00: trj3_s1-ug.stp #697: ignore points in singularity
1392 if ( mySurf->IsDegenerated ( points(ind), Precision::Confusion() ) )
1393 continue;
1394 OnBound = ( Abs ( Abs ( CurX - 0.5 * ( ul + uf ) ) - Up/2 ) <=
1395 Precision::PConfusion() );
1396 if ( ! start && Abs ( Abs ( CurX - PrevX ) - Up/2 ) <= 0.01*Up )
1397 break;
1398 start = Standard_False;
1399 PrevX = CurX;
1400 PrevOnBound = OnBound;
1401 }
1402 // if found, adjust seam part
2c26a53d 1403 if ( ind <= pnt2d.Length() ) {
7fd59977 1404 PrevX = ( myAdjustOverDegen ? uf : ul );
1405 Standard_Real dU = Up/2 + Precision::PConfusion();
1406 if ( PrevOnBound ) {
1407 pnt2d(ind-1).SetX ( PrevX );
1408 for ( Standard_Integer j=ind-2; j >0; j-- ) {
1409 Standard_Real CurX = pnt2d(j).X();
1410 while ( CurX < PrevX - dU ) pnt2d(j).SetX ( CurX += Up );
1411 while ( CurX > PrevX + dU ) pnt2d(j).SetX ( CurX -= Up );
1412 }
1413 }
1414 else if ( OnBound ) {
1415 pnt2d(ind).SetX ( PrevX );
2c26a53d 1416 for ( Standard_Integer j=ind+1; j <= pnt2d.Length(); j++ ) {
7fd59977 1417 Standard_Real CurX = pnt2d(j).X();
1418 while ( CurX < PrevX - dU ) pnt2d(j).SetX ( CurX += Up );
1419 while ( CurX > PrevX + dU ) pnt2d(j).SetX ( CurX -= Up );
1420 }
1421 }
1422 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE4);
1423 }
1424 }
1425 }
1426 else if ( mySurf->IsVClosed(myPreci) ) {//#78 rln 12.03.99 S4135
1427 mySurf->IsDegenerated ( gp_Pnt(0,0,0), myPreci ); // pour calculer les dgnr
1428 if ( mySurf->NbSingularities(myPreci) > 0 ) { //rln S4135
1429 // 1st, find gap point (degenerated pole)
1430 Standard_Real PrevY=0.;
1431 Standard_Integer OnBound=0, PrevOnBound=0;
1432 Standard_Integer ind; // svv #1
1433 Standard_Boolean start = Standard_True;
2c26a53d 1434 for ( ind=1; ind <= pnt2d.Length(); ind++ ) {
7fd59977 1435 Standard_Real CurY = pnt2d(ind).Y();
1436 // abv 16 Mar 00: trj3_s1-ug.stp #697: ignore points in singularity
1437 if ( mySurf->IsDegenerated ( points(ind), Precision::Confusion() ) )
1438 continue;
1439 OnBound = ( Abs ( Abs ( CurY - 0.5 * ( vl + vf ) ) - Vp/2 ) <=
1440 Precision::PConfusion() );
1441 if ( ! start && Abs ( Abs ( CurY - PrevY ) - Vp/2 ) <= 0.01*Vp )
1442 break;
1443 start = Standard_False;
1444 PrevY = CurY;
1445 PrevOnBound = OnBound;
1446 }
1447 // if found, adjust seam part
2c26a53d 1448 if ( ind <= pnt2d.Length() ) {
7fd59977 1449 PrevY = ( myAdjustOverDegen ? vf : vl );
1450 Standard_Real dV = Vp/2 + Precision::PConfusion();
1451 if ( PrevOnBound ) {
1452 pnt2d(ind-1).SetY ( PrevY );
1453 for ( Standard_Integer j=ind-2; j >0; j-- ) {
1454 Standard_Real CurY = pnt2d(j).Y();
1455 while ( CurY < PrevY - dV ) pnt2d(j).SetY ( CurY += Vp );
1456 while ( CurY > PrevY + dV ) pnt2d(j).SetY ( CurY -= Vp );
1457 }
1458 }
1459 else if ( OnBound ) {
1460 pnt2d(ind).SetY ( PrevY );
2c26a53d 1461 for ( Standard_Integer j=ind+1; j <= pnt2d.Length(); j++ ) {
7fd59977 1462 Standard_Real CurY = pnt2d(j).Y();
1463 while ( CurY < PrevY - dV ) pnt2d(j).SetY ( CurY += Vp );
1464 while ( CurY > PrevY + dV ) pnt2d(j).SetY ( CurY -= Vp );
1465 }
1466 }
1467 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE4);
1468 }
1469 }
1470 }
1471 }
1472
1473 //:q9: fill cashe
1474 myNbCashe = 2;
1475 if(ChangeCycle) { // msv 10.08.04: avoid using of uninitialised field
1476 //if(myCashe3d[0].Distance(points(1))>Precision::Confusion() &&
1477 // myCashe3d[1].Distance(points(1))>Precision::Confusion()) {
1478 myCashe3d[0] = points(1);
2c26a53d 1479 myCashe3d[1] = points.Last();
7fd59977 1480 myCashe2d[0] = pnt2d(1);
2c26a53d 1481 myCashe2d[1] = pnt2d.Last();
7fd59977 1482 }
1483 else {
1484 myCashe3d[1] = points(1);
2c26a53d 1485 myCashe3d[0] = points.Last();
7fd59977 1486 myCashe2d[1] = pnt2d(1);
2c26a53d 1487 myCashe2d[0] = pnt2d.Last();
7fd59977 1488 }
7fd59977 1489 return isDone;
1490}
1491
1492//=======================================================================
1493//function : ApproximatePCurve
1494//purpose :
1495//=======================================================================
1496
1497Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::ApproximatePCurve(const Standard_Integer /*nbrPnt*/,
1498 Handle(TColgp_HArray1OfPnt2d)& points2d,
1499 Handle(TColStd_HArray1OfReal)& params,
1500 const Handle(Geom_Curve)& /*orig*/) const
1501{
1502// Standard_Real resol = Min(mySurf->Adaptor3d()->VResolution(myPreci), mySurf->Adaptor3d()->UResolution(myPreci));
1503 Standard_Real theTolerance2d = myPreci; // (100*nbrPnt);//resol;
1504 Handle(Geom2d_Curve) C2d;
1505 try {
1506 OCC_CATCH_SIGNALS
1507 CheckPoints2d (points2d, params, theTolerance2d);
1508 Standard_Integer numberPnt = points2d->Length();
1509
1510 TColgp_Array1OfPnt points3d(1,numberPnt);
1511 gp_Pnt2d pnt2d;
1512 gp_Pnt pnt;
1513 Standard_Integer i; // svv #1
1514 for( i = 1; i <= numberPnt; i++) {
1515 pnt2d = points2d->Value(i);
1516 pnt.SetCoord(pnt2d.X(),pnt2d.Y(),0);
1517 points3d(i) = pnt;
1518 }
1519
1520 GeomAPI_PointsToBSpline appr(points3d, params->Array1(), 1, 10, GeomAbs_C1, theTolerance2d);
1521 Handle(Geom_BSplineCurve) crv3d = appr.Curve();
1522 Standard_Integer NbPoles = crv3d->NbPoles();
1523 TColgp_Array1OfPnt poles3d (1, NbPoles);
1524 TColgp_Array1OfPnt2d poles2d (1, NbPoles);
1525 crv3d->Poles(poles3d);
1526 for( i = 1; i <= NbPoles; i++) {
1527 pnt2d.SetCoord(poles3d(i).X(),poles3d(i).Y());
1528 poles2d(i) = pnt2d;
1529 }
1530 TColStd_Array1OfReal weights (1,NbPoles);
1531 TColStd_Array1OfInteger multiplicities (1,crv3d->NbKnots());
1532 TColStd_Array1OfReal knots(1,crv3d->NbKnots());
1533 crv3d->Knots(knots);
1534 crv3d->Weights(weights);
1535 crv3d->Multiplicities(multiplicities);
1536 C2d = new Geom2d_BSplineCurve ( poles2d, weights, knots, multiplicities, crv3d->Degree(), crv3d->IsPeriodic());
1537 return C2d;
1538 }
9775fa61 1539 catch(Standard_Failure const& anException) {
1540#ifdef OCCT_DEBUG
1541 //:s5
7fd59977 1542 // debug ...
1543 Standard_Integer nbp = params->Length();
1544 Standard_Integer nb2 = points2d->Length();
04232180 1545 std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::ApproximatePCurve(): Exception: ";
1546 anException.Print(std::cout);
1547 std::cout<<"Pb Geom2dAPI_Approximate, tol2d="<<theTolerance2d<<" NbParams="<<nbp<<" NbPnts="<<nb2<<std::endl;
7fd59977 1548// if (nb2 > nbp) nb2 = nbp;
1549// Standard_Real rbp,rb2; rbp = nbp; rb2 = nb2;
1550// // dbl.AddString ("NbP2d/NbParams puis X Y Param -> mini");
1551// dbl.AddReals (rb2,rbp);
1552// for (Standard_Integer i = 1; i <= nb2; i ++) {
1553// gp_XYZ quoi (points2d->Value(i).X(),points2d->Value(i).Y(),params->Value(i) );
1554// dbl.AddXYZ (quoi);
1555// }
1556#endif
9775fa61 1557 (void)anException;
7fd59977 1558 C2d.Nullify();
1559 }
1560 return C2d;
1561}
1562
1563//=======================================================================
1564//function : InterpolatePCurve
1565//purpose :
1566//=======================================================================
1567
1568 Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::InterpolatePCurve(const Standard_Integer nbrPnt,
1569 Handle(TColgp_HArray1OfPnt2d)& points2d,
1570 Handle(TColStd_HArray1OfReal)& params,
1571 const Handle(Geom_Curve)& /*orig*/) const
1572{
1573 Handle(Geom2d_Curve) C2d; // NULL si echec
1574 Standard_Real theTolerance2d = myPreci / (100 * nbrPnt);
1575 try {
1576 OCC_CATCH_SIGNALS
1577 // on verifie d abord s il n y a pas de points confondus
1578 // si besoin on retouche les valeurs ...
1579 CheckPoints2d (points2d, params, theTolerance2d);
1580 Geom2dAPI_Interpolate myInterPol2d (points2d, params,
1581 Standard_False, theTolerance2d);
1582 myInterPol2d.Perform();
1583 if (myInterPol2d.IsDone()) C2d = myInterPol2d.Curve();
1584 }
9775fa61 1585 catch(Standard_Failure const& anException) {
1586#ifdef OCCT_DEBUG
1587//:s5
7fd59977 1588// // debug ...
1589 Standard_Integer nbp = params->Length();
1590 Standard_Integer nb2 = points2d->Length();
04232180 1591 std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::InterpolatePCurve(): Exception: ";
1592 anException.Print(std::cout);
1593 std::cout<<"Pb Geom2dAPI_Interpolate, tol2d="<<theTolerance2d<<" NbParams="<<nbp<<" NbPnts="<<nb2<<std::endl;
7fd59977 1594// if (nb2 > nbp) nb2 = nbp;
1595// Standard_Real rbp,rb2; rbp = nbp; rb2 = nb2;
1596// // dbl.AddString ("NbP2d/NbParams puis X Y Param -> mini");
1597// dbl.AddReals (rb2,rbp);
1598// for (Standard_Integer i = 1; i <= nb2; i ++) {
1599// gp_XYZ quoi (points2d->Value(i).X(),points2d->Value(i).Y(),params->Value(i) );
1600// dbl.AddXYZ (quoi);
1601// }
1602#endif
9775fa61 1603 (void)anException;
7fd59977 1604 C2d.Nullify();
1605 }
1606 return C2d;
1607}
1608
1609//=======================================================================
1610//function : InterpolateCurve3d
1611//purpose :
1612//=======================================================================
1613
1614Handle(Geom_Curve) ShapeConstruct_ProjectCurveOnSurface::InterpolateCurve3d(const Standard_Integer,
1615 Handle(TColgp_HArray1OfPnt)& points,
1616 Handle(TColStd_HArray1OfReal)& params,
1617 const Handle(Geom_Curve)& /*orig*/) const
1618{
1619 Handle(Geom_Curve) C3d; // NULL si echec
1620 try {
1621 OCC_CATCH_SIGNALS
1622 Standard_Real Tol = myPreci;
1623 CheckPoints(points, params, Tol);
1624 GeomAPI_Interpolate myInterPol(points, params, Standard_False, Tol);
1625 myInterPol.Perform();
1626 if (myInterPol.IsDone()) C3d = myInterPol.Curve();
1627 }
9775fa61 1628 catch(Standard_Failure const& anException) {
1629#ifdef OCCT_DEBUG
1630 //:s5
04232180 1631 std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::InterpolateCurve3d(): Exception: ";
1632 anException.Print(std::cout); std::cout << std::endl;
7fd59977 1633#endif
9775fa61 1634 (void)anException;
1635 C3d.Nullify();
7fd59977 1636 }
1637 return C3d;
1638}
1639
2c26a53d 1640//============================================================================================
1641//function : CorrectExtremity
1642//purpose : corrects first or last 2d point of future curve
1643// in the case when it coincids with a singularity of surface
1644//============================================================================================
1645
1646void ShapeConstruct_ProjectCurveOnSurface::CorrectExtremity(const Handle(Geom_Curve)& theC3d,
1647 const TColStd_SequenceOfReal& theParams,
1648 TColgp_SequenceOfPnt2d& thePnt2d,
1649 const Standard_Boolean theIsFirstPoint,
1650 const gp_Pnt2d& thePointOnIsoLine,
1651 const Standard_Boolean theIsUiso)
1652{
1653 Standard_Integer NbPnt = thePnt2d.Length();
1654 Standard_Integer IndCoord = (theIsUiso)? 2 : 1;
1655 Standard_Real SingularityCoord = thePointOnIsoLine.Coord(3-IndCoord);
1656 gp_Pnt2d EndPoint = (theIsFirstPoint)? thePnt2d(1) : thePnt2d(NbPnt);
1657 Standard_Real FinishCoord = EndPoint.Coord(3-IndCoord); //the constant coord of isoline
1658
1659 gp_Dir2d aDir = (theIsUiso)? gp::DY2d() : gp::DX2d();
1660 gp_Lin2d anIsoLine(EndPoint, aDir);
1661 IntRes2d_Domain Dom1, Dom2;
1662
1663 Standard_Boolean IsPeriodic = (theIsUiso)?
1664 mySurf->Surface()->IsVPeriodic() : mySurf->Surface()->IsUPeriodic();
1665
1666 gp_Pnt2d FirstPointOfLine, SecondPointOfLine;
1667 Standard_Real FinishParam, FirstParam, SecondParam;
1668
1669 if (theIsFirstPoint)
1670 {
1671 FirstPointOfLine = thePnt2d(3);
1672 SecondPointOfLine = thePnt2d(2);
1673 FinishParam = theParams(1);
1674 FirstParam = theParams(3);
1675 SecondParam = theParams(2);
1676 }
1677 else //last point
1678 {
1679 FirstPointOfLine = thePnt2d(NbPnt-2);
1680 SecondPointOfLine = thePnt2d(NbPnt-1);
1681 FinishParam = theParams(NbPnt);
1682 FirstParam = theParams(NbPnt-2);
1683 SecondParam = theParams(NbPnt-1);
1684 }
1685
1686 if (SingularityCoord > FinishCoord &&
1687 SecondPointOfLine.Coord(3-IndCoord) > FinishCoord)
1688 return; //the curve passes through the singularity, do nothing
1689 if (SingularityCoord < FinishCoord &&
1690 SecondPointOfLine.Coord(3-IndCoord) < FinishCoord)
1691 return; //the curve passes through the singularity, do nothing
1692 //Check correctness of <EndPoint>
1693 {
1694 const Standard_Real aPrevDist = Abs(SecondPointOfLine.Coord(IndCoord) - FirstPointOfLine.Coord(IndCoord));
1695 const Standard_Real aCurDist = Abs(EndPoint.Coord(IndCoord) - SecondPointOfLine.Coord(IndCoord));
1696 if (aCurDist <= 2 * aPrevDist)
1697 return;
1698 }
1699
1700 gp_Pnt2d FinishPoint = (theIsUiso)? gp_Pnt2d(FinishCoord, SecondPointOfLine.Y()) :
1701 gp_Pnt2d(SecondPointOfLine.X(), FinishCoord); //first approximation of <FinishPoint>
1702
1703 for (;;)
1704 {
1705 if (Abs(SecondPointOfLine.Coord(3-IndCoord) - FinishCoord) <= 2*Precision::PConfusion())
1706 break;
1707
1708 gp_Vec2d aVec(FirstPointOfLine, SecondPointOfLine);
1709 Standard_Real aSqMagnitude = aVec.SquareMagnitude();
1710 if (aSqMagnitude <= 1.e-32)
1711 break;
1712 aDir.SetCoord(aVec.X(), aVec.Y());
1713
1714 gp_Lin2d aLine(FirstPointOfLine, aDir);
1715 IntCurve_IntConicConic Intersector(anIsoLine, Dom1,
1716 aLine, Dom2,
1717 1.e-10, 1.e-10);
1718 if (Intersector.IsDone() && !Intersector.IsEmpty())
1719 {
1720 IntRes2d_IntersectionPoint IntPoint = Intersector.Point(1);
1721 FinishPoint = IntPoint.Value();
1722 }
1723 else
1724 FinishPoint = (theIsUiso)? gp_Pnt2d(FinishCoord, SecondPointOfLine.Y()) :
1725 gp_Pnt2d(SecondPointOfLine.X(), FinishCoord);
1726
1727 gp_Pnt2d PrevPoint = FirstPointOfLine;
1728 FirstPointOfLine = SecondPointOfLine;
1729 FirstParam = SecondParam;
1730 SecondParam = (FirstParam + FinishParam)/2;
1731 if (Abs(SecondParam - FirstParam) <= 2*Precision::PConfusion())
1732 break;
1733 gp_Pnt aP3d;
1734 theC3d->D0(SecondParam, aP3d);
1735 SecondPointOfLine = mySurf->NextValueOfUV(FirstPointOfLine, aP3d,
1736 myPreci, Precision::Confusion());
1737 if (IsPeriodic)
1738 AdjustSecondPointToFirstPoint(FirstPointOfLine, SecondPointOfLine, mySurf->Surface());
1739
1740 //Check <SecondPointOfLine> to be enough close to <FirstPointOfLine>
1741 //because when a projected point is too close to singularity,
1742 //the non-constant coordinate becomes random.
1743 const Standard_Real aPrevDist = Abs(FirstPointOfLine.Coord(IndCoord) - PrevPoint.Coord(IndCoord));
1744 const Standard_Real aCurDist = Abs(SecondPointOfLine.Coord(IndCoord) - FirstPointOfLine.Coord(IndCoord));
1745 if (aCurDist > 2 * aPrevDist)
1746 break;
1747 }
1748
1749 if (theIsFirstPoint)
1750 thePnt2d(1) = FinishPoint;
1751 else
1752 thePnt2d(NbPnt) = FinishPoint;
1753}
1754
1755//============================================================================================
1756//function : InsertAdditionalPointOrAdjust
1757//purpose : If the current point is too far from the previous point
1758// (more than half-period of surface), it can happen in two cases:
1759// 1. Real current step on corresponding coordinate is small, all we need is adjust;
1760// 2. Current step on corresponding coordinate is really bigger than half-period of
1761// surface in this parametric direction, so we must add additional point to exclude
1762// such big intervals between points in 2d space.
1763//============================================================================================
1764
1765void ShapeConstruct_ProjectCurveOnSurface::
1766InsertAdditionalPointOrAdjust(Standard_Boolean& ToAdjust,
1767 const Standard_Integer theIndCoord,
1768 const Standard_Real Period,
1769 const Standard_Real TolOnPeriod,
1770 Standard_Real& CurCoord,
1771 const Standard_Real prevCoord,
1772 const Handle(Geom_Curve)& c3d,
1773 Standard_Integer& theIndex,
1774 TColgp_SequenceOfPnt& points,
1775 TColStd_SequenceOfReal& params,
1776 TColgp_SequenceOfPnt2d& pnt2d)
1777{
1778 Standard_Real CorrectedCurCoord = ElCLib::InPeriod(CurCoord,
1779 prevCoord - Period/2,
1780 prevCoord + Period/2);
1781 if (!ToAdjust)
1782 {
1783 Standard_Real CurPar = params(theIndex);
1784 Standard_Real PrevPar = params(theIndex-1);
1785 Standard_Real MidPar = (PrevPar + CurPar)/2;
1786 gp_Pnt MidP3d;
1787 c3d->D0(MidPar, MidP3d);
1788 gp_Pnt2d MidP2d = mySurf->ValueOfUV(MidP3d, myPreci);
1789 Standard_Real MidCoord = MidP2d.Coord(theIndCoord);
1790 MidCoord = ElCLib::InPeriod(MidCoord, prevCoord - Period/2, prevCoord + Period/2);
1791 Standard_Real FirstCoord = prevCoord, LastCoord = CorrectedCurCoord;
1792 if (LastCoord < FirstCoord)
1793 {Standard_Real tmp = FirstCoord; FirstCoord = LastCoord; LastCoord = tmp;}
1794 if (LastCoord - FirstCoord <= TolOnPeriod)
1795 ToAdjust = Standard_True;
1796 else if (FirstCoord <= MidCoord && MidCoord <= LastCoord)
1797 ToAdjust = Standard_True;
1798 else //add mid point
1799 {
1800 //Standard_Real RefU = prevX;
1801 Standard_Boolean Success = Standard_True;
1802 Standard_Real FirstT = PrevPar; //params(i-1)
1803 Standard_Real LastT = CurPar; //params(i)
1804 MidCoord = MidP2d.Coord(theIndCoord);
1805 while (Abs(MidCoord - prevCoord) >= Period/2 - TolOnPeriod ||
1806 Abs(CurCoord - MidCoord) >= Period/2 - TolOnPeriod)
1807 {
1808 if (MidPar - FirstT <= Precision::PConfusion() ||
1809 LastT - MidPar <= Precision::PConfusion())
1810 {
1811 Success = Standard_False;
1812 break; //wrong choice
1813 }
1814 if (Abs(MidCoord - prevCoord) >= Period/2 - TolOnPeriod)
1815 LastT = (FirstT + LastT)/2;
1816 else
1817 FirstT = (FirstT + LastT)/2;
1818 MidPar = (FirstT + LastT)/2;
1819 c3d->D0(MidPar, MidP3d);
1820 MidP2d = mySurf->ValueOfUV(MidP3d, myPreci);
1821 MidCoord = MidP2d.Coord(theIndCoord);
1822 }
1823 if (Success)
1824 {
1825 points.InsertBefore(theIndex, MidP3d);
1826 params.InsertBefore(theIndex, MidPar);
1827 pnt2d.InsertBefore(theIndex, MidP2d);
1828 theIndex++;
1829 }
1830 else
1831 ToAdjust = Standard_True;
1832 } //add mid point
1833 } //if (!ToAdjust)
1834 if (ToAdjust)
1835 {
1836 CurCoord = CorrectedCurCoord;
1837 pnt2d(theIndex).SetCoord (theIndCoord, CurCoord);
1838 }
1839}
1840
7fd59977 1841//=======================================================================
1842//function : CheckPoints
1843//purpose :
1844//=======================================================================
1845
1846 void ShapeConstruct_ProjectCurveOnSurface::CheckPoints(Handle(TColgp_HArray1OfPnt)& points,Handle(TColStd_HArray1OfReal)& params,Standard_Real& preci) const
1847{
1848 Standard_Integer firstElem = points->Lower();
1849 Standard_Integer lastElem = points->Upper();
1850 Standard_Integer i;
1851 Standard_Integer nbPntDropped = 0;
1852 Standard_Integer lastValid = firstElem; // indice of last undropped point
1853
1854 // will store 0 when the point is to be removed, 1 otherwise
1855 TColStd_Array1OfInteger tmpParam(firstElem, lastElem);
2c26a53d 1856 for (i = firstElem; i<=lastElem ; i++)
1857 tmpParam.SetValue(i,1);
7ae65f0d 1858 Standard_Real DistMin2 = RealLast();
7fd59977 1859 gp_Pnt Prev = points->Value (lastValid);
1860 gp_Pnt Curr;
1861 for (i = firstElem + 1; i <= lastElem ; i ++) {
1862 Curr = points->Value(i);
7ae65f0d 1863 Standard_Real CurDist2 = Prev.SquareDistance(Curr);
1864 if (CurDist2 < gp::Resolution()) { // test 0
7fd59977 1865 nbPntDropped ++;
1866 if ( i == lastElem ) tmpParam.SetValue(lastValid, 0); // last point kept
1867 else tmpParam.SetValue(i, 0); // current dropped, lastValid unchanged
1868 } else {
7ae65f0d 1869 if (CurDist2 < DistMin2)
1870 DistMin2 = CurDist2;
7fd59977 1871 // lastValid becomes the current (i.e. i)
1872 lastValid = i;
1873 Prev = Curr;
1874 }
1875 }
7ae65f0d 1876 if (DistMin2 < RealLast())
1877 preci = 0.9 * Sqrt (DistMin2); // preci est la distance min entre les points on la reduit un peu
1878 if (nbPntDropped == 0)
7fd59977 1879 return;
7ae65f0d 1880
0797d9d3 1881#ifdef OCCT_DEBUG
04232180 1882 std::cout << "Warning : removing 3d points for interpolation" << std::endl;
7fd59977 1883#endif
1884 // Build new HArrays
1885 Standard_Integer newLast = lastElem - nbPntDropped;
1886 if ((newLast - firstElem + 1) < 2) {
0797d9d3 1887#ifdef OCCT_DEBUG
04232180 1888 std::cout << "Too many degenerated points for 3D interpolation" << std::endl;
7fd59977 1889#endif
1890 return;
1891 }
1892 Handle(TColgp_HArray1OfPnt) newPnts =
1893 new TColgp_HArray1OfPnt(firstElem, newLast);
1894 Handle(TColStd_HArray1OfReal) newParams =
1895 new TColStd_HArray1OfReal(firstElem, newLast);
1896 Standard_Integer newCurr = 1;
1897 for (i = firstElem; i<= lastElem ; i++) {
1898 if (tmpParam.Value(i) == 1) {
1899 newPnts->SetValue(newCurr, points->Value(i));
1900 newParams->SetValue(newCurr, params->Value(i));
1901 newCurr ++;
1902 }
1903 }
1904 points = newPnts;
1905 params = newParams;
7fd59977 1906 // on la reduit un peu
1907}
1908
1909//=======================================================================
1910//function : CheckPoints2d
1911//purpose :
1912//=======================================================================
1913
1914 void ShapeConstruct_ProjectCurveOnSurface::CheckPoints2d(Handle(TColgp_HArray1OfPnt2d)& points,
1915 Handle(TColStd_HArray1OfReal)& params,
1916 Standard_Real& preci) const
1917{
1918 Standard_Integer firstElem = points->Lower();
1919 Standard_Integer lastElem = points->Upper();
1920 Standard_Integer i;
1921 Standard_Integer nbPntDropped = 0;
1922 Standard_Integer lastValid = firstElem; // indice of last undropped point
1923
1924 // will store 0 when the point is to be removed, 1 otherwise
1925 TColStd_Array1OfInteger tmpParam(firstElem, lastElem);
1926 for (i = firstElem; i<=lastElem ; i++) {
1927 tmpParam.SetValue(i,1);
1928 }
7ae65f0d 1929 Standard_Real DistMin2 = RealLast();
7fd59977 1930 gp_Pnt2d Prev = points->Value(lastValid);
1931 gp_Pnt2d Curr;
1932 for (i = firstElem + 1; i<=lastElem ; i++) {
1933 Curr = points->Value(i);
7ae65f0d 1934 Standard_Real CurDist2 = Prev.SquareDistance(Curr);
1935 if (CurDist2 < gp::Resolution()) { // test 0
7fd59977 1936 nbPntDropped ++;
1937 if ( i == lastElem ) tmpParam.SetValue(lastValid, 0); // last point kept
1938 else tmpParam.SetValue(i, 0); // current dropped, lastValid unchanged
1939 } else {
7ae65f0d 1940 if (CurDist2 < DistMin2)
1941 DistMin2 = CurDist2;
7fd59977 1942 // lastValid becomes the current (i.e. i)
1943 lastValid = i;
1944 Prev = Curr;
1945 }
1946 }
7ae65f0d 1947 if (DistMin2 < RealLast())
1948 preci = 0.9 * Sqrt (DistMin2);
1949 if (nbPntDropped == 0)
7fd59977 1950 return;
7ae65f0d 1951
0797d9d3 1952#ifdef OCCT_DEBUG
04232180 1953 std::cout << "Warning : removing 2d points for interpolation" << std::endl;
7fd59977 1954#endif
1955 // Build new HArrays
1956 Standard_Integer newLast = lastElem - nbPntDropped;
1957 if ((newLast - firstElem + 1) < 2) {
0797d9d3 1958#ifdef OCCT_DEBUG
04232180 1959 std::cout << "Too many degenerated points for 2D interpolation" << std::endl;
7fd59977 1960#endif
1961 //pdn 12.02.99 S4135 Creating pcurve with minimal length.
1962 tmpParam.SetValue(firstElem,1);
1963 tmpParam.SetValue(lastElem,1);
1964 gp_XY lastPnt = points->Value(lastElem).XY();
1965 lastPnt.Add(gp_XY(preci,preci));
1966 points->SetValue(lastElem,lastPnt);
1967 newLast = firstElem+1;
1968 //return;
1969 }
1970 Handle(TColgp_HArray1OfPnt2d) newPnts =
1971 new TColgp_HArray1OfPnt2d(firstElem, newLast);
1972 Handle(TColStd_HArray1OfReal) newParams =
1973 new TColStd_HArray1OfReal(firstElem, newLast);
1974 Standard_Integer newCurr = 1;
1975 for (i = firstElem; i <= lastElem ; i++) {
1976 if (tmpParam.Value(i) == 1) {
0797d9d3 1977#ifdef OCCT_DEBUG
04232180 1978 std::cout << "Point " << i << " : " << points->Value(i).X() << " " << points->Value(i).Y() << " at param " << params->Value(i) << std::endl;
7fd59977 1979#endif
1980 newPnts->SetValue(newCurr, points->Value(i));
1981 newParams->SetValue(newCurr, params->Value(i));
1982 newCurr ++;
1983 }
1984 else {
0797d9d3 1985#ifdef OCCT_DEBUG
04232180 1986 std::cout << "Removed " << i << " : " << points->Value(i).X() << " " << points->Value(i).Y() << " at param " << params->Value(i) << std::endl;
7fd59977 1987#endif
1988 }
1989 }
1990 points = newPnts;
1991 params = newParams;
7fd59977 1992}
1993
1994//=======================================================================
1995//function : IsAnIsoparametric
1996//purpose :
1997//=======================================================================
1998//:S4030: modified for optimization
1999//:p9 abv 11 Mar 99: PRO7226 #489490: find nearest boundary instead of first one
2000
2001 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::IsAnIsoparametric(const Standard_Integer nbrPnt,
2c26a53d 2002 const TColgp_SequenceOfPnt& points,
2003 const TColStd_SequenceOfReal& params,
7fd59977 2004 Standard_Boolean& isoTypeU,
2005 Standard_Boolean& p1OnIso,
2006 gp_Pnt2d& valueP1,
2007 Standard_Boolean& p2OnIso,
2008 gp_Pnt2d& valueP2,
2009 Standard_Boolean& isoPar2d3d,
2010 Handle(Geom_Curve)& cIso,
2011 Standard_Real& t1,
2012 Standard_Real& t2,
2013 TColStd_Array1OfReal& pout) const
2014{
2015 try { // RAJOUT
2016 OCC_CATCH_SIGNALS
2017
2018 Standard_Real prec = Precision::Confusion();//myPreci;
2019
2020 Standard_Boolean isoParam = Standard_False;
2021 isoPar2d3d = Standard_False;
2022
2023 Standard_Real U1, U2, V1, V2;
2024 mySurf->Bounds(U1, U2, V1, V2);
2025
2026 if ( mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
2027 Handle(Geom_RectangularTrimmedSurface) sTrim =
2028 Handle(Geom_RectangularTrimmedSurface)::DownCast(mySurf->Surface());
2029 sTrim->Bounds(U1, U2, V1, V2);
2030 }
2031
2032 gp_Pnt pt;
2033 Standard_Integer mpt[2]; mpt[0] = mpt[1] = 0;
d20d815b 2034 Standard_Real t, tpar[2] = { 0.0, 0.0 }, isoValue=0.;
7fd59977 2035 Standard_Real mindist2;
2036 Standard_Real mind2[2];
2037 mindist2 = mind2[0] = mind2[1] = 4*prec*prec;
2038
2039 p1OnIso = Standard_False;
2040 p2OnIso = Standard_False;
2041 const Bnd_Box* aBox = 0;
2042
2043 for (Standard_Integer j=1; (j<=4) /*&& !isoParam*/; j++) {
2044 Standard_Real isoVal=0.;
2045 Standard_Boolean isoU=Standard_False; //szv#4:S4163:12Mar99 `isoU` must be Standard_Boolean
2046 Handle(Geom_Curve) cI;
2047 Standard_Real tt1, tt2;
2048
2049 if (j == 1 ) {
2050 if (Precision::IsInfinite(U1)) continue;
2051 cI = mySurf->UIso(U1);
2052 isoU = Standard_True;
2053 isoVal = U1;
2054 aBox = & mySurf->GetBoxUF();
2055 }
2056 else if (j == 2) {
2057 if (Precision::IsInfinite(U2)) continue;
2058 cI = mySurf->UIso(U2);
2059 isoU = Standard_True;
2060 isoVal = U2;
2061 aBox = & mySurf->GetBoxUL();
2062 }
2063 else if (j == 3) {
2064 if (Precision::IsInfinite(V1)) continue;
2065 cI = mySurf->VIso(V1);
2066 isoU = Standard_False;
2067 isoVal = V1;
2068 aBox = & mySurf->GetBoxVF();
2069 }
2070 else if (j == 4) {
2071 if (Precision::IsInfinite(V2)) continue;
2072 cI = mySurf->VIso(V2);
2073 isoU = Standard_False;
2074 isoVal = V2;
2075 aBox = & mySurf->GetBoxVL();
2076 }
2077 if(cI.IsNull())
2078 continue;
2079
2080 if (isoU) { tt1 = V1; tt2 = V2; }
2081 else { tt1 = U1; tt2 = U2; }
2082
2083 gp_Pnt ext1, ext2;
2084 cI->D0(tt1, ext1);
2085 cI->D0(tt2, ext2);
2086
2087// PATCH CKY 9-JUL-1998 : protection contre singularite
2088 gp_Pnt extmi;
2089 cI->D0( (tt1+tt2)/2,extmi);
2090 if (ext1.IsEqual(ext2,prec) && ext1.IsEqual(extmi,prec)) continue;
2091
2092 Standard_Boolean PtEQext1 = Standard_False;
2093 Standard_Boolean PtEQext2 = Standard_False;
2094
eb1ebea4 2095 Standard_Real currd2[2], tp[2] = {0, 0};
7fd59977 2096 Standard_Integer mp[2];
2097
2098 for (Standard_Integer i=0; i<2; i++) {
2099 mp[i] = 0;
2100 Standard_Integer k = (i == 0 ? 1 : nbrPnt);
2101
2102 // si ext1 == ext2 => valueP1 == valueP2 => vect null plus tard
2103 currd2[i] = points(k).SquareDistance ( ext1 );
2104 if ( currd2[i] <= prec*prec && !PtEQext1) {
2105 mp[i] = 1;
2106 tp[i] = tt1;
2107 PtEQext1 = Standard_True;
2108 continue;
2109 }
2110
2111 currd2[i] = points(k).SquareDistance ( ext2 );
2112 if ( currd2[i] <= prec*prec && !PtEQext2) {
2113 mp[i] = 2;
2114 tp[i] = tt2;
2115 PtEQext2 = Standard_True;
2116 continue;
2117 }
2118
2119 // On evite de projecter sur un iso degenere
2120 // on doit egalement le faire pour l apex du cone
2121 if (mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface)) && !isoU) {
2122 continue;
2123 }
2124
2125 if(aBox->IsOut(points(k))) continue;
2126
2127 Standard_Real Cf = cI->FirstParameter();
2128 Standard_Real Cl = cI->LastParameter();
2129 if (Precision::IsInfinite(Cf)) Cf = -1000;
2130 if (Precision::IsInfinite(Cl)) Cl = +1000;
2131
2132 ShapeAnalysis_Curve sac;
2133 Standard_Real dist = sac.Project (cI,points(k),prec,pt,t,Cf,Cl);
2134 currd2[i] = dist * dist;
2135 if ((dist <= prec) && (t>= Cf) && (t<=Cl)) {
2136 mp[i] = 3;
2137 tp[i] = t;
2138 }
2139 }
2140
2141 //:e7 abv 21 Apr 98: ProSTEP TR8, r0501_pe #56679:
2142 // avoid possible null-length curves
2143 if ( mp[0] >0 && mp[1] >0 &&
2144 Abs ( tp[0] - tp[1] ) < Precision::PConfusion() ) continue;
2145
2146
2147 if (mp[0] > 0 &&
2148 ( ! p1OnIso || currd2[0] < mind2[0] ) ) {
2149 p1OnIso = Standard_True;
15b54261 2150 mind2[0] = currd2[0]; // LP2.stp #105899: FLT_INVALID_OPERATION on Windows 7 VC 9 Release mode on the whole file
7fd59977 2151 if (isoU) valueP1.SetCoord(isoVal, tp[0]);
2152 else valueP1.SetCoord(tp[0], isoVal);
2153 }
2154
2155 if (mp[1] > 0 &&
2156 ( ! p2OnIso || currd2[1] < mind2[1] ) ) {
2157 p2OnIso = Standard_True;
2158 mind2[1] = currd2[1];
2159 if (isoU) valueP2.SetCoord(isoVal, tp[1]);
2160 else valueP2.SetCoord(tp[1], isoVal);
2161 }
2162
2163 if ( mp[0] <=0 || mp[1] <=0 ) continue;
2164
2165 Standard_Real md2 = currd2[0] + currd2[1];
2166 if ( mindist2 <= md2 ) continue;
2167
2168 mindist2 = md2;
2169 mpt[0] = mp[0];
2170 mpt[1] = mp[1];
2171 tpar[0] = tp[0];
2172 tpar[1] = tp[1];
2173 isoTypeU = isoU;
2174 isoValue = isoVal;
2175 cIso = cI;
2176 t1 = tt1;
2177 t2 = tt2;
2178 }
2179
2180 // probablely it concerns an isoparametrics
2181 if ( mpt[0] >0 && mpt[1] >0 ) {
2182
2183 p1OnIso = p2OnIso = Standard_True;
2184 if (isoTypeU) {
2185 valueP1.SetCoord(isoValue, tpar[0]);
2186 valueP2.SetCoord(isoValue, tpar[1]);
2187 }
2188 else {
2189 valueP1.SetCoord(tpar[0], isoValue);
2190 valueP2.SetCoord(tpar[1], isoValue);
2191 }
2192
2193 if ( mpt[0] != 3 && mpt[1] != 3 ) {
2194 isoPar2d3d = Standard_True;
2195 for (Standard_Integer i=2; i < nbrPnt && isoPar2d3d; i++){
2196 if (tpar[1] > tpar[0]) t = params(i);
2197 else t = t1+t2-params(i);
2198 cIso->D0(t, pt);
2199 if (!points(i).IsEqual(pt, prec)) isoPar2d3d = Standard_False;
2200 }
2201 }
2202
2203 if (isoPar2d3d) isoParam = Standard_True;
2204 else {
2205 Standard_Real prevParam = tpar[0];
2206 Standard_Real Cf, Cl;
2207 Standard_Boolean isoByDistance = Standard_True;
2208 Cf = cIso->FirstParameter();
2209 Cl = cIso->LastParameter();
2210 if (Precision::IsInfinite(Cf)) Cf = -1000;
2211 if (Precision::IsInfinite(Cl)) Cl = +1000;
2212
2213 ShapeAnalysis_Curve sac;
2214 for (Standard_Integer i=2; i < nbrPnt && isoByDistance; i++) {
2215 Standard_Real dist = sac.NextProject (prevParam,cIso,points(i),
2216 prec,pt,t,Cf,Cl,
2217 Standard_False); //:j8 abv 10.12.98: TR10 r0501_db.stp #9423: avoid adjusting to ends
2218 prevParam = t;
2219 pout(i)=t;
2220 if( (dist > prec) || (t < Cf) || (t > Cl) )
2221 isoByDistance = Standard_False;
2222 }
2223 if (isoByDistance) isoParam = Standard_True;
2224 }
2225 }
2226/* if (!isoParam) { CKY 29-mai-1997 : garder tout ce qu on peut ?
2227 p1OnIso = Standard_False;
2228 p2OnIso = Standard_False;
2229 } */
2230 return isoParam;
2231 } // RAJOUT
9775fa61 2232 catch(Standard_Failure const& anException) {
0797d9d3 2233#ifdef OCCT_DEBUG
9775fa61 2234// pb : on affiche ce qu on peut
7fd59977 2235 for (Standard_Integer numpnt = 1; numpnt <= nbrPnt; numpnt ++) {
04232180 2236 std::cout<<"["<<numpnt<<"]param="<<params(numpnt)<<" point=("<<
2237 points(numpnt).X()<<" "<<points(numpnt).Y()<<" "<<points(numpnt).Z()<<")"<<std::endl;
7fd59977 2238 }
04232180 2239 std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::IsAnIsoparametric(): Exception: ";
2240 anException.Print(std::cout); std::cout << std::endl;
7fd59977 2241#endif
9775fa61 2242 (void)anException;
7fd59977 2243 return Standard_False;
2244 }
7fd59977 2245}
2246
2247/* S4135 : BestExtremum is commented after IsAnIsoparametric works with Precision::Confusion()
2248//=======================================================================
2249//function : BestExtremum
2250//purpose : auxiliaire prenant le meilleur extremum si ISO car doute possible
2251//=======================================================================
2252
2253 gp_Pnt2d ShapeConstruct_ProjectCurveOnSurface::BestExtremum(const gp_Pnt2d& P2iso,const gp_Pnt& P3ext,const gp_Pnt& P3next) const
2254{
2255// P2iso a ete calcule depuis P3ext sur une iso externe de la surface
2256// En principe bon mais circularite possible ... et IsU/VClosed faillible
2257// (si baillement 1e-4 ou 1e-5, on est dedans !). DONC
2258// 1/ on privilegie l iso mais a tout hasard on verifie si Surf meilleur
2259// 2/ si iso, attention a la circularite (cas limite)
2260
2261// NB : si isoParam, on suppose que P2iso est bon (car il y en a 2). A voir...
2262
2263// D abord, calcul p2ext depuis la surface. choix surface/iso
2264 return P2iso;
2265 Standard_Real prec = Precision::Confusion();//myPreci;
2266 gp_Pnt2d P2cal = mySurf->ValueOfUV(P3ext, prec);
2267 gp_Pnt P3cal = mySurf->Value (P2cal);
2268 Standard_Real dcal = P3ext.Distance (P3cal);
2269 Standard_Real dnxt = P3ext.Distance (P3next);
2270 if (dcal > dnxt) return P2iso; // en fait protection sur BUG (PRO8468)
2271
2272// On choisit entre P2iso et P2cal, le plus proche de P2next ... !!!
2273 gp_Pnt2d P2next = mySurf->ValueOfUV(P3next, prec);
2274 if (P2next.Distance(P2cal) < P2next.Distance(P2iso)) return P2cal;
2275 return P2iso;
2276}
2277*/
Open CASCADE Technology repository