src/Geom/Geom_BezierSurface.cxx

   1 // File:        Geom_BezierSurface.cxx
   2 // Created:     Tue Mar  9 19:52:37 1993
   3 // Author:      JCV
   4 //              <fid@phylox>
   5 // Copyright:   Matra Datavision 1993
   6
   7 //File Geom_BezierSurface.cxx  Septembre 1990
   8 //Passage en classe persistante - 23/01/91
   9 //Modif suite a la deuxieme revue de projet toolkit Geometry -23/01/91
  10 // pmn : 21/10/95 ; Correction de la methode segment (PRO5853)
  11 // pmn : 31-Dec-96; Bonne gestion des poids (bug PRO4622)
  12 // xab : 07-Jul-97; le cache est instable en degree 21
  13 //       a partir du degree 15 on ne l'utilise plus
  14 // RBD : 15/10/98 ; Le cache est desormais defini sur [-1,1] (pro15537).
  15 // pmn : 10/12/98 ; Update de la methode segment (suite a la modif de cache).
  16
  17 #define No_Standard_OutOfRange
  18 #define No_Standard_DimensionError
  19
  20 #include <Geom_BezierSurface.ixx>
  21 #include <Geom_BezierCurve.hxx>
  22
  23 #include <gp.hxx>
  24 #include <gp_XYZ.hxx>
  25 #include <PLib.hxx>
  26 #include <BSplCLib.hxx>
  27 #include <BSplSLib.hxx>
  28 #include <Precision.hxx>
  29
  30 #include <Standard_OutOfRange.hxx>
  31 #include <Standard_RangeError.hxx>
  32 #include <Standard_DimensionError.hxx>
  33 #include <Standard_ConstructionError.hxx>
  34 #include <TColStd_Array1OfInteger.hxx>
  35
  36 //=======================================================================
  37 //function : Rational
  38 //purpose  : check rationality of an array of weights
  39 //=======================================================================
  40
  41 static void Rational(const TColStd_Array2OfReal& Weights,
  42                      Standard_Boolean& Urational,
  43                      Standard_Boolean& Vrational)
  44 {
  45   Standard_Integer I,J;
  46   J = Weights.LowerCol ();
  47   Vrational = Standard_False;
  48   while (!Vrational && J <= Weights.UpperCol()) {
  49     I = Weights.LowerRow();
  50     while (!Vrational && I <= Weights.UpperRow() - 1) {
  51       Vrational = (Abs(Weights (I, J) - Weights (I+1, J))
  52                    > Epsilon (Abs(Weights (I, J))));
  53       I++;
  54     }
  55     J++;
  56   }
  57
  58   I = Weights.LowerRow ();
  59   Urational = Standard_False;
  60   while (!Urational && I <= Weights.UpperRow()) {
  61     J = Weights.LowerCol();
  62     while (!Urational && J <= Weights.UpperCol() - 1) {
  63       Urational = (Abs(Weights (I, J) - Weights (I, J+1))
  64                    > Epsilon (Abs(Weights (I, J))));
  65       J++;
  66     }
  67     I++;
  68   }
  69 }
  70
  71 //=======================================================================
  72 //function : AddPoleCol
  73 //purpose  : Internal use only.
  74 //=======================================================================
  75
  76 static void AddPoleCol
  77   (const TColgp_Array2OfPnt& Poles,
  78    const TColgp_Array1OfPnt& PoleCol,
  79    const Standard_Integer    AfterIndex,
  80          TColgp_Array2OfPnt& NewPoles)
  81 {
  82   Standard_Integer InsertIndex = AfterIndex + NewPoles.LowerCol();
  83   Standard_Integer Offset = NewPoles.LowerRow() - PoleCol.Lower();
  84   Standard_Integer ColIndex = NewPoles.LowerCol();
  85   Standard_Integer RowIndex;
  86   while (ColIndex < InsertIndex) {
  87     RowIndex = NewPoles.LowerRow();
  88     while (RowIndex <= NewPoles.UpperRow()){
  89       NewPoles (RowIndex, ColIndex) =  Poles (RowIndex, ColIndex);
  90       RowIndex++;
  91     }
  92     ColIndex++;
  93   }
  94   RowIndex = NewPoles.LowerRow();
  95   while (RowIndex <= NewPoles.UpperRow()){
  96     NewPoles (RowIndex, ColIndex) = PoleCol (RowIndex - Offset);
  97     RowIndex++;
  98   }
  99   ColIndex++;
 100   while (ColIndex <= NewPoles.UpperCol()) {
 101     RowIndex = NewPoles.LowerRow();
 102     while (RowIndex <= NewPoles.UpperRow()){
 103       NewPoles (RowIndex, ColIndex) =  Poles (RowIndex, ColIndex - 1);
 104       RowIndex++;
 105     }
 106     ColIndex++;
 107   }
 108 }
 109
 110 //=======================================================================
 111 //function : AddRatPoleCol
 112 //purpose  : Internal use only.
 113 //=======================================================================
 114
 115 static void AddRatPoleCol
 116   (const TColgp_Array2OfPnt&   Poles,
 117    const TColStd_Array2OfReal& Weights,
 118    const TColgp_Array1OfPnt&   PoleCol,
 119    const TColStd_Array1OfReal& PoleWeightCol,
 120    const Standard_Integer      AfterIndex,
 121          TColgp_Array2OfPnt&   NewPoles,
 122          TColStd_Array2OfReal& NewWeights)
 123 {
 124   Standard_Integer InsertIndex = AfterIndex + NewPoles.LowerCol();
 125   Standard_Integer OffsetPol   = NewPoles.LowerRow() - PoleCol.Lower();
 126   Standard_Integer OffsetW     = NewWeights.LowerRow() - PoleWeightCol.Lower();
 127
 128   Standard_Integer ColIndex = NewPoles.LowerCol();
 129   Standard_Integer RowIndex;
 130   while (ColIndex < InsertIndex) {
 131     RowIndex = NewPoles.LowerRow();
 132     while (RowIndex <= NewPoles.UpperRow()){
 133       NewPoles (RowIndex, ColIndex) = Poles (RowIndex, ColIndex);
 134       NewWeights (RowIndex, ColIndex)  =  Weights (RowIndex, ColIndex);
 135       RowIndex++;
 136     }
 137     ColIndex++;
 138   }
 139   RowIndex = NewPoles.LowerRow();
 140   while (RowIndex <= NewPoles.UpperRow()){
 141     NewPoles (RowIndex, ColIndex) = PoleCol (RowIndex - OffsetPol);
 142     NewWeights (RowIndex, ColIndex) =  PoleWeightCol (RowIndex - OffsetW);
 143     RowIndex++;
 144   }
 145   ColIndex++;
 146   while (ColIndex <= NewPoles.UpperCol()) {
 147     RowIndex = NewPoles.LowerRow();
 148     while (RowIndex <= NewPoles.UpperRow()){
 149       NewPoles (RowIndex, ColIndex) = Poles (RowIndex, ColIndex - 1);
 150       RowIndex++;
 151       NewWeights (RowIndex, ColIndex)  = Weights (RowIndex, ColIndex - 1);
 152     }
 153     ColIndex++;
 154   }
 155 }
 156
 157 //=======================================================================
 158 //function : AddPoleRow
 159 //purpose  : Internal use only.
 160 //=======================================================================
 161
 162 static void AddPoleRow
 163   (const TColgp_Array2OfPnt& Poles,
 164    const TColgp_Array1OfPnt& PoleRow,
 165    const Standard_Integer    AfterIndex,
 166          TColgp_Array2OfPnt& NewPoles)
 167 {
 168   Standard_Integer InsertIndex = AfterIndex + NewPoles.LowerRow();
 169   Standard_Integer Offset = NewPoles.LowerCol() - PoleRow.Lower();
 170   Standard_Integer RowIndex = NewPoles.LowerRow();
 171   Standard_Integer ColIndex;
 172   while (RowIndex < InsertIndex) {
 173     ColIndex = NewPoles.LowerCol();
 174     while (ColIndex <= NewPoles.UpperCol()){
 175       NewPoles (RowIndex, ColIndex) = Poles (RowIndex, ColIndex);
 176       ColIndex++;
 177     }
 178     RowIndex++;
 179   }
 180   ColIndex = NewPoles.LowerCol();
 181   while (ColIndex <= NewPoles.UpperCol()){
 182     NewPoles (RowIndex, ColIndex) = PoleRow (ColIndex - Offset);
 183     ColIndex++;
 184   }
 185   RowIndex++;
 186   while (RowIndex <= NewPoles.UpperRow()) {
 187     ColIndex = NewPoles.LowerCol();
 188     while (ColIndex <= NewPoles.UpperCol()){
 189       NewPoles (RowIndex, ColIndex) = Poles (RowIndex - 1, ColIndex);
 190       ColIndex++;
 191     }
 192     RowIndex++;
 193   }
 194 }
 195
 196 //=======================================================================
 197 //function : AddRatPoleRow
 198 //purpose  :
 199 //=======================================================================
 200
 201 static void AddRatPoleRow
 202   (const TColgp_Array2OfPnt&   Poles,
 203    const TColStd_Array2OfReal& Weights,
 204    const TColgp_Array1OfPnt&   PoleRow,
 205    const TColStd_Array1OfReal& PoleWeightRow,
 206    const Standard_Integer      AfterIndex,
 207          TColgp_Array2OfPnt&   NewPoles,
 208          TColStd_Array2OfReal& NewWeights)
 209 {
 210   Standard_Integer InsertIndex = AfterIndex + NewPoles.LowerRow();
 211   Standard_Integer OffsetPol = NewPoles.LowerCol() - PoleRow.Lower();
 212   Standard_Integer OffsetW = NewWeights.LowerCol() - PoleWeightRow.Lower();
 213
 214   Standard_Integer ColIndex;
 215   Standard_Integer RowIndex = NewPoles.LowerRow();
 216   while (RowIndex < InsertIndex) {
 217     ColIndex = NewPoles.LowerCol();
 218     while (ColIndex <= NewPoles.UpperCol()){
 219       NewPoles (RowIndex, ColIndex) =  Poles (RowIndex, ColIndex);
 220       NewWeights (RowIndex, ColIndex)  =  Weights (RowIndex, ColIndex);
 221       ColIndex++;
 222     }
 223     RowIndex++;
 224   }
 225   ColIndex = NewPoles.LowerCol();
 226   while (ColIndex <= NewPoles.UpperCol()){
 227     NewPoles (RowIndex, ColIndex) = PoleRow (ColIndex - OffsetPol);
 228     NewWeights (RowIndex, ColIndex) = PoleWeightRow (ColIndex - OffsetW);
 229     ColIndex++;
 230   }
 231   RowIndex++;
 232   while (RowIndex <= NewPoles.UpperRow()) {
 233     ColIndex = NewPoles.LowerCol();
 234     while (ColIndex <= NewPoles.UpperCol()){
 235       NewPoles (RowIndex, ColIndex) = Poles (RowIndex - 1, ColIndex);
 236       NewWeights (RowIndex, ColIndex)  = Weights (RowIndex - 1, ColIndex);
 237       ColIndex++;
 238     }
 239     RowIndex++;
 240   }
 241 }
 242
 243 //=======================================================================
 244 //function : DeletePoleCol
 245 //purpose  :
 246 //=======================================================================
 247
 248 static void DeletePoleCol
 249   (const TColgp_Array2OfPnt& Poles,
 250    const Standard_Integer    Index,
 251          TColgp_Array2OfPnt& NewPoles)
 252 {
 253   Standard_Integer Offset = 0;
 254   Standard_Integer RowIndex;
 255   Standard_Integer ColIndex = NewPoles.LowerCol();
 256   while (ColIndex <= NewPoles.UpperCol()) {
 257     RowIndex = NewPoles.LowerRow();
 258     if (ColIndex == Index)  Offset = 1;
 259     while (RowIndex <= NewPoles.UpperRow()){
 260       NewPoles (RowIndex, ColIndex) = Poles (RowIndex, ColIndex + Offset);
 261       RowIndex++;
 262     }
 263     ColIndex++;
 264   }
 265 }
 266
 267 //=======================================================================
 268 //function : DeleteRatPoleCol
 269 //purpose  :
 270 //=======================================================================
 271
 272 static void DeleteRatPoleCol
 273   (const TColgp_Array2OfPnt&   Poles,
 274    const TColStd_Array2OfReal& Weights,
 275    const Standard_Integer      Index,
 276          TColgp_Array2OfPnt&   NewPoles,
 277          TColStd_Array2OfReal& NewWeights)
 278 {
 279   Standard_Integer Offset = 0;
 280   Standard_Integer RowIndex;
 281   Standard_Integer ColIndex = NewPoles.LowerCol();
 282   while (ColIndex <= NewPoles.UpperCol()) {
 283     RowIndex = NewPoles.LowerRow();
 284     if (ColIndex == Index) Offset = 1;
 285     while (RowIndex <= NewPoles.UpperRow()){
 286       NewPoles (RowIndex, ColIndex) = Poles (RowIndex, ColIndex + Offset);
 287       NewWeights (RowIndex, ColIndex) = Weights (RowIndex, ColIndex+Offset);
 288       RowIndex++;
 289     }
 290     ColIndex++;
 291   }
 292 }
 293
 294 //=======================================================================
 295 //function : DeletePoleRow
 296 //purpose  :
 297 //=======================================================================
 298
 299 static void DeletePoleRow
 300   (const TColgp_Array2OfPnt& Poles,
 301    const Standard_Integer    Index,
 302          TColgp_Array2OfPnt& NewPoles)
 303 {
 304   Standard_Integer Offset = 0;
 305   Standard_Integer ColIndex;
 306   Standard_Integer RowIndex = NewPoles.LowerRow();
 307   while (RowIndex <= NewPoles.UpperRow()) {
 308     ColIndex = NewPoles.LowerCol();
 309     if (RowIndex == Index)  Offset = 1;
 310     while (ColIndex <= NewPoles.UpperCol()){
 311       NewPoles (RowIndex, ColIndex) = Poles (RowIndex + Offset, ColIndex);
 312       ColIndex++;
 313     }
 314     RowIndex++;
 315   }
 316 }
 317
 318 //=======================================================================
 319 //function : DeleteRatPoleRow
 320 //purpose  :
 321 //=======================================================================
 322
 323 static void DeleteRatPoleRow
 324   (const TColgp_Array2OfPnt&   Poles,
 325    const TColStd_Array2OfReal& Weights,
 326    const Standard_Integer      Index,
 327          TColgp_Array2OfPnt&   NewPoles,
 328          TColStd_Array2OfReal& NewWeights)
 329 {
 330   Standard_Integer Offset = 0;
 331   Standard_Integer ColIndex;
 332   Standard_Integer RowIndex = NewPoles.LowerRow();
 333   while (RowIndex <= NewPoles.UpperRow()) {
 334     ColIndex = NewPoles.LowerCol();
 335     if (RowIndex == Index)  Offset = 1;
 336     while (ColIndex <= NewPoles.UpperCol()){
 337       NewPoles (RowIndex, ColIndex) =  Poles (RowIndex +  Offset, ColIndex);
 338       NewWeights (RowIndex, ColIndex) = Weights (RowIndex+Offset, ColIndex);
 339       ColIndex++;
 340     }
 341     RowIndex++;
 342   }
 343 }
 344
 345 //=======================================================================
 346 //function : Geom_BezierSurface
 347 //purpose  :
 348 //=======================================================================
 349
 350 Geom_BezierSurface::Geom_BezierSurface
 351 (const TColgp_Array2OfPnt& SurfacePoles):
 352  ucacheparameter(0.),
 353  vcacheparameter(0.),
 354  ucachespanlenght(1.),
 355  vcachespanlenght(1.),
 356  validcache(0),
 357  maxderivinvok(Standard_False)
 358 {
 359   Standard_Integer NbUPoles = SurfacePoles.ColLength();
 360   Standard_Integer NbVPoles = SurfacePoles.RowLength();
 361   if (NbUPoles < 2 || NbUPoles > MaxDegree()+1 ||
 362       NbVPoles < 2 || NbVPoles > MaxDegree()+1) {
 363     Standard_ConstructionError::Raise();
 364   }
 365
 366   Handle(TColgp_HArray2OfPnt) npoles =
 367     new TColgp_HArray2OfPnt   (1, NbUPoles, 1, NbVPoles);
 368
 369   urational = 0;
 370   vrational = 0;
 371
 372   npoles->ChangeArray2() = SurfacePoles;
 373
 374   // Init non rational
 375   Init(npoles,
 376        Handle(TColStd_HArray2OfReal)());
 377 }
 378
 379 //=======================================================================
 380 //function : Geom_BezierSurface
 381 //purpose  :
 382 //=======================================================================
 383
 384 Geom_BezierSurface::Geom_BezierSurface
 385   (const TColgp_Array2OfPnt&   SurfacePoles,
 386    const TColStd_Array2OfReal& PoleWeights  ):
 387  ucacheparameter(0.),
 388  vcacheparameter(0.),
 389  ucachespanlenght(1.),
 390  vcachespanlenght(1.),
 391  validcache(0),
 392  maxderivinvok(Standard_False)
 393 {
 394   Standard_Integer NbUPoles = SurfacePoles.ColLength();
 395   Standard_Integer NbVPoles = SurfacePoles.RowLength();
 396   if (NbUPoles < 2 || NbUPoles > MaxDegree()+1 ||
 397       NbVPoles < 2 || NbVPoles > MaxDegree()+1 ||
 398       NbVPoles != PoleWeights.RowLength()      ||
 399       NbUPoles != PoleWeights.ColLength()      ) {
 400     Standard_ConstructionError::Raise();
 401   }
 402
 403   Standard_Integer Row = PoleWeights.LowerRow();
 404   Standard_Integer Col = PoleWeights.LowerCol();
 405   while (Col <= PoleWeights.UpperCol()) {
 406     Row = PoleWeights.LowerRow();
 407     while (Row <= PoleWeights.UpperRow()) {
 408       if (PoleWeights(Row, Col) <= gp::Resolution()) {
 409         Standard_ConstructionError::Raise();
 410       }
 411       Row++;
 412     }
 413     Col++;
 414   }
 415
 416   Handle(TColgp_HArray2OfPnt)
 417     npoles = new TColgp_HArray2OfPnt   (1, NbUPoles, 1, NbVPoles);
 418   npoles->ChangeArray2() = SurfacePoles;
 419
 420   Standard_Integer I, J;
 421   urational = Standard_False;
 422   vrational = Standard_False;
 423   J = PoleWeights.LowerCol ();
 424   while (!vrational && J <= PoleWeights.UpperCol()) {
 425     I = PoleWeights.LowerRow();
 426     while (!vrational && I <= PoleWeights.UpperRow() - 1) {
 427       vrational = (Abs(PoleWeights (I, J) - PoleWeights (I+1, J))
 428                    > Epsilon (Abs(PoleWeights (I, J))));
 429       I++;
 430     }
 431     J++;
 432   }
 433   I = PoleWeights.LowerRow ();
 434   while (!urational && I <= PoleWeights.UpperRow()) {
 435     J = PoleWeights.LowerCol();
 436     while (!urational && J <= PoleWeights.UpperCol() - 1) {
 437       urational = (Abs(PoleWeights (I, J) - PoleWeights (I, J+1))
 438                    > Epsilon (Abs(PoleWeights (I, J))));
 439       J++;
 440     }
 441     I++;
 442   }
 443
 444
 445   Handle(TColStd_HArray2OfReal) nweights;
 446   if (urational || vrational) {
 447     nweights = new TColStd_HArray2OfReal (1, NbUPoles, 1, NbVPoles);
 448     nweights->ChangeArray2() = PoleWeights;
 449   }
 450
 451   // Init
 452   Init(npoles,nweights);
 453 }
 454
 455 //=======================================================================
 456 //function : Geom_BezierSurface
 457 //purpose  :
 458 //=======================================================================
 459
 460 Geom_BezierSurface::Geom_BezierSurface
 461   (const Handle(TColgp_HArray2OfPnt)&   SurfacePoles,
 462    const Handle(TColgp_HArray2OfPnt)&   SurfaceCoefs,
 463    const Handle(TColStd_HArray2OfReal)& PoleWeights,
 464    const Handle(TColStd_HArray2OfReal)& CoefWeights,
 465    const Standard_Boolean               IsURational,
 466    const Standard_Boolean               IsVRational)
 467 :maxderivinvok(Standard_False)
 468 {
 469   urational = IsURational;
 470   vrational = IsVRational;
 471   ucachespanlenght = 1.;
 472   vcachespanlenght = 1.;
 473   validcache = 1;
 474   ucacheparameter =0.;
 475   vcacheparameter = 0.;
 476   Standard_Integer NbUPoles = SurfacePoles->ColLength();
 477   Standard_Integer NbVPoles = SurfacePoles->RowLength();
 478
 479   poles    = new TColgp_HArray2OfPnt   (1,NbUPoles,
 480                                         1,NbVPoles) ;
 481   poles->ChangeArray2() = SurfacePoles->Array2();
 482
 483   coeffs   = new TColgp_HArray2OfPnt   (1,SurfaceCoefs->ColLength(),
 484                                         1,SurfaceCoefs->RowLength()) ;
 485   coeffs->ChangeArray2() = SurfaceCoefs->Array2();
 486
 487   if ( urational || vrational) {
 488     weights  = new TColStd_HArray2OfReal (1,NbUPoles,1,NbVPoles);
 489     weights->ChangeArray2() = PoleWeights->Array2();
 490
 491     wcoeffs  = new TColStd_HArray2OfReal (1,SurfaceCoefs->ColLength(),
 492                                           1,SurfaceCoefs->RowLength()) ;
 493     wcoeffs->ChangeArray2() = CoefWeights->Array2();
 494   }
 495 }
 496
 497 //=======================================================================
 498 //function : MaxDegree
 499 //purpose  :
 500 //=======================================================================
 501
 502 Standard_Integer Geom_BezierSurface::MaxDegree ()
 503 {
 504   return BSplCLib::MaxDegree();
 505 }
 506
 507 //=======================================================================
 508 //function : ExchangeUV
 509 //purpose  :
 510 //=======================================================================
 511
 512 void Geom_BezierSurface::ExchangeUV ()
 513 {
 514   Standard_Integer LC = poles->LowerCol();
 515   Standard_Integer UC = poles->UpperCol();
 516   Standard_Integer LR = poles->LowerRow();
 517   Standard_Integer UR = poles->UpperRow();
 518
 519   Handle(TColgp_HArray2OfPnt) npoles =
 520     new TColgp_HArray2OfPnt (LC, UC, LR, UR);
 521   Handle(TColStd_HArray2OfReal) nweights =
 522     new TColStd_HArray2OfReal (LC, UC, LR, UR);
 523
 524   const TColgp_Array2OfPnt   & spoles   = poles->Array2();
 525   const TColStd_Array2OfReal & sweights = weights->Array2();
 526
 527   TColgp_Array2OfPnt&   snpoles   = npoles->ChangeArray2();
 528   TColStd_Array2OfReal& snweights = nweights->ChangeArray2();
 529
 530   Standard_Integer i, j;
 531   for (i = LC; i <= UC; i++) {
 532     for (j = LR; j <= UR; j++) {
 533       snpoles   (i,j) = spoles   (j,i);
 534       snweights (i,j) = sweights (j,i);
 535     }
 536   }
 537
 538   poles   = npoles;
 539   weights = nweights;
 540
 541   Standard_Boolean temp = urational;
 542   urational = vrational;
 543   vrational = temp;
 544   coeffs  = new TColgp_HArray2OfPnt   (LC, UC, LR, UR);
 545   wcoeffs = new TColStd_HArray2OfReal (LC, UC, LR, UR);
 546
 547   UpdateCoefficients();
 548 }
 549
 550 //=======================================================================
 551 //function : Increase
 552 //purpose  :
 553 //=======================================================================
 554
 555 void Geom_BezierSurface::Increase (const Standard_Integer UDeg,
 556                                    const Standard_Integer VDeg)
 557 {
 558   if (UDeg < UDegree() || UDeg > Geom_BezierSurface::MaxDegree() ||
 559       VDeg < VDegree() || VDeg > Geom_BezierSurface::MaxDegree() )  {
 560     Standard_ConstructionError::Raise();
 561   }
 562
 563   Standard_Integer oldUDeg = UDegree();
 564   Standard_Integer oldVDeg = VDegree();
 565   Standard_Integer IncUDeg = UDeg - oldUDeg;
 566   Standard_Integer IncVDeg = VDeg - oldVDeg;
 567   if (IncUDeg == 0 && IncVDeg == 0) return;
 568
 569   TColStd_Array1OfReal biduknots(1,2); biduknots(1) = 0.; biduknots(2) = 1.;
 570   TColStd_Array1OfInteger bidumults(1,2); bidumults.Init(UDegree() + 1);
 571   TColStd_Array1OfReal bidvknots(1,2); bidvknots(1) = 0.; bidvknots(2) = 1.;
 572   TColStd_Array1OfInteger bidvmults(1,2); bidvmults.Init(VDegree() + 1);
 573   Handle(TColgp_HArray2OfPnt) npoles;
 574   Handle(TColStd_HArray2OfReal) nweights;
 575
 576   if(IncUDeg > 0){
 577     npoles = new TColgp_HArray2OfPnt( 1, UDeg + 1, 1, oldVDeg + 1);
 578
 579     if ( urational || vrational) {
 580       nweights = new TColStd_HArray2OfReal( 1, UDeg + 1, 1, VDegree() + 1);
 581
 582       BSplSLib::IncreaseDegree(1, oldUDeg, UDeg, 0,
 583                                poles->Array2(),
 584                                weights->Array2(),
 585                                biduknots, bidumults,
 586                                npoles->ChangeArray2(),
 587                                nweights->ChangeArray2(),
 588                                biduknots, bidumults);
 589       weights = nweights;
 590     }
 591     else {
 592       BSplSLib::IncreaseDegree(1, oldUDeg, UDeg, 0,
 593                                poles->Array2(),
 594                                *((TColStd_Array2OfReal*) NULL),
 595                                biduknots, bidumults,
 596                                npoles->ChangeArray2(),
 597                                *((TColStd_Array2OfReal*) NULL),
 598                                biduknots, bidumults);
 599     }
 600     poles   = npoles;
 601   }
 602   if(IncVDeg > 0){
 603     npoles = new TColgp_HArray2OfPnt( 1, UDeg + 1, 1, VDeg + 1);
 604
 605     if ( urational || vrational) {
 606       nweights = new TColStd_HArray2OfReal( 1, UDeg + 1, 1, VDeg + 1);
 607
 608       BSplSLib::IncreaseDegree(0, oldVDeg, VDeg, 0,
 609                                poles->Array2(),
 610                                weights->Array2(),
 611                                bidvknots, bidvmults,
 612                                npoles->ChangeArray2(),
 613                                nweights->ChangeArray2(),
 614                                bidvknots, bidvmults);
 615       weights = nweights;
 616     }
 617     else {
 618       BSplSLib::IncreaseDegree(0, oldVDeg, VDeg, 0,
 619                                poles->Array2(),
 620                                *((TColStd_Array2OfReal*) NULL),
 621                                bidvknots, bidvmults,
 622                                npoles->ChangeArray2(),
 623                                *((TColStd_Array2OfReal*) NULL),
 624                                bidvknots, bidvmults);
 625
 626     }
 627     poles   = npoles;
 628   }
 629   Init(npoles,nweights);
 630 }
 631
 632 //=======================================================================
 633 //function : InsertPoleColAfter
 634 //purpose  :
 635 //=======================================================================
 636
 637 void Geom_BezierSurface::InsertPoleColAfter
 638   (const Standard_Integer    VIndex,
 639    const TColgp_Array1OfPnt& CPoles)
 640 {
 641   const TColgp_Array2OfPnt & Poles = poles->Array2();
 642   if (VIndex < 1 || VIndex > Poles.RowLength())  Standard_OutOfRange::Raise();
 643   if (CPoles.Length() != Poles.ColLength()) {
 644     Standard_ConstructionError::Raise();
 645   }
 646
 647   Handle(TColgp_HArray2OfPnt) npoles =
 648     new TColgp_HArray2OfPnt(1,poles->ColLength(),1,poles->RowLength()+1);
 649
 650   Handle(TColStd_HArray2OfReal) nweights;
 651
 652   if (urational || vrational) {
 653     nweights =
 654       new TColStd_HArray2OfReal(1,poles->ColLength(),1,poles->RowLength()+1);
 655
 656     TColStd_Array1OfReal CWeights(nweights->LowerRow(),nweights->UpperRow());
 657     CWeights.Init(1.);
 658
 659     AddRatPoleCol (poles->Array2(), weights->Array2(),
 660                    CPoles, CWeights, VIndex,
 661                    npoles->ChangeArray2(), nweights->ChangeArray2());
 662   }
 663   else {
 664     AddPoleCol (poles->Array2(),
 665                 CPoles, VIndex,
 666                 npoles->ChangeArray2());
 667   }
 668   poles   = npoles;
 669   weights = nweights;
 670   coeffs  = new TColgp_HArray2OfPnt(1,poles->ColLength(),
 671                                     1,poles->RowLength());
 672   wcoeffs = new TColStd_HArray2OfReal(1,poles->ColLength(),
 673                                       1,poles->RowLength());
 674   UpdateCoefficients();
 675 }
 676
 677 //=======================================================================
 678 //function : InsertPoleColAfter
 679 //purpose  :
 680 //=======================================================================
 681
 682 void Geom_BezierSurface::InsertPoleColAfter
 683   (const Standard_Integer      VIndex,
 684    const TColgp_Array1OfPnt&   CPoles,
 685    const TColStd_Array1OfReal& CPoleWeights)
 686 {
 687   const TColgp_Array2OfPnt & Poles = poles->Array2();
 688   if (VIndex < 1 || VIndex > Poles.RowLength())  Standard_OutOfRange::Raise();
 689   if (CPoles.Length() != Poles.ColLength() ||
 690       CPoleWeights.Length() != CPoles.Length()) {
 691     Standard_ConstructionError::Raise();
 692     }
 693   Standard_Integer Index = CPoleWeights.Lower();
 694   while (Index <= CPoleWeights.Upper()) {
 695     if (CPoleWeights (Index) <= gp::Resolution()) {
 696       Standard_ConstructionError::Raise();
 697       }
 698     Index++;
 699   }
 700
 701   Handle(TColgp_HArray2OfPnt) npoles =
 702     new TColgp_HArray2OfPnt(1,poles->ColLength(),1,poles->RowLength()+1);
 703
 704   Handle(TColStd_HArray2OfReal) nweights =
 705     new TColStd_HArray2OfReal(1,poles->ColLength(),1,poles->RowLength()+1);
 706
 707   AddRatPoleCol (poles->Array2(), weights->Array2(),
 708                  CPoles, CPoleWeights, VIndex,
 709                  npoles->ChangeArray2(), nweights->ChangeArray2());
 710
 711   poles   = npoles;
 712   weights = nweights;
 713   coeffs  = new TColgp_HArray2OfPnt(1,poles->ColLength(),
 714                                     1,poles->RowLength());
 715   wcoeffs = new TColStd_HArray2OfReal(1,poles->ColLength(),
 716                                       1,poles->RowLength());
 717
 718   Rational(weights->Array2(), urational, vrational);
 719
 720   UpdateCoefficients();
 721 }
 722
 723 //=======================================================================
 724 //function : InsertPoleColBefore
 725 //purpose  :
 726 //=======================================================================
 727
 728 void Geom_BezierSurface::InsertPoleColBefore (const Standard_Integer VIndex,
 729                                               const TColgp_Array1OfPnt& CPoles)
 730 {
 731   InsertPoleColAfter(VIndex - 1, CPoles);
 732 }
 733
 734 //=======================================================================
 735 //function : InsertPoleColBefore
 736 //purpose  :
 737 //=======================================================================
 738
 739 void Geom_BezierSurface::InsertPoleColBefore
 740   (const Standard_Integer      VIndex,
 741    const TColgp_Array1OfPnt&   CPoles,
 742    const TColStd_Array1OfReal& CPoleWeights)
 743 {
 744   InsertPoleColAfter( VIndex - 1, CPoles, CPoleWeights);
 745 }
 746
 747 //=======================================================================
 748 //function : InsertPoleRowAfter
 749 //purpose  :
 750 //=======================================================================
 751
 752 void Geom_BezierSurface::InsertPoleRowAfter (const Standard_Integer UIndex,
 753                                              const TColgp_Array1OfPnt& CPoles)
 754 {
 755   const TColgp_Array2OfPnt & Poles = poles->Array2();
 756   if (UIndex < 1 || UIndex > Poles.ColLength())  Standard_OutOfRange::Raise();
 757   if (CPoles.Length() != Poles.RowLength()) {
 758     Standard_ConstructionError::Raise();
 759   }
 760
 761   Handle(TColgp_HArray2OfPnt) npoles =
 762     new TColgp_HArray2OfPnt(1,poles->ColLength()+1,1,poles->RowLength());
 763
 764   Handle(TColStd_HArray2OfReal) nweights;
 765
 766   if (urational || vrational) {
 767     nweights =
 768       new TColStd_HArray2OfReal(1,poles->ColLength()+1,1,poles->RowLength());
 769
 770 //    TColStd_Array1OfReal CWeights(nweights->LowerCol(),nweights->UpperCol(),
 771 //                                1.0); ???????????
 772     TColStd_Array1OfReal CWeights(1.0,
 773                                   nweights->LowerCol(),nweights->UpperCol());
 774
 775     AddRatPoleRow (poles->Array2(), weights->Array2(),
 776                    CPoles, CWeights, UIndex,
 777                    npoles->ChangeArray2(), nweights->ChangeArray2());
 778   }
 779   else {
 780     AddPoleRow (poles->Array2(),
 781                 CPoles, UIndex,
 782                 npoles->ChangeArray2());
 783   }
 784   poles   = npoles;
 785   weights = nweights;
 786   coeffs  = new TColgp_HArray2OfPnt(1,poles->ColLength(),
 787                                     1,poles->RowLength());
 788   wcoeffs = new TColStd_HArray2OfReal(1,poles->ColLength(),
 789                                       1,poles->RowLength());
 790
 791   UpdateCoefficients();
 792 }
 793
 794 //=======================================================================
 795 //function : InsertPoleRowAfter
 796 //purpose  :
 797 //=======================================================================
 798
 799 void Geom_BezierSurface::InsertPoleRowAfter
 800   (const Standard_Integer      UIndex,
 801    const TColgp_Array1OfPnt&   CPoles,
 802    const TColStd_Array1OfReal& CPoleWeights)
 803 {
 804   const TColgp_Array2OfPnt & Poles = poles->Array2();
 805   if (UIndex < 1 || UIndex > Poles.ColLength())  Standard_OutOfRange::Raise();
 806   if (CPoles.Length() != Poles.RowLength() ||
 807       CPoleWeights.Length() != CPoles.Length()) {
 808     Standard_ConstructionError::Raise();
 809   }
 810   Standard_Integer Index = CPoleWeights.Lower();
 811   while (Index <= CPoleWeights.Upper()) {
 812     if (CPoleWeights(Index) <= gp::Resolution()) {
 813       Standard_ConstructionError::Raise();
 814     }
 815     Index++;
 816   }
 817
 818   Handle(TColgp_HArray2OfPnt) npoles =
 819     new TColgp_HArray2OfPnt(1,poles->ColLength()+1,1,poles->RowLength());
 820
 821   Handle(TColStd_HArray2OfReal) nweights =
 822     new TColStd_HArray2OfReal(1,poles->ColLength()+1,1,poles->RowLength());
 823
 824   AddRatPoleCol (poles->Array2(), weights->Array2(),
 825                  CPoles, CPoleWeights, UIndex,
 826                  npoles->ChangeArray2(), nweights->ChangeArray2());
 827
 828   poles   = npoles;
 829   weights = nweights;
 830   coeffs  = new TColgp_HArray2OfPnt(1,poles->ColLength(),
 831                                     1,poles->RowLength());
 832   wcoeffs = new TColStd_HArray2OfReal(1,poles->ColLength(),
 833                                       1,poles->RowLength());
 834
 835   Rational(weights->Array2(), urational, vrational);
 836
 837   UpdateCoefficients();
 838 }
 839
 840 //=======================================================================
 841 //function : InsertPoleRowBefore
 842 //purpose  :
 843 //=======================================================================
 844
 845 void Geom_BezierSurface::InsertPoleRowBefore (const Standard_Integer UIndex,
 846                                               const TColgp_Array1OfPnt& CPoles)
 847 {
 848   InsertPoleRowAfter( UIndex - 1, CPoles);
 849 }
 850
 851 //=======================================================================
 852 //function : InsertPoleRowBefore
 853 //purpose  :
 854 //=======================================================================
 855
 856 void Geom_BezierSurface::InsertPoleRowBefore
 857   (const Standard_Integer      UIndex,
 858    const TColgp_Array1OfPnt&   CPoles,
 859    const TColStd_Array1OfReal& CPoleWeights)
 860 {
 861   InsertPoleRowAfter( UIndex - 1, CPoles, CPoleWeights);
 862 }
 863
 864 //=======================================================================
 865 //function : RemovePoleCol
 866 //purpose  :
 867 //=======================================================================
 868
 869 void Geom_BezierSurface::RemovePoleCol (const Standard_Integer VIndex)
 870 {
 871   const TColgp_Array2OfPnt & Poles = poles->Array2();
 872   if (VIndex < 1 || VIndex > Poles.RowLength())  Standard_OutOfRange::Raise();
 873   if (Poles.RowLength() <= 2)             Standard_ConstructionError::Raise();
 874
 875   Handle(TColgp_HArray2OfPnt) npoles =
 876     new TColgp_HArray2OfPnt(1,poles->ColLength(),1,poles->RowLength()-1);
 877
 878   Handle(TColStd_HArray2OfReal) nweights;
 879
 880   if (urational || vrational) {
 881     nweights =
 882       new TColStd_HArray2OfReal(1,poles->ColLength(),1,poles->RowLength()-1);
 883
 884     DeleteRatPoleCol (poles->Array2(), weights->Array2(),
 885                       VIndex,
 886                       npoles->ChangeArray2(), nweights->ChangeArray2());
 887     // Mise a jour de la rationalite
 888     Rational(nweights->Array2(), urational, vrational);
 889   }
 890   else {
 891     DeletePoleCol (poles->Array2(),
 892                    VIndex,
 893                    npoles->ChangeArray2());
 894   }
 895   poles   = npoles;
 896   weights = nweights;
 897   coeffs  = new TColgp_HArray2OfPnt(1,poles->ColLength(),
 898                                     1,poles->RowLength());
 899   wcoeffs = new TColStd_HArray2OfReal(1,poles->ColLength(),
 900                                       1,poles->RowLength());
 901   UpdateCoefficients();
 902 }
 903
 904 //=======================================================================
 905 //function : RemovePoleRow
 906 //purpose  :
 907 //=======================================================================
 908
 909 void Geom_BezierSurface::RemovePoleRow (const Standard_Integer UIndex)
 910 {
 911   const TColgp_Array2OfPnt & Poles = poles->Array2();
 912   if (UIndex < 1 || UIndex > Poles.ColLength()) Standard_OutOfRange::Raise();
 913   if (Poles.ColLength() <= 2)            Standard_ConstructionError::Raise();
 914
 915   Handle(TColgp_HArray2OfPnt) npoles =
 916     new TColgp_HArray2OfPnt(1,poles->ColLength()-1,1,poles->RowLength());
 917
 918   Handle(TColStd_HArray2OfReal) nweights;
 919
 920   if (urational || vrational) {
 921     nweights =
 922       new TColStd_HArray2OfReal(1,poles->ColLength()-1,1,poles->RowLength());
 923
 924     DeleteRatPoleRow (poles->Array2(), weights->Array2(),
 925                       UIndex,
 926                       npoles->ChangeArray2(), nweights->ChangeArray2());
 927
 928     // Mise a jour de la rationalite
 929     Rational(nweights->Array2(), urational, vrational);
 930   }
 931   else {
 932     DeletePoleRow (poles->Array2(),
 933                    UIndex,
 934                    npoles->ChangeArray2());
 935   }
 936   poles   = npoles;
 937   weights = nweights;
 938   coeffs  = new TColgp_HArray2OfPnt(1,poles->ColLength(),
 939                                     1,poles->RowLength());
 940   wcoeffs = new TColStd_HArray2OfReal(1,poles->ColLength(),
 941                                       1,poles->RowLength());
 942   UpdateCoefficients();
 943 }
 944
 945 //=======================================================================
 946 //function : Segment
 947 //purpose  :
 948 //=======================================================================
 949
 950 void Geom_BezierSurface::Segment
 951   (const Standard_Real U1,
 952    const Standard_Real U2,
 953    const Standard_Real V1,
 954    const Standard_Real V2)
 955 {
 956   Standard_Boolean rat = (urational || vrational);
 957   Handle(TColgp_HArray2OfPnt)  Coefs;
 958   Handle(TColStd_HArray2OfReal) WCoefs;
 959
 960   if (validcache == 0) UpdateCoefficients(0., 0.);
 961
 962   // Attention si udeg <= vdeg u et v sont intervertis
 963   // dans les coeffs, il faut donc tout transposer.
 964   if(UDegree() <= VDegree()) {
 965     Standard_Integer ii, jj;
 966     Coefs = new  (TColgp_HArray2OfPnt)(1,UDegree()+1,1,VDegree()+1);
 967     if (rat) {
 968       WCoefs = new  (TColStd_HArray2OfReal)(1,UDegree()+1,1,VDegree()+1);
 969     }
 970     for (ii=1; ii<=UDegree()+1; ii++)
 971       for (jj=1; jj<=VDegree()+1; jj++) {
 972         Coefs->SetValue(ii, jj, coeffs->Value(jj,ii));
 973         if (rat)  WCoefs->SetValue(ii, jj, wcoeffs->Value(jj,ii));
 974       }
 975   }
 976   else {
 977    Coefs = coeffs;
 978    if (rat)  {WCoefs =  wcoeffs;}
 979  }
 980
 981 // Trim dans la base cannonique et Update des Poles et Coeffs
 982
 983 // PMN : tranfo sur les parametres
 984   Standard_Real ufirst = 2*(U1 - 0.5),
 985                 ulast =  2*(U2 - 0.5),
 986                 vfirst = 2*(V1 - 0.5),
 987                 vlast  = 2*(V2 - 0.5);
 988   if (rat) {
 989     PLib::UTrimming (ufirst, ulast, Coefs->ChangeArray2(),
 990                      WCoefs->ChangeArray2());
 991     PLib::VTrimming (vfirst, vlast, Coefs->ChangeArray2(),
 992                      WCoefs->ChangeArray2());
 993     PLib::CoefficientsPoles(Coefs->Array2(),
 994                             WCoefs->Array2(),
 995                             poles->ChangeArray2(),
 996                             weights->ChangeArray2());
 997   }
 998   else {
 999     PLib::UTrimming (ufirst, ulast, Coefs->ChangeArray2(),
1000                      *((TColStd_Array2OfReal*) NULL));
1001     PLib::VTrimming (vfirst, vlast, Coefs->ChangeArray2(),
1002                      *((TColStd_Array2OfReal*) NULL));
1003     PLib::CoefficientsPoles (Coefs->Array2(),
1004                              *((TColStd_Array2OfReal*) NULL),
1005                              poles->ChangeArray2(),
1006                              *((TColStd_Array2OfReal*) NULL));
1007   }
1008   UpdateCoefficients();
1009 }
1010
1011 //=======================================================================
1012 //function : SetPole
1013 //purpose  :
1014 //=======================================================================
1015
1016 void Geom_BezierSurface::SetPole
1017   (const Standard_Integer UIndex,
1018    const Standard_Integer VIndex,
1019    const gp_Pnt&          P)
1020 {
1021   TColgp_Array2OfPnt & Poles = poles->ChangeArray2();
1022   if (UIndex < 1                 ||
1023       UIndex > Poles.ColLength() ||
1024       VIndex < 1                 ||
1025       VIndex > Poles.RowLength() )    Standard_OutOfRange::Raise();
1026
1027   Poles (UIndex, VIndex) = P;
1028   UpdateCoefficients();
1029 }
1030
1031 //=======================================================================
1032 //function : SetPole
1033 //purpose  :
1034 //=======================================================================
1035
1036 void Geom_BezierSurface::SetPole
1037   (const Standard_Integer UIndex,
1038    const Standard_Integer VIndex,
1039    const gp_Pnt&          P,
1040    const Standard_Real    Weight)
1041 {
1042
1043   if (Weight <= gp::Resolution())
1044     Standard_ConstructionError::Raise("Geom_BezierSurface::SetPole");
1045   if (UIndex < 1                 ||
1046       UIndex > poles->ColLength() ||
1047       VIndex < 1                 ||
1048       VIndex > poles->RowLength())
1049     Standard_OutOfRange::Raise("Geom_BezierSurface::SetPole");
1050
1051   poles->SetValue(UIndex, VIndex, P);
1052
1053   SetWeight(UIndex, VIndex,  Weight); //L'update des coeff est fait la dedans
1054 }
1055
1056 //=======================================================================
1057 //function : SetPoleCol
1058 //purpose  :
1059 //=======================================================================
1060
1061 void Geom_BezierSurface::SetPoleCol
1062   (const Standard_Integer      VIndex,
1063    const TColgp_Array1OfPnt&   CPoles,
1064    const TColStd_Array1OfReal& CPoleWeights)
1065 {
1066   TColgp_Array2OfPnt & Poles = poles->ChangeArray2();
1067   if (VIndex < 1 || VIndex > Poles.RowLength()) Standard_OutOfRange::Raise();
1068
1069   if (CPoles.Lower() < 1                     ||
1070       CPoles.Lower() > Poles.ColLength()     ||
1071       CPoles.Upper() < 1                     ||
1072       CPoles.Upper() > Poles.ColLength()     ||
1073       CPoleWeights.Lower() != CPoles.Lower() ||
1074       CPoleWeights.Upper() != CPoles.Upper()) {
1075     Standard_ConstructionError::Raise();
1076   }
1077
1078   Standard_Integer I;
1079   for (I = CPoles.Lower();  I <= CPoles.Upper(); I++) {
1080     Poles (I, VIndex) = CPoles (I);
1081   }
1082   SetWeightCol(VIndex, CPoleWeights); //Avec l'update
1083 }
1084
1085 //=======================================================================
1086 //function : SetPoleCol
1087 //purpose  :
1088 //=======================================================================
1089
1090 void Geom_BezierSurface::SetPoleCol (const Standard_Integer      VIndex,
1091                                      const TColgp_Array1OfPnt&   CPoles)
1092 {
1093   TColgp_Array2OfPnt & Poles = poles->ChangeArray2();
1094   if (VIndex < 1 || VIndex > Poles.RowLength())  Standard_OutOfRange::Raise();
1095
1096   if (CPoles.Lower() < 1                 ||
1097       CPoles.Lower() > Poles.ColLength() ||
1098       CPoles.Upper() < 1                 ||
1099       CPoles.Upper() > Poles.ColLength()) {
1100     Standard_ConstructionError::Raise();
1101   }
1102   for (Standard_Integer I = CPoles.Lower(); I <= CPoles.Upper(); I++) {
1103     Poles (I, VIndex) = CPoles (I);
1104   }
1105
1106   UpdateCoefficients();
1107 }
1108
1109 //=======================================================================
1110 //function : SetPoleRow
1111 //purpose  :
1112 //=======================================================================
1113
1114 void Geom_BezierSurface::SetPoleRow (const Standard_Integer    UIndex,
1115                                      const TColgp_Array1OfPnt& CPoles)
1116 {
1117   TColgp_Array2OfPnt & Poles = poles->ChangeArray2();
1118   if (UIndex < 1 || UIndex > Poles.ColLength())  Standard_OutOfRange::Raise();
1119
1120   if (CPoles.Lower() < 1                 ||
1121       CPoles.Lower() > Poles.RowLength() ||
1122       CPoles.Upper() < 1                 ||
1123       CPoles.Upper() > Poles.RowLength())  Standard_ConstructionError::Raise();
1124
1125   for (Standard_Integer I = CPoles.Lower(); I <= CPoles.Upper(); I++) {
1126     Poles (UIndex, I) = CPoles (I);
1127   }
1128   UpdateCoefficients();
1129 }
1130
1131 //=======================================================================
1132 //function : SetPoleRow
1133 //purpose  :
1134 //=======================================================================
1135
1136 void Geom_BezierSurface::SetPoleRow
1137   (const Standard_Integer      UIndex,
1138    const TColgp_Array1OfPnt&   CPoles,
1139    const TColStd_Array1OfReal& CPoleWeights)
1140 {
1141   TColgp_Array2OfPnt & Poles = poles->ChangeArray2();
1142   if (UIndex < 1 || UIndex > Poles.ColLength())   Standard_OutOfRange::Raise();
1143
1144   if (CPoles.Lower() < 1                     ||
1145       CPoles.Lower() > Poles.RowLength()     ||
1146       CPoles.Upper() < 1                     ||
1147       CPoles.Upper() > Poles.RowLength()     ||
1148       CPoleWeights.Lower() != CPoles.Lower() ||
1149       CPoleWeights.Upper() != CPoles.Upper()) {
1150     Standard_ConstructionError::Raise();
1151   }
1152
1153   Standard_Integer I;
1154
1155   for (I = CPoles.Lower(); I <= CPoles.Upper(); I++) {
1156     Poles   (UIndex, I) = CPoles (I);
1157   }
1158
1159   SetWeightRow(UIndex, CPoleWeights); //Avec l'update
1160 }
1161
1162 //=======================================================================
1163 //function : SetWeight
1164 //purpose  :
1165 //=======================================================================
1166
1167 void Geom_BezierSurface::SetWeight (const Standard_Integer UIndex,
1168                                     const Standard_Integer VIndex,
1169                                     const Standard_Real    Weight)
1170 {
1171    // compute new rationality
1172   Standard_Boolean wasrat = (urational||vrational);
1173   if (!wasrat) {
1174     // a weight of 1. does not turn to rational
1175     if (Abs(Weight - 1.) <= gp::Resolution()) {
1176       UpdateCoefficients(); //Pour l'appel via SetPole
1177       return;
1178     }
1179
1180     // set weights of 1.
1181     weights = new TColStd_HArray2OfReal (1, poles->ColLength(),
1182                                          1, poles->RowLength(), 1.);
1183     wcoeffs = new TColStd_HArray2OfReal (1, poles->ColLength(),
1184                                          1, poles->RowLength());
1185   }
1186
1187   TColStd_Array2OfReal & Weights = weights->ChangeArray2();
1188   if (Weight <= gp::Resolution())
1189     Standard_ConstructionError::Raise("Geom_BezierSurface::SetWeight");
1190
1191   if (UIndex < 1                   ||
1192       UIndex > Weights.ColLength() ||
1193       VIndex < 1                   ||
1194       VIndex > Weights.RowLength())    Standard_OutOfRange::Raise();
1195
1196   if (Abs (Weight - Weights (UIndex, VIndex)) > gp::Resolution()) {
1197     Weights (UIndex, VIndex) = Weight;
1198     Rational(Weights, urational, vrational);
1199   }
1200
1201  // is it turning into non rational
1202   if (wasrat) {
1203     if (!(urational || vrational)) {
1204       weights.Nullify();
1205       wcoeffs.Nullify();
1206     }
1207   }
1208
1209   UpdateCoefficients(); //Dans tous cas :Attention a SetPoleCol !
1210 }
1211
1212 //=======================================================================
1213 //function : SetWeightCol
1214 //purpose  :
1215 //=======================================================================
1216
1217 void Geom_BezierSurface::SetWeightCol
1218   (const Standard_Integer      VIndex,
1219    const TColStd_Array1OfReal& CPoleWeights)
1220 {
1221   Standard_Integer I;
1222    // compute new rationality
1223   Standard_Boolean wasrat = (urational||vrational);
1224   if (!wasrat) {
1225     // set weights of 1.
1226     weights = new TColStd_HArray2OfReal (1, poles->ColLength(),
1227                                          1, poles->RowLength(), 1.);
1228     wcoeffs = new TColStd_HArray2OfReal (1, poles->ColLength(),
1229                                          1, poles->RowLength());
1230   }
1231
1232   TColStd_Array2OfReal & Weights = weights->ChangeArray2();
1233   if (VIndex < 1 || VIndex > Weights.RowLength()) Standard_OutOfRange::Raise();
1234
1235   if (CPoleWeights.Length() !=  Weights.ColLength())  {
1236     Standard_ConstructionError::Raise("Geom_BezierSurface::SetWeightCol");
1237   }
1238
1239   I = CPoleWeights.Lower();
1240   while (I <= CPoleWeights.Upper()) {
1241     if (CPoleWeights(I) <= gp::Resolution()) {
1242       Standard_ConstructionError::Raise();
1243     }
1244     Weights (I, VIndex) = CPoleWeights (I);
1245     I++;
1246   }
1247
1248  Rational(Weights, urational, vrational);
1249
1250  // is it turning into non rational
1251   if (wasrat) {
1252     if (!(urational || vrational)) {
1253       weights.Nullify();
1254       wcoeffs.Nullify();
1255     }
1256   }
1257
1258   UpdateCoefficients();
1259 }
1260
1261 //=======================================================================
1262 //function : SetWeightRow
1263 //purpose  :
1264 //=======================================================================
1265
1266 void Geom_BezierSurface::SetWeightRow
1267   (const Standard_Integer      UIndex,
1268    const TColStd_Array1OfReal& CPoleWeights)
1269 {
1270   Standard_Integer I;
1271    // compute new rationality
1272   Standard_Boolean wasrat = (urational||vrational);
1273   if (!wasrat) {
1274     // set weights of 1.
1275     weights = new TColStd_HArray2OfReal (1, poles->ColLength(),
1276                                          1, poles->RowLength(), 1.);
1277     wcoeffs = new TColStd_HArray2OfReal (1, poles->ColLength(),
1278                                          1, poles->RowLength());
1279   }
1280
1281   TColStd_Array2OfReal & Weights = weights->ChangeArray2();
1282   if (UIndex < 1 || UIndex > Weights.ColLength())
1283     Standard_OutOfRange::Raise("Geom_BezierSurface::SetWeightRow");
1284   if (CPoleWeights.Lower() < 1 ||
1285       CPoleWeights.Lower() > Weights.RowLength() ||
1286       CPoleWeights.Upper() < 1 ||
1287       CPoleWeights.Upper() > Weights.RowLength()  ) {
1288     Standard_ConstructionError::Raise("Geom_BezierSurface::SetWeightRow");
1289   }
1290
1291   I = CPoleWeights.Lower();
1292   while (I <= CPoleWeights.Upper()) {
1293     if (CPoleWeights(I) <= gp::Resolution())  {
1294       Standard_ConstructionError::Raise();
1295     }
1296     Weights (UIndex, I) = CPoleWeights (I);
1297     I++;
1298   }
1299
1300   Rational(Weights, urational, vrational);
1301
1302  // is it turning into non rational
1303   if (wasrat) {
1304     if (!(urational || vrational)) {
1305       weights.Nullify();
1306       wcoeffs.Nullify();
1307     }
1308   }
1309
1310   UpdateCoefficients();
1311 }
1312
1313 //=======================================================================
1314 //function : UReverse
1315 //purpose  :
1316 //=======================================================================
1317
1318 void Geom_BezierSurface::UReverse ()
1319 {
1320   gp_Pnt Pol;
1321   Standard_Integer Row,Col;
1322   TColgp_Array2OfPnt & Poles = poles->ChangeArray2();
1323   if (urational || vrational) {
1324     TColStd_Array2OfReal & Weights = weights->ChangeArray2();
1325     Standard_Real W;
1326     for (Col = 1; Col <= Poles.RowLength(); Col++) {
1327       for (Row = 1; Row <= IntegerPart (Poles.ColLength() / 2); Row++) {
1328         W = Weights (Row, Col);
1329         Weights (Row, Col) = Weights (Poles.ColLength()- Row + 1, Col);
1330         Weights (Poles.ColLength() - Row + 1, Col) = W;
1331         Pol = Poles (Row, Col);
1332         Poles (Row, Col) = Poles (Poles.ColLength() - Row + 1, Col);
1333         Poles (Poles.ColLength() - Row + 1, Col) = Pol;
1334       }
1335     }
1336   }
1337   else {
1338     for (Col = 1; Col <= Poles.RowLength(); Col++) {
1339       for (Row = 1; Row <= IntegerPart (Poles.ColLength() / 2); Row++) {
1340         Pol = Poles (Row, Col);
1341         Poles (Row, Col) = Poles (Poles.ColLength() - Row + 1, Col);
1342         Poles (Poles.ColLength() - Row + 1, Col) = Pol;
1343       }
1344     }
1345   }
1346   UpdateCoefficients();
1347 }
1348
1349 //=======================================================================
1350 //function : UReversedParameter
1351 //purpose  :
1352 //=======================================================================
1353
1354 Standard_Real Geom_BezierSurface::UReversedParameter
1355   ( const Standard_Real U) const
1356 {
1357   return ( 1. - U);
1358 }
1359
1360 //=======================================================================
1361 //function : VReverse
1362 //purpose  :
1363 //=======================================================================
1364
1365 void Geom_BezierSurface::VReverse ()
1366 {
1367   gp_Pnt Pol;
1368   Standard_Integer Row,Col;
1369   TColgp_Array2OfPnt & Poles = poles->ChangeArray2();
1370   if (urational || vrational) {
1371     TColStd_Array2OfReal & Weights = weights->ChangeArray2();
1372     Standard_Real W;
1373     for (Row = 1; Row <= Poles.ColLength(); Row++) {
1374       for (Col = 1; Col <= IntegerPart (Poles.RowLength()/2); Col++) {
1375         W = Weights (Row, Col);
1376         Weights (Row, Col) = Weights (Row, Poles.RowLength() - Col + 1);
1377         Weights (Row, Poles.RowLength() - Col + 1) = W;
1378         Pol = Poles (Row, Col);
1379         Poles (Row, Col) = Poles (Row, Poles.RowLength() - Col + 1);
1380         Poles (Row, Poles.RowLength() - Col + 1) = Pol;
1381       }
1382     }
1383   }
1384   else {
1385     for (Row = 1; Row <= Poles.ColLength(); Row++) {
1386       for (Col = 1; Col <= IntegerPart(Poles.RowLength()/2); Col++) {
1387         Pol = Poles (Row, Col);
1388         Poles (Row, Col)= Poles (Row, Poles.RowLength() - Col + 1);
1389         Poles (Row, Poles.RowLength() - Col + 1) = Pol;
1390       }
1391     }
1392   }
1393   UpdateCoefficients();
1394 }
1395
1396 //=======================================================================
1397 //function : VReversedParameter
1398 //purpose  :
1399 //=======================================================================
1400
1401 Standard_Real Geom_BezierSurface::VReversedParameter
1402   ( const Standard_Real V) const
1403 {
1404   return ( 1. - V);
1405 }
1406
1407 //=======================================================================
1408 //function : Bounds
1409 //purpose  :
1410 //=======================================================================
1411
1412 void Geom_BezierSurface::Bounds (Standard_Real& U1,
1413                                  Standard_Real& U2,
1414                                  Standard_Real& V1,
1415                                  Standard_Real& V2) const
1416 {
1417   U1 = 0.0;
1418   U2 = 1.0;
1419   V1 = 0.0;
1420   V2 = 1.0;
1421 }
1422
1423 //=======================================================================
1424 //function : Continuity
1425 //purpose  :
1426 //=======================================================================
1427
1428 GeomAbs_Shape Geom_BezierSurface::Continuity () const
1429 {
1430   return GeomAbs_CN;
1431 }
1432
1433 //=======================================================================
1434 //function : D0
1435 //purpose  :
1436 //=======================================================================
1437
1438 void Geom_BezierSurface::D0 (const Standard_Real U,
1439                              const Standard_Real V,
1440                                    gp_Pnt&       P ) const
1441 {
1442
1443   if (validcache == 1) {
1444     //
1445     // XAB : cet algorithme devient instable pour les hauts degres
1446     // RBD : Beaucoup moins maintenant avec le calcul d'un nouveau cache
1447     //       sur [-1,1].
1448     //
1449     Standard_Real uparameter_11 = (2*ucacheparameter + ucachespanlenght)/2,
1450                   uspanlenght_11 = ucachespanlenght/2,
1451                   vparameter_11 = (2*vcacheparameter + vcachespanlenght)/2,
1452                   vspanlenght_11 = vcachespanlenght/2 ;
1453     if (urational || vrational) {
1454       BSplSLib::CacheD0(U, V, UDegree(), VDegree(),
1455                         uparameter_11, vparameter_11,
1456                         uspanlenght_11, vspanlenght_11,
1457                         coeffs->Array2(),
1458                         wcoeffs->Array2(),
1459                         P);
1460     }
1461     else {
1462       BSplSLib::CacheD0(U, V, UDegree(), VDegree(),
1463                         uparameter_11, vparameter_11,
1464                         uspanlenght_11, vspanlenght_11,
1465                         coeffs->Array2(),
1466                         *((TColStd_Array2OfReal*) NULL),
1467                         P);
1468     }
1469   }
1470   else {
1471     Standard_Real array_u[2] ;
1472     Standard_Real array_v[2] ;
1473     Standard_Integer mult_u[2] ;
1474     Standard_Integer mult_v[2] ;
1475     TColStd_Array1OfReal biduknots(array_u[0],1,2); biduknots(1) = 0.; biduknots(2) = 1.;
1476     TColStd_Array1OfInteger bidumults(mult_u[0],1,2); bidumults.Init(UDegree() + 1);
1477     TColStd_Array1OfReal bidvknots(array_v[0],1,2); bidvknots(1) = 0.; bidvknots(2) = 1.;
1478     TColStd_Array1OfInteger bidvmults(mult_v[0],1,2); bidvmults.Init(VDegree() + 1);
1479     if (urational || vrational) {
1480         BSplSLib::D0(U, V, 1,1,poles->Array2(),
1481                      weights->Array2(),
1482                      biduknots,bidvknots,bidumults,bidvmults,
1483                      UDegree(),VDegree(),
1484                      urational,vrational,Standard_False,Standard_False,
1485                      P) ;
1486     }
1487     else {
1488
1489       BSplSLib::D0(U, V, 1,1,poles->Array2(),
1490                    *((TColStd_Array2OfReal*) NULL),
1491                    biduknots,bidvknots,bidumults,bidvmults,
1492                    UDegree(),VDegree(),
1493                    urational,vrational,Standard_False,Standard_False,
1494                    P) ;
1495     }
1496   }
1497 }
1498
1499 //=======================================================================
1500 //function : D1
1501 //purpose  :
1502 //=======================================================================
1503
1504 void Geom_BezierSurface::D1
1505   (const Standard_Real U,
1506    const Standard_Real V,
1507          gp_Pnt&       P,
1508          gp_Vec&       D1U,
1509          gp_Vec&       D1V ) const
1510 {
1511
1512   if (validcache == 1) {
1513     //
1514     // XAB : cet algorithme devient instable pour les hauts degres
1515     // RBD : Beaucoup moins maintenant avec le calcul d'un nouveau cache
1516     //       sur [-1,1].
1517     //
1518     Standard_Real uparameter_11 = (2*ucacheparameter + ucachespanlenght)/2,
1519                   uspanlenght_11 = ucachespanlenght/2,
1520                   vparameter_11 = (2*vcacheparameter + vcachespanlenght)/2,
1521                   vspanlenght_11 = vcachespanlenght/2 ;
1522    if (urational || vrational) {
1523      BSplSLib::CacheD1(U, V, UDegree(), VDegree(),
1524                        uparameter_11, vparameter_11,
1525                        uspanlenght_11, vspanlenght_11,
1526                        coeffs->Array2(),
1527                        wcoeffs->Array2(),
1528                        P, D1U, D1V);
1529    }
1530    else {
1531      BSplSLib::CacheD1(U, V, UDegree(), VDegree(),
1532                        uparameter_11, vparameter_11,
1533                        uspanlenght_11, vspanlenght_11,
1534                        coeffs->Array2(),
1535                        *((TColStd_Array2OfReal*) NULL),
1536                        P, D1U, D1V);
1537    }
1538   }
1539   else {
1540     Standard_Real array_u[2] ;
1541     Standard_Real array_v[2] ;
1542     Standard_Integer mult_u[2] ;
1543     Standard_Integer mult_v[2] ;
1544     TColStd_Array1OfReal biduknots(array_u[0],1,2); biduknots(1) = 0.; biduknots(2) = 1.;
1545     TColStd_Array1OfInteger bidumults(mult_u[0],1,2); bidumults.Init(UDegree() + 1);
1546     TColStd_Array1OfReal bidvknots(array_v[0],1,2); bidvknots(1) = 0.; bidvknots(2) = 1.;
1547     TColStd_Array1OfInteger bidvmults(mult_v[0],1,2); bidvmults.Init(VDegree() + 1);
1548     if (urational || vrational) {
1549       BSplSLib::D1(U, V, 1,1,poles->Array2(),
1550                    weights->Array2(),
1551                    biduknots,bidvknots,bidumults,bidvmults,
1552                    UDegree(),VDegree(),
1553                    urational,vrational,Standard_False,Standard_False,
1554                    P,D1U, D1V) ;
1555     }
1556     else {
1557       BSplSLib::D1(U, V, 1,1,poles->Array2(),
1558                    *((TColStd_Array2OfReal*) NULL),
1559                    biduknots,bidvknots,bidumults,bidvmults,
1560                    UDegree(),VDegree(),
1561                    urational,vrational,Standard_False,Standard_False,
1562                    P,D1U, D1V) ;
1563     }
1564   }
1565 }
1566
1567 //=======================================================================
1568 //function : D2
1569 //purpose  :
1570 //=======================================================================
1571
1572 void Geom_BezierSurface::D2
1573   (const Standard_Real U,
1574    const Standard_Real V,
1575          gp_Pnt&       P,
1576          gp_Vec&       D1U, gp_Vec& D1V,
1577          gp_Vec&       D2U, gp_Vec& D2V, gp_Vec& D2UV ) const
1578 {
1579
1580   if (validcache == 1) {
1581     //
1582     // XAB : cet algorithme devient instable pour les hauts degres
1583     // RBD : Beaucoup moins maintenant avec le calcul d'un nouveau cache
1584     //       sur [-1,1].
1585     //
1586     Standard_Real uparameter_11 = (2*ucacheparameter + ucachespanlenght)/2,
1587                   uspanlenght_11 = ucachespanlenght/2,
1588                   vparameter_11 = (2*vcacheparameter + vcachespanlenght)/2,
1589                   vspanlenght_11 = vcachespanlenght/2 ;
1590     if (urational || vrational) {
1591       //-- ATTENTION a l'ORDRE d'appel ds BSPLSLIB
1592       BSplSLib::CacheD2(U, V, UDegree(), VDegree(),
1593                         uparameter_11, vparameter_11,
1594                         uspanlenght_11, vspanlenght_11,
1595                         coeffs->Array2(),
1596                         wcoeffs->Array2(),
1597                         P, D1U, D1V, D2U, D2UV , D2V);
1598     }
1599     else {
1600       //-- ATTENTION a l'ORDRE d'appel ds BSPLSLIB
1601       BSplSLib::CacheD2(U, V, UDegree(), VDegree(),
1602                         uparameter_11, vparameter_11,
1603                         uspanlenght_11, vspanlenght_11,
1604                         coeffs->Array2(),
1605                         *((TColStd_Array2OfReal*) NULL),
1606                         P, D1U, D1V, D2U, D2UV , D2V);
1607     }
1608   }
1609   else {
1610     Standard_Real array_u[2] ;
1611     Standard_Real array_v[2] ;
1612     Standard_Integer mult_u[2] ;
1613     Standard_Integer mult_v[2] ;
1614     TColStd_Array1OfReal biduknots(array_u[0],1,2); biduknots(1) = 0.; biduknots(2) = 1.;
1615     TColStd_Array1OfInteger bidumults(mult_u[0],1,2); bidumults.Init(UDegree() + 1);
1616     TColStd_Array1OfReal bidvknots(array_v[0],1,2); bidvknots(1) = 0.; bidvknots(2) = 1.;
1617     TColStd_Array1OfInteger bidvmults(mult_v[0],1,2); bidvmults.Init(VDegree() + 1);
1618     if (urational || vrational) {
1619       //-- ATTENTION a l'ORDRE d'appel ds BSPLSLIB
1620       BSplSLib::D2(U, V, 1,1,poles->Array2(),
1621                    weights->Array2(),
1622                    biduknots,bidvknots,bidumults,bidvmults,
1623                    UDegree(),VDegree(),
1624                    urational,vrational,Standard_False,Standard_False,
1625                    P,D1U, D1V, D2U, D2V , D2UV) ;
1626     }
1627     else {
1628       //-- ATTENTION a l'ORDRE d'appel ds BSPLSLIB
1629       BSplSLib::D2(U, V, 1,1,poles->Array2(),
1630                    *((TColStd_Array2OfReal*) NULL),
1631                    biduknots,bidvknots,bidumults,bidvmults,
1632                    UDegree(),VDegree(),
1633                    urational,vrational,Standard_False,Standard_False,
1634                    P,D1U, D1V, D2U, D2V, D2UV ) ;
1635     }
1636   }
1637 }
1638
1639 //=======================================================================
1640 //function : D3
1641 //purpose  :
1642 //=======================================================================
1643
1644 void Geom_BezierSurface::D3
1645   (const Standard_Real U, const Standard_Real V,
1646    gp_Pnt& P,
1647    gp_Vec& D1U, gp_Vec& D1V,
1648    gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV,
1649    gp_Vec& D3U, gp_Vec& D3V, gp_Vec& D3UUV, gp_Vec& D3UVV) const
1650 {
1651   TColStd_Array1OfReal biduknots(1,2); biduknots(1) = 0.; biduknots(2) = 1.;
1652   TColStd_Array1OfInteger bidumults(1,2); bidumults.Init(UDegree() + 1);
1653   TColStd_Array1OfReal bidvknots(1,2); bidvknots(1) = 0.; bidvknots(2) = 1.;
1654   TColStd_Array1OfInteger bidvmults(1,2); bidvmults.Init(VDegree() + 1);
1655   if (urational || vrational) {
1656     BSplSLib::D3 (U, V, 0, 0, poles->Array2(),
1657                   weights->Array2(),
1658                   biduknots, bidvknots, bidumults, bidvmults,
1659                   UDegree(), VDegree(), urational, vrational, 0, 0,
1660                   P,
1661                   D1U, D1V,
1662                   D2U, D2V, D2UV,
1663                   D3U, D3V, D3UUV, D3UVV);
1664   }
1665   else {
1666     BSplSLib::D3 (U, V, 0, 0, poles->Array2(),
1667                   *((TColStd_Array2OfReal*) NULL),
1668                   biduknots, bidvknots, bidumults, bidvmults,
1669                   UDegree(), VDegree(), urational, vrational, 0, 0,
1670                   P,
1671                   D1U, D1V,
1672                   D2U, D2V, D2UV,
1673                   D3U, D3V, D3UUV, D3UVV);
1674   }
1675 }
1676
1677 //=======================================================================
1678 //function : DN
1679 //purpose  :
1680 //=======================================================================
1681
1682 gp_Vec Geom_BezierSurface::DN
1683   (const Standard_Real    U,
1684    const Standard_Real    V,
1685    const Standard_Integer Nu,
1686    const Standard_Integer Nv) const
1687 {
1688   Standard_RangeError_Raise_if (Nu + Nv < 1 || Nv < 0 || Nu <0, " ");
1689   gp_Vec Derivative;
1690   TColStd_Array1OfReal biduknots(1,2); biduknots(1) = 0.; biduknots(2) = 1.;
1691   TColStd_Array1OfInteger bidumults(1,2); bidumults.Init(UDegree() + 1);
1692   TColStd_Array1OfReal bidvknots(1,2); bidvknots(1) = 0.; bidvknots(2) = 1.;
1693   TColStd_Array1OfInteger bidvmults(1,2); bidvmults.Init(VDegree() + 1);
1694   if (urational || vrational) {
1695     BSplSLib::DN (U, V, Nu, Nv, 0, 0, poles->Array2(),
1696                   weights->Array2(),
1697                   biduknots, bidvknots, bidumults, bidvmults,
1698                   UDegree(), VDegree(), urational, vrational, 0, 0,
1699                   Derivative);
1700   }
1701   else {
1702     BSplSLib::DN (U, V, Nu, Nv, 0, 0, poles->Array2(),
1703                   *((TColStd_Array2OfReal*) NULL),
1704                   biduknots, bidvknots, bidumults, bidvmults,
1705                   UDegree(), VDegree(), urational, vrational, 0, 0,
1706                   Derivative);
1707   }
1708   return Derivative;
1709 }
1710
1711 //=======================================================================
1712 //function : NbUPoles
1713 //purpose  :
1714 //=======================================================================
1715
1716 Standard_Integer Geom_BezierSurface::NbUPoles () const
1717 {
1718   return poles->ColLength();
1719 }
1720
1721 //=======================================================================
1722 //function : NbVPoles
1723 //purpose  :
1724 //=======================================================================
1725
1726 Standard_Integer Geom_BezierSurface::NbVPoles () const
1727 {
1728   return poles->RowLength();
1729 }
1730
1731 //=======================================================================
1732 //function : Pole
1733 //purpose  :
1734 //=======================================================================
1735
1736 gp_Pnt Geom_BezierSurface::Pole (const Standard_Integer UIndex,
1737                                  const Standard_Integer VIndex) const
1738 {
1739   Standard_OutOfRange_Raise_if
1740     (UIndex < 1 || UIndex > poles->ColLength() ||
1741      VIndex < 1 || VIndex > poles->RowLength(), " ");
1742   return poles->Value (UIndex + poles->LowerRow() - 1,
1743                        VIndex + poles->LowerCol() - 1);
1744 }
1745
1746 //=======================================================================
1747 //function : Poles
1748 //purpose  :
1749 //=======================================================================
1750
1751 void Geom_BezierSurface::Poles (TColgp_Array2OfPnt& P) const
1752 {
1753   Standard_DimensionError_Raise_if
1754     (P.RowLength() != poles->RowLength() ||
1755      P.ColLength() != poles->ColLength(), " ");
1756   P = poles->Array2();
1757 }
1758
1759 //=======================================================================
1760 //function : UDegree
1761 //purpose  :
1762 //=======================================================================
1763
1764 Standard_Integer Geom_BezierSurface::UDegree () const
1765 {
1766   return poles->ColLength() - 1;
1767 }
1768
1769 //=======================================================================
1770 //function : UIso
1771 //purpose  :
1772 //=======================================================================
1773
1774 Handle(Geom_Curve) Geom_BezierSurface::UIso (const Standard_Real U) const
1775 {
1776   TColStd_Array1OfReal biduknots(1,2); biduknots(1) = 0.; biduknots(2) = 1.;
1777   TColStd_Array1OfInteger bidumults(1,2); bidumults.Init(UDegree() + 1);
1778
1779   Handle(Geom_BezierCurve) UIsoCurve;
1780   const TColgp_Array2OfPnt & Poles = poles->Array2();
1781   TColgp_Array1OfPnt VCurvePoles (Poles.LowerCol() , Poles.UpperCol());
1782   if (urational || vrational) {
1783     const TColStd_Array2OfReal & Weights = weights->Array2();
1784     TColStd_Array1OfReal VCurveWeights
1785       (Weights.LowerCol() , Weights.UpperCol());
1786     BSplSLib::Iso (U, 1, Poles,
1787                    Weights,
1788                    biduknots, bidumults,
1789                    UDegree(), 0, VCurvePoles, VCurveWeights);
1790     if (urational)
1791       UIsoCurve = new Geom_BezierCurve (VCurvePoles, VCurveWeights);
1792     else
1793       UIsoCurve = new Geom_BezierCurve (VCurvePoles);
1794   }
1795   else {
1796     BSplSLib::Iso (U, 1, Poles,
1797                    *((TColStd_Array2OfReal*) NULL),
1798                    biduknots, bidumults,
1799                    UDegree(), 0, VCurvePoles, PLib::NoWeights());
1800     UIsoCurve = new Geom_BezierCurve (VCurvePoles);
1801   }
1802   return UIsoCurve;
1803 }
1804
1805 //=======================================================================
1806 //function : VDegree
1807 //purpose  :
1808 //=======================================================================
1809
1810 Standard_Integer Geom_BezierSurface::VDegree () const
1811 {
1812   return poles->RowLength() - 1;
1813 }
1814
1815 //=======================================================================
1816 //function : VIso
1817 //purpose  :
1818 //=======================================================================
1819
1820 Handle(Geom_Curve) Geom_BezierSurface::VIso (const Standard_Real V) const
1821 {
1822   TColStd_Array1OfReal bidvknots(1,2); bidvknots(1) = 0.; bidvknots(2) = 1.;
1823   TColStd_Array1OfInteger bidvmults(1,2); bidvmults.Init(VDegree() + 1);
1824
1825   Handle(Geom_BezierCurve) VIsoCurve;
1826   const TColgp_Array2OfPnt & Poles = poles->Array2();
1827   TColgp_Array1OfPnt VCurvePoles (Poles.LowerRow() , Poles.UpperRow());
1828   if (vrational || urational) {
1829     const TColStd_Array2OfReal & Weights = weights->Array2();
1830     TColStd_Array1OfReal VCurveWeights
1831       (Weights.LowerRow() , Weights.UpperRow());
1832     BSplSLib::Iso (V, 0, Poles,
1833                    Weights,
1834                    bidvknots, bidvmults,
1835                    VDegree(), 0, VCurvePoles, VCurveWeights);
1836     if (vrational)
1837       VIsoCurve = new Geom_BezierCurve (VCurvePoles, VCurveWeights);
1838     else
1839       VIsoCurve = new Geom_BezierCurve (VCurvePoles);
1840   }
1841   else {
1842     BSplSLib::Iso (V, 0, Poles,
1843                    *((TColStd_Array2OfReal*) NULL),
1844                    bidvknots, bidvmults,
1845                    VDegree(), 0, VCurvePoles, PLib::NoWeights());
1846     VIsoCurve = new Geom_BezierCurve (VCurvePoles);
1847   }
1848   return VIsoCurve;
1849 }
1850
1851 //=======================================================================
1852 //function : Weight
1853 //purpose  :
1854 //=======================================================================
1855
1856 Standard_Real Geom_BezierSurface::Weight (const Standard_Integer UIndex,
1857                                           const Standard_Integer VIndex) const
1858 {
1859   Standard_OutOfRange_Raise_if (
1860            UIndex < 1 || UIndex > weights->ColLength() ||
1861            VIndex < 1 || VIndex > weights->RowLength(), " ");
1862
1863   if (urational || vrational)
1864     return weights->Value (UIndex, VIndex);
1865   else
1866     return 1;
1867 }
1868
1869 //=======================================================================
1870 //function : Weights
1871 //purpose  :
1872 //=======================================================================
1873
1874 void Geom_BezierSurface::Weights (TColStd_Array2OfReal& W ) const
1875 {
1876   Standard_DimensionError_Raise_if (
1877            W.RowLength() != weights->RowLength() ||
1878            W.ColLength() != weights->ColLength(), " " );
1879   if (urational || vrational)
1880     W = weights->Array2();
1881   else
1882     W.Init(1.);
1883 }
1884
1885 //=======================================================================
1886 //function : IsCNu
1887 //purpose  :
1888 //=======================================================================
1889
1890 Standard_Boolean Geom_BezierSurface::IsCNu (const Standard_Integer ) const
1891 {
1892   return Standard_True;
1893 }
1894
1895 //=======================================================================
1896 //function : IsCNv
1897 //purpose  :
1898 //=======================================================================
1899
1900 Standard_Boolean Geom_BezierSurface::IsCNv (const Standard_Integer ) const
1901 {
1902   return Standard_True;
1903 }
1904
1905 //=======================================================================
1906 //function : IsURational
1907 //purpose  :
1908 //=======================================================================
1909
1910 Standard_Boolean Geom_BezierSurface::IsURational () const
1911 {
1912   return urational;
1913 }
1914
1915 //=======================================================================
1916 //function : IsVRational
1917 //purpose  :
1918 //=======================================================================
1919
1920 Standard_Boolean Geom_BezierSurface::IsVRational () const
1921 {
1922   return vrational;
1923 }
1924
1925 //=======================================================================
1926 //function : Transform
1927 //purpose  :
1928 //=======================================================================
1929
1930 void Geom_BezierSurface::Transform (const gp_Trsf& T)
1931 {
1932   TColgp_Array2OfPnt & Poles = poles->ChangeArray2();
1933
1934   for (Standard_Integer I = 1; I <= Poles.ColLength(); I++) {
1935
1936     for (Standard_Integer J = 1; J <= Poles.RowLength(); J++) {
1937       Poles (I, J).Transform (T);
1938     }
1939   }
1940   UpdateCoefficients();
1941 }
1942
1943 //=======================================================================
1944 //function : IsUClosed
1945 //purpose  :
1946 //=======================================================================
1947
1948 Standard_Boolean Geom_BezierSurface::IsUClosed () const
1949 {
1950   const TColgp_Array2OfPnt & Poles = poles->Array2();
1951   Standard_Boolean Closed = Standard_True;
1952   Standard_Integer Lower  = Poles.LowerRow();
1953   Standard_Integer Upper  = Poles.UpperRow();
1954   Standard_Integer Length = Poles.RowLength();
1955   Standard_Integer j      = Poles.LowerCol();
1956
1957   while (Closed && j <= Length) {
1958     Closed = (Poles (Lower,j).Distance (Poles (Upper,j)) <= Precision::Confusion());
1959     j++;
1960   }
1961   return Closed;
1962 }
1963
1964 //=======================================================================
1965 //function : IsVClosed
1966 //purpose  :
1967 //=======================================================================
1968
1969 Standard_Boolean Geom_BezierSurface::IsVClosed () const
1970 {
1971   const TColgp_Array2OfPnt & Poles = poles->Array2();
1972   Standard_Boolean Closed = Standard_True;
1973   Standard_Integer Lower  = Poles.LowerCol();
1974   Standard_Integer Upper  = Poles.UpperCol();
1975   Standard_Integer Length = Poles.ColLength();
1976   Standard_Integer i      = Poles.LowerRow();
1977   while (Closed && i <= Length) {
1978     Closed = (Poles (i,Lower).Distance (Poles (i,Upper)) <= Precision::Confusion());
1979     i++;
1980   }
1981   return Closed;
1982 }
1983
1984 //=======================================================================
1985 //function : IsUPeriodic
1986 //purpose  :
1987 //=======================================================================
1988
1989 Standard_Boolean Geom_BezierSurface::IsUPeriodic () const
1990 {
1991   return Standard_False;
1992 }
1993
1994 //=======================================================================
1995 //function : IsVPeriodic
1996 //purpose  :
1997 //=======================================================================
1998
1999 Standard_Boolean Geom_BezierSurface::IsVPeriodic () const
2000 {
2001   return Standard_False;
2002 }
2003
2004 //=======================================================================
2005 //function : Resolution
2006 //purpose  :
2007 //=======================================================================
2008
2009 void Geom_BezierSurface::Resolution(const Standard_Real  Tolerance3D,
2010                                     Standard_Real&       UTolerance,
2011                                     Standard_Real&       VTolerance)
2012 {
2013   if(!maxderivinvok){
2014     TColStd_Array1OfReal biduknots(1,2); biduknots(1) = 0.; biduknots(2) = 1.;
2015     TColStd_Array1OfInteger bidumults(1,2); bidumults.Init(UDegree() + 1);
2016     TColStd_Array1OfReal bidvknots(1,2); bidvknots(1) = 0.; bidvknots(2) = 1.;
2017     TColStd_Array1OfInteger bidvmults(1,2); bidvmults.Init(VDegree() + 1);
2018     if(urational || vrational){
2019       BSplSLib::Resolution(poles->Array2(),
2020                            weights->Array2(),
2021                            biduknots,
2022                            bidvknots,
2023                            bidumults,
2024                            bidvmults,
2025                            UDegree(),
2026                            VDegree(),
2027                            urational,
2028                            vrational,
2029                            0,
2030                            0,
2031                            1.,
2032                            umaxderivinv,
2033                            vmaxderivinv) ;
2034     }
2035     else{
2036       BSplSLib::Resolution(poles->Array2(),
2037                            *((TColStd_Array2OfReal*) NULL),
2038                            biduknots,
2039                            bidvknots,
2040                            bidumults,
2041                            bidvmults,
2042                            UDegree(),
2043                            VDegree(),
2044                            urational,
2045                            vrational,
2046                            0,
2047                            0,
2048                            1.,
2049                            umaxderivinv,
2050                            vmaxderivinv) ;
2051     }
2052     maxderivinvok = 1;
2053   }
2054   UTolerance = Tolerance3D * umaxderivinv;
2055   VTolerance = Tolerance3D * vmaxderivinv;
2056 }
2057
2058 //=======================================================================
2059 //function : Copy
2060 //purpose  :
2061 //=======================================================================
2062
2063 Handle(Geom_Geometry) Geom_BezierSurface::Copy() const
2064 {
2065   Handle(Geom_BezierSurface) S = new Geom_BezierSurface
2066     (poles, coeffs, weights, wcoeffs, urational, vrational);
2067   return S;
2068 }
2069
2070 //=======================================================================
2071 //function : Init
2072 //purpose  :
2073 //=======================================================================
2074
2075 void Geom_BezierSurface::Init
2076   (const Handle(TColgp_HArray2OfPnt)&   Poles,
2077    const Handle(TColStd_HArray2OfReal)& Weights)
2078 {
2079   Standard_Integer NbUPoles = Poles->ColLength();
2080   Standard_Integer NbVPoles = Poles->RowLength();
2081
2082   Standard_Integer maxcls = Max(NbUPoles, NbVPoles);
2083   Standard_Integer mincls = Min(NbUPoles, NbVPoles);
2084
2085   // set fields
2086   poles   = Poles;
2087   coeffs  = new TColgp_HArray2OfPnt  (1,maxcls,1,mincls);
2088
2089   if (urational || vrational) {
2090     weights = Weights;
2091     wcoeffs = new TColStd_HArray2OfReal (1,maxcls,1,mincls);
2092   }
2093   else {
2094     weights.Nullify();
2095     wcoeffs.Nullify();
2096   }
2097
2098   UpdateCoefficients();
2099 }
2100
2101
2102 //=======================================================================
2103 //function : UpdateCoefficients
2104 //purpose  :
2105 //=======================================================================
2106
2107 void  Geom_BezierSurface::UpdateCoefficients(const Standard_Real ,
2108                                              const Standard_Real )
2109 {
2110   maxderivinvok = Standard_False;
2111   ucacheparameter = 0.;
2112   TColStd_Array1OfReal biduflatknots(BSplCLib::FlatBezierKnots(UDegree()),
2113                                      1, 2*(UDegree()+1));
2114   vcacheparameter = 0.;
2115   TColStd_Array1OfReal bidvflatknots(BSplCLib::FlatBezierKnots(VDegree()),
2116                                      1, 2*(VDegree()+1));
2117
2118   Standard_Real uparameter_11 = (2*ucacheparameter + ucachespanlenght)/2,
2119                 uspanlenght_11 = ucachespanlenght/2,
2120                 vparameter_11 = (2*vcacheparameter + vcachespanlenght)/2,
2121                 vspanlenght_11 = vcachespanlenght/2 ;
2122
2123   if ( urational || vrational ) {
2124     BSplSLib::BuildCache(uparameter_11,vparameter_11,
2125                          uspanlenght_11,vspanlenght_11,0,0,
2126                          UDegree(),VDegree(),0,0,
2127                          biduflatknots,bidvflatknots,
2128                          poles->Array2(),
2129                          weights->Array2(),
2130                          coeffs->ChangeArray2(),
2131                          wcoeffs->ChangeArray2());
2132   }
2133   else {
2134     BSplSLib::BuildCache(uparameter_11,vparameter_11,
2135                          uspanlenght_11,vspanlenght_11,0,0,
2136                          UDegree(),VDegree(),0,0,
2137                          biduflatknots,bidvflatknots,
2138                          poles->Array2(),
2139                          *((TColStd_Array2OfReal*) NULL),
2140                          coeffs->ChangeArray2(),
2141                          *((TColStd_Array2OfReal*) NULL));
2142   }
2143   validcache = 1;
2144 }
2145