src/Geom/Geom_BSplineCurve.cxx

   1 // Created on: 1993-03-09
   2 // Created by: JCV
   3 // Copyright (c) 1993-1999 Matra Datavision
   4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
   5 //
   6 // This file is part of Open CASCADE Technology software library.
   7 //
   8 // This library is free software; you can redistribute it and/or modify it under
   9 // the terms of the GNU Lesser General Public License version 2.1 as published
  10 // by the Free Software Foundation, with special exception defined in the file
  11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
  12 // distribution for complete text of the license and disclaimer of any warranty.
  13 //
  14 // Alternatively, this file may be used under the terms of Open CASCADE
  15 // commercial license or contractual agreement.
  16
  17 //Avril 1991 : constructeurs + methodes de lecture.
  18 //Mai 1991   : revue des specifs + debut de realisation des classes tool =>
  19 //             implementation des methodes Set et calcul du point courant.
  20 //Juillet 1991 : voir egalement File Geom_BSplineCurve_1.cxx
  21 //Juin    1992 : mise a plat des valeurs nodales - amelioration des
  22 //               performances sur calcul du point courant
  23
  24 //RLE Aug 1993  Remove Swaps, Remove typedefs, Update BSplCLib
  25 //              debug periodic, IncreaseDegree
  26 // 21-Mar-95 : xab implemented cache
  27 // 14-Mar-96 : xab implemented MovePointAndTangent
  28 // 13-Oct-96 : pmn Bug dans SetPeriodic (PRO6088) et Segment (PRO6250)
  29
  30 #define No_Standard_OutOfRange
  31
  32
  33 #include <BSplCLib.hxx>
  34 #include <BSplCLib_KnotDistribution.hxx>
  35 #include <BSplCLib_MultDistribution.hxx>
  36 #include <ElCLib.hxx>
  37 #include <Geom_BSplineCurve.hxx>
  38 #include <Geom_Geometry.hxx>
  39 #include <Geom_UndefinedDerivative.hxx>
  40 #include <gp.hxx>
  41 #include <gp_Pnt.hxx>
  42 #include <gp_Trsf.hxx>
  43 #include <gp_Vec.hxx>
  44 #include <Standard_ConstructionError.hxx>
  45 #include <Standard_DimensionError.hxx>
  46 #include <Standard_DomainError.hxx>
  47 #include <Standard_NoSuchObject.hxx>
  48 #include <Standard_NotImplemented.hxx>
  49 #include <Standard_OutOfRange.hxx>
  50 #include <Standard_RangeError.hxx>
  51 #include <Standard_Real.hxx>
  52 #include <Standard_Type.hxx>
  53
  54 IMPLEMENT_STANDARD_RTTIEXT(Geom_BSplineCurve,Geom_BoundedCurve)
  55
  56 //=======================================================================
  57 //function : CheckCurveData
  58 //purpose  : Internal use only
  59 //=======================================================================
  60 static void CheckCurveData
  61 (const TColgp_Array1OfPnt&         CPoles,
  62  const TColStd_Array1OfReal&       CKnots,
  63  const TColStd_Array1OfInteger&    CMults,
  64  const Standard_Integer            Degree,
  65  const Standard_Boolean            Periodic)
  66 {
  67   if (Degree < 1 || Degree > Geom_BSplineCurve::MaxDegree()) {
  68     throw Standard_ConstructionError("BSpline curve: invalid degree");
  69   }
  70
  71   if (CPoles.Length() < 2)                throw Standard_ConstructionError("BSpline curve: at least 2 poles required");
  72   if (CKnots.Length() != CMults.Length()) throw Standard_ConstructionError("BSpline curve: Knot and Mult array size mismatch");
  73
  74   for (Standard_Integer I = CKnots.Lower(); I < CKnots.Upper(); I++) {
  75     if (CKnots (I+1) - CKnots (I) <= Epsilon (Abs(CKnots (I)))) {
  76       throw Standard_ConstructionError("BSpline curve: Knots interval values too close");
  77     }
  78   }
  79
  80   if (CPoles.Length() != BSplCLib::NbPoles(Degree,Periodic,CMults))
  81     throw Standard_ConstructionError("BSpline curve: # Poles and degree mismatch");
  82 }
  83
  84 //! Check rationality of an array of weights
  85 static Standard_Boolean Rational (const TColStd_Array1OfReal& theWeights)
  86 {
  87   for (Standard_Integer i = theWeights.Lower(); i < theWeights.Upper(); i++)
  88   {
  89     if (Abs (theWeights[i] - theWeights[i + 1]) > gp::Resolution())
  90     {
  91       return Standard_True;
  92     }
  93   }
  94   return Standard_False;
  95 }
  96
  97 //=======================================================================
  98 //function : Copy
  99 //purpose  :
 100 //=======================================================================
 101
 102 Handle(Geom_Geometry) Geom_BSplineCurve::Copy() const
 103 {
 104   Handle(Geom_BSplineCurve) C;
 105   if (IsRational())
 106     C = new Geom_BSplineCurve(poles->Array1(),
 107                               weights->Array1(),
 108                               knots->Array1(),
 109                               mults->Array1(),
 110                               deg,periodic);
 111   else
 112     C = new Geom_BSplineCurve(poles->Array1(),
 113                               knots->Array1(),
 114                               mults->Array1(),
 115                               deg,periodic);
 116   return C;
 117 }
 118
 119 //=======================================================================
 120 //function : Geom_BSplineCurve
 121 //purpose  :
 122 //=======================================================================
 123
 124 Geom_BSplineCurve::Geom_BSplineCurve
 125 (const TColgp_Array1OfPnt&       Poles,
 126  const TColStd_Array1OfReal&     Knots,
 127  const TColStd_Array1OfInteger&  Mults,
 128  const Standard_Integer          Degree,
 129  const Standard_Boolean          Periodic) :
 130  rational(Standard_False),
 131  periodic(Periodic),
 132  deg(Degree),
 133  maxderivinvok(Standard_False)
 134 {
 135   // check
 136
 137   CheckCurveData(Poles,
 138                  Knots,
 139                  Mults,
 140                  Degree,
 141                  Periodic);
 142
 143   // copy arrays
 144
 145   poles =  new TColgp_HArray1OfPnt(1,Poles.Length());
 146   poles->ChangeArray1() = Poles;
 147
 148
 149   knots = new TColStd_HArray1OfReal(1,Knots.Length());
 150   knots->ChangeArray1() = Knots;
 151
 152   mults = new TColStd_HArray1OfInteger(1,Mults.Length());
 153   mults->ChangeArray1() = Mults;
 154
 155   UpdateKnots();
 156 }
 157
 158 //=======================================================================
 159 //function : Geom_BSplineCurve
 160 //purpose  :
 161 //=======================================================================
 162
 163 Geom_BSplineCurve::Geom_BSplineCurve
 164 (const TColgp_Array1OfPnt&      Poles,
 165  const TColStd_Array1OfReal&    Weights,
 166  const TColStd_Array1OfReal&    Knots,
 167  const TColStd_Array1OfInteger& Mults,
 168  const Standard_Integer         Degree,
 169  const Standard_Boolean         Periodic,
 170  const Standard_Boolean         CheckRational)  :
 171  rational(Standard_True),
 172  periodic(Periodic),
 173  deg(Degree),
 174  maxderivinvok(Standard_False)
 175
 176 {
 177
 178   // check
 179
 180   CheckCurveData(Poles,
 181                  Knots,
 182                  Mults,
 183                  Degree,
 184                  Periodic);
 185
 186   if (Weights.Length() != Poles.Length())
 187     throw Standard_ConstructionError("Geom_BSplineCurve: Weights and Poles array size mismatch");
 188
 189   Standard_Integer i;
 190   for (i = Weights.Lower(); i <= Weights.Upper(); i++) {
 191     if (Weights(i) <= gp::Resolution())
 192       throw Standard_ConstructionError("Geom_BSplineCurve: Weights values too small");
 193   }
 194
 195   // check really rational
 196   if (CheckRational)
 197     rational = Rational(Weights);
 198
 199   // copy arrays
 200
 201   poles =  new TColgp_HArray1OfPnt(1,Poles.Length());
 202   poles->ChangeArray1() = Poles;
 203   if (rational) {
 204     weights =  new TColStd_HArray1OfReal(1,Weights.Length());
 205     weights->ChangeArray1() = Weights;
 206   }
 207
 208   knots = new TColStd_HArray1OfReal(1,Knots.Length());
 209   knots->ChangeArray1() = Knots;
 210
 211   mults = new TColStd_HArray1OfInteger(1,Mults.Length());
 212   mults->ChangeArray1() = Mults;
 213
 214   UpdateKnots();
 215 }
 216
 217 //=======================================================================
 218 //function : MaxDegree
 219 //purpose  :
 220 //=======================================================================
 221
 222 Standard_Integer Geom_BSplineCurve::MaxDegree ()
 223 {
 224   return BSplCLib::MaxDegree();
 225 }
 226
 227 //=======================================================================
 228 //function : IncreaseDegree
 229 //purpose  :
 230 //=======================================================================
 231
 232 void Geom_BSplineCurve::IncreaseDegree  (const Standard_Integer Degree)
 233 {
 234   if (Degree == deg) return;
 235
 236   if (Degree < deg || Degree > Geom_BSplineCurve::MaxDegree()) {
 237     throw Standard_ConstructionError("BSpline curve: IncreaseDegree: bad degree value");
 238   }
 239   Standard_Integer FromK1 = FirstUKnotIndex ();
 240   Standard_Integer ToK2   = LastUKnotIndex  ();
 241
 242   Standard_Integer Step   = Degree - deg;
 243
 244   Handle(TColgp_HArray1OfPnt) npoles = new
 245     TColgp_HArray1OfPnt(1,poles->Length() + Step * (ToK2-FromK1));
 246
 247   Standard_Integer nbknots = BSplCLib::IncreaseDegreeCountKnots
 248     (deg,Degree,periodic,mults->Array1());
 249
 250   Handle(TColStd_HArray1OfReal) nknots =
 251     new TColStd_HArray1OfReal(1,nbknots);
 252
 253   Handle(TColStd_HArray1OfInteger) nmults =
 254     new TColStd_HArray1OfInteger(1,nbknots);
 255
 256   Handle(TColStd_HArray1OfReal) nweights;
 257   if (IsRational())
 258   {
 259     nweights = new TColStd_HArray1OfReal(1,npoles->Upper());
 260   }
 261   BSplCLib::IncreaseDegree (deg, Degree, periodic,
 262                             poles->Array1(), !nweights.IsNull() ? &weights->Array1() : BSplCLib::NoWeights(),
 263                             knots->Array1(), mults->Array1(),
 264                             npoles->ChangeArray1(), !nweights.IsNull() ? &nweights->ChangeArray1() : BSplCLib::NoWeights(),
 265                             nknots->ChangeArray1(),nmults->ChangeArray1());
 266   deg     = Degree;
 267   poles   = npoles;
 268   weights = nweights;
 269   knots   = nknots;
 270   mults   = nmults;
 271   UpdateKnots();
 272 }
 273
 274 //=======================================================================
 275 //function : IncreaseMultiplicity
 276 //purpose  :
 277 //=======================================================================
 278
 279 void Geom_BSplineCurve::IncreaseMultiplicity  (const Standard_Integer Index,
 280                                                const Standard_Integer M)
 281 {
 282   TColStd_Array1OfReal k(1,1);
 283   k(1) = knots->Value(Index);
 284   TColStd_Array1OfInteger m(1,1);
 285   m(1) = M - mults->Value(Index);
 286   InsertKnots(k,m,Epsilon(1.),Standard_True);
 287 }
 288
 289 //=======================================================================
 290 //function : IncreaseMultiplicity
 291 //purpose  :
 292 //=======================================================================
 293
 294 void Geom_BSplineCurve::IncreaseMultiplicity  (const Standard_Integer I1,
 295                                                const Standard_Integer I2,
 296                                                const Standard_Integer M)
 297 {
 298   Handle(TColStd_HArray1OfReal)  tk = knots;
 299   TColStd_Array1OfReal k((knots->Array1())(I1),I1,I2);
 300   TColStd_Array1OfInteger m(I1,I2);
 301   Standard_Integer i;
 302   for (i = I1; i <= I2; i++)
 303     m(i) = M - mults->Value(i);
 304   InsertKnots(k,m,Epsilon(1.),Standard_True);
 305 }
 306
 307 //=======================================================================
 308 //function : IncrementMultiplicity
 309 //purpose  :
 310 //=======================================================================
 311
 312 void Geom_BSplineCurve::IncrementMultiplicity
 313 (const Standard_Integer I1,
 314  const Standard_Integer I2,
 315  const Standard_Integer Step)
 316 {
 317   Handle(TColStd_HArray1OfReal) tk = knots;
 318   TColStd_Array1OfReal    k((knots->Array1())(I1),I1,I2);
 319   TColStd_Array1OfInteger m(I1,I2) ;
 320   m.Init(Step);
 321   InsertKnots(k,m,Epsilon(1.),Standard_True);
 322 }
 323
 324 //=======================================================================
 325 //function : InsertKnot
 326 //purpose  :
 327 //=======================================================================
 328
 329 void Geom_BSplineCurve::InsertKnot
 330 (const Standard_Real U,
 331  const Standard_Integer M,
 332  const Standard_Real ParametricTolerance,
 333  const Standard_Boolean Add)
 334 {
 335   TColStd_Array1OfReal k(1,1);
 336   k(1) = U;
 337   TColStd_Array1OfInteger m(1,1);
 338   m(1) = M;
 339   InsertKnots(k,m,ParametricTolerance,Add);
 340 }
 341
 342 //=======================================================================
 343 //function : InsertKnots
 344 //purpose  :
 345 //=======================================================================
 346
 347 void  Geom_BSplineCurve::InsertKnots(const TColStd_Array1OfReal& Knots,
 348                                      const TColStd_Array1OfInteger& Mults,
 349                                      const Standard_Real Epsilon,
 350                                      const Standard_Boolean Add)
 351 {
 352   // Check and compute new sizes
 353   Standard_Integer nbpoles,nbknots;
 354
 355   if (!BSplCLib::PrepareInsertKnots(deg,periodic,
 356                                     knots->Array1(),mults->Array1(),
 357                                     Knots,&Mults,nbpoles,nbknots,Epsilon,Add))
 358     throw Standard_ConstructionError("Geom_BSplineCurve::InsertKnots");
 359
 360   if (nbpoles == poles->Length()) return;
 361
 362   Handle(TColgp_HArray1OfPnt)      npoles = new TColgp_HArray1OfPnt(1,nbpoles);
 363   Handle(TColStd_HArray1OfReal)    nknots = knots;
 364   Handle(TColStd_HArray1OfInteger) nmults = mults;
 365
 366   if (nbknots != knots->Length()) {
 367     nknots = new TColStd_HArray1OfReal(1,nbknots);
 368     nmults = new TColStd_HArray1OfInteger(1,nbknots);
 369   }
 370
 371   Handle(TColStd_HArray1OfReal) nweights;
 372   if (rational)
 373   {
 374     nweights = new TColStd_HArray1OfReal(1,nbpoles);
 375   }
 376
 377   BSplCLib::InsertKnots (deg, periodic,
 378                          poles->Array1(), !nweights.IsNull() ? &weights->Array1() : BSplCLib::NoWeights(),
 379                          knots->Array1(), mults->Array1(),
 380                          Knots, &Mults,
 381                          npoles->ChangeArray1(), !nweights.IsNull() ? &nweights->ChangeArray1() : BSplCLib::NoWeights(),
 382                          nknots->ChangeArray1(), nmults->ChangeArray1(),
 383                          Epsilon, Add);
 384   weights = nweights;
 385   poles = npoles;
 386   knots = nknots;
 387   mults = nmults;
 388   UpdateKnots();
 389 }
 390
 391 //=======================================================================
 392 //function : RemoveKnot
 393 //purpose  :
 394 //=======================================================================
 395
 396 Standard_Boolean  Geom_BSplineCurve::RemoveKnot(const Standard_Integer Index,
 397                                                 const Standard_Integer M,
 398                                                 const Standard_Real Tolerance)
 399 {
 400   if (M < 0) return Standard_True;
 401
 402   Standard_Integer I1  = FirstUKnotIndex ();
 403   Standard_Integer I2  = LastUKnotIndex  ();
 404
 405   if ( !periodic && (Index <= I1 || Index >= I2) ) {
 406     throw Standard_OutOfRange("BSpline curve: RemoveKnot: index out of range");
 407   }
 408   else if ( periodic  && (Index < I1 || Index > I2)) {
 409     throw Standard_OutOfRange("BSpline curve: RemoveKnot: index out of range");
 410   }
 411
 412   const TColgp_Array1OfPnt   & oldpoles   = poles->Array1();
 413
 414   Standard_Integer step = mults->Value(Index) - M;
 415   if (step <= 0) return Standard_True;
 416
 417   Handle(TColgp_HArray1OfPnt) npoles =
 418     new TColgp_HArray1OfPnt(1,oldpoles.Length()-step);
 419
 420   Handle(TColStd_HArray1OfReal)    nknots  = knots;
 421   Handle(TColStd_HArray1OfInteger) nmults  = mults;
 422
 423   if (M == 0) {
 424     nknots = new TColStd_HArray1OfReal(1,knots->Length()-1);
 425     nmults = new TColStd_HArray1OfInteger(1,knots->Length()-1);
 426   }
 427
 428   Handle(TColStd_HArray1OfReal) nweights;
 429   if (IsRational())
 430   {
 431     nweights = new TColStd_HArray1OfReal(1,npoles->Length());
 432   }
 433
 434   if (!BSplCLib::RemoveKnot (Index, M, deg, periodic,
 435                              poles->Array1(), !nweights.IsNull() ? &weights->Array1() : BSplCLib::NoWeights(),
 436                              knots->Array1(),mults->Array1(),
 437                              npoles->ChangeArray1(), !nweights.IsNull() ? &nweights->ChangeArray1() : BSplCLib::NoWeights(),
 438                              nknots->ChangeArray1(),nmults->ChangeArray1(),
 439                              Tolerance))
 440   {
 441     return Standard_False;
 442   }
 443
 444   weights = nweights;
 445   poles = npoles;
 446   knots = nknots;
 447   mults = nmults;
 448
 449   UpdateKnots();
 450   maxderivinvok = 0;
 451   return Standard_True;
 452 }
 453
 454 //=======================================================================
 455 //function : Reverse
 456 //purpose  :
 457 //=======================================================================
 458
 459 void Geom_BSplineCurve::Reverse ()
 460 {
 461   BSplCLib::Reverse(knots->ChangeArray1());
 462   BSplCLib::Reverse(mults->ChangeArray1());
 463   Standard_Integer last;
 464   if (periodic)
 465     last = flatknots->Upper() - deg - 1;
 466   else
 467     last = poles->Upper();
 468   BSplCLib::Reverse(poles->ChangeArray1(),last);
 469   if (rational)
 470     BSplCLib::Reverse(weights->ChangeArray1(),last);
 471   UpdateKnots();
 472 }
 473
 474 //=======================================================================
 475 //function : ReversedParameter
 476 //purpose  :
 477 //=======================================================================
 478
 479 Standard_Real Geom_BSplineCurve::ReversedParameter
 480 (const Standard_Real U) const
 481 {
 482   return (FirstParameter() + LastParameter() - U);
 483 }
 484
 485 //=======================================================================
 486 //function : Segment
 487 //purpose  :
 488 //=======================================================================
 489
 490 void Geom_BSplineCurve::Segment(const Standard_Real U1,
 491                                 const Standard_Real U2,
 492                                 const Standard_Real theTolerance)
 493 {
 494   if (U2 < U1)
 495     throw Standard_DomainError("Geom_BSplineCurve::Segment");
 496
 497   Standard_Real NewU1, NewU2;
 498   Standard_Real U,DU=0,aDDU=0;
 499   Standard_Integer index;
 500   Standard_Boolean wasPeriodic = periodic;
 501
 502   TColStd_Array1OfReal    Knots(1,2);
 503   TColStd_Array1OfInteger Mults(1,2);
 504
 505   // define param ditance to keep (eap, Apr 18 2002, occ311)
 506   if (periodic) {
 507     Standard_Real Period = LastParameter() - FirstParameter();
 508     DU = U2 - U1;
 509     if (DU - Period > Precision::PConfusion())
 510       throw Standard_DomainError("Geom_BSplineCurve::Segment");
 511     if (DU > Period)
 512       DU = Period;
 513     aDDU = DU;
 514   }
 515
 516   index = 0;
 517   BSplCLib::LocateParameter(deg,knots->Array1(),mults->Array1(),
 518                             U1,periodic,knots->Lower(),knots->Upper(),
 519                             index,NewU1);
 520   index = 0;
 521   BSplCLib::LocateParameter(deg,knots->Array1(),mults->Array1(),
 522                             U2,periodic,knots->Lower(),knots->Upper(),
 523                             index,NewU2);
 524
 525   //-- DBB
 526   Standard_Real aNu2 = NewU2;
 527   //-- DBB
 528
 529
 530   Knots( 1) = Min( NewU1, NewU2);
 531   Knots( 2) = Max( NewU1, NewU2);
 532   Mults( 1) = Mults( 2) = deg;
 533
 534   Standard_Real AbsUMax = Max(Abs(NewU1),Abs(NewU2));
 535
 536 //  Modified by Sergey KHROMOV - Fri Apr 11 12:15:40 2003 Begin
 537   AbsUMax = Max(AbsUMax, Max(Abs(FirstParameter()),Abs(LastParameter())));
 538 //  Modified by Sergey KHROMOV - Fri Apr 11 12:15:40 2003 End
 539
 540   Standard_Real Eps = Max(Epsilon(AbsUMax), theTolerance);
 541
 542   InsertKnots( Knots, Mults, Eps);
 543
 544   if (periodic) { // set the origine at NewU1
 545     index = 0;
 546     BSplCLib::LocateParameter(deg,knots->Array1(),mults->Array1(),
 547                               U1,periodic,knots->Lower(),knots->Upper(),
 548                               index,U);
 549     // Test si l'insertion est Ok et decalage sinon.
 550     if ( Abs(knots->Value(index+1)-U) <= Eps) // <= pour etre homogene a InsertKnots
 551       index++;
 552     SetOrigin(index);
 553     SetNotPeriodic();
 554     NewU2 = NewU1 + DU;
 555   }
 556
 557   // compute index1 and index2 to set the new knots and mults
 558   Standard_Integer index1 = 0, index2 = 0;
 559   Standard_Integer FromU1 = knots->Lower();
 560   Standard_Integer ToU2   = knots->Upper();
 561   BSplCLib::LocateParameter(deg,knots->Array1(),mults->Array1(),
 562                             NewU1,periodic,FromU1,ToU2,index1,U);
 563   if ( Abs(knots->Value(index1+1)-U) <= Eps)
 564     index1++;
 565
 566   BSplCLib::LocateParameter(deg,knots->Array1(),mults->Array1(),
 567                             NewU2,periodic,FromU1,ToU2,index2,U);
 568   if ( Abs(knots->Value(index2+1)-U) <= Eps || index2 == index1)
 569     index2++;
 570
 571   Standard_Integer nbknots = index2 - index1 + 1;
 572
 573   Handle(TColStd_HArray1OfReal)
 574     nknots = new TColStd_HArray1OfReal(1,nbknots);
 575   Handle(TColStd_HArray1OfInteger)
 576     nmults = new TColStd_HArray1OfInteger(1,nbknots);
 577
 578   // to restore changed U1
 579   if (DU > 0) // if was periodic
 580     DU = NewU1 - U1;
 581
 582   Standard_Integer i , k = 1;
 583   for ( i = index1; i<= index2; i++) {
 584     nknots->SetValue(k, knots->Value(i) - DU);
 585     nmults->SetValue(k, mults->Value(i));
 586     k++;
 587   }
 588   nmults->SetValue(      1, deg + 1);
 589   nmults->SetValue(nbknots, deg + 1);
 590
 591
 592   // compute index1 and index2 to set the new poles and weights
 593   Standard_Integer pindex1
 594     = BSplCLib::PoleIndex(deg,index1,periodic,mults->Array1());
 595   Standard_Integer pindex2
 596     = BSplCLib::PoleIndex(deg,index2,periodic,mults->Array1());
 597
 598   pindex1++;
 599   pindex2 = Min( pindex2+1, poles->Length());
 600
 601   Standard_Integer nbpoles  = pindex2 - pindex1 + 1;
 602
 603   Handle(TColStd_HArray1OfReal)
 604     nweights = new TColStd_HArray1OfReal(1,nbpoles);
 605   Handle(TColgp_HArray1OfPnt)
 606     npoles = new TColgp_HArray1OfPnt(1,nbpoles);
 607
 608   k = 1;
 609   if ( rational) {
 610     nweights = new TColStd_HArray1OfReal( 1, nbpoles);
 611     for ( i = pindex1; i <= pindex2; i++) {
 612       npoles->SetValue(k, poles->Value(i));
 613       nweights->SetValue(k, weights->Value(i));
 614       k++;
 615     }
 616   }
 617   else {
 618     for ( i = pindex1; i <= pindex2; i++) {
 619       npoles->SetValue(k, poles->Value(i));
 620       k++;
 621     }
 622   }
 623
 624   //-- DBB
 625   if( wasPeriodic ) {
 626     nknots->ChangeValue(nknots->Lower()) = U1;
 627     if( aNu2 < U2 ) {
 628       nknots->ChangeValue(nknots->Upper()) = U1 + aDDU;
 629     }
 630   }
 631   //-- DBB
 632
 633   knots = nknots;
 634   mults = nmults;
 635   poles = npoles;
 636   if ( rational)
 637     weights = nweights;
 638
 639   maxderivinvok = 0;
 640   UpdateKnots();
 641 }
 642
 643 //=======================================================================
 644 //function : SetKnot
 645 //purpose  :
 646 //=======================================================================
 647
 648 void Geom_BSplineCurve::SetKnot
 649 (const Standard_Integer Index,
 650  const Standard_Real K)
 651 {
 652   if (Index < 1 || Index > knots->Length())     throw Standard_OutOfRange("BSpline curve: SetKnot: Index and #knots mismatch");
 653   Standard_Real DK = Abs(Epsilon (K));
 654   if (Index == 1) {
 655     if (K >= knots->Value(2) - DK) throw Standard_ConstructionError("BSpline curve: SetKnot: K out of range");
 656   }
 657   else if (Index == knots->Length()) {
 658     if (K <= knots->Value (knots->Length()-1) + DK)  {
 659       throw Standard_ConstructionError("BSpline curve: SetKnot: K out of range");
 660     }
 661   }
 662   else {
 663     if (K <= knots->Value(Index-1) + DK ||
 664         K >= knots->Value(Index+1) - DK ) {
 665       throw Standard_ConstructionError("BSpline curve: SetKnot: K out of range");
 666     }
 667   }
 668   if (K != knots->Value (Index)) {
 669     knots->SetValue (Index, K);
 670     maxderivinvok = 0;
 671     UpdateKnots();
 672   }
 673 }
 674
 675 //=======================================================================
 676 //function : SetKnots
 677 //purpose  :
 678 //=======================================================================
 679
 680 void Geom_BSplineCurve::SetKnots
 681 (const TColStd_Array1OfReal& K)
 682 {
 683   CheckCurveData(poles->Array1(),K,mults->Array1(),deg,periodic);
 684   knots->ChangeArray1() = K;
 685   maxderivinvok = 0;
 686   UpdateKnots();
 687 }
 688
 689 //=======================================================================
 690 //function : SetKnot
 691 //purpose  :
 692 //=======================================================================
 693
 694 void Geom_BSplineCurve::SetKnot
 695 (const Standard_Integer Index,
 696  const Standard_Real K,
 697  const Standard_Integer M)
 698 {
 699   IncreaseMultiplicity (Index, M);
 700   SetKnot (Index, K);
 701 }
 702
 703 //=======================================================================
 704 //function : SetPeriodic
 705 //purpose  :
 706 //=======================================================================
 707
 708 void Geom_BSplineCurve::SetPeriodic ()
 709 {
 710   Standard_Integer first = FirstUKnotIndex();
 711   Standard_Integer last  = LastUKnotIndex();
 712
 713   Handle(TColStd_HArray1OfReal) tk = knots;
 714   TColStd_Array1OfReal    cknots((knots->Array1())(first),first,last);
 715   knots = new TColStd_HArray1OfReal(1,cknots.Length());
 716   knots->ChangeArray1() = cknots;
 717
 718   Handle(TColStd_HArray1OfInteger) tm = mults;
 719   TColStd_Array1OfInteger cmults((mults->Array1())(first),first,last);
 720   cmults(first) = cmults(last) = Min(deg, Max( cmults(first), cmults(last)));
 721   mults = new TColStd_HArray1OfInteger(1,cmults.Length());
 722   mults->ChangeArray1() = cmults;
 723
 724   // compute new number of poles;
 725   Standard_Integer nbp = BSplCLib::NbPoles(deg,Standard_True,cmults);
 726
 727   Handle(TColgp_HArray1OfPnt) tp = poles;
 728   TColgp_Array1OfPnt cpoles((poles->Array1())(1),1,nbp);
 729   poles = new TColgp_HArray1OfPnt(1,nbp);
 730   poles->ChangeArray1() = cpoles;
 731
 732   if (rational) {
 733     Handle(TColStd_HArray1OfReal) tw = weights;
 734     TColStd_Array1OfReal cweights((weights->Array1())(1),1,nbp);
 735     weights = new TColStd_HArray1OfReal(1,nbp);
 736     weights->ChangeArray1() = cweights;
 737   }
 738
 739   periodic = Standard_True;
 740
 741   maxderivinvok = 0;
 742   UpdateKnots();
 743 }
 744
 745 //=======================================================================
 746 //function : SetOrigin
 747 //purpose  :
 748 //=======================================================================
 749
 750 void Geom_BSplineCurve::SetOrigin(const Standard_Integer Index)
 751 {
 752   if (!periodic)
 753     throw Standard_NoSuchObject("Geom_BSplineCurve::SetOrigin");
 754
 755   Standard_Integer i,k;
 756   Standard_Integer first = FirstUKnotIndex();
 757   Standard_Integer last  = LastUKnotIndex();
 758
 759   if ((Index < first) || (Index > last))
 760     throw Standard_DomainError("Geom_BSplineCurve::SetOrigin");
 761
 762   Standard_Integer nbknots = knots->Length();
 763   Standard_Integer nbpoles = poles->Length();
 764
 765   Handle(TColStd_HArray1OfReal) nknots =
 766     new TColStd_HArray1OfReal(1,nbknots);
 767   TColStd_Array1OfReal& newknots = nknots->ChangeArray1();
 768
 769   Handle(TColStd_HArray1OfInteger) nmults =
 770     new TColStd_HArray1OfInteger(1,nbknots);
 771   TColStd_Array1OfInteger& newmults = nmults->ChangeArray1();
 772
 773   // set the knots and mults
 774   Standard_Real period = knots->Value(last) - knots->Value(first);
 775   k = 1;
 776   for ( i = Index; i <= last ; i++) {
 777     newknots(k) = knots->Value(i);
 778     newmults(k) = mults->Value(i);
 779     k++;
 780   }
 781   for ( i = first+1; i <= Index; i++) {
 782     newknots(k) = knots->Value(i) + period;
 783     newmults(k) = mults->Value(i);
 784     k++;
 785   }
 786
 787   Standard_Integer index = 1;
 788   for (i = first+1; i <= Index; i++)
 789     index += mults->Value(i);
 790
 791   // set the poles and weights
 792   Handle(TColgp_HArray1OfPnt) npoles =
 793     new TColgp_HArray1OfPnt(1,nbpoles);
 794   Handle(TColStd_HArray1OfReal) nweights =
 795     new TColStd_HArray1OfReal(1,nbpoles);
 796   TColgp_Array1OfPnt   & newpoles   = npoles->ChangeArray1();
 797   TColStd_Array1OfReal & newweights = nweights->ChangeArray1();
 798   first = poles->Lower();
 799   last  = poles->Upper();
 800   if ( rational) {
 801     k = 1;
 802     for ( i = index; i <= last; i++) {
 803       newpoles(k)   = poles->Value(i);
 804       newweights(k) = weights->Value(i);
 805       k++;
 806     }
 807     for ( i = first; i < index; i++) {
 808       newpoles(k)   = poles->Value(i);
 809       newweights(k) = weights->Value(i);
 810       k++;
 811     }
 812   }
 813   else {
 814     k = 1;
 815     for ( i = index; i <= last; i++) {
 816       newpoles(k) = poles->Value(i);
 817       k++;
 818     }
 819     for ( i = first; i < index; i++) {
 820       newpoles(k) = poles->Value(i);
 821       k++;
 822     }
 823   }
 824
 825   poles = npoles;
 826   knots = nknots;
 827   mults = nmults;
 828   if (rational)
 829     weights = nweights;
 830   maxderivinvok = 0;
 831   UpdateKnots();
 832 }
 833
 834 //=======================================================================
 835 //function : SetOrigin
 836 //purpose  :
 837 //=======================================================================
 838
 839 void Geom_BSplineCurve::SetOrigin(const Standard_Real U,
 840                                   const Standard_Real Tol)
 841 {
 842   if (!periodic)
 843     throw Standard_NoSuchObject("Geom_BSplineCurve::SetOrigin");
 844   //U est il dans la period.
 845   Standard_Real uf = FirstParameter(), ul = LastParameter();
 846   Standard_Real u = U, period = ul - uf;
 847   while (Tol < (uf-u)) u += period;
 848   while (Tol > (ul-u)) u -= period;
 849
 850   if(Abs(U-u)>Tol) { //On reparametre la courbe
 851     Standard_Real delta = U-u;
 852     uf += delta;
 853     ul += delta;
 854     TColStd_Array1OfReal& kn = knots->ChangeArray1();
 855     Standard_Integer fk = kn.Lower(), lk = kn.Upper();
 856     for(Standard_Integer i = fk; i <= lk; i++){
 857       kn.ChangeValue(i) += delta;
 858     }
 859     UpdateKnots();
 860   }
 861   if(Abs(U-uf)<Tol) return;
 862
 863   TColStd_Array1OfReal& kn = knots->ChangeArray1();
 864   Standard_Integer fk = kn.Lower(), lk = kn.Upper(),ik=0;
 865   Standard_Real delta = RealLast();
 866   for(Standard_Integer i = fk; i<= lk; i++){
 867     Standard_Real dki = kn.Value(i)-U;
 868     if(Abs(dki)<Abs(delta)){
 869       ik = i;
 870       delta = dki;
 871     }
 872   }
 873   if(Abs(delta)>Tol){
 874     InsertKnot(U);
 875     if(delta < 0.) ik++;
 876   }
 877   SetOrigin(ik);
 878 }
 879
 880 //=======================================================================
 881 //function : SetNotPeriodic
 882 //purpose  :
 883 //=======================================================================
 884
 885 void Geom_BSplineCurve::SetNotPeriodic ()
 886 {
 887   if ( periodic) {
 888     Standard_Integer NbKnots, NbPoles;
 889     BSplCLib::PrepareUnperiodize( deg, mults->Array1(),NbKnots,NbPoles);
 890
 891     Handle(TColgp_HArray1OfPnt) npoles
 892       = new TColgp_HArray1OfPnt(1,NbPoles);
 893
 894     Handle(TColStd_HArray1OfReal) nknots
 895       = new TColStd_HArray1OfReal(1,NbKnots);
 896
 897     Handle(TColStd_HArray1OfInteger) nmults
 898       = new TColStd_HArray1OfInteger(1,NbKnots);
 899
 900     Handle(TColStd_HArray1OfReal) nweights;
 901     if (IsRational())
 902     {
 903       nweights = new TColStd_HArray1OfReal(1,NbPoles);
 904     }
 905
 906     BSplCLib::Unperiodize (deg,
 907                            mults->Array1(), knots->Array1(), poles->Array1(),
 908                            !nweights.IsNull() ? &weights->Array1() : BSplCLib::NoWeights(),
 909                            nmults->ChangeArray1(), nknots->ChangeArray1(), npoles->ChangeArray1(),
 910                            !nweights.IsNull() ? &nweights->ChangeArray1() : BSplCLib::NoWeights());
 911     poles   = npoles;
 912     weights = nweights;
 913     mults   = nmults;
 914     knots   = nknots;
 915     periodic = Standard_False;
 916
 917     maxderivinvok = 0;
 918     UpdateKnots();
 919   }
 920 }
 921
 922 //=======================================================================
 923 //function : SetPole
 924 //purpose  :
 925 //=======================================================================
 926
 927 void Geom_BSplineCurve::SetPole
 928 (const Standard_Integer Index,
 929  const gp_Pnt& P)
 930 {
 931   if (Index < 1 || Index > poles->Length()) throw Standard_OutOfRange("BSpline curve: SetPole: index and #pole mismatch");
 932   poles->SetValue (Index, P);
 933   maxderivinvok = 0;
 934 }
 935
 936 //=======================================================================
 937 //function : SetPole
 938 //purpose  :
 939 //=======================================================================
 940
 941 void Geom_BSplineCurve::SetPole
 942 (const Standard_Integer Index,
 943  const gp_Pnt& P,
 944  const Standard_Real W)
 945 {
 946   SetPole(Index,P);
 947   SetWeight(Index,W);
 948 }
 949
 950 //=======================================================================
 951 //function : SetWeight
 952 //purpose  :
 953 //=======================================================================
 954
 955 void Geom_BSplineCurve::SetWeight
 956 (const Standard_Integer Index,
 957  const Standard_Real W)
 958 {
 959   if (Index < 1 || Index > poles->Length())   throw Standard_OutOfRange("BSpline curve: SetWeight: Index and #pole mismatch");
 960
 961   if (W <= gp::Resolution ())     throw Standard_ConstructionError("BSpline curve: SetWeight: Weight too small");
 962
 963
 964   Standard_Boolean rat = IsRational() || (Abs(W - 1.) > gp::Resolution());
 965
 966   if ( rat) {
 967     if (rat && !IsRational()) {
 968       weights = new TColStd_HArray1OfReal(1,poles->Length());
 969       weights->Init(1.);
 970     }
 971
 972     weights->SetValue (Index, W);
 973
 974     if (IsRational()) {
 975       rat = Rational(weights->Array1());
 976       if (!rat) weights.Nullify();
 977     }
 978
 979     rational = !weights.IsNull();
 980   }
 981   maxderivinvok = 0;
 982 }
 983
 984 //=======================================================================
 985 //function : MovePoint
 986 //purpose  :
 987 //=======================================================================
 988
 989 void Geom_BSplineCurve::MovePoint(const Standard_Real U,
 990                                   const gp_Pnt& P,
 991                                   const Standard_Integer Index1,
 992                                   const Standard_Integer Index2,
 993                                   Standard_Integer& FirstModifiedPole,
 994                                   Standard_Integer& LastmodifiedPole)
 995 {
 996   if (Index1 < 1 || Index1 > poles->Length() ||
 997       Index2 < 1 || Index2 > poles->Length() || Index1 > Index2) {
 998     throw Standard_OutOfRange("BSpline curve: MovePoint: Index and #pole mismatch");
 999   }
1000   TColgp_Array1OfPnt npoles(1, poles->Length());
1001   gp_Pnt P0;
1002   D0(U, P0);
1003   gp_Vec Displ(P0, P);
1004   BSplCLib::MovePoint (U, Displ, Index1, Index2, deg, poles->Array1(),
1005                        rational ? &weights->Array1() : BSplCLib::NoWeights(), flatknots->Array1(),
1006                        FirstModifiedPole, LastmodifiedPole, npoles);
1007   if (FirstModifiedPole) {
1008     poles->ChangeArray1() = npoles;
1009     maxderivinvok = 0;
1010   }
1011 }
1012
1013 //=======================================================================
1014 //function : MovePointAndTangent
1015 //purpose  :
1016 //=======================================================================
1017
1018 void Geom_BSplineCurve::MovePointAndTangent(const Standard_Real    U,
1019                                             const gp_Pnt&          P,
1020                                             const gp_Vec&          Tangent,
1021                                             const Standard_Real    Tolerance,
1022                                             const Standard_Integer StartingCondition,
1023                                             const Standard_Integer EndingCondition,
1024                                             Standard_Integer&      ErrorStatus)
1025 {
1026   Standard_Integer ii ;
1027   if (IsPeriodic()) {
1028     //
1029     // for the time being do not deal with periodic curves
1030     //
1031     SetNotPeriodic() ;
1032   }
1033   TColgp_Array1OfPnt new_poles(1, poles->Length());
1034   gp_Pnt P0;
1035
1036
1037   gp_Vec delta_derivative;
1038   D1(U, P0,
1039      delta_derivative) ;
1040   gp_Vec delta(P0, P);
1041   for (ii = 1 ; ii <= 3 ; ii++) {
1042     delta_derivative.SetCoord(ii, Tangent.Coord(ii)-delta_derivative.Coord(ii));
1043   }
1044   BSplCLib::MovePointAndTangent(U,
1045                                 delta,
1046                                 delta_derivative,
1047                                 Tolerance,
1048                                 deg,
1049                                 StartingCondition,
1050                                 EndingCondition,
1051                                 poles->Array1(),
1052                                 rational ? &weights->Array1() : BSplCLib::NoWeights(),
1053                                 flatknots->Array1(),
1054                                 new_poles,
1055                                 ErrorStatus) ;
1056   if (!ErrorStatus) {
1057     poles->ChangeArray1() = new_poles;
1058     maxderivinvok = 0;
1059   }
1060 }
1061
1062 //=======================================================================
1063 //function : UpdateKnots
1064 //purpose  :
1065 //=======================================================================
1066
1067 void Geom_BSplineCurve::UpdateKnots()
1068 {
1069   rational = !weights.IsNull();
1070
1071   Standard_Integer MaxKnotMult = 0;
1072   BSplCLib::KnotAnalysis(deg,
1073                          periodic,
1074                          knots->Array1(),
1075                          mults->Array1(),
1076                          knotSet, MaxKnotMult);
1077
1078   if (knotSet == GeomAbs_Uniform && !periodic)  {
1079     flatknots = knots;
1080   }
1081   else {
1082     flatknots = new TColStd_HArray1OfReal
1083       (1, BSplCLib::KnotSequenceLength(mults->Array1(),deg,periodic));
1084
1085     BSplCLib::KnotSequence(knots->Array1(),
1086                            mults->Array1(),
1087                            deg,periodic,
1088                            flatknots->ChangeArray1());
1089   }
1090
1091   if (MaxKnotMult == 0)  smooth = GeomAbs_CN;
1092   else {
1093     switch (deg - MaxKnotMult) {
1094     case 0 :   smooth = GeomAbs_C0;   break;
1095     case 1 :   smooth = GeomAbs_C1;   break;
1096     case 2 :   smooth = GeomAbs_C2;   break;
1097     case 3 :   smooth = GeomAbs_C3;   break;
1098       default :  smooth = GeomAbs_C3;   break;
1099     }
1100   }
1101 }
1102
1103 //=======================================================================
1104 //function : Normalizes the parameters if the curve is periodic
1105 //purpose  : that is compute the cache so that it is valid
1106 //=======================================================================
1107
1108 void Geom_BSplineCurve::PeriodicNormalization(Standard_Real&  Parameter) const
1109 {
1110   Standard_Real Period ;
1111
1112   if (periodic) {
1113     Period = flatknots->Value(flatknots->Upper() - deg) - flatknots->Value (deg + 1) ;
1114     while (Parameter > flatknots->Value(flatknots->Upper()-deg)) {
1115       Parameter -= Period ;
1116     }
1117     while (Parameter < flatknots->Value((deg + 1))) {
1118       Parameter +=  Period ;
1119     }
1120   }
1121 }
1122