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