1 // Copyright (c) 1996-1999 Matra Datavision
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
4 // This file is part of Open CASCADE Technology software library.
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
16 #include <Adaptor3d_Curve.hxx>
17 #include <Bnd_Box.hxx>
19 #include <BndLib_Add3dCurve.hxx>
21 #include <Geom_BezierCurve.hxx>
22 #include <Geom_BSplineCurve.hxx>
23 #include <GeomAbs_CurveType.hxx>
24 #include <GeomAdaptor_Curve.hxx>
25 #include <gp_Circ.hxx>
26 #include <gp_Elips.hxx>
27 #include <gp_Hypr.hxx>
29 #include <gp_Parab.hxx>
31 #include <Precision.hxx>
32 #include <TColgp_Array1OfPnt.hxx>
33 #include <TColStd_Array1OfInteger.hxx>
34 #include <TColStd_Array1OfReal.hxx>
35 #include <math_MultipleVarFunction.hxx>
36 #include <math_Function.hxx>
37 #include <math_PSO.hxx>
38 #include <math_BrentMinimum.hxx>
40 static Standard_Integer NbSamples(const Adaptor3d_Curve& C,
41 const Standard_Real Umin,
42 const Standard_Real Umax);
44 static Standard_Real AdjustExtr(const Adaptor3d_Curve& C,
45 const Standard_Real UMin,
46 const Standard_Real UMax,
47 const Standard_Real Extr0,
48 const Standard_Integer CoordIndx,
49 const Standard_Real Tol,
50 const Standard_Boolean IsMin);
53 //=======================================================================
54 //function : reduceSplineBox
55 //purpose : This method intended to reduce box in case of
56 // bezier and bspline curve.
57 //=======================================================================
58 static void reduceSplineBox(const Adaptor3d_Curve& theCurve,
59 const Bnd_Box& theOrigBox,
60 Bnd_Box & theReducedBox)
62 // Guaranteed bounding box based on poles of bspline.
64 Standard_Real aPolesXMin, aPolesYMin, aPolesZMin,
65 aPolesXMax, aPolesYMax, aPolesZMax;
67 if (theCurve.GetType() == GeomAbs_BSplineCurve)
69 Handle(Geom_BSplineCurve) aC = theCurve.BSpline();
70 const TColgp_Array1OfPnt& aPoles = aC->Poles();
72 for(Standard_Integer anIdx = aPoles.Lower();
73 anIdx <= aPoles.Upper();
76 aPolesBox.Add(aPoles.Value(anIdx));
79 if (theCurve.GetType() == GeomAbs_BezierCurve)
81 Handle(Geom_BezierCurve) aC = theCurve.Bezier();
82 const TColgp_Array1OfPnt& aPoles = aC->Poles();
84 for(Standard_Integer anIdx = aPoles.Lower();
85 anIdx <= aPoles.Upper();
88 aPolesBox.Add(aPoles.Value(anIdx));
92 aPolesBox.Get(aPolesXMin, aPolesYMin, aPolesZMin,
93 aPolesXMax, aPolesYMax, aPolesZMax);
95 Standard_Real x, y, z, X, Y, Z;
96 theOrigBox.Get(x, y, z, X, Y, Z);
114 theReducedBox.Update(x, y, z, X, Y, Z);
117 //=======================================================================
120 //=======================================================================
121 void BndLib_Add3dCurve::Add( const Adaptor3d_Curve& C,
122 const Standard_Real Tol,
125 BndLib_Add3dCurve::Add(C,
132 static Standard_Real FillBox(Bnd_Box& B, const Adaptor3d_Curve& C,
133 const Standard_Real first, const Standard_Real last,
134 const Standard_Integer N)
137 C.D0(first,P1); B.Add(P1);
138 Standard_Real p = first, dp = last-first, tol= 0.;
139 if(Abs(dp) > Precision::PConfusion()){
142 for(i = 1; i <= N; i++){
143 p += dp; C.D0(p,P2); B.Add(P2);
144 p += dp; C.D0(p,P3); B.Add(P3);
145 gp_Pnt Pc((P1.XYZ()+P3.XYZ())/2.0);
146 tol = Max(tol,Pc.Distance(P2));
150 C.D0(first,P1); B.Add(P1);
151 C.D0(last,P3); B.Add(P3);
157 //=======================================================================
160 //=======================================================================
162 void BndLib_Add3dCurve::Add( const Adaptor3d_Curve& C,
163 const Standard_Real U1,
164 const Standard_Real U2,
165 const Standard_Real Tol,
168 static Standard_Real weakness = 1.5; //OCC566(apo)
169 Standard_Real tol = 0.0;
170 switch (C.GetType()) {
174 BndLib::Add(C.Line(),U1,U2,Tol,B);
179 BndLib::Add(C.Circle(),U1,U2,Tol,B);
182 case GeomAbs_Ellipse:
184 BndLib::Add(C.Ellipse(),U1,U2,Tol,B);
187 case GeomAbs_Hyperbola:
189 BndLib::Add(C.Hyperbola(),U1,U2,Tol,B);
192 case GeomAbs_Parabola:
194 BndLib::Add(C.Parabola(),U1,U2,Tol,B);
197 case GeomAbs_BezierCurve:
199 Handle(Geom_BezierCurve) Bz = C.Bezier();
200 Standard_Integer N = Bz->Degree();
201 GeomAdaptor_Curve GACurve(Bz);
203 tol = FillBox(B1,GACurve,U1,U2,N);
204 B1.Enlarge(weakness*tol);
205 reduceSplineBox(C, B1, B);
209 case GeomAbs_BSplineCurve:
211 Handle(Geom_BSplineCurve) Bs = C.BSpline();
212 if(Abs(Bs->FirstParameter() - U1) > Precision::Parametric(Tol)||
213 Abs(Bs->LastParameter() - U2) > Precision::Parametric(Tol)) {
215 Handle(Geom_Geometry) G = Bs->Copy();
216 Handle(Geom_BSplineCurve) Bsaux (Handle(Geom_BSplineCurve)::DownCast (G));
217 Standard_Real u1 = U1, u2 = U2;
218 //// modified by jgv, 24.10.01 for BUC61031 ////
219 if (Bsaux->IsPeriodic())
220 ElCLib::AdjustPeriodic( Bsaux->FirstParameter(), Bsaux->LastParameter(), Precision::PConfusion(), u1, u2 );
222 ////////////////////////////////////////////////
223 // modified by NIZHNY-EAP Fri Dec 3 14:29:14 1999 ___BEGIN___
224 // To avoid exeption in Segment
225 if(Bsaux->FirstParameter() > U1) u1 = Bsaux->FirstParameter();
226 if(Bsaux->LastParameter() < U2 ) u2 = Bsaux->LastParameter();
227 // modified by NIZHNY-EAP Fri Dec 3 14:29:18 1999 ___END___
229 Bsaux->Segment(u1, u2);
234 Standard_Integer k, k1 = Bs->FirstUKnotIndex(), k2 = Bs->LastUKnotIndex(),
235 N = Bs->Degree(), NbKnots = Bs->NbKnots();
236 TColStd_Array1OfReal Knots(1,NbKnots);
238 GeomAdaptor_Curve GACurve(Bs);
239 Standard_Real first = Knots(k1), last;
240 for(k = k1 + 1; k <= k2; k++){
242 tol = Max(FillBox(B1,GACurve,first,last,N), tol);
247 B1.Enlarge(weakness*tol);
248 reduceSplineBox(C, B1, B);
257 static Standard_Integer N = 33;
258 tol = FillBox(B1,C,U1,U2,N);
259 B1.Enlarge(weakness*tol);
260 Standard_Real x, y, z, X, Y, Z;
261 B1.Get(x, y, z, X, Y, Z);
262 B.Update(x, y, z, X, Y, Z);
268 //=======================================================================
269 //function : AddOptimal
271 //=======================================================================
273 void BndLib_Add3dCurve::AddOptimal( const Adaptor3d_Curve& C,
274 const Standard_Real Tol,
277 BndLib_Add3dCurve::AddOptimal(C,
283 //=======================================================================
284 //function : AddOptimal
286 //=======================================================================
288 void BndLib_Add3dCurve::AddOptimal( const Adaptor3d_Curve& C,
289 const Standard_Real U1,
290 const Standard_Real U2,
291 const Standard_Real Tol,
294 switch (C.GetType()) {
298 BndLib::Add(C.Line(),U1,U2,Tol,B);
303 BndLib::Add(C.Circle(),U1,U2,Tol,B);
306 case GeomAbs_Ellipse:
308 BndLib::Add(C.Ellipse(),U1,U2,Tol,B);
311 case GeomAbs_Hyperbola:
313 BndLib::Add(C.Hyperbola(),U1,U2,Tol,B);
316 case GeomAbs_Parabola:
318 BndLib::Add(C.Parabola(),U1,U2,Tol,B);
323 AddGenCurv(C, U1, U2, Tol, B);
328 //=======================================================================
329 //function : AddGenCurv
331 //=======================================================================
332 void BndLib_Add3dCurve::AddGenCurv(const Adaptor3d_Curve& C,
333 const Standard_Real UMin,
334 const Standard_Real UMax,
335 const Standard_Real Tol,
338 Standard_Integer Nu = NbSamples(C, UMin, UMax);
340 Standard_Real CoordMin[3] = {RealLast(), RealLast(), RealLast()};
341 Standard_Real CoordMax[3] = {-RealLast(), -RealLast(), -RealLast()};
342 Standard_Real DeflMax[3] = {-RealLast(), -RealLast(), -RealLast()};
345 Standard_Integer i, k;
346 Standard_Real du = (UMax-UMin)/(Nu-1), du2 = du / 2.;
347 NCollection_Array1<gp_XYZ> aPnts(1, Nu);
349 for (i = 1, u = UMin; i <= Nu; i++, u += du)
354 for(k = 0; k < 3; ++k)
356 if(CoordMin[k] > P.Coord(k+1))
358 CoordMin[k] = P.Coord(k+1);
360 if(CoordMax[k] < P.Coord(k+1))
362 CoordMax[k] = P.Coord(k+1);
368 gp_XYZ aPm = 0.5 * (aPnts(i-1) + aPnts(i));
370 gp_XYZ aD = (P.XYZ() - aPm);
371 for(k = 0; k < 3; ++k)
373 if(CoordMin[k] > P.Coord(k+1))
375 CoordMin[k] = P.Coord(k+1);
377 if(CoordMax[k] < P.Coord(k+1))
379 CoordMax[k] = P.Coord(k+1);
381 Standard_Real d = Abs(aD.Coord(k+1));
391 Standard_Real eps = Max(Tol, Precision::Confusion());
392 for(k = 0; k < 3; ++k)
394 Standard_Real d = DeflMax[k];
399 Standard_Real CMin = CoordMin[k];
400 Standard_Real CMax = CoordMax[k];
401 for(i = 1; i <= Nu; ++i)
403 if(aPnts(i).Coord(k+1) - CMin < d)
405 Standard_Real umin, umax;
406 umin = UMin + Max(0, i-2) * du;
407 umax = UMin + Min(Nu-1, i) * du;
408 Standard_Real cmin = AdjustExtr(C, umin, umax,
409 CMin, k + 1, eps, Standard_True);
415 else if(CMax - aPnts(i).Coord(k+1) < d)
417 Standard_Real umin, umax;
418 umin = UMin + Max(0, i-2) * du;
419 umax = UMin + Min(Nu-1, i) * du;
420 Standard_Real cmax = AdjustExtr(C, umin, umax,
421 CMax, k + 1, eps, Standard_False);
432 B.Add(gp_Pnt(CoordMin[0], CoordMin[1], CoordMin[2]));
433 B.Add(gp_Pnt(CoordMax[0], CoordMax[1], CoordMax[2]));
437 class CurvMaxMinCoordMVar : public math_MultipleVarFunction
440 CurvMaxMinCoordMVar(const Adaptor3d_Curve& theCurve,
441 const Standard_Real UMin,
442 const Standard_Real UMax,
443 const Standard_Integer CoordIndx,
444 const Standard_Real Sign)
448 myCoordIndx(CoordIndx),
453 Standard_Boolean Value (const math_Vector& X,
456 if (!CheckInputData(X(1)))
458 return Standard_False;
460 gp_Pnt aP = myCurve.Value(X(1));
462 F = mySign * aP.Coord(myCoordIndx);
464 return Standard_True;
469 Standard_Integer NbVariables() const
475 CurvMaxMinCoordMVar & operator = (const CurvMaxMinCoordMVar & theOther);
477 Standard_Boolean CheckInputData(Standard_Real theParam)
479 if (theParam < myUMin ||
481 return Standard_False;
482 return Standard_True;
485 const Adaptor3d_Curve& myCurve;
486 Standard_Real myUMin;
487 Standard_Real myUMax;
488 Standard_Integer myCoordIndx;
489 Standard_Real mySign;
492 class CurvMaxMinCoord : public math_Function
495 CurvMaxMinCoord(const Adaptor3d_Curve& theCurve,
496 const Standard_Real UMin,
497 const Standard_Real UMax,
498 const Standard_Integer CoordIndx,
499 const Standard_Real Sign)
503 myCoordIndx(CoordIndx),
508 Standard_Boolean Value (const Standard_Real X,
511 if (!CheckInputData(X))
513 return Standard_False;
515 gp_Pnt aP = myCurve.Value(X);
517 F = mySign * aP.Coord(myCoordIndx);
519 return Standard_True;
523 CurvMaxMinCoord & operator = (const CurvMaxMinCoord & theOther);
525 Standard_Boolean CheckInputData(Standard_Real theParam)
527 if (theParam < myUMin ||
529 return Standard_False;
530 return Standard_True;
533 const Adaptor3d_Curve& myCurve;
534 Standard_Real myUMin;
535 Standard_Real myUMax;
536 Standard_Integer myCoordIndx;
537 Standard_Real mySign;
540 //=======================================================================
541 //function : AdjustExtr
543 //=======================================================================
545 Standard_Real AdjustExtr(const Adaptor3d_Curve& C,
546 const Standard_Real UMin,
547 const Standard_Real UMax,
548 const Standard_Real Extr0,
549 const Standard_Integer CoordIndx,
550 const Standard_Real Tol,
551 const Standard_Boolean IsMin)
553 Standard_Real aSign = IsMin ? 1.:-1.;
554 Standard_Real extr = aSign * Extr0;
556 Standard_Real uTol = Max(C.Resolution(Tol), Precision::PConfusion());
557 Standard_Real Du = (C.LastParameter() - C.FirstParameter());
559 Standard_Real reltol = uTol / Max(Abs(UMin), Abs(UMax));
560 if(UMax - UMin < 0.01 * Du)
563 math_BrentMinimum anOptLoc(reltol, 100, uTol);
564 CurvMaxMinCoord aFunc(C, UMin, UMax, CoordIndx, aSign);
565 anOptLoc.Perform(aFunc, UMin, (UMin+UMax)/2., UMax);
566 if(anOptLoc.IsDone())
568 extr = anOptLoc.Minimum();
573 Standard_Integer aNbParticles = Max(8, RealToInt(32 * (UMax - UMin) / Du));
574 Standard_Real maxstep = (UMax - UMin) / (aNbParticles + 1);
576 math_Vector aLowBorder(1,1);
577 math_Vector aUppBorder(1,1);
578 math_Vector aSteps(1,1);
579 aLowBorder(1) = UMin;
580 aUppBorder(1) = UMax;
581 aSteps(1) = Min(0.1 * Du, maxstep);
583 CurvMaxMinCoordMVar aFunc(C, UMin, UMax, CoordIndx, aSign);
584 math_PSO aFinder(&aFunc, aLowBorder, aUppBorder, aSteps, aNbParticles);
585 aFinder.Perform(aSteps, extr, aT);
587 math_BrentMinimum anOptLoc(reltol, 100, uTol);
588 CurvMaxMinCoord aFunc1(C, UMin, UMax, CoordIndx, aSign);
589 anOptLoc.Perform(aFunc1, Max(aT(1) - aSteps(1), UMin), aT(1), Min(aT(1) + aSteps(1), UMax));
591 if(anOptLoc.IsDone())
593 extr = anOptLoc.Minimum();
600 //=======================================================================
601 //function : NbSamples
603 //=======================================================================
605 Standard_Integer NbSamples(const Adaptor3d_Curve& C,
606 const Standard_Real Umin,
607 const Standard_Real Umax)
610 GeomAbs_CurveType Type = C.GetType();
612 case GeomAbs_BezierCurve:
615 //By default parametric range of Bezier curv is [0, 1]
616 Standard_Real du = Umax - Umin;
619 N = RealToInt(du*N) + 1;
624 case GeomAbs_BSplineCurve:
626 const Handle(Geom_BSplineCurve)& BC = C.BSpline();
627 N = 2 * (BC->Degree() + 1)*(BC->NbKnots() -1);
628 Standard_Real umin = BC->FirstParameter(),
629 umax = BC->LastParameter();
630 Standard_Real du = (Umax - Umin) / (umax - umin);
633 N = RealToInt(du*N) + 1;