0025706: SIGSEGV after making existing BSplineCurve rational
[occt.git] / src / Geom / Geom_BSplineSurface.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 // 14-Mar-96 : xab  portage hp
18 // pmn : 28-Jun-96 Distinction entre la continuite en U et V (bug PRO4625)
19 // pmn : 07-Jan-97 Centralisation des verif rational (PRO6834)
20 //       et ajout des InvalideCache() dans les SetWeight*(PRO6833)
21 // RBD : 15-10-98 ; Le cache est maintenant calcule sur [-1,1] (pro15537).
22 // jct : 19-01-99 ; permutation de urational et vrational dans Rational.
23 #define No_Standard_OutOfRange
24
25 #include <Geom_BSplineSurface.ixx>
26
27 #include <gp.hxx>
28 #include <BSplCLib.hxx>
29 #include <BSplSLib.hxx>
30
31 #include <Standard_ConstructionError.hxx>
32 #include <Standard_NotImplemented.hxx>
33 #include <Standard_OutOfRange.hxx>
34 #include <Precision.hxx>
35
36 //=======================================================================
37 //function : CheckSurfaceData
38 //purpose  : Internal use only.
39 //=======================================================================
40
41 static void CheckSurfaceData
42 (const TColgp_Array2OfPnt&      SPoles,
43  const TColStd_Array1OfReal&    SUKnots,
44  const TColStd_Array1OfReal&    SVKnots,
45  const TColStd_Array1OfInteger& SUMults,
46  const TColStd_Array1OfInteger& SVMults,
47  const Standard_Integer         UDegree,
48  const Standard_Integer         VDegree,
49  const Standard_Boolean         UPeriodic,
50  const Standard_Boolean         VPeriodic)
51 {
52   if (UDegree < 1 || UDegree > Geom_BSplineSurface::MaxDegree () || 
53       VDegree < 1 || VDegree > Geom_BSplineSurface::MaxDegree ()) {
54     Standard_ConstructionError::Raise("Geom_BSplineSurface");
55   }
56   if (SPoles.ColLength () < 2 || SPoles.RowLength () < 2) {
57     Standard_ConstructionError::Raise("Geom_BSplineSurface");
58   }
59
60   if (SUKnots.Length() != SUMults.Length() ||
61       SVKnots.Length() != SVMults.Length()) {
62     Standard_ConstructionError::Raise("Geom_BSplineSurface");
63   }
64
65   Standard_Integer i;
66   for (i = SUKnots.Lower(); i < SUKnots.Upper(); i++) {
67     if (SUKnots(i+1) - SUKnots(i) <= Epsilon(Abs(SUKnots(i)))) {
68       Standard_ConstructionError::Raise("Geom_BSplineSurface");
69     }
70   }
71
72   for (i = SVKnots.Lower(); i < SVKnots.Upper(); i++) {
73     if (SVKnots(i+1) - SVKnots(i) <= Epsilon(Abs(SVKnots(i)))) {
74       Standard_ConstructionError::Raise("Geom_BSplineSurface");
75     }
76   }
77   
78   if (SPoles.ColLength() != BSplCLib::NbPoles(UDegree,UPeriodic,SUMults))
79     Standard_ConstructionError::Raise("Geom_BSplineSurface");
80
81   if (SPoles.RowLength() != BSplCLib::NbPoles(VDegree,VPeriodic,SVMults))
82     Standard_ConstructionError::Raise("Geom_BSplineSurface");
83 }
84
85 //=======================================================================
86 //function : Rational
87 //purpose  : Internal use only.
88 //=======================================================================
89
90 static void Rational(const TColStd_Array2OfReal& Weights,
91                      Standard_Boolean& Urational,
92                      Standard_Boolean& Vrational)
93 {
94   Standard_Integer I,J;
95   J = Weights.LowerCol ();
96   Vrational = Standard_False;
97   while (!Vrational && J <= Weights.UpperCol()) {
98     I = Weights.LowerRow();
99     while (!Vrational && I <= Weights.UpperRow() - 1) {
100       Vrational = (Abs(Weights (I, J) - Weights (I+1, J)) 
101                    > Epsilon (Abs(Weights (I, J))));
102       I++;
103     }
104     J++;
105   }
106
107   I = Weights.LowerRow ();
108   Urational = Standard_False;
109   while (!Urational && I <= Weights.UpperRow()) {
110     J = Weights.LowerCol();
111     while (!Urational && J <= Weights.UpperCol() - 1) {
112       Urational = (Abs(Weights (I, J) - Weights (I, J+1))
113                    > Epsilon (Abs(Weights (I, J))));
114       J++;
115     }
116     I++;
117   }
118 }
119
120 //=======================================================================
121 //function : Copy
122 //purpose  : 
123 //=======================================================================
124
125 Handle(Geom_Geometry) Geom_BSplineSurface::Copy () const
126 {
127   Handle(Geom_BSplineSurface) S;
128   if (urational || vrational) 
129     S = new Geom_BSplineSurface (poles->Array2() , weights->Array2(), 
130                                  uknots->Array1(), vknots->Array1(), 
131                                  umults->Array1(), vmults->Array1(), 
132                                  udeg     , vdeg, 
133                                  uperiodic, vperiodic);
134   else
135     S = new Geom_BSplineSurface (poles->Array2(),
136                                  uknots->Array1(), vknots->Array1(), 
137                                  umults->Array1(), vmults->Array1(), 
138                                  udeg     , vdeg, 
139                                  uperiodic, vperiodic);
140   return S;
141 }
142
143 //=======================================================================
144 //function : Geom_BSplineSurface
145 //purpose  : 
146 //=======================================================================
147
148 Geom_BSplineSurface::Geom_BSplineSurface
149 (const TColgp_Array2OfPnt&      Poles, 
150  const TColStd_Array1OfReal&    UKnots, 
151  const TColStd_Array1OfReal&    VKnots,
152  const TColStd_Array1OfInteger& UMults, 
153  const TColStd_Array1OfInteger& VMults,
154  const Standard_Integer         UDegree, 
155  const Standard_Integer         VDegree,
156  const Standard_Boolean         UPeriodic,
157  const Standard_Boolean         VPeriodic
158  ) :
159  urational(Standard_False),
160  vrational(Standard_False),
161  uperiodic(UPeriodic),
162  vperiodic(VPeriodic),
163  udeg(UDegree),
164  vdeg(VDegree),
165  maxderivinvok(0)
166
167 {
168   Standard_Integer MinDegree,
169   MaxDegree ;
170
171   // check
172   
173   CheckSurfaceData(Poles,
174                    UKnots   , VKnots,
175                    UMults   , VMults,
176                    UDegree  , VDegree,
177                    UPeriodic, VPeriodic);
178
179   // copy arrays
180
181   poles   = new TColgp_HArray2OfPnt(1,Poles.ColLength(),
182                                     1,Poles.RowLength());
183   poles->ChangeArray2() = Poles;
184
185   weights = new TColStd_HArray2OfReal (1,Poles.ColLength(),
186                                        1,Poles.RowLength(), 1.0);
187
188   uknots  = new TColStd_HArray1OfReal    (1,UKnots.Length());
189   uknots->ChangeArray1() = UKnots;
190
191   umults  = new TColStd_HArray1OfInteger (1,UMults.Length());
192   umults->ChangeArray1() = UMults;
193
194   vknots  = new TColStd_HArray1OfReal    (1,VKnots.Length());
195   vknots->ChangeArray1() = VKnots;
196
197   vmults  = new TColStd_HArray1OfInteger (1,VMults.Length());
198   vmults->ChangeArray1() = VMults;
199   MinDegree = Min(udeg,vdeg) ;
200   MaxDegree = Max(udeg,vdeg) ;
201   cachepoles = new TColgp_HArray2OfPnt(1,MaxDegree + 1,
202                                        1,MinDegree + 1) ;
203   
204   cacheweights.Nullify() ;
205   ucacheparameter = 0.0e0 ;
206   vcacheparameter = 0.0e0 ;
207   ucachespanlenght = 1.0e0 ;
208   vcachespanlenght = 1.0e0 ;
209   ucachespanindex = 0 ;
210   vcachespanindex = 0 ;
211   validcache = 0 ;
212
213   UpdateUKnots();
214   UpdateVKnots();
215 }
216
217 //=======================================================================
218 //function : Geom_BSplineSurface
219 //purpose  : 
220 //=======================================================================
221
222 Geom_BSplineSurface::Geom_BSplineSurface
223 (const TColgp_Array2OfPnt&      Poles,
224  const TColStd_Array2OfReal&    Weights,
225  const TColStd_Array1OfReal&    UKnots,
226  const TColStd_Array1OfReal&    VKnots,
227  const TColStd_Array1OfInteger& UMults, 
228  const TColStd_Array1OfInteger& VMults,
229  const Standard_Integer         UDegree,
230  const Standard_Integer         VDegree,
231  const Standard_Boolean         UPeriodic,
232  const Standard_Boolean         VPeriodic) :
233  urational(Standard_False),
234  vrational(Standard_False),
235  uperiodic(UPeriodic),
236  vperiodic(VPeriodic),
237  udeg(UDegree),
238  vdeg(VDegree),
239  maxderivinvok(0)
240 {
241   Standard_Integer MinDegree,
242   MaxDegree ;
243   // check weights
244
245   if (Weights.ColLength() != Poles.ColLength())
246     Standard_ConstructionError::Raise("Geom_BSplineSurface");
247
248   if (Weights.RowLength() != Poles.RowLength())
249     Standard_ConstructionError::Raise("Geom_BSplineSurface");
250
251   Standard_Integer i,j;
252   for (i = Weights.LowerRow(); i <= Weights.UpperRow(); i++) {
253     for (j = Weights.LowerCol(); j <= Weights.UpperCol(); j++) {
254       if (Weights(i,j) <= gp::Resolution())  
255         Standard_ConstructionError::Raise("Geom_BSplineSurface");
256     }
257   }
258   
259   // check really rational
260   
261   Rational(Weights, urational, vrational);
262
263   // check
264   
265   CheckSurfaceData(Poles,
266                    UKnots   , VKnots,
267                    UMults   , VMults,
268                    UDegree  , VDegree,
269                    UPeriodic, VPeriodic);
270
271   // copy arrays
272
273   poles   = new TColgp_HArray2OfPnt(1,Poles.ColLength(),
274                                     1,Poles.RowLength());
275   poles->ChangeArray2() = Poles;
276
277   weights = new TColStd_HArray2OfReal (1,Poles.ColLength(),
278                                        1,Poles.RowLength());
279   weights->ChangeArray2() = Weights;
280
281   uknots  = new TColStd_HArray1OfReal    (1,UKnots.Length());
282   uknots->ChangeArray1() = UKnots;
283
284   umults  = new TColStd_HArray1OfInteger (1,UMults.Length());
285   umults->ChangeArray1() = UMults;
286
287   vknots  = new TColStd_HArray1OfReal    (1,VKnots.Length());
288   vknots->ChangeArray1() = VKnots;
289
290   vmults  = new TColStd_HArray1OfInteger (1,VMults.Length());
291   vmults->ChangeArray1() = VMults;
292   MinDegree = Min(udeg,vdeg) ;
293   MaxDegree = Max(udeg,vdeg) ;
294   cachepoles = new TColgp_HArray2OfPnt(1,MaxDegree + 1,
295                                        1,MinDegree + 1) ;
296   if (urational || vrational) {
297     cacheweights = new TColStd_HArray2OfReal (1,MaxDegree + 1,
298                                               1,MinDegree + 1);
299   }
300   ucacheparameter = 0.0e0 ;
301   vcacheparameter = 0.0e0 ;
302   ucachespanlenght = 1.0e0 ;
303   vcachespanlenght = 1.0e0 ;
304   ucachespanindex = 0 ;
305   vcachespanindex = 0 ;
306   validcache = 0 ;
307
308   UpdateUKnots();
309   UpdateVKnots();
310 }
311
312 //=======================================================================
313 //function : ExchangeUV
314 //purpose  : 
315 //=======================================================================
316
317 void Geom_BSplineSurface::ExchangeUV ()
318 {
319   Standard_Integer LC = poles->LowerCol();
320   Standard_Integer UC = poles->UpperCol();
321   Standard_Integer LR = poles->LowerRow();
322   Standard_Integer UR = poles->UpperRow();
323
324   Handle(TColgp_HArray2OfPnt)   npoles   = 
325     new TColgp_HArray2OfPnt  (LC, UC, LR, UR);
326   Handle(TColStd_HArray2OfReal) nweights = 
327     new TColStd_HArray2OfReal (LC, UC, LR, UR);
328
329   const TColgp_Array2OfPnt   & spoles   = poles->Array2();
330   const TColStd_Array2OfReal & sweights = weights->Array2();
331   
332   TColgp_Array2OfPnt&   snpoles   = npoles->ChangeArray2();
333   TColStd_Array2OfReal& snweights = nweights->ChangeArray2();
334   
335   Standard_Integer i, j;
336   for (i = LC; i <= UC; i++) {
337     for (j = LR; j <= UR; j++) {
338       snpoles   (i,j) = spoles   (j,i);
339       snweights (i,j) = sweights (j,i);
340     }
341   }
342
343   poles   = npoles;
344   weights = nweights;
345
346   Standard_Boolean temp = urational;
347   urational = vrational;
348   vrational = temp;
349
350   temp = uperiodic;
351   uperiodic = vperiodic;
352   vperiodic = temp;
353
354   Standard_Integer tempdeg = udeg;
355   udeg = vdeg;
356   vdeg = tempdeg;
357   
358
359   Handle(TColStd_HArray1OfReal) tempknots = uknots;
360   uknots = vknots;
361   vknots = tempknots;
362
363   Handle(TColStd_HArray1OfInteger) tempmults = umults;
364   umults = vmults;
365   vmults = tempmults;
366
367   UpdateUKnots();
368   UpdateVKnots();
369 }
370
371 //=======================================================================
372 //function : IncreaseDegree
373 //purpose  : 
374 //=======================================================================
375
376 void Geom_BSplineSurface::IncreaseDegree (const Standard_Integer UDegree,
377                                           const Standard_Integer VDegree)
378
379   if (UDegree != udeg) {
380     if ( UDegree < udeg || UDegree > Geom_BSplineSurface::MaxDegree())
381       Standard_ConstructionError::Raise();
382     
383     Standard_Integer FromK1 = FirstUKnotIndex();
384     Standard_Integer ToK2   = LastUKnotIndex();
385
386     Standard_Integer Step   = UDegree - udeg;
387
388     Handle(TColgp_HArray2OfPnt) npoles = new
389       TColgp_HArray2OfPnt( 1, poles->ColLength() + Step * (ToK2 - FromK1),
390                           1, poles->RowLength());
391
392     Standard_Integer nbknots = BSplCLib::IncreaseDegreeCountKnots
393       (udeg,UDegree,uperiodic,umults->Array1());
394
395     Handle(TColStd_HArray1OfReal) nknots = 
396       new TColStd_HArray1OfReal(1,nbknots);
397     
398     Handle(TColStd_HArray1OfInteger) nmults = 
399       new TColStd_HArray1OfInteger(1,nbknots);
400
401     Handle(TColStd_HArray2OfReal) nweights 
402       = new TColStd_HArray2OfReal(1,npoles->ColLength(),
403                                   1,npoles->RowLength(), 1.);
404
405     if (urational || vrational) {
406       
407       BSplSLib::IncreaseDegree
408         (Standard_True, udeg, UDegree, uperiodic,
409          poles->Array2(),weights->Array2(),
410          uknots->Array1(),umults->Array1(),
411          npoles->ChangeArray2(),nweights->ChangeArray2(),
412          nknots->ChangeArray1(),nmults->ChangeArray1());
413     }
414     else {
415
416       BSplSLib::IncreaseDegree
417         (Standard_True, udeg, UDegree, uperiodic,
418          poles->Array2(),BSplSLib::NoWeights(),
419          uknots->Array1(),umults->Array1(),
420          npoles->ChangeArray2(),*((TColStd_Array2OfReal*) NULL),
421          nknots->ChangeArray1(),nmults->ChangeArray1());
422     }
423     udeg    = UDegree;
424     poles   = npoles;
425     weights = nweights;
426     uknots  = nknots;
427     umults  = nmults;
428     UpdateUKnots();
429   }
430
431   if (VDegree != vdeg) {
432     if ( VDegree < vdeg || VDegree > Geom_BSplineSurface::MaxDegree())
433       Standard_ConstructionError::Raise();
434     
435     Standard_Integer FromK1 = FirstVKnotIndex();
436     Standard_Integer ToK2   = LastVKnotIndex();
437
438     Standard_Integer Step   = VDegree - vdeg;
439
440     Handle(TColgp_HArray2OfPnt) npoles = new
441       TColgp_HArray2OfPnt( 1, poles->ColLength(),
442                           1, poles->RowLength() + Step * (ToK2 - FromK1));
443
444     Standard_Integer nbknots = BSplCLib::IncreaseDegreeCountKnots
445       (vdeg,VDegree,vperiodic,vmults->Array1());
446
447     Handle(TColStd_HArray1OfReal) nknots = 
448       new TColStd_HArray1OfReal(1,nbknots);
449     
450     Handle(TColStd_HArray1OfInteger) nmults = 
451       new TColStd_HArray1OfInteger(1,nbknots);
452
453     Handle(TColStd_HArray2OfReal) nweights
454       = new TColStd_HArray2OfReal(1,npoles->ColLength(),
455                                   1,npoles->RowLength(), 1.);
456
457     if (urational || vrational) {
458       
459       BSplSLib::IncreaseDegree
460         (Standard_False, vdeg, VDegree, vperiodic,
461          poles->Array2(),weights->Array2(),
462          vknots->Array1(),vmults->Array1(),
463          npoles->ChangeArray2(),nweights->ChangeArray2(),
464          nknots->ChangeArray1(),nmults->ChangeArray1());
465     }
466     else {
467
468       BSplSLib::IncreaseDegree
469         (Standard_False, vdeg, VDegree, vperiodic,
470          poles->Array2(),BSplSLib::NoWeights(),
471          vknots->Array1(),vmults->Array1(),
472          npoles->ChangeArray2(),*((TColStd_Array2OfReal*) NULL),
473          nknots->ChangeArray1(),nmults->ChangeArray1());
474     }
475     vdeg    = VDegree;
476     poles   = npoles;
477     weights = nweights;
478     vknots  = nknots;
479     vmults  = nmults;
480     UpdateVKnots();
481   }
482 }
483
484 //=======================================================================
485 //function : IncreaseUMultiplicity
486 //purpose  : 
487 //=======================================================================
488
489 void Geom_BSplineSurface::IncreaseUMultiplicity
490 (const Standard_Integer UIndex, 
491  const Standard_Integer M)
492 {
493   TColStd_Array1OfReal k(1,1);
494   k(1) = uknots->Value(UIndex);
495   TColStd_Array1OfInteger m(1,1);
496   m(1) = M - umults->Value(UIndex);
497   InsertUKnots(k,m,Epsilon(1.),Standard_True);
498 }
499
500 //=======================================================================
501 //function : IncreaseUMultiplicity
502 //purpose  : 
503 //=======================================================================
504
505 void Geom_BSplineSurface::IncreaseUMultiplicity
506 (const Standard_Integer FromI1, 
507  const Standard_Integer ToI2,
508  const Standard_Integer M)
509 {
510   Handle(TColStd_HArray1OfReal) tk = uknots;
511   TColStd_Array1OfReal k((uknots->Array1())(FromI1),FromI1,ToI2);
512   TColStd_Array1OfInteger m(FromI1, ToI2);
513   for (Standard_Integer i = FromI1; i <= ToI2; i++) 
514     m(i) = M - umults->Value(i);
515   InsertUKnots(k,m,Epsilon(1.),Standard_True);
516 }
517
518 //=======================================================================
519 //function : IncreaseVMultiplicity
520 //purpose  : 
521 //=======================================================================
522
523 void Geom_BSplineSurface::IncreaseVMultiplicity
524 (const Standard_Integer VIndex, 
525  const Standard_Integer M)
526 {
527   TColStd_Array1OfReal k(1,1);
528   k(1) = vknots->Value(VIndex);
529   TColStd_Array1OfInteger m(1,1);
530   m(1) = M - vmults->Value(VIndex);
531   InsertVKnots(k,m,Epsilon(1.),Standard_True);
532 }
533
534 //=======================================================================
535 //function : IncreaseVMultiplicity
536 //purpose  : 
537 //=======================================================================
538
539 void Geom_BSplineSurface::IncreaseVMultiplicity
540 (const Standard_Integer FromI1,
541  const Standard_Integer ToI2,
542  const Standard_Integer M)
543 {
544   Handle(TColStd_HArray1OfReal) tk = vknots;
545   TColStd_Array1OfReal k((vknots->Array1())(FromI1),FromI1,ToI2);
546   TColStd_Array1OfInteger m(FromI1,ToI2);
547   for (Standard_Integer i = FromI1; i <= ToI2; i++)
548     m(i) = M - vmults->Value(i);
549   InsertVKnots(k,m,Epsilon(1.),Standard_True);
550 }
551
552 //=======================================================================
553 //function : Segment
554 //purpose  : 
555 //=======================================================================
556
557 void Geom_BSplineSurface::Segment(const Standard_Real U1, 
558                                   const Standard_Real U2,
559                                   const Standard_Real V1,
560                                   const Standard_Real V2) 
561 {
562
563   Standard_DomainError_Raise_if ( (U2 < U1) || (V2 < V1),
564                                  "Geom_BSplineCurve::Segment");
565   Standard_Real deltaU = Max(Abs(U2),Abs(U1));
566   Standard_Real EpsU = Epsilon(deltaU);
567   deltaU = U2 - U1;
568   
569   Standard_Real deltaV = Max(Abs(V2),Abs(V1));
570   Standard_Real EpsV = Epsilon(deltaV);
571   deltaV = V2 - V1;
572
573   Standard_Real NewU1, NewU2, NewV1, NewV2;
574   Standard_Real U,V;
575   Standard_Integer indexU, indexV;
576
577   // inserting the UKnots
578   TColStd_Array1OfReal    UKnots(1,2);
579   TColStd_Array1OfInteger UMults(1,2);
580
581   indexU = 0;
582   BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
583                             U1,uperiodic,uknots->Lower(),uknots->Upper(),
584                             indexU,NewU1);
585   indexU = 0;
586   BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
587                             U2,uperiodic,uknots->Lower(),uknots->Upper(),
588                             indexU,NewU2);
589   UKnots( 1) = Min( NewU1, NewU2);
590   UKnots( 2) = Max( NewU1, NewU2);
591   UMults( 1) = UMults( 2) = udeg;
592   InsertUKnots( UKnots, UMults, EpsU);
593
594   // Inserting the VKnots
595   TColStd_Array1OfReal    VKnots(1,2);
596   TColStd_Array1OfInteger VMults(1,2);
597
598   indexV = 0;
599   BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
600                             V1,vperiodic,vknots->Lower(),vknots->Upper(),
601                             indexV,NewV1);
602   indexV = 0;
603   BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
604                             V2,vperiodic,vknots->Lower(),vknots->Upper(),
605                             indexV,NewV2);
606   VKnots( 1) = Min( NewV1, NewV2);
607   VKnots( 2) = Max( NewV1, NewV2);
608   VMults( 1) = VMults( 2) = vdeg;
609   InsertVKnots( VKnots, VMults, EpsV);
610
611
612   if (uperiodic) { // set the origine at NewU1
613     Standard_Integer index = 0;
614     BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
615                               U1,uperiodic,uknots->Lower(),uknots->Upper(),
616                               index,U);
617     if ( Abs(uknots->Value(index+1)-U) <= EpsU)
618       index++;
619     SetUOrigin(index);
620     SetUNotPeriodic();
621   }
622
623   // compute index1 and index2 to set the new knots and mults 
624   Standard_Integer index1U = 0, index2U = 0;
625   Standard_Integer FromU1 = uknots->Lower();
626   Standard_Integer ToU2   = uknots->Upper();
627   BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
628                             NewU1,uperiodic,FromU1,ToU2,index1U,U);
629   BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
630                             NewU1 + deltaU,uperiodic,FromU1,ToU2,index2U,U);
631   if ( Abs(uknots->Value(index2U+1)-U) <= EpsU)
632     index2U++;
633   
634   Standard_Integer nbuknots = index2U - index1U + 1;
635
636   Handle(TColStd_HArray1OfReal) 
637     nuknots = new TColStd_HArray1OfReal(1,nbuknots);
638   Handle(TColStd_HArray1OfInteger) 
639     numults = new TColStd_HArray1OfInteger(1,nbuknots);
640
641   Standard_Integer i , k = 1;
642   for ( i = index1U; i<= index2U; i++) {
643     nuknots->SetValue(k, uknots->Value(i));
644     numults->SetValue(k, umults->Value(i));
645     k++;
646   }
647   numults->SetValue(       1, udeg + 1);
648   numults->SetValue(nbuknots, udeg + 1);
649
650
651   if (vperiodic) { // set the origine at NewV1
652     Standard_Integer index = 0;
653     BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
654                               V1,vperiodic,vknots->Lower(),vknots->Upper(),
655                               index,V);
656     if ( Abs(vknots->Value(index+1)-V) <= EpsV)
657       index++;
658     SetVOrigin(index);
659     SetVNotPeriodic();
660   }
661
662   // compute index1 and index2 to set the new knots and mults 
663   Standard_Integer index1V = 0, index2V = 0;
664   Standard_Integer FromV1 = vknots->Lower();
665   Standard_Integer ToV2   = vknots->Upper();
666   BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
667                             NewV1,vperiodic,FromV1,ToV2,index1V,V);
668   BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
669                             NewV1 + deltaV,vperiodic,FromV1,ToV2,index2V,V);
670   if ( Abs(vknots->Value(index2V+1)-V) <= EpsV)
671     index2V++;
672   
673   Standard_Integer nbvknots = index2V - index1V + 1;
674
675   Handle(TColStd_HArray1OfReal) 
676     nvknots = new TColStd_HArray1OfReal(1,nbvknots);
677   Handle(TColStd_HArray1OfInteger) 
678     nvmults = new TColStd_HArray1OfInteger(1,nbvknots);
679
680   k = 1;
681   for ( i = index1V; i<= index2V; i++) {
682     nvknots->SetValue(k, vknots->Value(i));
683     nvmults->SetValue(k, vmults->Value(i));
684     k++;
685   }
686   nvmults->SetValue(       1, vdeg + 1);
687   nvmults->SetValue(nbvknots, vdeg + 1);
688
689
690   // compute index1 and index2 to set the new poles and weights
691   Standard_Integer pindex1U 
692     = BSplCLib::PoleIndex(udeg,index1U,uperiodic,umults->Array1());
693   Standard_Integer pindex2U 
694     = BSplCLib::PoleIndex(udeg,index2U,uperiodic,umults->Array1());
695
696   pindex1U++;
697   pindex2U = Min( pindex2U+1, poles->ColLength());
698
699   Standard_Integer nbupoles  = pindex2U - pindex1U + 1;
700
701   // compute index1 and index2 to set the new poles and weights
702   Standard_Integer pindex1V 
703     = BSplCLib::PoleIndex(vdeg,index1V,vperiodic,vmults->Array1());
704   Standard_Integer pindex2V 
705     = BSplCLib::PoleIndex(vdeg,index2V,vperiodic,vmults->Array1());
706
707   pindex1V++;
708   pindex2V = Min( pindex2V+1, poles->RowLength());
709
710   Standard_Integer nbvpoles  = pindex2V - pindex1V + 1;
711
712
713   Handle(TColStd_HArray2OfReal) nweights; 
714
715   Handle(TColgp_HArray2OfPnt)
716     npoles = new TColgp_HArray2OfPnt(1,nbupoles,1,nbvpoles);
717
718   k = 1;
719   Standard_Integer j, l;
720   if ( urational || vrational) {
721     nweights = new TColStd_HArray2OfReal( 1,nbupoles,1,nbvpoles);
722     for ( i = pindex1U; i <= pindex2U; i++) {
723       l = 1;
724       for ( j = pindex1V; j <= pindex2V; j++) {
725         npoles->SetValue(k,l, poles->Value(i,j));
726         nweights->SetValue(k,l, weights->Value(i,j));
727         l++;
728       }
729       k++;
730     }
731   }
732   else {
733     for ( i = pindex1U; i <= pindex2U; i++) {
734       l = 1;
735       for ( j = pindex1V; j <= pindex2V; j++) {
736         npoles->SetValue(k,l, poles->Value(i,j));
737         l++;
738       }
739       k++;
740     }
741   }
742
743   uknots = nuknots;
744   umults = numults;
745   vknots = nvknots;
746   vmults = nvmults;
747   poles = npoles;
748   if ( urational || vrational) 
749     weights = nweights;
750   else 
751     weights = new TColStd_HArray2OfReal (1,poles->ColLength(),
752                                          1,poles->RowLength(), 1.0);
753
754   maxderivinvok = 0;
755   UpdateUKnots();
756   UpdateVKnots();
757
758 }
759
760 //=======================================================================
761 //function : CheckAndSegment
762 //purpose  : 
763 //=======================================================================
764
765 void Geom_BSplineSurface::CheckAndSegment(const Standard_Real U1, 
766                                           const Standard_Real U2,
767                                           const Standard_Real V1,
768                                           const Standard_Real V2) 
769 {
770
771   Standard_DomainError_Raise_if ( (U2 < U1) || (V2 < V1),
772                                  "Geom_BSplineCurve::Segment");
773   Standard_Real deltaU = Max(Abs(U2),Abs(U1));
774   Standard_Real EpsU = Epsilon(deltaU);
775   deltaU = U2 - U1;
776   
777   Standard_Real deltaV = Max(Abs(V2),Abs(V1));
778   Standard_Real EpsV = Epsilon(deltaV);
779   deltaV = V2 - V1;
780
781   Standard_Real NewU1, NewU2, NewV1, NewV2;
782   Standard_Real U,V;
783   Standard_Integer indexU, indexV;
784
785   Standard_Boolean segment_in_U = Standard_True;
786   Standard_Boolean segment_in_V = Standard_True;
787   segment_in_U = ( Abs(U1 - uknots->Value(uknots->Lower())) > EpsU )
788                         || ( Abs(U2 - uknots->Value(uknots->Upper())) > EpsU );
789   segment_in_V = ( Abs(V1 - vknots->Value(vknots->Lower())) > EpsV )
790                         || ( Abs(V2 - vknots->Value(vknots->Upper())) > EpsV );
791
792   indexU = 0;
793   BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
794                             U1,uperiodic,uknots->Lower(),uknots->Upper(),
795                             indexU,NewU1);
796   indexU = 0;
797   BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
798                             U2,uperiodic,uknots->Lower(),uknots->Upper(),
799                             indexU,NewU2);
800   if (segment_in_U) {
801     // inserting the UKnots
802     TColStd_Array1OfReal    UKnots(1,2);
803     TColStd_Array1OfInteger UMults(1,2);
804     UKnots( 1) = Min( NewU1, NewU2);
805     UKnots( 2) = Max( NewU1, NewU2);
806     UMults( 1) = UMults( 2) = udeg;
807     
808     InsertUKnots( UKnots, UMults, EpsU);
809   }
810
811   indexV = 0;
812   BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
813                             V1,vperiodic,vknots->Lower(),vknots->Upper(),
814                             indexV,NewV1);
815   indexV = 0;
816   BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
817                             V2,vperiodic,vknots->Lower(),vknots->Upper(),
818                             indexV,NewV2);
819   if (segment_in_V) {
820     // Inserting the VKnots
821     TColStd_Array1OfReal    VKnots(1,2);
822     TColStd_Array1OfInteger VMults(1,2);
823     
824     VKnots( 1) = Min( NewV1, NewV2);
825     VKnots( 2) = Max( NewV1, NewV2);
826     VMults( 1) = VMults( 2) = vdeg;
827     InsertVKnots( VKnots, VMults, EpsV);
828   }
829
830   if (uperiodic && segment_in_U) { // set the origine at NewU1
831     Standard_Integer index = 0;
832     BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
833                               U1,uperiodic,uknots->Lower(),uknots->Upper(),
834                               index,U);
835     if ( Abs(uknots->Value(index+1)-U) <= EpsU)
836       index++;
837     SetUOrigin(index);
838     SetUNotPeriodic();
839   }
840
841   // compute index1 and index2 to set the new knots and mults 
842   Standard_Integer index1U = 0, index2U = 0;
843   Standard_Integer FromU1 = uknots->Lower();
844   Standard_Integer ToU2   = uknots->Upper();
845   BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
846                             NewU1,uperiodic,FromU1,ToU2,index1U,U);
847   BSplCLib::LocateParameter(udeg,uknots->Array1(),umults->Array1(),
848                             NewU1 + deltaU,uperiodic,FromU1,ToU2,index2U,U);
849   if ( Abs(uknots->Value(index2U+1)-U) <= EpsU)
850     index2U++;
851   
852   Standard_Integer nbuknots = index2U - index1U + 1;
853
854   Handle(TColStd_HArray1OfReal) 
855     nuknots = new TColStd_HArray1OfReal(1,nbuknots);
856   Handle(TColStd_HArray1OfInteger) 
857     numults = new TColStd_HArray1OfInteger(1,nbuknots);
858
859   Standard_Integer i , k = 1;
860   for ( i = index1U; i<= index2U; i++) {
861     nuknots->SetValue(k, uknots->Value(i));
862     numults->SetValue(k, umults->Value(i));
863     k++;
864   }
865   if (segment_in_U) {
866     numults->SetValue(       1, udeg + 1);
867     numults->SetValue(nbuknots, udeg + 1);
868   }
869
870   if (vperiodic&& segment_in_V) { // set the origine at NewV1
871     Standard_Integer index = 0;
872     BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
873                               V1,vperiodic,vknots->Lower(),vknots->Upper(),
874                               index,V);
875     if ( Abs(vknots->Value(index+1)-V) <= EpsV)
876       index++;
877     SetVOrigin(index);
878     SetVNotPeriodic();
879   }
880
881   // compute index1 and index2 to set the new knots and mults 
882   Standard_Integer index1V = 0, index2V = 0;
883   Standard_Integer FromV1 = vknots->Lower();
884   Standard_Integer ToV2   = vknots->Upper();
885   BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
886                             NewV1,vperiodic,FromV1,ToV2,index1V,V);
887   BSplCLib::LocateParameter(vdeg,vknots->Array1(),vmults->Array1(),
888                             NewV1 + deltaV,vperiodic,FromV1,ToV2,index2V,V);
889   if ( Abs(vknots->Value(index2V+1)-V) <= EpsV)
890     index2V++;
891   
892   Standard_Integer nbvknots = index2V - index1V + 1;
893
894   Handle(TColStd_HArray1OfReal) 
895     nvknots = new TColStd_HArray1OfReal(1,nbvknots);
896   Handle(TColStd_HArray1OfInteger) 
897     nvmults = new TColStd_HArray1OfInteger(1,nbvknots);
898
899   k = 1;
900   for ( i = index1V; i<= index2V; i++) {
901     nvknots->SetValue(k, vknots->Value(i));
902     nvmults->SetValue(k, vmults->Value(i));
903     k++;
904   }
905   if (segment_in_V) {
906     nvmults->SetValue(       1, vdeg + 1);
907     nvmults->SetValue(nbvknots, vdeg + 1);
908   }
909
910   // compute index1 and index2 to set the new poles and weights
911   Standard_Integer pindex1U 
912     = BSplCLib::PoleIndex(udeg,index1U,uperiodic,umults->Array1());
913   Standard_Integer pindex2U 
914     = BSplCLib::PoleIndex(udeg,index2U,uperiodic,umults->Array1());
915
916   pindex1U++;
917   pindex2U = Min( pindex2U+1, poles->ColLength());
918
919   Standard_Integer nbupoles  = pindex2U - pindex1U + 1;
920
921   // compute index1 and index2 to set the new poles and weights
922   Standard_Integer pindex1V 
923     = BSplCLib::PoleIndex(vdeg,index1V,vperiodic,vmults->Array1());
924   Standard_Integer pindex2V 
925     = BSplCLib::PoleIndex(vdeg,index2V,vperiodic,vmults->Array1());
926
927   pindex1V++;
928   pindex2V = Min( pindex2V+1, poles->RowLength());
929
930   Standard_Integer nbvpoles  = pindex2V - pindex1V + 1;
931
932
933   Handle(TColStd_HArray2OfReal) nweights; 
934
935   Handle(TColgp_HArray2OfPnt)
936     npoles = new TColgp_HArray2OfPnt(1,nbupoles,1,nbvpoles);
937
938   k = 1;
939   Standard_Integer j, l;
940   if ( urational || vrational) {
941     nweights = new TColStd_HArray2OfReal( 1,nbupoles,1,nbvpoles);
942     for ( i = pindex1U; i <= pindex2U; i++) {
943       l = 1;
944       for ( j = pindex1V; j <= pindex2V; j++) {
945         npoles->SetValue(k,l, poles->Value(i,j));
946         nweights->SetValue(k,l, weights->Value(i,j));
947         l++;
948       }
949       k++;
950     }
951   }
952   else {
953     for ( i = pindex1U; i <= pindex2U; i++) {
954       l = 1;
955       for ( j = pindex1V; j <= pindex2V; j++) {
956         npoles->SetValue(k,l, poles->Value(i,j));
957         l++;
958       }
959       k++;
960     }
961   }
962
963   uknots = nuknots;
964   umults = numults;
965   vknots = nvknots;
966   vmults = nvmults;
967   poles = npoles;
968   if ( urational || vrational) 
969     weights = nweights;
970   else 
971     weights = new TColStd_HArray2OfReal (1,poles->ColLength(),
972                                          1,poles->RowLength(), 1.0);
973
974   maxderivinvok = 0;
975   UpdateUKnots();
976   UpdateVKnots();
977
978 }
979
980 //=======================================================================
981 //function : SetUKnot
982 //purpose  : 
983 //=======================================================================
984
985 void Geom_BSplineSurface::SetUKnot
986 (const Standard_Integer UIndex,
987  const Standard_Real    K      )
988 {
989   if (UIndex < 1 || UIndex > uknots->Length()) Standard_OutOfRange::Raise();
990
991   Standard_Integer NewIndex = UIndex;
992   Standard_Real DU = Abs(Epsilon (K));
993   if (UIndex == 1) {
994     if (K >= uknots->Value (2) - DU) Standard_ConstructionError::Raise();
995   }
996   else if (UIndex == uknots->Length()) {
997     if (K <= uknots->Value (uknots->Length()-1) + DU)  {
998       Standard_ConstructionError::Raise();
999     }
1000   }
1001   else {
1002     if (K <= uknots->Value (NewIndex-1) + DU || 
1003         K >= uknots->Value (NewIndex+1) - DU ) { 
1004       Standard_ConstructionError::Raise();
1005     } 
1006   }
1007   
1008   if (K != uknots->Value (NewIndex)) {
1009     uknots->SetValue (NewIndex, K);
1010     maxderivinvok = 0;
1011     UpdateUKnots();
1012   }
1013 }
1014
1015 //=======================================================================
1016 //function : SetUKnots
1017 //purpose  : 
1018 //=======================================================================
1019
1020 void Geom_BSplineSurface::SetUKnots (const TColStd_Array1OfReal& UK) {
1021
1022   Standard_Integer Lower = UK.Lower();
1023   Standard_Integer Upper = UK.Upper();
1024   if (Lower < 1 || Lower > uknots->Length() ||
1025       Upper < 1 || Upper > uknots->Length() ) {
1026     Standard_OutOfRange::Raise();
1027   }
1028   if (Lower > 1) {
1029     if (Abs (UK (Lower) - uknots->Value (Lower-1)) <= gp::Resolution()) {
1030       Standard_ConstructionError::Raise();
1031     }
1032   }
1033   if (Upper < uknots->Length ()) {
1034     if (Abs (UK (Upper) - uknots->Value (Upper+1)) <= gp::Resolution()) {
1035       Standard_ConstructionError::Raise();
1036     }
1037   }
1038   Standard_Real K1 = UK (Lower);
1039   for (Standard_Integer i = Lower; i <= Upper; i++) {
1040     uknots->SetValue (i, UK(i));
1041     if (i != Lower) {
1042       if (Abs (UK(i) - K1) <= gp::Resolution()) {
1043         Standard_ConstructionError::Raise();
1044       }
1045       K1 = UK (i);
1046     }
1047   }
1048
1049   maxderivinvok = 0;
1050   UpdateUKnots();
1051 }
1052
1053 //=======================================================================
1054 //function : SetUKnot
1055 //purpose  : 
1056 //=======================================================================
1057
1058 void Geom_BSplineSurface::SetUKnot
1059 (const Standard_Integer UIndex,
1060  const Standard_Real    K,
1061  const Standard_Integer M)
1062 {
1063   IncreaseUMultiplicity (UIndex, M);
1064   SetUKnot (UIndex, K);
1065 }
1066
1067 //=======================================================================
1068 //function : SetVKnot
1069 //purpose  : 
1070 //=======================================================================
1071
1072 void Geom_BSplineSurface::SetVKnot
1073 (const Standard_Integer VIndex,
1074  const Standard_Real    K)
1075 {
1076   if (VIndex < 1 || VIndex > vknots->Length())  Standard_OutOfRange::Raise();
1077   Standard_Integer NewIndex = VIndex + vknots->Lower() - 1;
1078   Standard_Real DV = Abs(Epsilon (K));
1079   if (VIndex == 1) {
1080     if (K >=  vknots->Value (2) - DV) {
1081       Standard_ConstructionError::Raise();
1082     }
1083   }
1084   else if (VIndex == vknots->Length()) {
1085     if (K <= vknots->Value (vknots->Length()-1) + DV)  {
1086       Standard_ConstructionError::Raise();
1087     }
1088   }
1089   else {
1090     if (K <= vknots->Value (NewIndex-1) + DV || 
1091         K >= vknots->Value (NewIndex+1) - DV ) { 
1092       Standard_ConstructionError::Raise();
1093     } 
1094   }
1095   
1096   maxderivinvok = 0;
1097   UpdateVKnots();
1098 }
1099
1100 //=======================================================================
1101 //function : SetVKnots
1102 //purpose  : 
1103 //=======================================================================
1104
1105 void Geom_BSplineSurface::SetVKnots (const TColStd_Array1OfReal& VK) {
1106
1107   Standard_Integer Lower = VK.Lower();
1108   Standard_Integer Upper = VK.Upper();
1109   if (Lower < 1 || Lower > vknots->Length() ||
1110       Upper < 1 || Upper > vknots->Length() ) {
1111     Standard_OutOfRange::Raise();
1112   }
1113   if (Lower > 1) {
1114     if (Abs (VK (Lower) - vknots->Value (Lower-1)) <= gp::Resolution()) {
1115       Standard_ConstructionError::Raise();
1116     }
1117   }
1118   if (Upper < vknots->Length ()) {
1119     if (Abs (VK (Upper) - vknots->Value (Upper+1)) <= gp::Resolution()) {
1120       Standard_ConstructionError::Raise();
1121     }
1122   }
1123   Standard_Real K1 = VK (Lower);
1124   for (Standard_Integer i = Lower; i <= Upper; i++) {
1125     vknots->SetValue (i, VK(i));
1126     if (i != Lower) {
1127       if (Abs (VK(i) - K1) <= gp::Resolution()) {
1128         Standard_ConstructionError::Raise();
1129       }
1130       K1 = VK (i);
1131     }
1132   }
1133
1134   maxderivinvok = 0;
1135   UpdateVKnots();
1136 }
1137
1138 //=======================================================================
1139 //function : SetVKnot
1140 //purpose  : 
1141 //=======================================================================
1142
1143 void Geom_BSplineSurface::SetVKnot
1144 (const Standard_Integer VIndex,
1145  const Standard_Real    K,
1146  const Standard_Integer M)
1147 {
1148   IncreaseVMultiplicity (VIndex, M);
1149   SetVKnot (VIndex, K);
1150 }
1151
1152 //=======================================================================
1153 //function : InsertUKnot
1154 //purpose  : 
1155 //=======================================================================
1156
1157 void Geom_BSplineSurface::InsertUKnot
1158 (const Standard_Real    U,
1159  const Standard_Integer M,
1160  const Standard_Real    ParametricTolerance,
1161  const Standard_Boolean Add)
1162 {
1163   TColStd_Array1OfReal k(1,1);
1164   k(1) = U;
1165   TColStd_Array1OfInteger m(1,1);
1166   m(1) = M;
1167   InsertUKnots(k,m,ParametricTolerance,Add);
1168 }
1169
1170 //=======================================================================
1171 //function : InsertVKnot
1172 //purpose  : 
1173 //=======================================================================
1174
1175 void Geom_BSplineSurface::InsertVKnot
1176 (const Standard_Real    V,
1177  const Standard_Integer M,
1178  const Standard_Real    ParametricTolerance,
1179  const Standard_Boolean Add)
1180 {
1181   TColStd_Array1OfReal k(1,1);
1182   k(1) = V;
1183   TColStd_Array1OfInteger m(1,1);
1184   m(1) = M;
1185   InsertVKnots(k,m,ParametricTolerance,Add);
1186 }
1187
1188 //=======================================================================
1189 //function : IncrementUMultiplicity
1190 //purpose  : 
1191 //=======================================================================
1192
1193 void  Geom_BSplineSurface::IncrementUMultiplicity
1194 (const Standard_Integer FromI1,
1195  const Standard_Integer ToI2,
1196  const Standard_Integer Step)
1197 {
1198   Handle(TColStd_HArray1OfReal) tk = uknots;
1199   TColStd_Array1OfReal k( (uknots->Array1())(FromI1), FromI1, ToI2);
1200   TColStd_Array1OfInteger m( FromI1, ToI2) ;
1201   m.Init(Step);
1202   InsertUKnots( k, m, Epsilon(1.));
1203 }
1204
1205 //=======================================================================
1206 //function : IncrementVMultiplicity
1207 //purpose  : 
1208 //=======================================================================
1209
1210 void  Geom_BSplineSurface::IncrementVMultiplicity
1211 (const Standard_Integer FromI1,
1212  const Standard_Integer ToI2,
1213  const Standard_Integer Step)
1214 {
1215   Handle(TColStd_HArray1OfReal) tk = vknots;
1216   TColStd_Array1OfReal k( (vknots->Array1())(FromI1), FromI1, ToI2);
1217
1218   TColStd_Array1OfInteger m( FromI1, ToI2) ;
1219   m.Init(Step);
1220   
1221   InsertVKnots( k, m, Epsilon(1.));
1222 }
1223
1224 //=======================================================================
1225 //function : UpdateUKnots
1226 //purpose  : 
1227 //=======================================================================
1228
1229 void Geom_BSplineSurface::UpdateUKnots()
1230 {
1231
1232   Standard_Integer MaxKnotMult = 0;
1233   BSplCLib::KnotAnalysis (udeg, uperiodic,
1234                 uknots->Array1(), 
1235                 umults->Array1(), 
1236                 uknotSet, MaxKnotMult);
1237   
1238   if (uknotSet == GeomAbs_Uniform && !uperiodic)  {
1239     ufknots = uknots;
1240   }
1241   else {
1242     ufknots = new TColStd_HArray1OfReal 
1243       (1, BSplCLib::KnotSequenceLength(umults->Array1(),udeg,uperiodic));
1244
1245     BSplCLib::KnotSequence (uknots->Array1(), 
1246                             umults->Array1(),
1247                             udeg,uperiodic,
1248                             ufknots->ChangeArray1());
1249   }
1250   
1251   if (MaxKnotMult == 0)  Usmooth = GeomAbs_CN;
1252   else {
1253     switch (udeg - MaxKnotMult) {
1254     case 0 :   Usmooth = GeomAbs_C0;   break;
1255     case 1 :   Usmooth = GeomAbs_C1;   break;
1256     case 2 :   Usmooth = GeomAbs_C2;   break;
1257     case 3 :   Usmooth = GeomAbs_C3;   break;
1258       default :  Usmooth = GeomAbs_C3;   break;
1259     }
1260   }
1261
1262   InvalidateCache() ;
1263 }
1264
1265 //=======================================================================
1266 //function : UpdateVKnots
1267 //purpose  : 
1268 //=======================================================================
1269
1270 void Geom_BSplineSurface::UpdateVKnots()
1271 {
1272   Standard_Integer MaxKnotMult = 0;
1273   BSplCLib::KnotAnalysis (vdeg, vperiodic,
1274                 vknots->Array1(), 
1275                 vmults->Array1(), 
1276                 vknotSet, MaxKnotMult);
1277   
1278   if (vknotSet == GeomAbs_Uniform && !vperiodic)  {
1279     vfknots = vknots;
1280   }
1281   else {
1282     vfknots = new TColStd_HArray1OfReal 
1283       (1, BSplCLib::KnotSequenceLength(vmults->Array1(),vdeg,vperiodic));
1284
1285     BSplCLib::KnotSequence (vknots->Array1(), 
1286                             vmults->Array1(),
1287                             vdeg,vperiodic,
1288                             vfknots->ChangeArray1());
1289   }
1290   
1291   if (MaxKnotMult == 0)  Vsmooth = GeomAbs_CN;
1292   else {
1293     switch (vdeg - MaxKnotMult) {
1294     case 0 :   Vsmooth = GeomAbs_C0;   break;
1295     case 1 :   Vsmooth = GeomAbs_C1;   break;
1296     case 2 :   Vsmooth = GeomAbs_C2;   break;
1297     case 3 :   Vsmooth = GeomAbs_C3;   break;
1298       default :  Vsmooth = GeomAbs_C3;   break;
1299     }
1300   }
1301   InvalidateCache() ;
1302 }
1303
1304 //=======================================================================
1305 //function :  InvalidateCache
1306 //purpose  : Invalidates the Cache of the surface
1307 //=======================================================================
1308
1309 void Geom_BSplineSurface::InvalidateCache()
1310 {
1311   validcache = 0 ;
1312 }
1313
1314 //=======================================================================
1315 //function : Normalizes the parameters if the curve is periodic
1316 //purpose  : that is compute the cache so that it is valid
1317 //=======================================================================
1318
1319 void Geom_BSplineSurface::PeriodicNormalization
1320 (Standard_Real&  Uparameter, 
1321  Standard_Real&  Vparameter) const 
1322 {
1323   Standard_Real Period, aMaxVal, aMinVal;
1324   
1325   if (uperiodic) {
1326     aMaxVal = ufknots->Value(ufknots->Upper() - udeg);
1327     aMinVal = ufknots->Value (udeg + 1);
1328     Standard_Real eps = Abs(Epsilon(Uparameter));
1329     Period =  aMaxVal - aMinVal;
1330
1331     if(Period <= eps) 
1332       Standard_OutOfRange::Raise("Geom_BSplineSurface::PeriodicNormalization: Uparameter is too great number");
1333
1334     Standard_Boolean isLess, isGreater;
1335     isLess = aMinVal - Uparameter > 0;
1336     isGreater = Uparameter - aMaxVal > 0;
1337     if (isLess || isGreater) {
1338       Standard_Real aDPar, aNbPer;
1339       aDPar = (isLess) ? (aMaxVal - Uparameter) : (aMinVal - Uparameter);
1340       modf(aDPar / Period, &aNbPer);
1341       Uparameter += aNbPer * Period;
1342     }
1343   }
1344   if (vperiodic) {
1345     aMaxVal = vfknots->Value(vfknots->Upper() - vdeg);
1346     aMinVal = vfknots->Value (vdeg + 1);
1347     Standard_Real eps = Abs(Epsilon(Vparameter));
1348     Period = aMaxVal - aMinVal;
1349
1350     if(Period <= eps) 
1351       Standard_OutOfRange::Raise("Geom_BSplineSurface::PeriodicNormalization: Vparameter is too great number");
1352
1353     Standard_Boolean isLess, isGreater;
1354     isLess = aMinVal - Vparameter > 0;
1355     isGreater = Vparameter - aMaxVal > 0;
1356     if (isLess || isGreater) {
1357       Standard_Real aDPar, aNbPer;
1358       aDPar = (isLess) ? (aMaxVal - Vparameter) : (aMinVal - Vparameter);
1359       modf(aDPar / Period, &aNbPer);
1360       Vparameter += aNbPer * Period;
1361     }
1362   }
1363 }
1364
1365 //=======================================================================
1366 //function : ValidateCache
1367 //purpose  : function that validates the cache of the surface
1368 //=======================================================================
1369
1370 void Geom_BSplineSurface::ValidateCache(const Standard_Real  Uparameter,
1371                                         const Standard_Real  Vparameter) 
1372 {
1373   Standard_Real NewParameter  ;
1374   Standard_Integer LocalIndex = 0  ;
1375   Standard_Integer MinDegree,
1376   MaxDegree ;
1377   //
1378   // check if the degree did not change
1379   //
1380
1381   MinDegree = Min(udeg,vdeg) ;
1382   MaxDegree = Max(udeg,vdeg) ;
1383   if (cachepoles->ColLength() < MaxDegree + 1 ||
1384       cachepoles->RowLength() < MinDegree + 1) {
1385     cachepoles = new TColgp_HArray2OfPnt(1,MaxDegree + 1,
1386                                          1,MinDegree + 1);
1387   }
1388   //
1389   // Verif + poussee pour les poids
1390   //
1391   if (urational || vrational) {
1392     if (cacheweights.IsNull()) {
1393       cacheweights  = new TColStd_HArray2OfReal(1,MaxDegree + 1,
1394                                                 1,MinDegree + 1);
1395     }
1396     else {
1397       if (cacheweights->ColLength() < MaxDegree + 1 ||
1398           cacheweights->RowLength() < MinDegree + 1) {
1399         cacheweights  = new TColStd_HArray2OfReal(1,MaxDegree + 1,
1400                                                   1,MinDegree + 1);
1401       }
1402     }
1403   }
1404   else if (!cacheweights.IsNull())
1405     cacheweights.Nullify();
1406
1407   BSplCLib::LocateParameter(udeg,
1408                             (ufknots->Array1()),
1409                             (BSplCLib::NoMults()),
1410                             Uparameter,
1411                             uperiodic,
1412                             LocalIndex,
1413                             NewParameter);
1414   ucachespanindex = LocalIndex ;
1415   if (Uparameter == ufknots->Value(LocalIndex + 1)) {
1416     
1417     LocalIndex += 1 ;
1418     ucacheparameter = ufknots->Value(LocalIndex) ;
1419     if (LocalIndex == ufknots->Upper() - udeg) {
1420       //
1421       // for the last span if the parameter is outside of 
1422       // the domain of the curve than use the last knot
1423       // and normalize with the last span Still set the
1424       // cachespanindex to flatknots->Upper() - deg so that
1425       // the IsCacheValid will know for sure we are extending
1426       // the Bspline 
1427       //
1428       
1429       ucachespanlenght = ufknots->Value(LocalIndex - 1) - ucacheparameter ;
1430     }
1431     else {
1432         ucachespanlenght = ufknots->Value(LocalIndex + 1) - ucacheparameter ;
1433       }
1434   }
1435   else {
1436       ucacheparameter = ufknots->Value(LocalIndex) ;
1437       ucachespanlenght = ufknots->Value(LocalIndex + 1) - ucacheparameter ;
1438     }
1439
1440   LocalIndex = 0 ;   
1441   BSplCLib::LocateParameter(vdeg,
1442                             (vfknots->Array1()),
1443                             (BSplCLib::NoMults()),
1444                             Vparameter,
1445                             vperiodic,
1446                             LocalIndex,
1447                             NewParameter);
1448   vcachespanindex = LocalIndex ;
1449   if (Vparameter == vfknots->Value(LocalIndex + 1)) {
1450     LocalIndex += 1 ;
1451     vcacheparameter = vfknots->Value(LocalIndex) ;
1452     if (LocalIndex == vfknots->Upper() - vdeg) {
1453       //
1454       // for the last span if the parameter is outside of 
1455       // the domain of the curve than use the last knot
1456       // and normalize with the last span Still set the
1457       // cachespanindex to flatknots->Upper() - deg so that
1458       // the IsCacheValid will know for sure we are extending
1459       // the Bspline 
1460       //
1461       
1462       vcachespanlenght = vfknots->Value(LocalIndex - 1) - vcacheparameter ;
1463     }
1464     else {
1465       vcachespanlenght = vfknots->Value(LocalIndex + 1) - vcacheparameter ;
1466     }
1467   }
1468   else {
1469     vcacheparameter = vfknots->Value(LocalIndex) ;
1470     vcachespanlenght = vfknots->Value(LocalIndex + 1) - vcacheparameter ;
1471   }
1472   
1473   Standard_Real uparameter_11 = (2*ucacheparameter + ucachespanlenght)/2,
1474                 uspanlenght_11 = ucachespanlenght/2,
1475                 vparameter_11 = (2*vcacheparameter + vcachespanlenght)/2,
1476                 vspanlenght_11 = vcachespanlenght/2 ;
1477   if (urational || vrational) {
1478     BSplSLib::BuildCache(uparameter_11,
1479                          vparameter_11,
1480                          uspanlenght_11,
1481                          vspanlenght_11,
1482                          uperiodic,
1483                          vperiodic,
1484                          udeg,
1485                          vdeg,
1486                          ucachespanindex,
1487                          vcachespanindex,
1488                          (ufknots->Array1()),
1489                          (vfknots->Array1()),
1490                          poles->Array2(),
1491                          weights->Array2(),
1492                          cachepoles->ChangeArray2(),
1493                          cacheweights->ChangeArray2()) ;
1494   }
1495   else {
1496     BSplSLib::BuildCache(uparameter_11,
1497                          vparameter_11,
1498                          uspanlenght_11,
1499                          vspanlenght_11,
1500                          uperiodic,
1501                          vperiodic,
1502                          udeg,
1503                          vdeg,
1504                          ucachespanindex,
1505                          vcachespanindex,
1506                          (ufknots->Array1()),
1507                          (vfknots->Array1()),
1508                          poles->Array2(),
1509                          *((TColStd_Array2OfReal*) NULL),
1510                          cachepoles->ChangeArray2(),
1511                          *((TColStd_Array2OfReal*) NULL)) ;
1512   }
1513   validcache = 1 ;
1514 }
1515
1516 //=======================================================================
1517 //function : IsCacheValid
1518 //purpose  : function that checks for the validity of the cache of the
1519 //           surface
1520 //=======================================================================
1521 Standard_Boolean Geom_BSplineSurface::IsCacheValid
1522 (const Standard_Real  U,
1523  const Standard_Real  V)  const 
1524 {
1525   //Roman Lygin 26.12.08, see comments in Geom_BSplineCurve::IsCacheValid()
1526   Standard_Real aDeltaU = U - ucacheparameter;
1527   Standard_Real aDeltaV = V - vcacheparameter;
1528
1529   return ( validcache &&
1530       (aDeltaU >= 0.0e0) &&
1531       ((aDeltaU < ucachespanlenght) || (ucachespanindex == ufknots->Upper() - udeg)) &&
1532       (aDeltaV >= 0.0e0) &&
1533       ((aDeltaV < vcachespanlenght) || (vcachespanindex == vfknots->Upper() - vdeg)) );
1534 }
1535
1536 //=======================================================================
1537 //function : SetWeight
1538 //purpose  : 
1539 //=======================================================================
1540
1541 void Geom_BSplineSurface::SetWeight (const Standard_Integer UIndex,
1542                                      const Standard_Integer VIndex,
1543                                      const Standard_Real    Weight)
1544 {
1545   if (Weight <= gp::Resolution())  Standard_ConstructionError::Raise(); 
1546   TColStd_Array2OfReal & Weights = weights->ChangeArray2();
1547   if (UIndex < 1 || UIndex > Weights.ColLength() ||
1548       VIndex < 1 || VIndex > Weights.RowLength() ) {
1549     Standard_OutOfRange::Raise();
1550   }
1551   Weights (UIndex+Weights.LowerRow()-1, VIndex+Weights.LowerCol()-1) = Weight;
1552   Rational(Weights, urational, vrational);
1553   InvalidateCache();
1554 }
1555
1556 //=======================================================================
1557 //function : SetWeightCol
1558 //purpose  : 
1559 //=======================================================================
1560
1561 void Geom_BSplineSurface::SetWeightCol
1562 (const Standard_Integer       VIndex, 
1563  const TColStd_Array1OfReal&  CPoleWeights)
1564 {
1565   TColStd_Array2OfReal & Weights = weights->ChangeArray2();   
1566   if (VIndex < 1 || VIndex > Weights.RowLength()) {
1567     Standard_OutOfRange::Raise();
1568   }
1569   if (CPoleWeights.Lower() < 1 || 
1570       CPoleWeights.Lower() > Weights.ColLength() ||
1571       CPoleWeights.Upper() < 1 ||
1572       CPoleWeights.Upper() > Weights.ColLength()  ) {
1573     Standard_ConstructionError::Raise();
1574   }
1575   Standard_Integer I = CPoleWeights.Lower();
1576   while (I <= CPoleWeights.Upper()) {
1577     if (CPoleWeights(I) <= gp::Resolution()) { 
1578       Standard_ConstructionError::Raise();
1579     }
1580     Weights (I+Weights.LowerRow()-1, VIndex+Weights.LowerCol()-1) = 
1581       CPoleWeights (I);
1582     I++;
1583   }
1584   // Verifie si c'est rationnel
1585   Rational(Weights, urational, vrational);
1586
1587   InvalidateCache();
1588 }
1589
1590 //=======================================================================
1591 //function : SetWeightRow
1592 //purpose  : 
1593 //=======================================================================
1594
1595 void Geom_BSplineSurface::SetWeightRow
1596 (const Standard_Integer       UIndex, 
1597  const TColStd_Array1OfReal&  CPoleWeights)
1598 {
1599   TColStd_Array2OfReal & Weights = weights->ChangeArray2();   
1600   if (UIndex < 1 || UIndex > Weights.ColLength()) {
1601     Standard_OutOfRange::Raise();
1602   }
1603   if (CPoleWeights.Lower() < 1 ||
1604       CPoleWeights.Lower() > Weights.RowLength() ||
1605       CPoleWeights.Upper() < 1 ||
1606       CPoleWeights.Upper() > Weights.RowLength()  ) {
1607     
1608     Standard_ConstructionError::Raise();
1609   }
1610   Standard_Integer I = CPoleWeights.Lower();
1611
1612   while (I <= CPoleWeights.Upper()) {
1613     if (CPoleWeights(I)<=gp::Resolution()) {
1614       Standard_ConstructionError::Raise();
1615     }
1616     Weights (UIndex+Weights.LowerRow()-1, I+Weights.LowerCol()-1) = 
1617       CPoleWeights (I);
1618     I++;
1619   }
1620   // Verifie si c'est rationnel
1621   Rational(Weights, urational, vrational);
1622   InvalidateCache();
1623 }
1624