1 // Copyright (c) 1999-2014 OPEN CASCADE SAS
3 // This file is part of Open CASCADE Technology software library.
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
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.
11 // Alternatively, this file may be used under the terms of Open CASCADE
12 // commercial license or contractual agreement.
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)
29 //:s5 abv 22.04.99 Adding debug printouts in catch {} blocks
30 //#1 svv 11.01.00 Porting on DEC
32 #include <Approx_CurveOnSurface.hxx>
33 #include <Geom2d_BSplineCurve.hxx>
34 #include <Geom2d_Circle.hxx>
35 #include <Geom2d_Curve.hxx>
36 #include <Geom2d_Ellipse.hxx>
37 #include <Geom2d_Hyperbola.hxx>
38 #include <Geom2d_Line.hxx>
39 #include <Geom2d_Parabola.hxx>
40 #include <Geom2d_TrimmedCurve.hxx>
41 #include <Geom2dAdaptor.hxx>
42 #include <Geom2dAPI_Interpolate.hxx>
43 #include <Geom_BezierSurface.hxx>
44 #include <Geom_BoundedCurve.hxx>
45 #include <Geom_BSplineCurve.hxx>
46 #include <Geom_BSplineSurface.hxx>
47 #include <Geom_Curve.hxx>
48 #include <Geom_OffsetSurface.hxx>
49 #include <Geom_Plane.hxx>
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>
55 #include <GeomAdaptor_HCurve.hxx>
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>
66 #include <ShapeAnalysis.hxx>
67 #include <ShapeAnalysis_Curve.hxx>
68 #include <ShapeAnalysis_Surface.hxx>
69 #include <ShapeConstruct_ProjectCurveOnSurface.hxx>
70 #include <ShapeExtend.hxx>
71 #include <Standard_ErrorHandler.hxx>
72 #include <Standard_Failure.hxx>
73 #include <Standard_Type.hxx>
74 #include <TColgp_Array1OfPnt.hxx>
75 #include <TColStd_Array1OfInteger.hxx>
78 IMPLEMENT_STANDARD_RTTIEXT(ShapeConstruct_ProjectCurveOnSurface,MMgt_TShared)
82 //=======================================================================
83 //function : ShapeConstruct_ProjectCurveOnSurface
85 //=======================================================================
87 ShapeConstruct_ProjectCurveOnSurface::ShapeConstruct_ProjectCurveOnSurface()
89 myPreci = Precision::Confusion();
90 myBuild = Standard_False;
91 myAdjustOverDegen = 1; //:c0 //szv#4:S4163:12Mar99 was boolean
95 //=======================================================================
98 //=======================================================================
100 void ShapeConstruct_ProjectCurveOnSurface::Init(const Handle(Geom_Surface)& surf,const Standard_Real preci)
102 Init (new ShapeAnalysis_Surface (surf), preci);
105 //=======================================================================
108 //=======================================================================
110 void ShapeConstruct_ProjectCurveOnSurface::Init(const Handle(ShapeAnalysis_Surface)& surf,const Standard_Real preci)
113 SetPrecision (preci);
116 //=======================================================================
117 //function : SetSurface
119 //=======================================================================
121 void ShapeConstruct_ProjectCurveOnSurface::SetSurface(const Handle(Geom_Surface)& surf)
123 SetSurface (new ShapeAnalysis_Surface (surf));
126 //=======================================================================
127 //function : SetSurface
129 //=======================================================================
131 void ShapeConstruct_ProjectCurveOnSurface::SetSurface(const Handle(ShapeAnalysis_Surface)& surf)
133 if ( mySurf == surf ) return;
138 //=======================================================================
139 //function : SetPrecision
141 //=======================================================================
143 void ShapeConstruct_ProjectCurveOnSurface::SetPrecision(const Standard_Real preci)
148 //=======================================================================
149 //function : BuildCurveMode
151 //=======================================================================
153 Standard_Boolean& ShapeConstruct_ProjectCurveOnSurface::BuildCurveMode()
158 //=======================================================================
159 //function : AdjustOverDegenMode
161 //=======================================================================
164 //szv#4:S4163:12Mar99 was Boolean
165 Standard_Integer& ShapeConstruct_ProjectCurveOnSurface::AdjustOverDegenMode()
167 return myAdjustOverDegen;
171 //=======================================================================
172 //function : NbSurfIntervals
173 //purpose : work-around of bug in standard method
174 // GeomAdaptor_Surface->NbIntervals() (PRO16346)
175 //=======================================================================
177 static Standard_Integer NbSurfIntervals(const Handle(GeomAdaptor_HSurface)& GAS, const GeomAbs_Shape cont)
179 Standard_Integer NbU = 0;
180 if (GAS->GetType() == GeomAbs_SurfaceOfExtrusion) {
181 // extract the surface
182 Handle(Geom_SurfaceOfLinearExtrusion) surf = Handle(Geom_SurfaceOfLinearExtrusion)::DownCast(GAS->ChangeSurface().Surface());
183 // build a 3d adaptor curve
184 GeomAdaptor_Curve Adaptor3dCurve(surf->BasisCurve(), GAS->FirstUParameter(), GAS->LastUParameter());
185 if (Adaptor3dCurve.GetType() == GeomAbs_BSplineCurve)
186 NbU = Adaptor3dCurve.NbIntervals(cont);
189 NbU = GAS->NbUIntervals(cont);
190 return NbU * (GAS->NbVIntervals(cont));
193 //=======================================================================
196 //=======================================================================
198 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::Status (const ShapeExtend_Status theStatus) const
200 return ShapeExtend::DecodeStatus (myStatus, theStatus);
203 //=======================================================================
206 //=======================================================================
207 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::Perform (Handle(Geom_Curve)& c3d,
208 const Standard_Real First,
209 const Standard_Real Last,
210 Handle(Geom2d_Curve)& c2d,
212 const Standard_Integer,
213 const Standard_Integer)
215 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_OK);
216 //Standard_Boolean OK = Standard_True; //szv#4:S4163:12Mar99 not needed
218 if (mySurf.IsNull()) {
220 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_FAIL1);
221 return Standard_False;
223 // Projection Analytique
224 Handle(Geom_Curve) crv3dtrim = c3d;
225 if ( ! c3d->IsKind(STANDARD_TYPE(Geom_BoundedCurve)) )
226 crv3dtrim = new Geom_TrimmedCurve ( c3d, First, Last );
227 c2d = ProjectAnalytic ( crv3dtrim );
229 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE1);
230 return Standard_True;
233 // Projection par approximation
235 // discretize the 3d curve
237 Standard_Integer nbrPnt;
239 // $$$$ :92 abv 28 Jan 98 see PRO10107, big BSplineCurve C0
240 Standard_Integer nbPini = NCONTROL; // as in BRepCheck_Edge (RLN/Nijni)
241 // 20; // number of points for interpolation, should be "parametric dependent"
243 //:92 abv 28 Jan 98: if curve is BSpline with many intervals,
244 // increase number of points to provide at least Degree()+1 points per interval
245 Handle(Geom_BSplineCurve) bspl;
246 if ( c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve)) ) {
247 Handle(Geom_TrimmedCurve) ctrim = Handle(Geom_TrimmedCurve)::DownCast(c3d);
248 bspl = Handle(Geom_BSplineCurve)::DownCast ( ctrim->BasisCurve() );
250 else bspl = Handle(Geom_BSplineCurve)::DownCast ( c3d );
251 if ( ! bspl.IsNull() ) {
252 Standard_Integer nint = 0;
253 for ( Standard_Integer i=1; i < bspl->NbKnots(); i++ )
254 if ( bspl->Knot(i+1) > First && bspl->Knot(i) < Last ) nint++;
255 Standard_Integer minPnt = nint * ( bspl->Degree() + 1 );
256 while ( nbPini < minPnt ) nbPini += NCONTROL - 1;
258 if ( nbPini > NCONTROL )
259 cout << "Warning: number of points for projecting is " << nbPini << endl;
263 // $$$$ end :92 (big BSplineCurve C0)
265 // this number should be "parametric dependent"
266 TColgp_Array1OfPnt points(1, nbPini);
267 TColStd_Array1OfReal params(1, nbPini);
268 NCollection_Sequence<Standard_Real> aKnotCoeffs;
270 Standard_Integer iPnt;
272 // In case of bspline compute parametrization speed on each
273 // knot interval inside [aFirstParam, aLastParam].
274 // If quotient = (MaxSpeed / MinSpeed) >= aMaxQuotientCoeff then
275 // use PerformByProjLib algorithm.
278 Standard_Real aFirstParam = First; // First parameter of current interval.
279 Standard_Real aLastParam = Last; // Last parameter of current interval.
281 // First index computation.
282 Standard_Integer anIdx = 1;
283 for(; anIdx <= bspl->NbKnots() && aFirstParam < Last; anIdx++)
285 if(bspl->Knot(anIdx) > First)
291 GeomAdaptor_Curve aC3DAdaptor(c3d);
292 Standard_Real aMinParSpeed = Precision::Infinite(); // Minimal parameterization speed.
293 for(; anIdx <= bspl->NbKnots() && aFirstParam < Last; anIdx++)
295 // Fill current knot interval.
296 aLastParam = Min(Last, bspl->Knot(anIdx));
297 Standard_Integer aNbIntPnts = NCONTROL;
298 // Number of inner points is adapted according to the length of the interval
299 // to avoid a lot of calculations on small range of parameters.
302 const Standard_Real aLenThres = 1.e-2;
303 const Standard_Real aLenRatio =
304 (aLastParam - aFirstParam) / (bspl->Knot(anIdx) - bspl->Knot(anIdx - 1));
305 if (aLenRatio < aLenThres)
307 aNbIntPnts = Standard_Integer(aLenRatio / aLenThres * aNbIntPnts);
312 Standard_Real aStep = (aLastParam - aFirstParam) / (aNbIntPnts - 1);
313 Standard_Integer anIntIdx;
315 // Start filling from first point.
316 aC3DAdaptor.D0(aFirstParam, p3d1);
318 Standard_Real aLength3d = 0.0;
319 for(anIntIdx = 1; anIntIdx < aNbIntPnts; anIntIdx++)
321 Standard_Real aParam = aFirstParam + aStep * anIntIdx;
322 aC3DAdaptor.D0 (aParam, p3d2);
323 const Standard_Real aDist = p3d2.Distance(p3d1);
328 aMinParSpeed = Min(aMinParSpeed, aDist / aStep);
330 const Standard_Real aCoeff = aLength3d / (aLastParam - aFirstParam);
331 if (Abs(aCoeff) > gp::Resolution())
332 aKnotCoeffs.Append(aCoeff);
333 aFirstParam = aLastParam;
336 Standard_Real anEvenlyCoeff = 0;
337 if (aKnotCoeffs.Size() > 0)
339 anEvenlyCoeff = *std::max_element(aKnotCoeffs.begin(), aKnotCoeffs.end()) /
340 *std::min_element(aKnotCoeffs.begin(), aKnotCoeffs.end());
343 const Standard_Real aMaxQuotientCoeff = 1500.0;
344 if (anEvenlyCoeff > aMaxQuotientCoeff &&
345 aMinParSpeed > Precision::Confusion() )
347 PerformByProjLib(c3d, First, Last, c2d);
348 // PerformByProjLib fail detection:
351 return Status (ShapeExtend_DONE);
356 Standard_Real deltaT, t;
357 deltaT = (Last - First) / (nbPini-1);
359 for (iPnt = 1; iPnt <= nbPini; iPnt ++)
361 if (iPnt == 1) t = First;
362 else if (iPnt == nbPini) t = Last;
363 else t = First + (iPnt - 1) * deltaT;
370 // CALCUL par approximation
371 TColgp_Array1OfPnt2d pnt2d(1, nbrPnt);
372 ApproxPCurve (nbrPnt,points,params,pnt2d,c2d); //szv#4:S4163:12Mar99 OK not needed
374 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE2);
375 return Standard_True;
376 }// cas particulier d iso
378 // INTERPOLATION du resultat
381 Handle(TColgp_HArray1OfPnt) thePnts = new TColgp_HArray1OfPnt (1, nbPini);
382 Handle(TColStd_HArray1OfReal) theParams = new TColStd_HArray1OfReal(1, nbPini);
383 for (iPnt = 1; iPnt <= nbPini ; iPnt ++) {
384 thePnts->SetValue(iPnt, points(iPnt));
385 theParams->SetValue(iPnt, params(iPnt));
388 Handle(Geom_Curve) newc3d = InterpolateCurve3d (nbPini,thePnts,theParams, c3d);
389 if ( newc3d.IsNull() ) myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_FAIL2);
391 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE3);
396 Handle(TColgp_HArray1OfPnt2d) thePnts2d = new TColgp_HArray1OfPnt2d(1, nbPini);
397 Handle(TColStd_HArray1OfReal) theParams2d = new TColStd_HArray1OfReal(1, nbPini);
398 for (iPnt = 1; iPnt <= nbPini ; iPnt ++) {
399 theParams2d->SetValue(iPnt, params(iPnt));
400 thePnts2d->SetValue(iPnt, pnt2d(iPnt));
403 c2d = InterpolatePCurve (nbPini, thePnts2d, theParams2d, c3d);
404 // c2d = ApproximatePCurve (nbPini, thePnts2d, theParams2d, c3d);
405 // Faut-il aussi reprendre la C3D ?
406 myStatus |= ShapeExtend::EncodeStatus (c2d.IsNull() ? ShapeExtend_FAIL1 : ShapeExtend_DONE2);
407 return Status (ShapeExtend_DONE);
410 //=======================================================================
411 //function : PerformByProjLib
413 //=======================================================================
415 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::PerformByProjLib(Handle(Geom_Curve)& c3d,
416 const Standard_Real First,
417 const Standard_Real Last,
418 Handle(Geom2d_Curve)& c2d,
419 const GeomAbs_Shape /*continuity*/,
420 const Standard_Integer /*maxdeg */,
421 const Standard_Integer /*nbinterval */)
423 //Standard_Boolean OK = Standard_True; //szv#4:S4163:12Mar99 unused
425 if (mySurf.IsNull()) {
426 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL1);
427 return Standard_False;
433 Handle(GeomAdaptor_HSurface) GAS = mySurf->Adaptor3d();
434 Handle(GeomAdaptor_HCurve) GAC = new GeomAdaptor_HCurve (c3d,First,Last);
435 ProjLib_ProjectedCurve Projector(GAS, GAC);
437 switch (Projector.GetType())
440 c2d = new Geom2d_Line(Projector.Line());
442 case GeomAbs_Circle :
443 c2d = new Geom2d_Circle(Projector.Circle());
445 case GeomAbs_Ellipse :
446 c2d = new Geom2d_Ellipse(Projector.Ellipse());
448 case GeomAbs_Parabola :
449 c2d = new Geom2d_Parabola(Projector.Parabola());
451 case GeomAbs_Hyperbola :
452 c2d = new Geom2d_Hyperbola(Projector.Hyperbola());
454 case GeomAbs_BSplineCurve :
455 c2d = Projector.BSpline();
458 // Not possible, handling added to avoid gcc warning.
464 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL2);
465 return Standard_False;
469 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_DONE1);
470 return Standard_True;
474 catch(Standard_Failure const& anException) {
476 cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::PerformByProjLib(): Exception: ";
477 anException.Print(cout); cout << endl;
480 myStatus = ShapeExtend::EncodeStatus (ShapeExtend_FAIL3);
483 return Standard_False;
486 //=======================================================================
487 //function : PerformAdvanced
489 //=======================================================================
491 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::PerformAdvanced (Handle(Geom_Curve)& c3d,
492 const Standard_Real First,
493 const Standard_Real Last,
494 Handle(Geom2d_Curve)& c2d)
496 Standard_Boolean hasResult = Standard_False;
497 Standard_Integer nbintervals;
499 Standard_Boolean isStandard = (mySurf->Adaptor3d()->GetType() != GeomAbs_Cylinder);
500 // && (mySurf->Adaptor3d()->GetType() != GeomAbs_SurfaceOfRevolution);
502 if (isStandard) isStandard = !mySurf->HasSingularities(myPreci);
504 Handle(GeomAdaptor_HSurface) GAS = mySurf->Adaptor3d();
505 Handle(GeomAdaptor_HCurve) GAC = new GeomAdaptor_HCurve (c3d,First,Last);
506 nbintervals = NbSurfIntervals(GAS, GeomAbs_C1);//+GAC->NbIntervals(GeomAbs_C3);
507 isStandard = (nbintervals < 2);
510 hasResult = PerformByProjLib(c3d, First, Last, c2d);
512 if (!hasResult) hasResult = Perform (c3d, First, Last, c2d);
516 //=======================================================================
517 //function : ProjectAnalytic
519 //=======================================================================
521 Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::ProjectAnalytic(const Handle(Geom_Curve)& c3d) const
523 Handle(Geom2d_Curve) result;
525 //:k1 abv 16 Dec 98: limit analytic cases by Plane surfaces only
526 // This is necessary for K4L since it fails on other surfaces
527 // when general method GeomProjLib::Curve2d() is used
528 // Projection is done as in BRep_Tool and BRepCheck_Edge
529 Handle(Geom_Surface) surf = mySurf->Surface();
530 Handle(Geom_Plane) Plane = Handle(Geom_Plane)::DownCast ( surf );
531 if ( Plane.IsNull() ) {
532 Handle(Geom_RectangularTrimmedSurface) RTS =
533 Handle(Geom_RectangularTrimmedSurface)::DownCast ( surf );
534 if ( ! RTS.IsNull() ) Plane = Handle(Geom_Plane)::DownCast ( RTS->BasisSurface() );
536 Handle(Geom_OffsetSurface) OS =
537 Handle(Geom_OffsetSurface)::DownCast ( surf );
539 Plane = Handle(Geom_Plane)::DownCast ( OS->BasisSurface() );
542 if ( ! Plane.IsNull() ) {
543 Handle(Geom_Curve) ProjOnPlane =
544 GeomProjLib::ProjectOnPlane (c3d, Plane,
545 Plane->Position().Direction(), Standard_True);
546 Handle(GeomAdaptor_HCurve) HC = new GeomAdaptor_HCurve ( ProjOnPlane );
547 ProjLib_ProjectedCurve Proj ( mySurf->Adaptor3d(), HC );
549 result = Geom2dAdaptor::MakeCurve(Proj);
550 if ( result.IsNull() ) return result;
551 if ( result->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve)) ) {
552 Handle(Geom2d_TrimmedCurve) TC = Handle(Geom2d_TrimmedCurve)::DownCast ( result );
553 result = TC->BasisCurve();
562 //! Fix possible period jump and handle walking period parameter.
563 static Standard_Boolean fixPeriodictyTroubles(gp_Pnt2d *thePnt, // pointer to gp_Pnt2d[4] beginning
564 Standard_Integer theIdx, // Index of objective coord: 1 ~ X, 2 ~ Y
565 Standard_Real thePeriod, // Period on objective coord
566 Standard_Integer theSavedPoint, // Point number to choose period
567 Standard_Real theSavedParam) // Param from cashe to choose period
569 Standard_Real aSavedParam;
570 Standard_Integer aSavedPoint;
571 Standard_Real aMinParam = 0.0, aMaxParam = thePeriod;
572 if (theSavedPoint < 0) {
573 // normalize to first period by default
574 aSavedParam = 0.5 * thePeriod;
578 aSavedParam = theSavedParam;
579 aSavedPoint = theSavedPoint;
580 while (aMinParam > aSavedParam) {
581 aMinParam -= thePeriod;
582 aMaxParam -= thePeriod;
584 while (aMaxParam < aSavedParam) {
585 aMinParam += thePeriod;
586 aMaxParam += thePeriod;
590 Standard_Real aFixIsoParam = aMinParam;
591 Standard_Boolean isIsoLine = Standard_False;
592 if (aMaxParam - aSavedParam < Precision::PConfusion() ||
593 aSavedParam - aMinParam < Precision::PConfusion()) {
594 aFixIsoParam = aSavedParam;
595 isIsoLine = Standard_True;
597 // normalize all coordinates to [aMinParam, aMaxParam)
598 for (Standard_Integer i = 0; i < 4; i++) {
599 Standard_Real aParam = thePnt[i].Coord(theIdx);
600 Standard_Real aShift = ShapeAnalysis::AdjustToPeriod(aParam, aMinParam, aMaxParam);
602 // Walk over period coord -> not walking on another isoline in parameter space.
604 if (aMaxParam - aParam < Precision::PConfusion() || aParam - aMinParam < Precision::PConfusion())
605 aParam = aFixIsoParam;
608 if (aMaxParam - aParam < Precision::PConfusion())
610 if (aParam - aMinParam < Precision::PConfusion())
614 thePnt[i].SetCoord(theIdx, aParam);
617 // find possible period jump and increasing or decreasing coordinates vector
618 Standard_Boolean isJump = Standard_False;
619 Standard_Real aPrevDiff = 0.0;
620 Standard_Real aSumDiff = 1.0;
621 for (Standard_Integer i = 0; i < 3; i++) {
622 Standard_Real aDiff = thePnt[i + 1].Coord(theIdx) - thePnt[i].Coord(theIdx);
623 if (aDiff < -Precision::PConfusion()) {
626 //if first derivative changes its sign then period jump may exists in this place
627 if (aDiff * aPrevDiff < -Precision::PConfusion()) {
628 isJump = Standard_True;
634 return Standard_False;
636 if (aSumDiff > 0) { // decreasing sequence (parameters decrease twice(--) and one period jump(+))
637 for (Standard_Integer i = aSavedPoint; i > 0; i--)
638 if (thePnt[i].Coord(theIdx) > thePnt[i - 1].Coord(theIdx)) {
639 thePnt[i - 1].SetCoord(theIdx, thePnt[i - 1].Coord(theIdx) + thePeriod);
641 for (Standard_Integer i = aSavedPoint; i < 3; i++)
642 if (thePnt[i].Coord(theIdx) < thePnt[i + 1].Coord(theIdx)) {
643 thePnt[i + 1].SetCoord(theIdx, thePnt[i + 1].Coord(theIdx) - thePeriod);
646 else {// increasing sequence (parameters increase twice(++) and one period jump(-))
647 for (Standard_Integer i = aSavedPoint; i > 0; i--)
648 if (thePnt[i].Coord(theIdx) < thePnt[i - 1].Coord(theIdx)) {
649 thePnt[i - 1].SetCoord(theIdx, thePnt[i - 1].Coord(theIdx) - thePeriod);
651 for (Standard_Integer i = aSavedPoint; i < 3; i++)
652 if (thePnt[i].Coord(theIdx) > thePnt[i + 1].Coord(theIdx)) {
653 thePnt[i + 1].SetCoord(theIdx, thePnt[i + 1].Coord(theIdx) + thePeriod);
657 // Do not return false, because for nonlinear 2d curves vector of parameters
658 // may change its first derivative and shifted parameters will be broken for this case.
659 return Standard_True;
662 //=======================================================================
665 //=======================================================================
667 Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::getLine(
668 const TColgp_Array1OfPnt& thepoints,
669 const TColStd_Array1OfReal& theparams,
670 TColgp_Array1OfPnt2d& thePnt2ds,
671 Standard_Real theTol,
672 Standard_Boolean &isRecompute,
673 Standard_Boolean &isFromCashe) const
675 Standard_Integer nb = thepoints.Length();
677 aP[0] = thepoints(1);
678 aP[1] = thepoints(2);
679 aP[2] = thepoints(nb - 1);
680 aP[3] = thepoints(nb);
682 Standard_Integer i = 0;
684 Standard_Real aTol2 = theTol * theTol;
685 Standard_Boolean isPeriodicU = mySurf->Surface()->IsUPeriodic();
686 Standard_Boolean isPeriodicV = mySurf->Surface()->IsVPeriodic();
689 // Protection against bad "tolerance" shapes.
692 theTol = Precision::Confusion();
693 aTol2 = theTol * theTol;
695 if (aTol2 < Precision::SquareConfusion())
696 aTol2 = Precision::SquareConfusion();
697 Standard_Real anOldTol2 = aTol2;
698 // auxiliary variables to choose period for connection with previous 2dcurve (if exist)
699 Standard_Integer aSavedPointNum = -1;
700 gp_Pnt2d aSavedPoint;
702 // project first and last points
706 for (j = 0; j < myNbCashe; ++j)
708 if ( myCashe3d[j].SquareDistance (aP[i] ) < aTol2)
710 aP2d[i] = mySurf->NextValueOfUV (myCashe2d[j], aP[i], theTol, theTol);
712 aSavedPoint = myCashe2d[j];
714 isFromCashe = Standard_True;
721 aP2d[i] = mySurf->ValueOfUV(aP[i], theTol);
724 Standard_Real aDist = mySurf->Gap();
725 Standard_Real aCurDist = aDist * aDist;
726 if (aTol2 < aDist * aDist)
732 if ( isPeriodicU || isPeriodicV )
734 // Compute second and last but one c2d points.
735 for(i = 1; i < 3; i++)
738 for (j = 0; j < myNbCashe; ++j)
740 if ( myCashe3d[j].SquareDistance (aP[i] ) < aTol2)
742 aP2d[i] = mySurf->NextValueOfUV (myCashe2d[j], aP[i], theTol, theTol);
744 aSavedPoint = myCashe2d[j];
751 aP2d[i] = mySurf->ValueOfUV(aP[i], theTol);
754 Standard_Real aDist = mySurf->Gap();
755 Standard_Real aCurDist = aDist * aDist;
756 if (aTol2 < aDist * aDist)
764 isRecompute = fixPeriodictyTroubles(&aP2d[0], 1 /* X Coord */, mySurf->Surface()->UPeriod(), aSavedPointNum, aSavedPoint.X());
769 isRecompute = fixPeriodictyTroubles(&aP2d[0], 2 /* Y Coord */, mySurf->Surface()->VPeriod(), aSavedPointNum, aSavedPoint.Y());
773 if (isRecompute && mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface))) {
774 // Do not try to make line, because in this case may be very special case when 3d curve
775 // go over the pole of, e.g., sphere, and partly lies along seam
776 // (see ApproxPCurve() for more information).
780 thePnt2ds.SetValue(1, aP2d[0]);
781 thePnt2ds.SetValue(nb, aP2d[3]);
783 // Restore old tolerance in 2d space to avoid big gap cases.
785 // Check that straight line in 2d with parameterisation as in 3d will fit
786 // fit 3d curve at all points.
787 Standard_Real dPar = theparams(nb) - theparams(1);
788 if ( Abs(dPar) < Precision::PConfusion() )
790 gp_Vec2d aVec0 (aP2d[0], aP2d[3]);
791 gp_Vec2d aVec = aVec0 / dPar;
792 Handle(Geom_Surface) aSurf = mySurf->Surface();
793 Standard_Boolean isNormalCheck = aSurf->IsCNu(1) && aSurf->IsCNv(1);
795 for(i = 1; i <= nb; i++)
797 gp_XY aCurPoint = aP2d[0].XY() + aVec.XY() * (theparams(i) - theparams(1));
799 gp_Vec aNormalVec, aDu, aDv;
800 aSurf->D1(aCurPoint.X(), aCurPoint.Y(), aCurP, aDu, aDv);
801 aNormalVec = aDu ^ aDv;
802 if (aNormalVec.SquareMagnitude() < Precision::SquareConfusion()) {
803 isNormalCheck = Standard_False;
806 gp_Lin aNormalLine(aCurP, gp_Dir(aNormalVec));
807 Standard_Real aDist = aNormalLine.Distance(thepoints(i));
812 if (!isNormalCheck) {
813 Standard_Real aFirstPointDist = mySurf->Surface()->Value(aP2d[0].X(), aP2d[0].Y()).
814 SquareDistance(thepoints(1));
815 aTol2 = Max(aTol2, aTol2 * 2 * aFirstPointDist);
816 for(i = 2; i < nb; i++)
818 gp_XY aCurPoint = aP2d[0].XY() + aVec.XY() * (theparams(i) - theparams(1));
820 aSurf->D0(aCurPoint.X(), aCurPoint.Y(), aCurP);
821 Standard_Real aDist1 = aCurP.SquareDistance(thepoints(i));
823 if(Abs (aFirstPointDist - aDist1) > aTol2)
828 // check if pcurve can be represented by Geom2d_Line (parameterised by length)
829 Standard_Real aLLength = aVec0.Magnitude();
830 if ( Abs (aLLength - dPar) <= Precision::PConfusion() )
832 gp_XY aDirL = aVec0.XY() / aLLength;
833 gp_Pnt2d aPL (aP2d[0].XY() - theparams(1) * aDirL);
834 return new Geom2d_Line (aPL, gp_Dir2d(aDirL));
837 // create straight bspline
838 TColgp_Array1OfPnt2d aPoles(1, 2);
842 TColStd_Array1OfReal aKnots(1,2);
843 aKnots(1) = theparams(1);
844 aKnots(2) = theparams(theparams.Length());
846 TColStd_Array1OfInteger aMults(1,2);
849 Standard_Integer aDegree = 1;
850 Handle(Geom2d_BSplineCurve) abspl2d =
851 new Geom2d_BSplineCurve (aPoles, aKnots, aMults, aDegree);
855 //=======================================================================
856 //function : ApproxPCurve
858 //=======================================================================
860 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::ApproxPCurve(const Standard_Integer nbrPnt,
861 const TColgp_Array1OfPnt& points,
862 const TColStd_Array1OfReal& params,
863 TColgp_Array1OfPnt2d& pnt2d,
864 Handle(Geom2d_Curve)& c2d)
866 // for performance, first try to handle typical case when pcurve is straight
867 Standard_Boolean isRecompute = Standard_False;
868 Standard_Boolean isFromCasheLine = Standard_False;
869 c2d = getLine(points, params, pnt2d, myPreci, isRecompute, isFromCasheLine);
873 Standard_Boolean ChangeCycle = Standard_False;
874 if(myNbCashe>0 && myCashe3d[0].Distance(points(1)) > myCashe3d[0].Distance(points(nbrPnt)) &&
875 myCashe3d[0].Distance(points(nbrPnt))<Precision::Confusion())
876 ChangeCycle = Standard_True;
879 myCashe3d[0] = points(1);
880 myCashe3d[1] = points(nbrPnt);
881 myCashe2d[0] = pnt2d(1);
882 myCashe2d[1] = pnt2d(nbrPnt);
885 myCashe3d[1] = points(1);
886 myCashe3d[0] = points(nbrPnt);
887 myCashe2d[1] = pnt2d(1);
888 myCashe2d[0] = pnt2d(nbrPnt);
890 return Standard_True;
892 Standard_Boolean isDone = Standard_True;
893 // test if the curve 3d is a boundary of the surface
894 // (only for Bezier or BSpline surface)
896 Standard_Boolean isoParam, isoPar2d3d, isoTypeU, p1OnIso, p2OnIso, isoclosed;
897 gp_Pnt2d valueP1, valueP2;
898 Handle(Geom_Curve) cIso;
899 Standard_Real t1, t2;
901 Handle(Standard_Type) sType = mySurf->Surface()->DynamicType();
902 Standard_Boolean isAnalytic = Standard_True;
903 if (sType == STANDARD_TYPE(Geom_BezierSurface) || sType == STANDARD_TYPE(Geom_BSplineSurface)) isAnalytic = Standard_False;
904 Standard_Real uf, ul, vf, vl;
905 mySurf->Surface()->Bounds(uf, ul, vf, vl);
906 isoclosed = Standard_False;
907 TColStd_Array1OfReal pout(1, nbrPnt);
909 isoParam = IsAnIsoparametric(nbrPnt, points, params,
910 isoTypeU, p1OnIso, valueP1, p2OnIso, valueP2,
911 isoPar2d3d, cIso, t1, t2, pout);
913 // projection of the points on surfaces
918 Standard_Real isoValue=0., isoPar1=0., isoPar2=0., tPar=0., tdeb,tfin;
919 Standard_Real Cf, Cl, parf, parl; //szv#4:S4163:12Mar99 dist not needed
921 // Le calcul part-il dans le bon sens, c-a-d deb et fin dans le bon ordre ?
922 // Si uclosed et iso en V, attention isoPar1 ET/OU 2 peut toucher la fermeture
925 isoValue = valueP1.X();
926 isoPar1 = valueP1.Y();
927 isoPar2 = valueP2.Y();
928 isoclosed = mySurf->IsVClosed(myPreci);//#78 rln 12.03.99 S4135
929 parf = vf; parl = vl;
932 isoValue = valueP1.Y();
933 isoPar1 = valueP1.X();
934 isoPar2 = valueP2.X();
935 isoclosed = mySurf->IsUClosed(myPreci);//#78 rln 12.03.99 S4135
936 parf = uf; parl = ul;
938 if (!isoPar2d3d && !isAnalytic) {
939 Cf = cIso->FirstParameter();
940 Cl = cIso->LastParameter();
941 if (Precision::IsInfinite(Cf)) Cf = -1000;
942 if (Precision::IsInfinite(Cl)) Cl = +1000;
943 //pdn S4030 optimizing and fix isopar case on PRO41323
945 // dist = ShapeAnalysis_Curve().Project (cIso,points(2),myPreci,pt,tdeb,Cf,Cl);
946 // Chacun des par1 ou par2 est-il sur un bord. Attention first/last : recaler
947 if (isoclosed && (isoPar1 == parf || isoPar1 == parl)) {
948 if (Abs(tdeb-parf) < Abs(tdeb-parl)) isoPar1 = parf;
950 if (isoTypeU) valueP1.SetY (isoPar1);
951 else valueP1.SetX (isoPar1);
953 if (isoclosed && (isoPar2 == parf || isoPar2 == parl)) {
954 //pdn S4030 optimizing and fix isopar case on PRO41323
955 tfin = pout(nbrPnt-1);
956 //dist = ShapeAnalysis_Curve().Project (cIso,points(nbrPnt-1),myPreci,pt,tfin,Cf,Cl);
957 if (Abs(tfin-parf) < Abs(tfin-parl)) isoPar2 = parf;
959 if (isoTypeU) valueP2.SetY (isoPar2);
960 else valueP2.SetX (isoPar2);
963 // Interversion Par1/Par2 (ne veut que si les 2 sont sur les bords ...)
964 // Est-ce encore necessaire apres ce qui vient d etre fait ?
966 // PTV 05.02.02 fix for translation face from 12_hp_mouse (PARASOLID) face 24008
967 // if curve is periodic do not change the points
968 // skl change "if" for pout(nbrPnt-1) 19.11.2003
970 if( (Abs(tdeb-isoPar1)>Abs(tdeb-isoPar2)) &&
971 (Abs(pout(nbrPnt-1)-isoPar2)>Abs(pout(nbrPnt-1)-isoPar1)) ) {
972 gp_Pnt2d valueTmp = valueP1;
973 valueP1 = valueP2; valueP2 = valueTmp;
975 isoValue = valueP1.X();
976 isoPar1 = valueP1.Y();
977 isoPar2 = valueP2.Y();
980 isoValue = valueP1.Y();
981 isoPar1 = valueP1.X();
982 isoPar2 = valueP2.X();
984 // Fin calcul sens de courbe iso
986 } // end of fix check 05.02.02
990 // Si pas isoParam, on a quand meme du p1OnIso/p2OnIso possible ... !!!
991 // (utile pour detromper bug de projection). Mais detromper aussi circularite
993 //if (p1OnIso) valueP1 =
994 //BestExtremum (valueP1,points(1),points(2));
995 //if (p2OnIso) valueP2 =
996 //BestExtremum (valueP2,points(nbrPnt),points(nbrPnt-1));
999 Standard_Real gap = myPreci; //:q1
1000 Standard_Boolean ChangeCycle = Standard_False; //skl for OCC3430
1001 // auxiliaruy variables to shift 2dcurve, according to previous
1002 Standard_Boolean isFromCashe = Standard_False;
1003 gp_Pnt2d aSavedPoint;
1004 if( myNbCashe>0 && myCashe3d[0].Distance(points(1))>myCashe3d[0].Distance(points(nbrPnt)) )
1005 //if(myCashe3d[0].Distance(points(nbrPnt))<myPreci)
1006 if(myCashe3d[0].Distance(points(nbrPnt))<Precision::Confusion())
1007 ChangeCycle = Standard_True;
1008 //for( i = 1; i <= nbrPnt; i ++) {
1009 for(Standard_Integer ii=1; ii<=nbrPnt; ii++) {
1010 if(ChangeCycle) //skl for OCC3430
1018 if (isoPar2 > isoPar1) tPar = params(i);
1019 else tPar = t1 + t2 - params(i);
1020 } else if (!isAnalytic) {
1021 // projection to iso
1022 if (i==1) tPar = isoPar1;
1023 else if (i==nbrPnt) tPar = isoPar2;
1026 //:S4030 ShapeAnalysis_Curve().Project (cIso,p3d,myPreci,pt,tPar,Cf,Cl); //szv#4:S4163:12Mar99 `dist=` not needed
1030 if (!isoPar2d3d && isAnalytic) {
1031 if (i == 1) p2d = valueP1;
1032 else if (i == nbrPnt) p2d = valueP2;
1034 p2d = mySurf->NextValueOfUV(p2d,p3d, myPreci, //%12 pdn 15.02.99 optimizing
1035 Precision::Confusion()+1000*gap); //:q1
1036 gap = mySurf->Gap();
1039 if(isoTypeU) { p2d.SetX(isoValue); p2d.SetY(tPar); }
1040 else { p2d.SetX(tPar); p2d.SetY(isoValue); }
1045 if ( (i == 1) && p1OnIso) p2d = valueP1;
1046 else if( (i == nbrPnt) && p2OnIso) p2d = valueP2;
1047 else {// general case (not an iso) mais attention aux singularites !
1048 // first and last points are already computed by getLine()
1049 if ( (i == 1 || i == nbrPnt))
1054 gap = mySurf->Gap();
1056 isFromCashe = isFromCasheLine;
1063 //:q9 abv 23 Mar 99: use cashe as 1st approach
1064 Standard_Integer j; // svv #1
1065 for (j = 0; j < myNbCashe; ++j)
1067 if ( myCashe3d[j].SquareDistance ( p3d ) < myPreci*myPreci )
1069 p2d = mySurf->NextValueOfUV (myCashe2d[j], p3d, myPreci, Precision::Confusion()+gap);
1072 isFromCashe = Standard_True;
1073 aSavedPoint = myCashe2d[j];
1080 p2d = mySurf->ValueOfUV(p3d, myPreci);
1085 p2d = mySurf->NextValueOfUV (p2d, p3d, myPreci, //:S4030: optimizing
1086 Precision::Confusion()+1000*gap); //:q1
1088 gap = mySurf->Gap();
1094 p2d.SetXY ( 2. * p2d.XY() - pnt2d(i+1).XY() );
1096 p2d.SetXY ( 2. * p2d.XY() - pnt2d(i-1).XY() );
1100 //pdn %12 11.02.99 PRO9234 entity 15402
1102 mySurf->ProjectDegenerated(nbrPnt,points,pnt2d,myPreci,Standard_True);
1103 mySurf->ProjectDegenerated(nbrPnt,points,pnt2d,myPreci,Standard_False);
1106 // attention aux singularites ... (hors cas iso qui les traite deja)
1109 // if (mySurf->ProjectDegenerated (points(1),myPreci,pnt2d (2),p2d))
1111 // p2d = pnt2d (nbrPnt);
1112 // if (mySurf->ProjectDegenerated (points(nbrPnt),myPreci,pnt2d (nbrPnt-1),p2d))
1113 // pnt2d (nbrPnt) = p2d;
1116 // Si la surface est UCLosed et VClosed, on recadre les points
1117 // algo un peu complique, on retarde l implementation
1118 Standard_Real Up = ul - uf;
1119 Standard_Real Vp = vl - vf;
1120 Standard_Real dist2d;
1122 if (mySurf->IsUClosed(myPreci) && mySurf->IsVClosed(myPreci)) {//#78 rln 12.03.99 S4135
1123 cout << "WARNING : Recadrage incertain sur U & VClosed" << endl;
1126 // Si la surface est UCLosed, on recadre les points
1127 if (mySurf->IsUClosed(myPreci)) {//#78 rln 12.03.99 S4135
1128 // Premier point dans le domain [uf, ul]
1129 Standard_Real prevX, firstX = pnt2d (1).X();
1131 // do not shift 2dcurve, if it connects to previous
1132 while (firstX < uf) { firstX += Up; pnt2d (1).SetX(firstX); }
1133 while (firstX > ul) { firstX -= Up; pnt2d (1).SetX(firstX); }
1135 // shift first point, according to cashe
1136 if (mySurf->Surface()->IsUPeriodic() && isFromCashe) {
1137 Standard_Real aMinParam = uf, aMaxParam = ul;
1138 while (aMinParam > aSavedPoint.X()) {
1142 while (aMaxParam < aSavedPoint.X()) {
1146 Standard_Real aShift = ShapeAnalysis::AdjustToPeriod(firstX, aMinParam, aMaxParam);
1148 pnt2d(1).SetX(firstX);
1152 //:97 abv 1 Feb 98: treat case when curve is whole out of surface bounds
1153 Standard_Real minX = firstX, maxX = firstX;
1155 // On decalle toujours le suivant
1156 for (i = 2; i <= nbrPnt; i++) {
1157 // dist2d = pnt2d (i-1).Distance(pnt2d (i));
1158 Standard_Real CurX = pnt2d (i).X();
1159 dist2d = Abs (CurX - prevX);
1160 if (dist2d > ( Up / 2) ) {
1161 if (CurX > prevX + Up/2) {
1162 while (CurX > prevX + Up/2) { CurX -= Up; pnt2d (i).SetX (CurX); }
1163 } else if (CurX < prevX - Up/2) {
1164 while (CurX < prevX - Up/2) { CurX += Up; pnt2d (i).SetX (CurX); }
1169 if ( minX > CurX ) minX = CurX; //:97
1170 else if ( maxX < CurX ) maxX = CurX; //:97
1175 // do not shift 2dcurve, if it connects to previous
1176 Standard_Real midX = 0.5 * ( minX + maxX );
1177 Standard_Real shiftX=0.;
1178 if ( midX > ul ) shiftX = -Up;
1179 else if ( midX < uf ) shiftX = Up;
1181 for ( i=1; i <= nbrPnt; i++ ) pnt2d(i).SetX ( pnt2d(i).X() + shiftX );
1184 // Si la surface est VCLosed, on recadre les points
1185 // Same code as UClosed : optimisation souhaitable !!
1186 // CKY : d abord un code IDENTIQUE A UClosed; PUIS le special Seam ...
1187 // Si la surface est UCLosed, on recadre les points
1189 //#69 rln 01.03.99 S4135 bm2_sd_t4-A.stp entity 30
1190 //#78 rln 12.03.99 S4135
1191 if (mySurf->IsVClosed(myPreci) || mySurf->Surface()->IsKind (STANDARD_TYPE (Geom_SphericalSurface))) {
1192 // Premier point dans le domain [vf, vl]
1193 Standard_Real prevY, firstY = pnt2d (1).Y();
1195 // do not shift 2dcurve, if it connects to previous
1196 while (firstY < vf) { firstY += Vp; pnt2d (1).SetY(firstY); }
1197 while (firstY > vl) { firstY -= Vp; pnt2d (1).SetY(firstY); }
1199 // shift first point, according to cashe
1200 if (mySurf->Surface()->IsVPeriodic() && isFromCashe) {
1201 Standard_Real aMinParam = vf, aMaxParam = vl;
1202 while (aMinParam > aSavedPoint.Y()) {
1206 while (aMaxParam < aSavedPoint.Y()) {
1210 Standard_Real aShift = ShapeAnalysis::AdjustToPeriod(firstY, aMinParam, aMaxParam);
1212 pnt2d(1).SetY(firstY);
1216 //:97 abv 1 Feb 98: treat case when curve is whole out of surface bounds
1217 Standard_Real minY = firstY, maxY = firstY;
1219 // On decalle toujours le suivant
1220 for (i = 2; i <= nbrPnt; i ++) {
1221 // dist2d = pnt2d (i-1).Distance(pnt2d (i));
1222 Standard_Real CurY = pnt2d (i).Y();
1223 dist2d = Abs (CurY - prevY);
1224 if (dist2d > ( Vp / 2) ) {
1225 if (CurY > prevY + Vp/2) {
1226 while (CurY > prevY + Vp/2) { CurY -= Vp; pnt2d (i).SetY (CurY); }
1227 } else if (CurY < prevY - Vp/2) {
1228 while (CurY < prevY - Vp/2) { CurY += Vp; pnt2d (i).SetY (CurY); }
1232 if ( minY > CurY ) minY = CurY; //:97
1233 else if ( maxY < CurY ) maxY = CurY; //:97
1238 // do not shift 2dcurve, if it connects to previous
1239 Standard_Real midY = 0.5 * ( minY + maxY );
1240 Standard_Real shiftY=0.;
1241 if ( midY > vl ) shiftY = -Vp;
1242 else if ( midY < vf ) shiftY = Vp;
1244 for ( i=1; i <= nbrPnt; i++ ) pnt2d(i).SetY ( pnt2d(i).Y() + shiftY );
1248 //#69 rln 01.03.99 S4135 bm2_sd_t4-A.stp entity 30
1249 //#78 rln 12.03.99 S4135
1250 if (mySurf->IsVClosed(myPreci) || mySurf->Surface()->IsKind (STANDARD_TYPE (Geom_SphericalSurface))) {
1251 for (i = 2; i <= nbrPnt; i++) {
1252 //#1 rln 11/02/98 ca_exhaust.stp entity #9869 dist2d = pnt2d (i-1).Distance(pnt2d (i));
1253 dist2d = Abs (pnt2d(i).Y() - pnt2d(i - 1).Y());
1254 if (dist2d > ( Vp / 2) ) {
1255 // ATTENTION : il faut regarder ou le decalage se fait.
1256 // si plusieurs points sont decalles, il faut plusieurs passes
1257 // pour obtenir un resultat correct.
1258 // NOT YET IMPLEMENTED
1260 // one of those point is incorrectly placed
1261 // i.e on the wrong side of the "seam"
1262 // on prend le point le plus pres des bords vf ou vl
1263 Standard_Boolean prevOnFirst = Standard_False;
1264 Standard_Boolean prevOnLast = Standard_False;
1265 Standard_Boolean currOnFirst = Standard_False;
1266 Standard_Boolean currOnLast = Standard_False;
1268 // .X ? plutot .Y , non ?
1269 Standard_Real distPrevVF = Abs(pnt2d (i-1).Y() - vf);
1270 Standard_Real distPrevVL = Abs(pnt2d (i-1).Y() - vl);
1271 Standard_Real distCurrVF = Abs(pnt2d (i).Y() - vf);
1272 Standard_Real distCurrVL = Abs(pnt2d (i).Y() - vl);
1274 Standard_Real theMin = distPrevVF;
1275 prevOnFirst = Standard_True;
1276 if (distPrevVL < theMin) {
1277 theMin = distPrevVL;
1278 prevOnFirst = Standard_False;
1279 prevOnLast = Standard_True;
1281 if (distCurrVF < theMin) {
1282 theMin = distCurrVF;
1283 prevOnFirst = Standard_False;
1284 prevOnLast = Standard_False;
1285 currOnFirst = Standard_True;
1287 if (distCurrVL < theMin) {
1288 theMin = distCurrVL;
1289 prevOnFirst = Standard_False;
1290 prevOnLast = Standard_False;
1291 currOnFirst = Standard_False;
1292 currOnLast = Standard_True;
1294 // Modifs RLN/Nijni 3-DEC-1997
1296 // on decalle le point (i-1) en V Last
1297 gp_Pnt2d newPrev(pnt2d (i-1).X(), vf); // instead of vl RLN/Nijni
1298 pnt2d (i-1) = newPrev;
1300 else if (prevOnLast) {
1301 // on decalle le point (i-1) en V first
1302 gp_Pnt2d newPrev(pnt2d (i-1).X(), vl); // instead of vf RLN/Nijni
1303 pnt2d (i-1) = newPrev;
1305 else if (currOnFirst) {
1306 // on decalle le point (i) en V Last
1307 gp_Pnt2d newCurr(pnt2d (i).X(),vf); // instead of vl RLN/Nijni
1308 pnt2d (i) = newCurr;
1310 else if (currOnLast) {
1311 // on decalle le point (i) en V First
1312 gp_Pnt2d newCurr(pnt2d (i).X(), vl); // instead of vf RLN/Nijni
1313 pnt2d (i) = newCurr;
1317 dist2d = pnt2d (i-1).Distance(pnt2d (i));
1318 if (dist2d > ( Vp / 2) ) {
1319 cout << "Echec dans le recadrage" << endl;
1326 //:c0 abv 20 Feb 98: treat very special case when 3d curve
1327 // go over the pole of, e.g., sphere, and partly lies along seam.
1328 // 2d representation of such a curve should consist of 3 parts - one on
1329 // regular part of surface (interior), one part along degenerated boundary
1330 // and one along seam.
1331 // Since it cannot be adjusted later by arranging pcurves (curve is single),
1332 // to fix it it is nesessary to have a possibility of adjusting seam
1333 // part of such curve either to left or right boundary of surface.
1334 // Test is performed only if flag AdjustOverDegen is not -1.
1335 // If AdjustOverDegen is True, seam part of curve is adjusted to
1336 // the left, and if False - to the right parametric boundary
1337 // If treated case is detected, flag DONE4 is set to status
1338 // NOTE: currently, precision is Precision::PConfusion() since it
1339 // is enough on encountered example
1340 // (ug_turbine-A.stp from ProSTEP Benchmark #3, entities ##2470 & 5680)
1341 // (r1001_ac.stp from Test Rally #10, face #35027 and others)
1342 if ( myAdjustOverDegen != -1 ) {
1343 if ( mySurf->IsUClosed(myPreci) ) {//#78 rln 12.03.99 S4135
1344 mySurf->IsDegenerated ( gp_Pnt(0,0,0), myPreci ); // pour calculer les dgnr
1345 if ( mySurf->NbSingularities(myPreci) > 0 ) { //rln S4135
1346 // 1st, find gap point (degenerated pole)
1347 Standard_Real PrevX=0.;
1348 Standard_Integer OnBound=0, PrevOnBound=0;
1349 Standard_Integer ind; // svv #1
1350 Standard_Boolean start = Standard_True;
1351 for ( ind=1; ind <= nbrPnt; ind++ ) {
1352 Standard_Real CurX = pnt2d(ind).X();
1353 // abv 16 Mar 00: trj3_s1-ug.stp #697: ignore points in singularity
1354 if ( mySurf->IsDegenerated ( points(ind), Precision::Confusion() ) )
1356 OnBound = ( Abs ( Abs ( CurX - 0.5 * ( ul + uf ) ) - Up/2 ) <=
1357 Precision::PConfusion() );
1358 if ( ! start && Abs ( Abs ( CurX - PrevX ) - Up/2 ) <= 0.01*Up )
1360 start = Standard_False;
1362 PrevOnBound = OnBound;
1364 // if found, adjust seam part
1365 if ( ind <= nbrPnt ) {
1366 PrevX = ( myAdjustOverDegen ? uf : ul );
1367 Standard_Real dU = Up/2 + Precision::PConfusion();
1368 if ( PrevOnBound ) {
1369 pnt2d(ind-1).SetX ( PrevX );
1370 for ( Standard_Integer j=ind-2; j >0; j-- ) {
1371 Standard_Real CurX = pnt2d(j).X();
1372 while ( CurX < PrevX - dU ) pnt2d(j).SetX ( CurX += Up );
1373 while ( CurX > PrevX + dU ) pnt2d(j).SetX ( CurX -= Up );
1376 else if ( OnBound ) {
1377 pnt2d(ind).SetX ( PrevX );
1378 for ( Standard_Integer j=ind+1; j <= nbrPnt; j++ ) {
1379 Standard_Real CurX = pnt2d(j).X();
1380 while ( CurX < PrevX - dU ) pnt2d(j).SetX ( CurX += Up );
1381 while ( CurX > PrevX + dU ) pnt2d(j).SetX ( CurX -= Up );
1384 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE4);
1388 else if ( mySurf->IsVClosed(myPreci) ) {//#78 rln 12.03.99 S4135
1389 mySurf->IsDegenerated ( gp_Pnt(0,0,0), myPreci ); // pour calculer les dgnr
1390 if ( mySurf->NbSingularities(myPreci) > 0 ) { //rln S4135
1391 // 1st, find gap point (degenerated pole)
1392 Standard_Real PrevY=0.;
1393 Standard_Integer OnBound=0, PrevOnBound=0;
1394 Standard_Integer ind; // svv #1
1395 Standard_Boolean start = Standard_True;
1396 for ( ind=1; ind <= nbrPnt; ind++ ) {
1397 Standard_Real CurY = pnt2d(ind).Y();
1398 // abv 16 Mar 00: trj3_s1-ug.stp #697: ignore points in singularity
1399 if ( mySurf->IsDegenerated ( points(ind), Precision::Confusion() ) )
1401 OnBound = ( Abs ( Abs ( CurY - 0.5 * ( vl + vf ) ) - Vp/2 ) <=
1402 Precision::PConfusion() );
1403 if ( ! start && Abs ( Abs ( CurY - PrevY ) - Vp/2 ) <= 0.01*Vp )
1405 start = Standard_False;
1407 PrevOnBound = OnBound;
1409 // if found, adjust seam part
1410 if ( ind <= nbrPnt ) {
1411 PrevY = ( myAdjustOverDegen ? vf : vl );
1412 Standard_Real dV = Vp/2 + Precision::PConfusion();
1413 if ( PrevOnBound ) {
1414 pnt2d(ind-1).SetY ( PrevY );
1415 for ( Standard_Integer j=ind-2; j >0; j-- ) {
1416 Standard_Real CurY = pnt2d(j).Y();
1417 while ( CurY < PrevY - dV ) pnt2d(j).SetY ( CurY += Vp );
1418 while ( CurY > PrevY + dV ) pnt2d(j).SetY ( CurY -= Vp );
1421 else if ( OnBound ) {
1422 pnt2d(ind).SetY ( PrevY );
1423 for ( Standard_Integer j=ind+1; j <= nbrPnt; j++ ) {
1424 Standard_Real CurY = pnt2d(j).Y();
1425 while ( CurY < PrevY - dV ) pnt2d(j).SetY ( CurY += Vp );
1426 while ( CurY > PrevY + dV ) pnt2d(j).SetY ( CurY -= Vp );
1429 myStatus |= ShapeExtend::EncodeStatus (ShapeExtend_DONE4);
1437 if(ChangeCycle) { // msv 10.08.04: avoid using of uninitialised field
1438 //if(myCashe3d[0].Distance(points(1))>Precision::Confusion() &&
1439 // myCashe3d[1].Distance(points(1))>Precision::Confusion()) {
1440 myCashe3d[0] = points(1);
1441 myCashe3d[1] = points(nbrPnt);
1442 myCashe2d[0] = pnt2d(1);
1443 myCashe2d[1] = pnt2d(nbrPnt);
1446 myCashe3d[1] = points(1);
1447 myCashe3d[0] = points(nbrPnt);
1448 myCashe2d[1] = pnt2d(1);
1449 myCashe2d[0] = pnt2d(nbrPnt);
1454 //=======================================================================
1455 //function : ApproximatePCurve
1457 //=======================================================================
1459 Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::ApproximatePCurve(const Standard_Integer /*nbrPnt*/,
1460 Handle(TColgp_HArray1OfPnt2d)& points2d,
1461 Handle(TColStd_HArray1OfReal)& params,
1462 const Handle(Geom_Curve)& /*orig*/) const
1464 // Standard_Real resol = Min(mySurf->Adaptor3d()->VResolution(myPreci), mySurf->Adaptor3d()->UResolution(myPreci));
1465 Standard_Real theTolerance2d = myPreci; // (100*nbrPnt);//resol;
1466 Handle(Geom2d_Curve) C2d;
1469 CheckPoints2d (points2d, params, theTolerance2d);
1470 Standard_Integer numberPnt = points2d->Length();
1472 TColgp_Array1OfPnt points3d(1,numberPnt);
1475 Standard_Integer i; // svv #1
1476 for( i = 1; i <= numberPnt; i++) {
1477 pnt2d = points2d->Value(i);
1478 pnt.SetCoord(pnt2d.X(),pnt2d.Y(),0);
1482 GeomAPI_PointsToBSpline appr(points3d, params->Array1(), 1, 10, GeomAbs_C1, theTolerance2d);
1483 Handle(Geom_BSplineCurve) crv3d = appr.Curve();
1484 Standard_Integer NbPoles = crv3d->NbPoles();
1485 TColgp_Array1OfPnt poles3d (1, NbPoles);
1486 TColgp_Array1OfPnt2d poles2d (1, NbPoles);
1487 crv3d->Poles(poles3d);
1488 for( i = 1; i <= NbPoles; i++) {
1489 pnt2d.SetCoord(poles3d(i).X(),poles3d(i).Y());
1492 TColStd_Array1OfReal weights (1,NbPoles);
1493 TColStd_Array1OfInteger multiplicities (1,crv3d->NbKnots());
1494 TColStd_Array1OfReal knots(1,crv3d->NbKnots());
1495 crv3d->Knots(knots);
1496 crv3d->Weights(weights);
1497 crv3d->Multiplicities(multiplicities);
1498 C2d = new Geom2d_BSplineCurve ( poles2d, weights, knots, multiplicities, crv3d->Degree(), crv3d->IsPeriodic());
1501 catch(Standard_Failure const& anException) {
1505 Standard_Integer nbp = params->Length();
1506 Standard_Integer nb2 = points2d->Length();
1507 cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::ApproximatePCurve(): Exception: ";
1508 anException.Print(cout);
1509 cout<<"Pb Geom2dAPI_Approximate, tol2d="<<theTolerance2d<<" NbParams="<<nbp<<" NbPnts="<<nb2<<endl;
1510 // if (nb2 > nbp) nb2 = nbp;
1511 // Standard_Real rbp,rb2; rbp = nbp; rb2 = nb2;
1512 // // dbl.AddString ("NbP2d/NbParams puis X Y Param -> mini");
1513 // dbl.AddReals (rb2,rbp);
1514 // for (Standard_Integer i = 1; i <= nb2; i ++) {
1515 // gp_XYZ quoi (points2d->Value(i).X(),points2d->Value(i).Y(),params->Value(i) );
1516 // dbl.AddXYZ (quoi);
1525 //=======================================================================
1526 //function : InterpolatePCurve
1528 //=======================================================================
1530 Handle(Geom2d_Curve) ShapeConstruct_ProjectCurveOnSurface::InterpolatePCurve(const Standard_Integer nbrPnt,
1531 Handle(TColgp_HArray1OfPnt2d)& points2d,
1532 Handle(TColStd_HArray1OfReal)& params,
1533 const Handle(Geom_Curve)& /*orig*/) const
1535 Handle(Geom2d_Curve) C2d; // NULL si echec
1536 Standard_Real theTolerance2d = myPreci / (100 * nbrPnt);
1539 // on verifie d abord s il n y a pas de points confondus
1540 // si besoin on retouche les valeurs ...
1541 CheckPoints2d (points2d, params, theTolerance2d);
1542 Geom2dAPI_Interpolate myInterPol2d (points2d, params,
1543 Standard_False, theTolerance2d);
1544 myInterPol2d.Perform();
1545 if (myInterPol2d.IsDone()) C2d = myInterPol2d.Curve();
1547 catch(Standard_Failure const& anException) {
1551 Standard_Integer nbp = params->Length();
1552 Standard_Integer nb2 = points2d->Length();
1553 cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::InterpolatePCurve(): Exception: ";
1554 anException.Print(cout);
1555 cout<<"Pb Geom2dAPI_Interpolate, tol2d="<<theTolerance2d<<" NbParams="<<nbp<<" NbPnts="<<nb2<<endl;
1556 // if (nb2 > nbp) nb2 = nbp;
1557 // Standard_Real rbp,rb2; rbp = nbp; rb2 = nb2;
1558 // // dbl.AddString ("NbP2d/NbParams puis X Y Param -> mini");
1559 // dbl.AddReals (rb2,rbp);
1560 // for (Standard_Integer i = 1; i <= nb2; i ++) {
1561 // gp_XYZ quoi (points2d->Value(i).X(),points2d->Value(i).Y(),params->Value(i) );
1562 // dbl.AddXYZ (quoi);
1571 //=======================================================================
1572 //function : InterpolateCurve3d
1574 //=======================================================================
1576 Handle(Geom_Curve) ShapeConstruct_ProjectCurveOnSurface::InterpolateCurve3d(const Standard_Integer,
1577 Handle(TColgp_HArray1OfPnt)& points,
1578 Handle(TColStd_HArray1OfReal)& params,
1579 const Handle(Geom_Curve)& /*orig*/) const
1581 Handle(Geom_Curve) C3d; // NULL si echec
1584 Standard_Real Tol = myPreci;
1585 CheckPoints(points, params, Tol);
1586 GeomAPI_Interpolate myInterPol(points, params, Standard_False, Tol);
1587 myInterPol.Perform();
1588 if (myInterPol.IsDone()) C3d = myInterPol.Curve();
1590 catch(Standard_Failure const& anException) {
1593 cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::InterpolateCurve3d(): Exception: ";
1594 anException.Print(cout); cout << endl;
1602 //=======================================================================
1603 //function : CheckPoints
1605 //=======================================================================
1607 void ShapeConstruct_ProjectCurveOnSurface::CheckPoints(Handle(TColgp_HArray1OfPnt)& points,Handle(TColStd_HArray1OfReal)& params,Standard_Real& preci) const
1609 Standard_Integer firstElem = points->Lower();
1610 Standard_Integer lastElem = points->Upper();
1612 Standard_Integer nbPntDropped = 0;
1613 Standard_Integer lastValid = firstElem; // indice of last undropped point
1615 // will store 0 when the point is to be removed, 1 otherwise
1616 TColStd_Array1OfInteger tmpParam(firstElem, lastElem);
1617 for (i = firstElem; i<=lastElem ; i++) tmpParam.SetValue(i,1);
1618 Standard_Real DistMin2 = RealLast();
1619 gp_Pnt Prev = points->Value (lastValid);
1621 for (i = firstElem + 1; i <= lastElem ; i ++) {
1622 Curr = points->Value(i);
1623 Standard_Real CurDist2 = Prev.SquareDistance(Curr);
1624 if (CurDist2 < gp::Resolution()) { // test 0
1626 if ( i == lastElem ) tmpParam.SetValue(lastValid, 0); // last point kept
1627 else tmpParam.SetValue(i, 0); // current dropped, lastValid unchanged
1629 if (CurDist2 < DistMin2)
1630 DistMin2 = CurDist2;
1631 // lastValid becomes the current (i.e. i)
1636 if (DistMin2 < RealLast())
1637 preci = 0.9 * Sqrt (DistMin2); // preci est la distance min entre les points on la reduit un peu
1638 if (nbPntDropped == 0)
1642 cout << "Warning : removing 3d points for interpolation" << endl;
1644 // Build new HArrays
1645 Standard_Integer newLast = lastElem - nbPntDropped;
1646 if ((newLast - firstElem + 1) < 2) {
1648 cout << "Too many degenerated points for 3D interpolation" << endl;
1652 Handle(TColgp_HArray1OfPnt) newPnts =
1653 new TColgp_HArray1OfPnt(firstElem, newLast);
1654 Handle(TColStd_HArray1OfReal) newParams =
1655 new TColStd_HArray1OfReal(firstElem, newLast);
1656 Standard_Integer newCurr = 1;
1657 for (i = firstElem; i<= lastElem ; i++) {
1658 if (tmpParam.Value(i) == 1) {
1659 newPnts->SetValue(newCurr, points->Value(i));
1660 newParams->SetValue(newCurr, params->Value(i));
1666 // on la reduit un peu
1669 //=======================================================================
1670 //function : CheckPoints2d
1672 //=======================================================================
1674 void ShapeConstruct_ProjectCurveOnSurface::CheckPoints2d(Handle(TColgp_HArray1OfPnt2d)& points,
1675 Handle(TColStd_HArray1OfReal)& params,
1676 Standard_Real& preci) const
1678 Standard_Integer firstElem = points->Lower();
1679 Standard_Integer lastElem = points->Upper();
1681 Standard_Integer nbPntDropped = 0;
1682 Standard_Integer lastValid = firstElem; // indice of last undropped point
1684 // will store 0 when the point is to be removed, 1 otherwise
1685 TColStd_Array1OfInteger tmpParam(firstElem, lastElem);
1686 for (i = firstElem; i<=lastElem ; i++) {
1687 tmpParam.SetValue(i,1);
1689 Standard_Real DistMin2 = RealLast();
1690 gp_Pnt2d Prev = points->Value(lastValid);
1692 for (i = firstElem + 1; i<=lastElem ; i++) {
1693 Curr = points->Value(i);
1694 Standard_Real CurDist2 = Prev.SquareDistance(Curr);
1695 if (CurDist2 < gp::Resolution()) { // test 0
1697 if ( i == lastElem ) tmpParam.SetValue(lastValid, 0); // last point kept
1698 else tmpParam.SetValue(i, 0); // current dropped, lastValid unchanged
1700 if (CurDist2 < DistMin2)
1701 DistMin2 = CurDist2;
1702 // lastValid becomes the current (i.e. i)
1707 if (DistMin2 < RealLast())
1708 preci = 0.9 * Sqrt (DistMin2);
1709 if (nbPntDropped == 0)
1713 cout << "Warning : removing 2d points for interpolation" << endl;
1715 // Build new HArrays
1716 Standard_Integer newLast = lastElem - nbPntDropped;
1717 if ((newLast - firstElem + 1) < 2) {
1719 cout << "Too many degenerated points for 2D interpolation" << endl;
1721 //pdn 12.02.99 S4135 Creating pcurve with minimal length.
1722 tmpParam.SetValue(firstElem,1);
1723 tmpParam.SetValue(lastElem,1);
1724 gp_XY lastPnt = points->Value(lastElem).XY();
1725 lastPnt.Add(gp_XY(preci,preci));
1726 points->SetValue(lastElem,lastPnt);
1727 newLast = firstElem+1;
1730 Handle(TColgp_HArray1OfPnt2d) newPnts =
1731 new TColgp_HArray1OfPnt2d(firstElem, newLast);
1732 Handle(TColStd_HArray1OfReal) newParams =
1733 new TColStd_HArray1OfReal(firstElem, newLast);
1734 Standard_Integer newCurr = 1;
1735 for (i = firstElem; i <= lastElem ; i++) {
1736 if (tmpParam.Value(i) == 1) {
1738 cout << "Point " << i << " : " << points->Value(i).X() << " " << points->Value(i).Y() << " at param " << params->Value(i) << endl;
1740 newPnts->SetValue(newCurr, points->Value(i));
1741 newParams->SetValue(newCurr, params->Value(i));
1746 cout << "Removed " << i << " : " << points->Value(i).X() << " " << points->Value(i).Y() << " at param " << params->Value(i) << endl;
1754 //=======================================================================
1755 //function : IsAnIsoparametric
1757 //=======================================================================
1758 //:S4030: modified for optimization
1759 //:p9 abv 11 Mar 99: PRO7226 #489490: find nearest boundary instead of first one
1761 Standard_Boolean ShapeConstruct_ProjectCurveOnSurface::IsAnIsoparametric(const Standard_Integer nbrPnt,
1762 const TColgp_Array1OfPnt& points,
1763 const TColStd_Array1OfReal& params,
1764 Standard_Boolean& isoTypeU,
1765 Standard_Boolean& p1OnIso,
1767 Standard_Boolean& p2OnIso,
1769 Standard_Boolean& isoPar2d3d,
1770 Handle(Geom_Curve)& cIso,
1773 TColStd_Array1OfReal& pout) const
1778 Standard_Real prec = Precision::Confusion();//myPreci;
1780 Standard_Boolean isoParam = Standard_False;
1781 isoPar2d3d = Standard_False;
1783 Standard_Real U1, U2, V1, V2;
1784 mySurf->Bounds(U1, U2, V1, V2);
1786 if ( mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
1787 Handle(Geom_RectangularTrimmedSurface) sTrim =
1788 Handle(Geom_RectangularTrimmedSurface)::DownCast(mySurf->Surface());
1789 sTrim->Bounds(U1, U2, V1, V2);
1793 Standard_Integer mpt[2]; mpt[0] = mpt[1] = 0;
1794 Standard_Real t, tpar[2] = { 0.0, 0.0 }, isoValue=0.;
1795 Standard_Real mindist2;
1796 Standard_Real mind2[2];
1797 mindist2 = mind2[0] = mind2[1] = 4*prec*prec;
1799 p1OnIso = Standard_False;
1800 p2OnIso = Standard_False;
1801 const Bnd_Box* aBox = 0;
1803 for (Standard_Integer j=1; (j<=4) /*&& !isoParam*/; j++) {
1804 Standard_Real isoVal=0.;
1805 Standard_Boolean isoU=Standard_False; //szv#4:S4163:12Mar99 `isoU` must be Standard_Boolean
1806 Handle(Geom_Curve) cI;
1807 Standard_Real tt1, tt2;
1810 if (Precision::IsInfinite(U1)) continue;
1811 cI = mySurf->UIso(U1);
1812 isoU = Standard_True;
1814 aBox = & mySurf->GetBoxUF();
1817 if (Precision::IsInfinite(U2)) continue;
1818 cI = mySurf->UIso(U2);
1819 isoU = Standard_True;
1821 aBox = & mySurf->GetBoxUL();
1824 if (Precision::IsInfinite(V1)) continue;
1825 cI = mySurf->VIso(V1);
1826 isoU = Standard_False;
1828 aBox = & mySurf->GetBoxVF();
1831 if (Precision::IsInfinite(V2)) continue;
1832 cI = mySurf->VIso(V2);
1833 isoU = Standard_False;
1835 aBox = & mySurf->GetBoxVL();
1840 if (isoU) { tt1 = V1; tt2 = V2; }
1841 else { tt1 = U1; tt2 = U2; }
1847 // PATCH CKY 9-JUL-1998 : protection contre singularite
1849 cI->D0( (tt1+tt2)/2,extmi);
1850 if (ext1.IsEqual(ext2,prec) && ext1.IsEqual(extmi,prec)) continue;
1852 Standard_Boolean PtEQext1 = Standard_False;
1853 Standard_Boolean PtEQext2 = Standard_False;
1855 Standard_Real currd2[2], tp[2] = {0, 0};
1856 Standard_Integer mp[2];
1858 for (Standard_Integer i=0; i<2; i++) {
1860 Standard_Integer k = (i == 0 ? 1 : nbrPnt);
1862 // si ext1 == ext2 => valueP1 == valueP2 => vect null plus tard
1863 currd2[i] = points(k).SquareDistance ( ext1 );
1864 if ( currd2[i] <= prec*prec && !PtEQext1) {
1867 PtEQext1 = Standard_True;
1871 currd2[i] = points(k).SquareDistance ( ext2 );
1872 if ( currd2[i] <= prec*prec && !PtEQext2) {
1875 PtEQext2 = Standard_True;
1879 // On evite de projecter sur un iso degenere
1880 // on doit egalement le faire pour l apex du cone
1881 if (mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface)) && !isoU) {
1885 if(aBox->IsOut(points(k))) continue;
1887 Standard_Real Cf = cI->FirstParameter();
1888 Standard_Real Cl = cI->LastParameter();
1889 if (Precision::IsInfinite(Cf)) Cf = -1000;
1890 if (Precision::IsInfinite(Cl)) Cl = +1000;
1892 ShapeAnalysis_Curve sac;
1893 Standard_Real dist = sac.Project (cI,points(k),prec,pt,t,Cf,Cl);
1894 currd2[i] = dist * dist;
1895 if ((dist <= prec) && (t>= Cf) && (t<=Cl)) {
1901 //:e7 abv 21 Apr 98: ProSTEP TR8, r0501_pe #56679:
1902 // avoid possible null-length curves
1903 if ( mp[0] >0 && mp[1] >0 &&
1904 Abs ( tp[0] - tp[1] ) < Precision::PConfusion() ) continue;
1908 ( ! p1OnIso || currd2[0] < mind2[0] ) ) {
1909 p1OnIso = Standard_True;
1910 mind2[0] = currd2[0]; // LP2.stp #105899: FLT_INVALID_OPERATION on Windows 7 VC 9 Release mode on the whole file
1911 if (isoU) valueP1.SetCoord(isoVal, tp[0]);
1912 else valueP1.SetCoord(tp[0], isoVal);
1916 ( ! p2OnIso || currd2[1] < mind2[1] ) ) {
1917 p2OnIso = Standard_True;
1918 mind2[1] = currd2[1];
1919 if (isoU) valueP2.SetCoord(isoVal, tp[1]);
1920 else valueP2.SetCoord(tp[1], isoVal);
1923 if ( mp[0] <=0 || mp[1] <=0 ) continue;
1925 Standard_Real md2 = currd2[0] + currd2[1];
1926 if ( mindist2 <= md2 ) continue;
1940 // probablely it concerns an isoparametrics
1941 if ( mpt[0] >0 && mpt[1] >0 ) {
1943 p1OnIso = p2OnIso = Standard_True;
1945 valueP1.SetCoord(isoValue, tpar[0]);
1946 valueP2.SetCoord(isoValue, tpar[1]);
1949 valueP1.SetCoord(tpar[0], isoValue);
1950 valueP2.SetCoord(tpar[1], isoValue);
1953 if ( mpt[0] != 3 && mpt[1] != 3 ) {
1954 isoPar2d3d = Standard_True;
1955 for (Standard_Integer i=2; i < nbrPnt && isoPar2d3d; i++){
1956 if (tpar[1] > tpar[0]) t = params(i);
1957 else t = t1+t2-params(i);
1959 if (!points(i).IsEqual(pt, prec)) isoPar2d3d = Standard_False;
1963 if (isoPar2d3d) isoParam = Standard_True;
1965 Standard_Real prevParam = tpar[0];
1966 Standard_Real Cf, Cl;
1967 Standard_Boolean isoByDistance = Standard_True;
1968 Cf = cIso->FirstParameter();
1969 Cl = cIso->LastParameter();
1970 if (Precision::IsInfinite(Cf)) Cf = -1000;
1971 if (Precision::IsInfinite(Cl)) Cl = +1000;
1973 ShapeAnalysis_Curve sac;
1974 for (Standard_Integer i=2; i < nbrPnt && isoByDistance; i++) {
1975 Standard_Real dist = sac.NextProject (prevParam,cIso,points(i),
1977 Standard_False); //:j8 abv 10.12.98: TR10 r0501_db.stp #9423: avoid adjusting to ends
1980 if( (dist > prec) || (t < Cf) || (t > Cl) )
1981 isoByDistance = Standard_False;
1983 if (isoByDistance) isoParam = Standard_True;
1986 /* if (!isoParam) { CKY 29-mai-1997 : garder tout ce qu on peut ?
1987 p1OnIso = Standard_False;
1988 p2OnIso = Standard_False;
1992 catch(Standard_Failure const& anException) {
1994 // pb : on affiche ce qu on peut
1995 for (Standard_Integer numpnt = 1; numpnt <= nbrPnt; numpnt ++) {
1996 cout<<"["<<numpnt<<"]param="<<params(numpnt)<<" point=("<<
1997 points(numpnt).X()<<" "<<points(numpnt).Y()<<" "<<points(numpnt).Z()<<")"<<endl;
1999 cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::IsAnIsoparametric(): Exception: ";
2000 anException.Print(cout); cout << endl;
2003 return Standard_False;
2007 /* S4135 : BestExtremum is commented after IsAnIsoparametric works with Precision::Confusion()
2008 //=======================================================================
2009 //function : BestExtremum
2010 //purpose : auxiliaire prenant le meilleur extremum si ISO car doute possible
2011 //=======================================================================
2013 gp_Pnt2d ShapeConstruct_ProjectCurveOnSurface::BestExtremum(const gp_Pnt2d& P2iso,const gp_Pnt& P3ext,const gp_Pnt& P3next) const
2015 // P2iso a ete calcule depuis P3ext sur une iso externe de la surface
2016 // En principe bon mais circularite possible ... et IsU/VClosed faillible
2017 // (si baillement 1e-4 ou 1e-5, on est dedans !). DONC
2018 // 1/ on privilegie l iso mais a tout hasard on verifie si Surf meilleur
2019 // 2/ si iso, attention a la circularite (cas limite)
2021 // NB : si isoParam, on suppose que P2iso est bon (car il y en a 2). A voir...
2023 // D abord, calcul p2ext depuis la surface. choix surface/iso
2025 Standard_Real prec = Precision::Confusion();//myPreci;
2026 gp_Pnt2d P2cal = mySurf->ValueOfUV(P3ext, prec);
2027 gp_Pnt P3cal = mySurf->Value (P2cal);
2028 Standard_Real dcal = P3ext.Distance (P3cal);
2029 Standard_Real dnxt = P3ext.Distance (P3next);
2030 if (dcal > dnxt) return P2iso; // en fait protection sur BUG (PRO8468)
2032 // On choisit entre P2iso et P2cal, le plus proche de P2next ... !!!
2033 gp_Pnt2d P2next = mySurf->ValueOfUV(P3next, prec);
2034 if (P2next.Distance(P2cal) < P2next.Distance(P2iso)) return P2cal;