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