src/BndLib/BndLib_Add2dCurve.cxx

   1 // Copyright (c) 1996-1999 Matra Datavision
   2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
   3 //
   4 // This file is part of Open CASCADE Technology software library.
   5 //
   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.
  11 //
  12 // Alternatively, this file may be used under the terms of Open CASCADE
  13 // commercial license or contractual agreement.
  14
  15
  16 #include <Adaptor2d_Curve2d.hxx>
  17 #include <Bnd_Box2d.hxx>
  18 #include <BndLib_Add2dCurve.hxx>
  19 #include <Geom2d_BezierCurve.hxx>
  20 #include <Geom2d_BSplineCurve.hxx>
  21 #include <Geom2d_Circle.hxx>
  22 #include <Geom2d_Conic.hxx>
  23 #include <Geom2d_Curve.hxx>
  24 #include <Geom2d_Ellipse.hxx>
  25 #include <Geom2d_Geometry.hxx>
  26 #include <Geom2d_Hyperbola.hxx>
  27 #include <Geom2d_Line.hxx>
  28 #include <Geom2d_OffsetCurve.hxx>
  29 #include <Geom2d_Parabola.hxx>
  30 #include <Geom2d_TrimmedCurve.hxx>
  31 #include <Geom2dAdaptor_Curve.hxx>
  32 #include <gp.hxx>
  33 #include <Precision.hxx>
  34 #include <Standard_Type.hxx>
  35 #include <math_MultipleVarFunction.hxx>
  36 #include <math_Function.hxx>
  37 #include <math_BrentMinimum.hxx>
  38 #include <math_PSO.hxx>
  39
  40 //=======================================================================
  41 //function : BndLib_Box2dCurve
  42 //purpose  :
  43 //=======================================================================
  44 class BndLib_Box2dCurve  {
  45  public:
  46   BndLib_Box2dCurve();
  47
  48   virtual ~BndLib_Box2dCurve();
  49
  50   void SetCurve(const Handle(Geom2d_Curve)& aC);
  51
  52   const Handle(Geom2d_Curve)& Curve() const;
  53
  54   void SetRange(const Standard_Real aT1,
  55                 const Standard_Real aT2);
  56
  57   void Range(Standard_Real& aT1,
  58              Standard_Real& aT2) const;
  59
  60   const Bnd_Box2d& Box() const;
  61
  62   void Perform();
  63
  64   void PerformOptimal(const Standard_Real Tol);
  65
  66   void Clear();
  67
  68   Standard_Integer ErrorStatus() const;
  69   //
  70   //-----------------------------
  71  protected:
  72   void CheckData();
  73   void GetInfoBase();
  74   void PerformLineConic();
  75   void PerformBezier();
  76   void PerformBSpline();
  77   void PerformOther();
  78   void D0(const Standard_Real, gp_Pnt2d&);
  79   //
  80   void Compute(const Handle(Geom2d_Conic)&,
  81                const GeomAbs_CurveType,
  82                const Standard_Real,
  83                const Standard_Real,
  84                Bnd_Box2d& aBox2D);
  85   //
  86   static
  87     Standard_Integer Compute(const Handle(Geom2d_Conic)&,
  88                              const GeomAbs_CurveType,
  89                              Standard_Real *);
  90   static
  91     Standard_Boolean IsTypeBase(const Handle(Geom2d_Curve)& ,
  92                               GeomAbs_CurveType& );
  93   static
  94     Standard_Real AdjustToPeriod(const Standard_Real ,
  95                                  const Standard_Real );
  96   //
  97   void PerformOnePoint();
  98   //
  99   void PerformGenCurv(const Standard_Real Tol = Precision::PConfusion());
 100   //
 101   Standard_Integer NbSamples();
 102   //
 103   Standard_Real AdjustExtr(const Standard_Real UMin,
 104                            const Standard_Real UMax,
 105                            const Standard_Real Extr0,
 106                            const Standard_Integer CoordIndx,
 107                            const Standard_Real Tol,
 108                            const Standard_Boolean IsMin);
 109   //-----------------------------
 110  protected:
 111   Handle(Geom2d_Curve) myCurve;
 112   Bnd_Box2d myBox;
 113   Standard_Integer myErrorStatus;
 114   Handle(Geom2d_Curve) myCurveBase;
 115   Standard_Real myOffsetBase;
 116   Standard_Boolean myOffsetFlag;
 117   Standard_Real myT1;
 118   Standard_Real myT2;
 119   GeomAbs_CurveType myTypeBase;
 120 };
 121 //
 122 class Curv2dMaxMinCoordMVar : public math_MultipleVarFunction
 123 {
 124 public:
 125   Curv2dMaxMinCoordMVar(const Handle(Geom2d_Curve)& theCurve,
 126                         const Standard_Real UMin,
 127                         const Standard_Real UMax,
 128                         const Standard_Integer CoordIndx,
 129                         const Standard_Real Sign)
 130 : myCurve(theCurve),
 131   myUMin(UMin),
 132   myUMax(UMax),
 133   myCoordIndx(CoordIndx),
 134   mySign(Sign)
 135   {
 136   }
 137
 138   Standard_Boolean Value (const math_Vector& X,
 139                                 Standard_Real& F)
 140   {
 141     if (!CheckInputData(X(1)))
 142     {
 143       return Standard_False;
 144     }
 145     gp_Pnt2d aP = myCurve->Value(X(1));
 146
 147     F = mySign * aP.Coord(myCoordIndx);
 148
 149     return Standard_True;
 150   }
 151
 152
 153
 154   Standard_Integer NbVariables() const
 155   {
 156     return 1;
 157   }
 158
 159 private:
 160   Curv2dMaxMinCoordMVar & operator = (const Curv2dMaxMinCoordMVar & theOther);
 161
 162   Standard_Boolean CheckInputData(Standard_Real theParam)
 163   {
 164     if (theParam < myUMin ||
 165         theParam > myUMax)
 166       return Standard_False;
 167     return Standard_True;
 168   }
 169
 170   const Handle(Geom2d_Curve)& myCurve;
 171   Standard_Real myUMin;
 172   Standard_Real myUMax;
 173   Standard_Integer myCoordIndx;
 174   Standard_Real mySign;
 175 };
 176 //
 177 class Curv2dMaxMinCoord : public math_Function
 178 {
 179 public:
 180   Curv2dMaxMinCoord(const Handle(Geom2d_Curve)& theCurve,
 181                      const Standard_Real UMin,
 182                      const Standard_Real UMax,
 183                      const Standard_Integer CoordIndx,
 184                      const Standard_Real Sign)
 185 : myCurve(theCurve),
 186   myUMin(UMin),
 187   myUMax(UMax),
 188   myCoordIndx(CoordIndx),
 189   mySign(Sign)
 190   {
 191   }
 192
 193   Standard_Boolean Value (const Standard_Real X,
 194                                 Standard_Real& F)
 195   {
 196     if (!CheckInputData(X))
 197     {
 198       return Standard_False;
 199     }
 200     gp_Pnt2d aP = myCurve->Value(X);
 201
 202     F = mySign * aP.Coord(myCoordIndx);
 203
 204     return Standard_True;
 205   }
 206
 207 private:
 208   Curv2dMaxMinCoord & operator = (const Curv2dMaxMinCoord & theOther);
 209
 210   Standard_Boolean CheckInputData(Standard_Real theParam)
 211   {
 212     if (theParam < myUMin ||
 213         theParam > myUMax)
 214       return Standard_False;
 215     return Standard_True;
 216   }
 217
 218   const Handle(Geom2d_Curve)& myCurve;
 219   Standard_Real myUMin;
 220   Standard_Real myUMax;
 221   Standard_Integer myCoordIndx;
 222   Standard_Real mySign;
 223 };
 224
 225 //=======================================================================
 226 //function :
 227 //purpose  :
 228 //=======================================================================
 229 BndLib_Box2dCurve::BndLib_Box2dCurve()
 230 {
 231   Clear();
 232 }
 233 //=======================================================================
 234 //function : ~
 235 //purpose  :
 236 //=======================================================================
 237 BndLib_Box2dCurve::~BndLib_Box2dCurve()
 238 {
 239 }
 240 //=======================================================================
 241 //function : Clear
 242 //purpose  :
 243 //=======================================================================
 244 void BndLib_Box2dCurve::Clear()
 245 {
 246   myBox.SetVoid();
 247   //
 248   myErrorStatus=-1;
 249   myTypeBase=GeomAbs_OtherCurve;
 250   myOffsetBase=0.;
 251   myOffsetFlag=Standard_False;
 252 }
 253 //=======================================================================
 254 //function : SetCurve
 255 //purpose  :
 256 //=======================================================================
 257 void BndLib_Box2dCurve::SetCurve(const Handle(Geom2d_Curve)& aC2D)
 258 {
 259   myCurve=aC2D;
 260 }
 261 //=======================================================================
 262 //function : Curve
 263 //purpose  :
 264 //=======================================================================
 265 const Handle(Geom2d_Curve)& BndLib_Box2dCurve::Curve()const
 266 {
 267   return myCurve;
 268 }
 269 //=======================================================================
 270 //function : SetRange
 271 //purpose  :
 272 //=======================================================================
 273 void BndLib_Box2dCurve::SetRange(const Standard_Real aT1,
 274                                  const Standard_Real aT2)
 275 {
 276   myT1=aT1;
 277   myT2=aT2;
 278 }
 279 //=======================================================================
 280 //function : tRange
 281 //purpose  :
 282 //=======================================================================
 283 void BndLib_Box2dCurve::Range(Standard_Real& aT1,
 284                               Standard_Real& aT2) const
 285 {
 286   aT1=myT1;
 287   aT2=myT2;
 288 }
 289 //=======================================================================
 290 //function : ErrorStatus
 291 //purpose  :
 292 //=======================================================================
 293 Standard_Integer BndLib_Box2dCurve::ErrorStatus()const
 294 {
 295   return myErrorStatus;
 296 }
 297 //=======================================================================
 298 //function : Box
 299 //purpose  :
 300 //=======================================================================
 301 const Bnd_Box2d& BndLib_Box2dCurve::Box()const
 302 {
 303   return myBox;
 304 }
 305 //=======================================================================
 306 //function : CheckData
 307 //purpose  :
 308 //=======================================================================
 309 void BndLib_Box2dCurve::CheckData()
 310 {
 311   myErrorStatus=0;
 312   //
 313   if(myCurve.IsNull()) {
 314     myErrorStatus=10;
 315     return;
 316   }
 317   //
 318   if(myT1>myT2) {
 319     myErrorStatus=12; // invalid range
 320     return;
 321   }
 322 }
 323 //=======================================================================
 324 //function : Perform
 325 //purpose  :
 326 //=======================================================================
 327 void BndLib_Box2dCurve::Perform()
 328 {
 329   Clear();
 330   //
 331   myErrorStatus=0;
 332   //
 333   CheckData();
 334   if(myErrorStatus) {
 335     return;
 336   }
 337   //
 338   if (myT1==myT2) {
 339     PerformOnePoint();
 340     return;
 341   }
 342   //
 343   GetInfoBase();
 344   if(myErrorStatus) {
 345     return;
 346   }
 347   //
 348   if (myTypeBase==GeomAbs_Line ||
 349       myTypeBase==GeomAbs_Circle ||
 350       myTypeBase==GeomAbs_Ellipse ||
 351       myTypeBase==GeomAbs_Parabola ||
 352       myTypeBase==GeomAbs_Hyperbola) { // LineConic
 353     PerformLineConic();
 354   }
 355   else if (myTypeBase==GeomAbs_BezierCurve) { // Bezier
 356     PerformBezier();
 357   }
 358   else if (myTypeBase==GeomAbs_BSplineCurve) { //B-Spline
 359     PerformBSpline();
 360   }
 361   else {
 362     myErrorStatus=11; // unknown type base
 363   }
 364 }
 365 //=======================================================================
 366 //function : PerformOptimal
 367 //purpose  :
 368 //=======================================================================
 369 void BndLib_Box2dCurve::PerformOptimal(const Standard_Real Tol)
 370 {
 371   Clear();
 372   myErrorStatus=0;
 373   CheckData();
 374
 375   if(myErrorStatus) {
 376     return;
 377   }
 378
 379   if (myT1==myT2) {
 380     PerformOnePoint();
 381     return;
 382   }
 383
 384   GetInfoBase();
 385   if(myErrorStatus) {
 386     return;
 387   }
 388
 389   if (myTypeBase==GeomAbs_Line ||
 390       myTypeBase==GeomAbs_Circle ||
 391       myTypeBase==GeomAbs_Ellipse ||
 392       myTypeBase==GeomAbs_Parabola ||
 393       myTypeBase==GeomAbs_Hyperbola) { // LineConic
 394     PerformLineConic();
 395   }
 396   else {
 397     PerformGenCurv(Tol);
 398   }
 399 }
 400 //=======================================================================
 401 //function : PerformOnePoint
 402 //purpose  :
 403 //=======================================================================
 404 void BndLib_Box2dCurve::PerformOnePoint()
 405 {
 406   gp_Pnt2d aP2D;
 407   //
 408   myCurve->D0(myT1, aP2D);
 409   myBox.Add(aP2D);
 410 }
 411 //=======================================================================
 412 //function : PerformBezier
 413 //purpose  :
 414 //=======================================================================
 415 void BndLib_Box2dCurve::PerformBezier()
 416 {
 417   if (myOffsetFlag) {
 418     PerformOther();
 419     return;
 420   }
 421   //
 422   Standard_Integer i, aNbPoles;
 423   Standard_Real aT1, aT2, aTb[2];
 424   gp_Pnt2d aP2D;
 425   Handle(Geom2d_Geometry) aG;
 426   Handle(Geom2d_BezierCurve) aCBz, aCBzSeg;
 427   //
 428   myErrorStatus=0;
 429   Bnd_Box2d& aBox2D=myBox;
 430   //
 431   aCBz=Handle(Geom2d_BezierCurve)::DownCast(myCurveBase);
 432   aT1=aCBz->FirstParameter();
 433   aT2=aCBz->LastParameter();
 434   //
 435   aTb[0]=myT1;
 436   if (aTb[0]<aT1) {
 437     aTb[0]=aT1;
 438   }
 439   //
 440   aTb[1]=myT2;
 441   if (aTb[1]>aT2) {
 442     aTb[1]=aT2;
 443   }
 444   //
 445   if (!(aT1==aTb[0] && aT2==aTb[1])) {
 446     aG=aCBz->Copy();
 447     //
 448     aCBzSeg=Handle(Geom2d_BezierCurve)::DownCast(aG);
 449     aCBzSeg->Segment(aTb[0], aTb[1]);
 450     aCBz=aCBzSeg;
 451   }
 452   //
 453   aNbPoles=aCBz->NbPoles();
 454   for (i=1; i<=aNbPoles; ++i) {
 455     aP2D=aCBz->Pole(i);
 456     aBox2D.Add(aP2D);
 457   }
 458 }
 459 //=======================================================================
 460 //function : PerformBSpline
 461 //purpose  :
 462 //=======================================================================
 463 void BndLib_Box2dCurve::PerformBSpline()
 464 {
 465   if (myOffsetFlag) {
 466     PerformOther();
 467     return;
 468   }
 469   //
 470   Standard_Integer i, aNbPoles;
 471   Standard_Real  aT1, aT2, aTb[2];
 472   gp_Pnt2d aP2D;
 473   Handle(Geom2d_Geometry) aG;
 474   Handle(Geom2d_BSplineCurve) aCBS, aCBSs;
 475   //
 476   myErrorStatus=0;
 477   Bnd_Box2d& aBox2D=myBox;
 478   //
 479   aCBS=Handle(Geom2d_BSplineCurve)::DownCast(myCurveBase);
 480   aT1=aCBS->FirstParameter();
 481   aT2=aCBS->LastParameter();
 482   //
 483   aTb[0]=myT1;
 484   if (aTb[0]<aT1) {
 485     aTb[0]=aT1;
 486   }
 487   aTb[1]=myT2;
 488   if (aTb[1]>aT2) {
 489     aTb[1]=aT2;
 490   }
 491
 492   if(aTb[1] < aTb[0])
 493   {
 494     aTb[0]=aT1;
 495     aTb[1]=aT2;
 496   }
 497
 498   //
 499   const Standard_Real eps = Precision::PConfusion();
 500   if (fabs(aT1-aTb[0]) > eps || fabs(aT2-aTb[1]) > eps) {
 501     aG=aCBS->Copy();
 502     //
 503     aCBSs=Handle(Geom2d_BSplineCurve)::DownCast(aG);
 504     aCBSs->Segment(aTb[0], aTb[1]);
 505     aCBS=aCBSs;
 506   }
 507   //
 508   aNbPoles=aCBS->NbPoles();
 509   for (i=1; i<=aNbPoles; ++i) {
 510     aP2D=aCBS->Pole(i);
 511     aBox2D.Add(aP2D);
 512   }
 513 }
 514 //=======================================================================
 515 //function : PerformOther
 516 //purpose  :
 517 //=======================================================================
 518 void BndLib_Box2dCurve::PerformOther()
 519 {
 520   Standard_Integer j, aNb;
 521   Standard_Real aT, dT;
 522   gp_Pnt2d aP2D;
 523   //
 524   aNb=33;
 525   dT=(myT2-myT1)/(aNb-1);
 526   //
 527   aT=myT1;
 528   for (j=0; j<aNb; ++j) {
 529     aT=j*dT;
 530     myCurve->D0(aT, aP2D);
 531     myBox.Add(aP2D);
 532   }
 533   myCurve->D0(myT2, aP2D);
 534   myBox.Add(aP2D);
 535 }
 536 //=======================================================================
 537 //function : NbSamples
 538 //purpose  :
 539 //=======================================================================
 540 Standard_Integer BndLib_Box2dCurve::NbSamples()
 541 {
 542   Standard_Integer N;
 543   switch (myTypeBase) {
 544   case GeomAbs_BezierCurve:
 545     {
 546       Handle(Geom2d_BezierCurve) aCBz=Handle(Geom2d_BezierCurve)::DownCast(myCurveBase);
 547       N = aCBz->NbPoles();
 548       //By default parametric range of Bezier curv is [0, 1]
 549       Standard_Real du = myT2 - myT1;
 550       if(du < .9)
 551       {
 552         N = RealToInt(du*N) + 1;
 553         N = Max(N, 5);
 554       }
 555       break;
 556     }
 557   case GeomAbs_BSplineCurve:
 558     {
 559       Handle(Geom2d_BSplineCurve) aCBS=Handle(Geom2d_BSplineCurve)::DownCast(myCurveBase);
 560       N = (aCBS->Degree() + 1)*(aCBS->NbKnots() -1);
 561       Standard_Real umin = aCBS->FirstParameter(),
 562                     umax = aCBS->LastParameter();
 563       Standard_Real du = (myT2 - myT1) / (umax - umin);
 564       if(du < .9)
 565       {
 566         N = RealToInt(du*N) + 1;
 567         N = Max(N, 5);
 568       }
 569       break;
 570     }
 571   default:
 572     N = 17;
 573   }
 574   return Min (23,N);
 575 }
 576 //=======================================================================
 577 //function : AdjustExtr
 578 //purpose  :
 579 //=======================================================================
 580 Standard_Real BndLib_Box2dCurve::AdjustExtr(const Standard_Real UMin,
 581                                             const Standard_Real UMax,
 582                                             const Standard_Real Extr0,
 583                                             const Standard_Integer CoordIndx,
 584                                             const Standard_Real Tol,
 585                                             const Standard_Boolean IsMin)
 586 {
 587   Standard_Real aSign = IsMin ? 1.:-1.;
 588   Standard_Real extr = aSign * Extr0;
 589   //
 590   Standard_Real Du = (myCurve->LastParameter() - myCurve->FirstParameter());
 591   //
 592   Geom2dAdaptor_Curve aGAC(myCurve);
 593   Standard_Real UTol = Max(aGAC.Resolution(Tol), Precision::PConfusion());
 594   Standard_Real reltol = UTol / Max(Abs(UMin), Abs(UMax));
 595   if(UMax - UMin < 0.01 * Du)
 596   {
 597     //It is suggested that function has one extremum on small interval
 598     math_BrentMinimum anOptLoc(reltol, 100, UTol);
 599     Curv2dMaxMinCoord aFunc(myCurve, UMin, UMax, CoordIndx, aSign);
 600     anOptLoc.Perform(aFunc, UMin, (UMin+UMax)/2., UMax);
 601     if(anOptLoc.IsDone())
 602     {
 603       extr = anOptLoc.Minimum();
 604       return aSign * extr;
 605     }
 606   }
 607   //
 608   Standard_Integer aNbParticles = Max(8, RealToInt(32 * (UMax - UMin) / Du));
 609   Standard_Real maxstep = (UMax - UMin) / (aNbParticles + 1);
 610   math_Vector aT(1,1);
 611   math_Vector aLowBorder(1,1);
 612   math_Vector aUppBorder(1,1);
 613   math_Vector aSteps(1,1);
 614   aLowBorder(1) = UMin;
 615   aUppBorder(1) = UMax;
 616   aSteps(1) = Min(0.1 * Du, maxstep);
 617
 618   Curv2dMaxMinCoordMVar aFunc(myCurve, UMin, UMax, CoordIndx, aSign);
 619   math_PSO aFinder(&aFunc, aLowBorder, aUppBorder, aSteps, aNbParticles);
 620   aFinder.Perform(aSteps, extr, aT);
 621   //
 622   math_BrentMinimum anOptLoc(reltol, 100, UTol);
 623   Curv2dMaxMinCoord aFunc1(myCurve, UMin, UMax, CoordIndx, aSign);
 624   anOptLoc.Perform(aFunc1, Max(aT(1) - aSteps(1), UMin), aT(1), Min(aT(1) + aSteps(1), UMax));
 625
 626   if(anOptLoc.IsDone())
 627   {
 628     extr = anOptLoc.Minimum();
 629     return aSign * extr;
 630   }
 631
 632   return aSign * extr;
 633 }
 634
 635 //=======================================================================
 636 //function : PerformGenCurv
 637 //purpose  :
 638 //=======================================================================
 639 void BndLib_Box2dCurve::PerformGenCurv(const Standard_Real Tol)
 640 {
 641   //
 642   Standard_Integer Nu = NbSamples();
 643   //
 644   Standard_Real CoordMin[2] = {RealLast(), RealLast()};
 645   Standard_Real CoordMax[2] = {-RealLast(), -RealLast()};
 646   Standard_Real DeflMax[2] = {-RealLast(), -RealLast()};
 647   //
 648   gp_Pnt2d P;
 649   Standard_Integer i, k;
 650   Standard_Real du = (myT2 - myT1)/(Nu-1), du2 = du / 2.;
 651   NCollection_Array1<gp_XY> aPnts(1, Nu);
 652   Standard_Real u;
 653   for (i = 1, u = myT1; i <= Nu; i++, u += du)
 654   {
 655     D0(u,P);
 656     aPnts(i) = P.XY();
 657     //
 658     for(k = 0; k < 2; ++k)
 659     {
 660       if(CoordMin[k] > P.Coord(k+1))
 661       {
 662         CoordMin[k] = P.Coord(k+1);
 663       }
 664       if(CoordMax[k] < P.Coord(k+1))
 665       {
 666         CoordMax[k] = P.Coord(k+1);
 667       }
 668     }
 669     //
 670     if(i > 1)
 671     {
 672       gp_XY aPm = 0.5 * (aPnts(i-1) + aPnts(i));
 673       D0(u - du2, P);
 674       gp_XY aD = (P.XY() - aPm);
 675       for(k = 0; k < 2; ++k)
 676       {
 677         if(CoordMin[k] > P.Coord(k+1))
 678         {
 679           CoordMin[k] = P.Coord(k+1);
 680         }
 681         if(CoordMax[k] < P.Coord(k+1))
 682         {
 683           CoordMax[k] = P.Coord(k+1);
 684         }
 685         Standard_Real d = Abs(aD.Coord(k+1));
 686         if(DeflMax[k] < d)
 687         {
 688           DeflMax[k] = d;
 689         }
 690       }
 691     }
 692   }
 693   //
 694   //Adjusting minmax
 695   for(k = 0; k < 2; ++k)
 696   {
 697     Standard_Real d = DeflMax[k];
 698     if(d <= Tol)
 699     {
 700       continue;
 701     }
 702     Standard_Real CMin = CoordMin[k];
 703     Standard_Real CMax = CoordMax[k];
 704     for(i = 1; i <= Nu; ++i)
 705     {
 706       if(aPnts(i).Coord(k+1) - CMin < d)
 707       {
 708         Standard_Real tmin, tmax;
 709         tmin = myT1 + Max(0, i-2) * du;
 710         tmax = myT1 + Min(Nu-1, i) * du;
 711         Standard_Real cmin = AdjustExtr(tmin, tmax,
 712                                         CMin, k + 1, Tol, Standard_True);
 713         if(cmin < CMin)
 714         {
 715           CMin = cmin;
 716         }
 717       }
 718       else if(CMax - aPnts(i).Coord(k+1) < d)
 719       {
 720         Standard_Real tmin, tmax;
 721         tmin = myT1 + Max(0, i-2) * du;
 722         tmax = myT1 + Min(Nu-1, i) * du;
 723         Standard_Real cmax = AdjustExtr(tmin, tmax,
 724                                         CMax, k + 1, Tol, Standard_False);
 725         if(cmax > CMax)
 726         {
 727           CMax = cmax;
 728         }
 729       }
 730     }
 731     CoordMin[k] = CMin;
 732     CoordMax[k] = CMax;
 733   }
 734
 735   myBox.Add(gp_Pnt2d(CoordMin[0], CoordMin[1]));
 736   myBox.Add(gp_Pnt2d(CoordMax[0], CoordMax[1]));
 737   myBox.Enlarge(Tol);
 738 }
 739 //=======================================================================
 740 //function : D0
 741 //purpose  :
 742 //=======================================================================
 743 void BndLib_Box2dCurve::D0(const Standard_Real aU,
 744                             gp_Pnt2d& aP2D)
 745 {
 746   gp_Vec2d aV1;
 747   //
 748   myCurveBase->D1(aU, aP2D, aV1);
 749   //
 750   if (myOffsetFlag) {
 751     Standard_Integer aIndex, aMaxDegree;
 752     Standard_Real aA, aB, aR, aRes;
 753     //
 754     aMaxDegree=9;
 755     aIndex = 2;
 756     aRes=gp::Resolution();
 757     //
 758     while (aV1.Magnitude() <= aRes && aIndex <= aMaxDegree) {
 759       aV1=myCurveBase->DN(aU, aIndex);
 760       ++aIndex;
 761     }
 762     //
 763     aA=aV1.Y();
 764     aB=-aV1.X();
 765     aR=sqrt(aA*aA+aB*aB);
 766     if(aR<=aRes) {
 767       myErrorStatus=13;
 768       return;
 769     }
 770     //
 771     aR=myOffsetBase/aR;
 772     aA=aA*aR;
 773     aB=aB*aR;
 774     aP2D.SetCoord(aP2D.X()+aA, aP2D.Y()+aB);
 775   }
 776   //
 777 }
 778 //=======================================================================
 779 //function : GetInfoBase
 780 //purpose  :
 781 //=======================================================================
 782 void BndLib_Box2dCurve::GetInfoBase()
 783 {
 784   Standard_Boolean bIsTypeBase;
 785   Standard_Integer  iTrimmed, iOffset;
 786   GeomAbs_CurveType aTypeB;
 787   Handle(Geom2d_Curve) aC2DB;
 788   Handle(Geom2d_TrimmedCurve) aCT2D;
 789   Handle(Geom2d_OffsetCurve) aCF2D;
 790   //
 791   myErrorStatus=0;
 792   myTypeBase=GeomAbs_OtherCurve;
 793   myOffsetBase=0;
 794   //
 795   aC2DB=myCurve;
 796   bIsTypeBase=IsTypeBase(aC2DB, aTypeB);
 797   if (bIsTypeBase) {
 798     myTypeBase=aTypeB;
 799     myCurveBase=myCurve;
 800     return;
 801   }
 802   //
 803   while(!bIsTypeBase) {
 804     iTrimmed=0;
 805     iOffset=0;
 806     aCT2D=Handle(Geom2d_TrimmedCurve)::DownCast(aC2DB);
 807     if (!aCT2D.IsNull()) {
 808       aC2DB=aCT2D->BasisCurve();
 809       ++iTrimmed;
 810     }
 811     //
 812     aCF2D=Handle(Geom2d_OffsetCurve)::DownCast(aC2DB);
 813     if (!aCF2D.IsNull()) {
 814       Standard_Real aOffset;
 815       //
 816       aOffset=aCF2D->Offset();
 817       myOffsetBase=myOffsetBase+aOffset;
 818       myOffsetFlag=Standard_True;
 819       //
 820       aC2DB=aCF2D->BasisCurve();
 821       ++iOffset;
 822     }
 823     //
 824     if (!(iTrimmed || iOffset)) {
 825       break;
 826     }
 827     //
 828     bIsTypeBase=IsTypeBase(aC2DB, aTypeB);
 829     if (bIsTypeBase) {
 830       myTypeBase=aTypeB;
 831       myCurveBase=aC2DB;
 832       return;
 833     }
 834   }
 835   //
 836   myErrorStatus=11; // unknown type base
 837 }
 838 //=======================================================================
 839 //function : IsTypeBase
 840 //purpose  :
 841 //=======================================================================
 842 Standard_Boolean BndLib_Box2dCurve::IsTypeBase
 843   (const Handle(Geom2d_Curve)& aC2D,
 844    GeomAbs_CurveType& aTypeB)
 845 {
 846   Standard_Boolean bRet;
 847   Handle(Standard_Type) aType;
 848   //
 849   bRet=Standard_True;
 850   //
 851   aType=aC2D->DynamicType();
 852   if (aType==STANDARD_TYPE(Geom2d_Line)) {
 853     aTypeB=GeomAbs_Line;
 854   }
 855   else if (aType==STANDARD_TYPE(Geom2d_Circle)) {
 856     aTypeB=GeomAbs_Circle;
 857   }
 858   else if (aType==STANDARD_TYPE(Geom2d_Ellipse)) {
 859     aTypeB=GeomAbs_Ellipse;
 860   }
 861   else if (aType==STANDARD_TYPE(Geom2d_Parabola)) {
 862     aTypeB=GeomAbs_Parabola;
 863   }
 864   else if (aType==STANDARD_TYPE(Geom2d_Hyperbola)) {
 865     aTypeB=GeomAbs_Hyperbola;
 866   }
 867   else if (aType==STANDARD_TYPE(Geom2d_BezierCurve)) {
 868     aTypeB=GeomAbs_BezierCurve;
 869   }
 870   else if (aType==STANDARD_TYPE(Geom2d_BSplineCurve)) {
 871     aTypeB=GeomAbs_BSplineCurve;
 872   }
 873   else {
 874     aTypeB=GeomAbs_OtherCurve;
 875     bRet=!bRet;
 876   }
 877   return bRet;
 878 }
 879 //=======================================================================
 880 //function : PerformLineConic
 881 //purpose  :
 882 //=======================================================================
 883 void BndLib_Box2dCurve::PerformLineConic()
 884 {
 885   Standard_Integer i, iInf[2];
 886   Standard_Real  aTb[2];
 887   gp_Pnt2d aP2D;
 888   //
 889   myErrorStatus=0;
 890   //
 891   Bnd_Box2d& aBox2D=myBox;
 892   //
 893   iInf[0]=0;
 894   iInf[1]=0;
 895   aTb[0]=myT1;
 896   aTb[1]=myT2;
 897   //
 898   for (i=0; i<2; ++i) {
 899     if (Precision::IsNegativeInfinite(aTb[i])) {
 900       D0(aTb[i], aP2D);
 901       aBox2D.Add(aP2D);
 902       ++iInf[0];
 903     }
 904     else if (Precision::IsPositiveInfinite(aTb[i])) {
 905       D0(aTb[i], aP2D);
 906       aBox2D.Add(aP2D);
 907       ++iInf[1];
 908     }
 909     else {
 910       D0(aTb[i], aP2D);
 911       aBox2D.Add(aP2D);
 912     }
 913   }
 914   //
 915   if (myTypeBase==GeomAbs_Line) {
 916     return;
 917   }
 918   //
 919   if (iInf[0] && iInf[1]) {
 920     return;
 921   }
 922   //-------------
 923   Handle(Geom2d_Conic) aConic2D;
 924   //
 925   aConic2D=Handle(Geom2d_Conic)::DownCast(myCurveBase);
 926   Compute(aConic2D, myTypeBase, aTb[0], aTb[1], aBox2D);
 927
 928 }
 929 //=======================================================================
 930 //function : Compute
 931 //purpose  :
 932 //=======================================================================
 933 void BndLib_Box2dCurve::Compute(const Handle(Geom2d_Conic)& aConic2D,
 934                                 const GeomAbs_CurveType aType,
 935                                 const Standard_Real aT1,
 936                                 const Standard_Real aT2,
 937                                 Bnd_Box2d& aBox2D)
 938 {
 939   Standard_Integer i, aNbT;
 940   Standard_Real pT[10], aT, aTwoPI, dT, aEps;
 941   gp_Pnt2d aP2D;
 942   //
 943   aNbT=Compute(aConic2D, aType, pT);
 944   //
 945   if (aType==GeomAbs_Parabola ||  aType==GeomAbs_Hyperbola) {
 946     for (i=0; i<aNbT; ++i) {
 947       aT=pT[i];
 948       if (aT>aT1 && aT<aT2) {
 949         D0(aT, aP2D);
 950         aBox2D.Add(aP2D);
 951       }
 952     }
 953     return;
 954   }
 955   //
 956   //aType==GeomAbs_Circle ||  aType==GeomAbs_Ellipse
 957   aEps=1.e-14;
 958   aTwoPI=2.*M_PI;
 959   dT=aT2-aT1;
 960   //
 961   Standard_Real aT1z = AdjustToPeriod (aT1, aTwoPI);
 962   if (fabs(aT1z)<aEps) {
 963     aT1z=0.;
 964   }
 965   //
 966   Standard_Real aT2z = aT1z + dT;
 967   if (fabs(aT2z-aTwoPI)<aEps) {
 968     aT2z=aTwoPI;
 969   }
 970   //
 971   for (i=0; i<aNbT; ++i) {
 972     aT = pT[i];
 973     // adjust aT to range [aT1z, aT1z + 2*PI]; note that pT[i] and aT1z
 974     // are adjusted to range [0, 2*PI], but aT2z can be greater than 2*PI
 975     aT = (aT < aT1z ? aT + aTwoPI : aT);
 976     if (aT <= aT2z) {
 977       D0(aT, aP2D);
 978       aBox2D.Add(aP2D);
 979     }
 980   }
 981 }
 982
 983 //=======================================================================
 984 //function : Compute
 985 //purpose  :
 986 //=======================================================================
 987 Standard_Integer BndLib_Box2dCurve::Compute
 988   (const Handle(Geom2d_Conic)& aConic2D,
 989    const GeomAbs_CurveType aType,
 990    Standard_Real *pT)
 991 {
 992   Standard_Integer iRet, i, j;
 993   Standard_Real aCosBt, aSinBt, aCosGm, aSinGm;
 994   Standard_Real aLx, aLy;
 995   //
 996   iRet=0;
 997   //
 998   const gp_Ax22d& aPos=aConic2D->Position();
 999   const gp_XY& aXDir=aPos.XDirection().XY();
1000   const gp_XY& aYDir=aPos.YDirection().XY();
1001   //
1002   aCosBt=aXDir.X();
1003   aSinBt=aXDir.Y();
1004   aCosGm=aYDir.X();
1005   aSinGm=aYDir.Y();
1006   //
1007   if (aType==GeomAbs_Circle ||  aType==GeomAbs_Ellipse) {
1008     Standard_Real aR1 = 0.0, aR2 = 0.0, aTwoPI = M_PI+M_PI;
1009     Standard_Real aA11 = 0.0, aA12 = 0.0, aA21 = 0.0, aA22 = 0.0;
1010     Standard_Real aBx = 0.0, aBy = 0.0, aB = 0.0, aCosFi = 0.0, aSinFi = 0.0, aFi = 0.0;
1011     //
1012     if(aType==GeomAbs_Ellipse) {
1013       Handle(Geom2d_Ellipse) aEL2D;
1014       //
1015       aEL2D=Handle(Geom2d_Ellipse)::DownCast(aConic2D);
1016       aR1=aEL2D->MajorRadius();
1017       aR2=aEL2D->MinorRadius();
1018     }
1019     else if(aType==GeomAbs_Circle) {
1020       Handle(Geom2d_Circle) aCR2D;
1021       //
1022       aCR2D=Handle(Geom2d_Circle)::DownCast(aConic2D);
1023       aR1=aCR2D->Radius();
1024       aR2=aR1;
1025     }
1026     //
1027     aA11=-aR1*aCosBt;
1028     aA12= aR2*aCosGm;
1029     aA21=-aR1*aSinBt;
1030     aA22= aR2*aSinGm;
1031     //
1032     for (i=0; i<2; ++i) {
1033       aLx=(!i) ? 0. : 1.;
1034       aLy=(!i) ? 1. : 0.;
1035       aBx=aLx*aA21-aLy*aA11;
1036       aBy=aLx*aA22-aLy*aA12;
1037       aB=sqrt(aBx*aBx+aBy*aBy);
1038       //
1039       aCosFi=aBx/aB;
1040       aSinFi=aBy/aB;
1041       //
1042       aFi=acos(aCosFi);
1043       if (aSinFi<0.) {
1044         aFi=aTwoPI-aFi;
1045       }
1046       //
1047       j=2*i;
1048       pT[j]=aTwoPI-aFi;
1049       pT[j]=AdjustToPeriod(pT[j], aTwoPI);
1050       //
1051       pT[j+1]=M_PI-aFi;
1052       pT[j+1]=AdjustToPeriod(pT[j+1], aTwoPI);
1053     }
1054     iRet=4;
1055   }//if (aType==GeomAbs_Ellipse) {
1056   //
1057   else if (aType==GeomAbs_Parabola) {
1058     Standard_Real aFc, aEps;
1059     Standard_Real aA1, aA2;
1060     Handle(Geom2d_Parabola) aPR2D;
1061     //
1062     aEps=1.e-12;
1063     //
1064     aPR2D=Handle(Geom2d_Parabola)::DownCast(aConic2D);
1065     aFc=aPR2D->Focal();
1066     //
1067     j=0;
1068     for (i=0; i<2; i++) {
1069       aLx=(!i) ? 0. : 1.;
1070       aLy=(!i) ? 1. : 0.;
1071       //
1072       aA2=aLx*aSinBt-aLy*aCosBt;
1073       if (fabs(aA2)<aEps) {
1074         continue;
1075       }
1076       //
1077       aA1=aLy*aCosGm-aLx*aSinGm;
1078       //
1079       pT[j]=2.*aFc*aA1/aA2;
1080       ++j;
1081     }
1082     iRet=j;
1083   }// else if (aType==GeomAbs_Parabola) {
1084   //
1085   else if (aType==GeomAbs_Hyperbola) {
1086     Standard_Integer k;
1087     Standard_Real aR1, aR2;
1088     Standard_Real aEps, aB1, aB2, aB12, aB22, aZ, aD;
1089     Handle(Geom2d_Hyperbola) aHP2D;
1090     //
1091     aEps=1.e-12;
1092     //
1093     aHP2D=Handle(Geom2d_Hyperbola)::DownCast(aConic2D);
1094     aR1=aHP2D->MajorRadius();
1095     aR2=aHP2D->MinorRadius();
1096     //
1097     j=0;
1098     for (i=0; i<2; i++) {
1099       aLx=(!i) ? 0. : 1.;
1100       aLy=(!i) ? 1. : 0.;
1101       //
1102       aB1=aR1*(aLx*aSinBt-aLy*aCosBt);
1103       aB2=aR2*(aLx*aSinGm-aLy*aCosGm);
1104       //
1105       if (fabs(aB1)<aEps) {
1106         continue;
1107       }
1108       //
1109       if (fabs(aB2)<aEps) {
1110         pT[j]=0.;
1111         ++j;
1112       }
1113       else {
1114         aB12=aB1*aB1;
1115         aB22=aB2*aB2;
1116         if (!(aB12>aB22)) {
1117           continue;
1118         }
1119         //
1120         aD=sqrt(aB12-aB22);
1121         //-------------
1122         for (k=-1; k<2; k+=2) {
1123           aZ=(aB1+k*aD)/aB2;
1124           if (fabs(aZ)<1.) {
1125             pT[j]=-log((1.+aZ)/(1.-aZ));
1126             ++j;
1127           }
1128         }
1129       }
1130     }
1131     iRet=j;
1132   }// else if (aType==GeomAbs_Hyperbola) {
1133   //
1134   return iRet;
1135 }
1136 //=======================================================================
1137 //function : AdjustToPeriod
1138 //purpose  :
1139 //=======================================================================
1140 Standard_Real BndLib_Box2dCurve::AdjustToPeriod(const Standard_Real aT,
1141                                                 const Standard_Real aPeriod)
1142 {
1143   Standard_Integer k;
1144   Standard_Real aTRet;
1145   //
1146   aTRet=aT;
1147   if (aT<0.) {
1148     k=1+(Standard_Integer)(-aT/aPeriod);
1149     aTRet=aT+k*aPeriod;
1150   }
1151   else if (aT>aPeriod) {
1152     k=(Standard_Integer)(aT/aPeriod);
1153     aTRet=aT-k*aPeriod;
1154   }
1155   if (aTRet==aPeriod) {
1156     aTRet=0.;
1157   }
1158   return aTRet;
1159 }
1160 //
1161 // myErrorStatus:
1162 //
1163 // -1 - object is just initialized
1164 // 10 - myCurve is Null
1165 // 12 - invalid range myT1 >  myT2l
1166 // 11 - unknown type of base curve
1167 // 13 - offset curve can not be computed
1168 //NMTTest
1169
1170 //=======================================================================
1171 //function : Add
1172 //purpose  :
1173 //=======================================================================
1174 void BndLib_Add2dCurve::Add(const Adaptor2d_Curve2d& aC,
1175                              const Standard_Real aTol,
1176                              Bnd_Box2d& aBox2D)
1177 {
1178   BndLib_Add2dCurve::Add(aC,
1179                           aC.FirstParameter(),
1180                           aC.LastParameter (),
1181                           aTol,
1182                           aBox2D);
1183 }
1184 //=======================================================================
1185 //function : Add
1186 //purpose  :
1187 //=======================================================================
1188 void BndLib_Add2dCurve::Add(const Adaptor2d_Curve2d& aC,
1189                              const Standard_Real aU1,
1190                              const Standard_Real aU2,
1191                              const Standard_Real aTol,
1192                              Bnd_Box2d& aBox2D)
1193 {
1194   Adaptor2d_Curve2d *pC=(Adaptor2d_Curve2d *)&aC;
1195   Geom2dAdaptor_Curve *pA=dynamic_cast<Geom2dAdaptor_Curve*>(pC);
1196   if (!pA) {
1197     Standard_Real U, DU;
1198     Standard_Integer N, j;
1199     gp_Pnt2d P;
1200     N = 33;
1201     U  = aU1;
1202     DU = (aU2-aU1)/(N-1);
1203     for (j=1; j<N; j++) {
1204       aC.D0(U,P);
1205       U+=DU;
1206       aBox2D.Add(P);
1207     }
1208     aC.D0(aU2,P);
1209     aBox2D.Add(P);
1210     aBox2D.Enlarge(aTol);
1211     return;
1212   }
1213   //
1214   const Handle(Geom2d_Curve)& aC2D=pA->Curve();
1215   //
1216   BndLib_Add2dCurve::Add(aC2D, aU1, aU2, aTol, aBox2D);
1217 }
1218 //=======================================================================
1219 //function : Add
1220 //purpose  :
1221 //=======================================================================
1222 void BndLib_Add2dCurve::Add(const Handle(Geom2d_Curve)& aC2D,
1223                              const Standard_Real aTol,
1224                              Bnd_Box2d& aBox2D)
1225 {
1226   Standard_Real aT1, aT2;
1227   //
1228   aT1=aC2D->FirstParameter();
1229   aT2=aC2D->LastParameter();
1230   //
1231   BndLib_Add2dCurve::Add(aC2D, aT1, aT2, aTol, aBox2D);
1232 }
1233
1234 //=======================================================================
1235 //function : Add
1236 //purpose  :
1237 //=======================================================================
1238 void BndLib_Add2dCurve::Add(const Handle(Geom2d_Curve)& aC2D,
1239                              const Standard_Real aT1,
1240                              const Standard_Real aT2,
1241                              const Standard_Real aTol,
1242                              Bnd_Box2d& aBox2D)
1243 {
1244   BndLib_Box2dCurve aBC;
1245   //
1246   aBC.SetCurve(aC2D);
1247   aBC.SetRange(aT1, aT2);
1248   //
1249   aBC.Perform();
1250   //
1251   const Bnd_Box2d& aBoxC=aBC.Box();
1252   aBox2D.Add(aBoxC);
1253   aBox2D.Enlarge(aTol);
1254 }
1255 //=======================================================================
1256 //function : AddOptimal
1257 //purpose  :
1258 //=======================================================================
1259 void BndLib_Add2dCurve::AddOptimal(const Handle(Geom2d_Curve)& aC2D,
1260                                    const Standard_Real aT1,
1261                                    const Standard_Real aT2,
1262                                    const Standard_Real aTol,
1263                                    Bnd_Box2d& aBox2D)
1264 {
1265   BndLib_Box2dCurve aBC;
1266   //
1267   aBC.SetCurve(aC2D);
1268   aBC.SetRange(aT1, aT2);
1269   //
1270   aBC.PerformOptimal(aTol);
1271   //
1272   const Bnd_Box2d& aBoxC=aBC.Box();
1273   aBox2D.Add(aBoxC);
1274   aBox2D.Enlarge(aTol);
1275 }