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