0024589: Draw Harness - disable floating-point exceptions by default
[occt.git] / src / Convert / Convert_ConicToBSplineCurve.cxx
1 // Copyright (c) 1995-1999 Matra Datavision
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
3 //
4 // This file is part of Open CASCADE Technology software library.
5 //
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
11 //
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
14
15 //JCV 16/10/91
16
17 #define No_Standard_OutOfRange
18
19 #include <Convert_ConicToBSplineCurve.ixx>
20 #include <Convert_CosAndSinEvalFunction.hxx>
21 #include <Convert_PolynomialCosAndSin.hxx>
22 #include <TColStd_Array1OfReal.hxx>
23 #include <TColStd_Array1OfInteger.hxx>
24 #include <TColgp_Array1OfPnt2d.hxx>
25 #include <TColgp_Array1OfPnt.hxx>
26 #include <TColStd_HArray1OfReal.hxx>
27 #include <TColStd_Array1OfReal.hxx>
28 #include <TColStd_HArray1OfInteger.hxx>
29 #include <TColgp_HArray1OfPnt2d.hxx>
30
31 #include <PLib.hxx>
32 #include <BSplCLib.hxx>
33 #include <Precision.hxx>
34 #include <gp.hxx>
35 #include <Standard_OutOfRange.hxx>
36 #include <Standard_ConstructionError.hxx>
37
38 //=======================================================================
39 //function : Convert_ConicToBSplineCurve
40 //purpose  : 
41 //=======================================================================
42
43 Convert_ConicToBSplineCurve::Convert_ConicToBSplineCurve 
44   (const Standard_Integer NbPoles,
45    const Standard_Integer NbKnots,
46    const Standard_Integer Degree  ) 
47 : degree (Degree)    ,  nbPoles (NbPoles)   ,  nbKnots (NbKnots)
48   
49
50   if (NbPoles >= 2) {
51     poles    = new TColgp_HArray1OfPnt2d (1, NbPoles)  ;
52    
53     weights  = new TColStd_HArray1OfReal (1, NbPoles)  ;
54   }
55   if (NbKnots >= 2) {
56     knots    = new TColStd_HArray1OfReal (1, NbKnots)  ;
57     mults    = new TColStd_HArray1OfInteger(1,NbKnots) ;
58   }
59 }
60
61
62 //=======================================================================
63 //function : Degree
64 //purpose  : 
65 //=======================================================================
66
67 Standard_Integer Convert_ConicToBSplineCurve::Degree () const 
68 {
69   return degree;
70 }
71
72 //=======================================================================
73 //function : NbPoles
74 //purpose  : 
75 //=======================================================================
76
77 Standard_Integer Convert_ConicToBSplineCurve::NbPoles () const 
78 {
79  return nbPoles; 
80 }
81
82 //=======================================================================
83 //function : NbKnots
84 //purpose  : 
85 //=======================================================================
86
87 Standard_Integer Convert_ConicToBSplineCurve::NbKnots () const 
88 {
89   return nbKnots;
90 }
91
92 //=======================================================================
93 //function : IsPeriodic
94 //purpose  : 
95 //=======================================================================
96
97 Standard_Boolean Convert_ConicToBSplineCurve::IsPeriodic() const 
98 {
99   return isperiodic;
100 }
101
102 //=======================================================================
103 //function : Pole
104 //purpose  : 
105 //=======================================================================
106
107 gp_Pnt2d Convert_ConicToBSplineCurve::Pole 
108   (const Standard_Integer Index) const 
109 {
110   if (Index < 1 || Index > nbPoles)
111     Standard_OutOfRange::Raise(" ");
112   return poles->Value (Index);
113 }
114
115
116 //=======================================================================
117 //function : Weight
118 //purpose  : 
119 //=======================================================================
120
121 Standard_Real Convert_ConicToBSplineCurve::Weight 
122   (const Standard_Integer Index) const 
123 {
124   if (Index < 1 || Index > nbPoles)
125     Standard_OutOfRange::Raise(" ");
126   return weights->Value (Index);
127 }
128
129
130 //=======================================================================
131 //function : Knot
132 //purpose  : 
133 //=======================================================================
134
135 Standard_Real Convert_ConicToBSplineCurve::Knot 
136   (const Standard_Integer Index) const 
137 {
138   if (Index < 1 || Index > nbKnots)
139     Standard_OutOfRange::Raise(" ");
140   return knots->Value (Index);
141 }
142
143
144 //=======================================================================
145 //function : Multiplicity
146 //purpose  : 
147 //=======================================================================
148
149 Standard_Integer Convert_ConicToBSplineCurve::Multiplicity
150   (const Standard_Integer Index) const 
151 {
152   if (Index < 1 || Index > nbKnots)
153     Standard_OutOfRange::Raise(" ");
154   return mults->Value (Index);
155 }
156 //=======================================================================
157 //function : CosAndSinRationalC1
158 //purpose  : evaluates U(t) and V(t) such that
159 //                      2      2
160 //                     U   -  V 
161 //   cos (theta(t)) = ----------
162 //                      2      2
163 //                     U   +  V 
164 //
165
166 //                      2 * U*V
167 //   sin (theta(t)) = ----------
168 //                      2      2
169 //                     U   +  V 
170 //                                                    2     2
171 //  such that the derivative at the domain bounds of U   + V   is 0.0e0 
172 //  with is helpfull when having to make a C1 BSpline  by merging two
173 //  BSpline toghether
174 //=======================================================================
175
176 void CosAndSinRationalC1(Standard_Real Parameter,
177                           const Standard_Integer         EvalDegree,
178                           const TColgp_Array1OfPnt2d&    EvalPoles,
179                           const TColStd_Array1OfReal&    EvalKnots,
180                           const TColStd_Array1OfInteger& EvalMults,
181                           Standard_Real Result[2]) 
182 {
183  gp_Pnt2d a_point ;
184  BSplCLib::D0(Parameter,
185               0,
186               EvalDegree,
187               Standard_False,
188               EvalPoles,
189               BSplCLib::NoWeights(),
190               EvalKnots,
191               EvalMults,
192               a_point) ;
193  Result[0] = a_point.Coord(1) ;
194  Result[1] = a_point.Coord(2) ;
195 }
196       
197
198 //=======================================================================
199 //function : CosAndSinQuasiAngular
200 //purpose  : evaluates U(t) and V(t) such that
201 //                      2      2
202 //                     U   -  V 
203 //   cos (theta(t)) = ----------
204 //                      2      2
205 //                     U   +  V 
206 //
207
208 //                      2 * U*V
209 //   sin (theta(t)) = ----------
210 //                      2      2
211 //                     U   +  V 
212 //=======================================================================
213
214 void  CosAndSinQuasiAngular(Standard_Real  Parameter,
215                             const Standard_Integer         EvalDegree,
216                             const TColgp_Array1OfPnt2d&    EvalPoles,
217 //                          const TColStd_Array1OfReal&    EvalKnots,
218                             const TColStd_Array1OfReal&    ,
219 //                          const TColStd_Array1OfInteger& EvalMults,
220                             const TColStd_Array1OfInteger& ,
221                             Standard_Real  Result[2])
222 {
223   Standard_Real 
224   param,
225   *coeff ;
226  
227    coeff = (Standard_Real *) &EvalPoles(EvalPoles.Lower()) ;
228 //
229 //   rational_function_coeff represent a rational approximation
230 //   of U ---> cotan( PI * U /2) between [0 1] 
231 //   rational_function_coeff[i][0] is the denominator
232 //   rational_function_coeff[i][1] is the numerator
233 //   
234    param = Parameter * 0.5e0 ;
235   PLib::NoDerivativeEvalPolynomial (param,
236                                     EvalDegree,
237                                     2,
238                                     EvalDegree << 1,
239                                     coeff[0],
240                                     Result[0]) ;
241 }
242    
243 //=======================================================================
244 //function : function that build the Bspline Representation of 
245 // an algorithmic description of the function cos and sin
246 //purpose  : 
247 //=======================================================================
248 void AlgorithmicCosAndSin(Standard_Integer               Degree,
249                           const TColStd_Array1OfReal&    FlatKnots,
250                           const Standard_Integer         EvalDegree,
251                           const TColgp_Array1OfPnt2d&    EvalPoles,
252                           const TColStd_Array1OfReal&    EvalKnots,
253                           const TColStd_Array1OfInteger& EvalMults,
254                           Convert_CosAndSinEvalFunction  Evaluator,
255                           TColStd_Array1OfReal&          CosNumerator,
256                           TColStd_Array1OfReal&          SinNumerator,
257                           TColStd_Array1OfReal&          Denominator) 
258 {
259    Standard_Integer order,
260    num_poles,
261    pivot_index_problem,
262    ii;
263  
264    Standard_Real  result[2],
265    inverse ;
266
267    order = Degree + 1 ;
268    num_poles = FlatKnots.Length() - order ;
269
270    if (num_poles != CosNumerator.Length() ||
271        num_poles != SinNumerator.Length() ||
272        num_poles != Denominator.Length() ) {
273       Standard_ConstructionError::Raise();
274       }
275    TColStd_Array1OfReal      parameters(1,num_poles)  ;
276    TColgp_Array1OfPnt        poles_array(1,num_poles) ;
277    TColStd_Array1OfInteger   contact_order_array(1,num_poles) ;  
278    BSplCLib::BuildSchoenbergPoints(Degree,
279                                    FlatKnots,
280                                    parameters) ;
281    for (ii = parameters.Lower() ; ii <= parameters.Upper() ; ii++) {
282      Evaluator(parameters(ii),
283                EvalDegree,
284                EvalPoles,
285                EvalKnots,
286                EvalMults,
287                result) ;
288      contact_order_array(ii) = 0 ;
289      
290      poles_array(ii).SetCoord(1,
291                               (result[1]*result[1] - result[0]*result[0]));
292      poles_array(ii).SetCoord(2,
293                               2.0e0  * result[1]* result[0]) ;
294      poles_array(ii).SetCoord(3,
295                               result[1]*result[1] + result[0] * result[0]) ;
296    }
297     BSplCLib::Interpolate(Degree,
298                           FlatKnots,
299                           parameters,
300                           contact_order_array,
301                           poles_array,
302                           pivot_index_problem) ;
303     for (ii = 1 ; ii <= num_poles ; ii++) {
304       inverse = 1.0e0 / poles_array(ii).Coord(3) ;
305       CosNumerator(ii) = poles_array(ii).Coord(1) * inverse ;
306       SinNumerator(ii) = poles_array(ii).Coord(2) * inverse ;
307       Denominator(ii)  = poles_array(ii).Coord(3) ;
308     }
309  }
310
311 //=======================================================================
312 //function : BuildCosAndSin
313 //purpose  : 
314 //=======================================================================
315
316 void Convert_ConicToBSplineCurve::BuildCosAndSin(
317           const Convert_ParameterisationType     Parameterisation,
318           const Standard_Real                    UFirst,
319           const Standard_Real                    ULast,
320           Handle(TColStd_HArray1OfReal)&         CosNumeratorPtr,
321           Handle(TColStd_HArray1OfReal)&         SinNumeratorPtr,
322           Handle(TColStd_HArray1OfReal)&         DenominatorPtr,
323           Standard_Integer&                      Degree,
324           Handle(TColStd_HArray1OfReal)&         KnotsPtr,
325           Handle(TColStd_HArray1OfInteger)&      MultsPtr)  const 
326 {
327   Standard_Real delta = ULast - UFirst,
328   direct,
329   inverse,
330   value1,
331   value2,
332   cos_beta,
333   sin_beta,
334   alpha=0,
335   alpha_2,
336   alpha_4,
337   tan_alpha_2,
338   beta,
339   p_param,
340   q_param,
341   param ;
342
343   Standard_Integer num_poles = 0,
344   ii,
345   jj,
346   num_knots = 1,
347   num_spans = 1,
348   num_flat_knots,
349   num_temp_knots,
350   temp_degree = 0,
351   tgt_theta_flag,
352   num_temp_poles,
353   order  = 0;
354
355   Convert_CosAndSinEvalFunction *EvaluatorPtr=NULL ;
356
357   tgt_theta_flag = 0 ;
358
359
360   switch (Parameterisation) {
361   case Convert_TgtThetaOver2: 
362     num_spans =
363       (Standard_Integer)IntegerPart( 1.2 * delta / M_PI) + 1;
364     
365     tgt_theta_flag = 1 ;
366     break ;
367   case Convert_TgtThetaOver2_1:
368     num_spans = 1 ;
369     if (delta > 0.9999 * M_PI) {
370       Standard_ConstructionError::Raise() ; 
371       }
372     tgt_theta_flag = 1 ;
373     break ;
374   case Convert_TgtThetaOver2_2:
375     num_spans = 2 ;
376     if (delta > 1.9999 * M_PI) {
377       Standard_ConstructionError::Raise() ;
378       }
379     tgt_theta_flag = 1 ;
380     break ;
381   
382   case Convert_TgtThetaOver2_3:
383     num_spans = 3 ;
384     tgt_theta_flag = 1 ;
385     break ;
386   case Convert_TgtThetaOver2_4:
387     num_spans = 4 ;
388     tgt_theta_flag = 1 ;
389     break ; 
390   case Convert_QuasiAngular:
391     num_poles = 7 ;
392     Degree    = 6 ;
393     num_spans = 1 ;
394     num_knots = 2 ;
395     order = Degree + 1 ;
396     break ;
397    case Convert_RationalC1:
398     Degree    = 4 ;
399     order = Degree + 1 ; 
400     num_poles = 8 ;
401     num_knots = 3 ;
402     num_spans = 2 ;
403     break ;
404    case Convert_Polynomial:
405     Degree    = 7 ;
406     num_poles = 8 ;
407     num_knots = 2 ;
408     num_spans = 1 ;
409     break ;
410   default:
411     break ;
412   }
413   if (tgt_theta_flag) {
414     alpha = delta / ( 2.0e0 * num_spans) ;    
415     Degree = 2 ;
416     num_poles = 2 * num_spans + 1;
417   }  
418   
419   CosNumeratorPtr = 
420     new TColStd_HArray1OfReal(1,num_poles) ;
421   SinNumeratorPtr =
422     new TColStd_HArray1OfReal(1,num_poles) ;      
423   DenominatorPtr =
424     new TColStd_HArray1OfReal(1,num_poles) ;      
425   KnotsPtr = 
426     new TColStd_HArray1OfReal(1,num_spans+1) ;
427   MultsPtr =
428     new TColStd_HArray1OfInteger(1,num_spans+1) ;
429   if (tgt_theta_flag) {
430
431     param = UFirst ;
432     CosNumeratorPtr->SetValue(1,Cos(UFirst)) ;
433     SinNumeratorPtr->SetValue(1,Sin(UFirst)) ;
434     DenominatorPtr ->SetValue(1,1.0e0) ;
435     KnotsPtr->SetValue(1,param) ;
436     MultsPtr->SetValue(1,Degree + 1) ;
437     direct = Cos(alpha) ;
438     inverse = 1.0e0 / direct ;
439     for (ii = 1 ; ii <= num_spans ; ii++ ) {
440       CosNumeratorPtr->SetValue(2 * ii, inverse * Cos(param + alpha)) ;
441       SinNumeratorPtr->SetValue(2 * ii, inverse * Sin(param + alpha)) ;
442       DenominatorPtr->SetValue(2 * ii,  direct) ; 
443       CosNumeratorPtr->SetValue(2 * ii + 1, Cos(param + 2 * alpha)) ;
444       SinNumeratorPtr->SetValue(2 * ii + 1, Sin(param + 2 * alpha)) ;
445       DenominatorPtr->SetValue(2 * ii + 1,  1.0e0) ;
446       KnotsPtr->SetValue(ii + 1, param + 2 * alpha) ;
447       MultsPtr->SetValue(ii + 1, 2) ;
448       param += 2 * alpha ;
449     }
450     MultsPtr->SetValue(num_spans + 1, Degree + 1) ;
451   }
452   else if (Parameterisation != Convert_Polynomial) {
453     alpha = ULast - UFirst ;
454     alpha *= 0.5e0 ;
455     beta = ULast + UFirst ;
456     beta *= 0.5e0 ;
457     cos_beta = Cos(beta) ;
458     sin_beta = Sin(beta) ;
459     num_flat_knots = num_poles + order ;
460
461     num_temp_poles = 4 ;
462     num_temp_knots = 3 ;
463     TColStd_Array1OfReal   flat_knots(1, num_flat_knots) ;
464     
465
466     TColgp_Array1OfPnt2d    temp_poles(1,num_temp_poles)  ;
467     TColStd_Array1OfReal    temp_knots(1,num_temp_knots)  ;
468     TColStd_Array1OfInteger temp_mults(1,num_temp_knots) ;
469
470     for (ii = 1 ; ii <= order ; ii++) {
471       flat_knots(ii) = -alpha ;
472       flat_knots(ii + num_poles) = alpha ;
473     }
474     KnotsPtr->SetValue(1,UFirst) ;
475     KnotsPtr->SetValue(num_knots, ULast) ;
476     MultsPtr->SetValue(1,order) ;
477     MultsPtr->SetValue(num_knots,order) ;
478     
479     switch (Parameterisation) {
480     case Convert_QuasiAngular:
481 //
482 //    we code here in temp_poles(xx).Coord(1) the following function V(t) 
483 //   and in temp_poles(xx).Coord(2) the function U(t) 
484 //                     3
485 //       V(t) = t + c t
486 //                     2
487 //       U(t) = 1 + b t
488 //            1
489 //       c = ---  + b   = q_param
490 //            3 
491 //                          3
492 //                     gamma
493 //            gamma +  ------  - tang gamma
494 //                      3
495 //       b =------------------------------    = p_param
496 //                 2 
497 //            gamma  (tang gamma - gamma) 
498 //
499 //     with gamma = alpha / 2
500 //
501 //
502      
503       alpha_2 = alpha * 0.5e0 ;
504       p_param = - 1.0e0 / (alpha_2 * alpha_2) ;
505      
506       if (alpha_2 <  M_PI * 0.5e0) 
507       {
508         if (alpha_2 < 1.0e-7)
509         {
510           // Fixed degenerate case, when obtain 0 / 0 uncertainty.
511           // According to Taylor aprroximation:
512           // b (gamma) = -6.0 / 15.0 + o(gamma^2)
513           p_param = -6.0 / 15.0;
514         }
515         else
516         {
517           tan_alpha_2 = Tan(alpha_2) ;
518           value1 = 3.0e0 * (tan_alpha_2 - alpha_2) ;
519           value1 = alpha_2 / value1 ;
520           p_param += value1 ;
521         }
522       }
523       q_param = (1.0e0 / 3.0e0)  + p_param ;
524       
525       
526       temp_degree = 3 ;
527       temp_poles(1).SetCoord(1,0.0e0);
528       temp_poles(2).SetCoord(1,1.0e0);
529       temp_poles(3).SetCoord(1,0.0e0) ;
530       temp_poles(4).SetCoord(1,q_param) ;
531
532       temp_poles(1).SetCoord(2, 1.0e0) ;
533       temp_poles(2).SetCoord(2, 0.0e0) ;
534       temp_poles(3).SetCoord(2, p_param) ;
535       temp_poles(4).SetCoord(2, 0.0e0);
536       EvaluatorPtr = &CosAndSinQuasiAngular ;
537       break ;
538     case  Convert_RationalC1:
539       for (ii = order + 1 ; ii <= num_poles ; ii++) {
540         flat_knots(ii) = 0.0e0 ;
541       }
542       KnotsPtr->SetValue(2,UFirst + alpha) ;
543       MultsPtr->SetValue(2,Degree -1) ;
544       temp_degree = 2 ;   
545       alpha_2 = alpha * 0.5e0 ;
546       alpha_4 = alpha * 0.25e0 ;
547       tan_alpha_2 = Tan(alpha_2) ;
548       jj = 1 ;
549       for (ii = 1 ; ii <= 2 ; ii++) {
550         temp_poles(1+ ii).SetCoord(2,1.0e0 + alpha_4 * tan_alpha_2) ;
551         temp_poles(jj).SetCoord(2,1.e0) ;
552         jj += 3 ;
553       }
554       temp_poles(1).SetCoord(1,-tan_alpha_2) ;  
555       temp_poles(2).SetCoord(1,alpha_4 - tan_alpha_2) ; 
556       temp_poles(3).SetCoord(1,-alpha_4 + tan_alpha_2) ; 
557       temp_poles(4).SetCoord(1,tan_alpha_2) ; 
558       temp_knots(1) = -alpha ;
559       temp_knots(2) = 0.0e0 ;
560       temp_knots(3) = alpha ;
561       temp_mults(1) = temp_degree + 1;
562       temp_mults(2) = 1 ;
563       temp_mults(3) = temp_degree + 1;
564      
565       EvaluatorPtr = &CosAndSinRationalC1 ;
566       break ;
567     default: 
568       break ;
569     }
570     AlgorithmicCosAndSin(Degree,
571                          flat_knots,
572                          temp_degree,
573                          temp_poles,
574                          temp_knots,
575                          temp_mults,
576                          *EvaluatorPtr,
577                          CosNumeratorPtr->ChangeArray1(),
578                          SinNumeratorPtr->ChangeArray1(),
579                          DenominatorPtr->ChangeArray1()) ;
580
581     for (ii = 1 ; ii <= num_poles ; ii++) {
582        value1 = cos_beta * CosNumeratorPtr->Value(ii) -
583          sin_beta * SinNumeratorPtr->Value(ii) ;
584        value2 = sin_beta * CosNumeratorPtr->Value(ii) +
585          cos_beta * SinNumeratorPtr->Value(ii) ;
586        CosNumeratorPtr->SetValue(ii,value1) ;
587        SinNumeratorPtr->SetValue(ii,value2) ;
588      }
589   }
590   else { // Convert_Polynomial
591
592     KnotsPtr->SetValue(1, 0.) ;
593     KnotsPtr->SetValue(num_knots, 1.);
594     MultsPtr->SetValue(1, num_poles);
595     MultsPtr->SetValue(num_knots, num_poles);
596     
597     BuildPolynomialCosAndSin(UFirst,ULast,num_poles,
598                              CosNumeratorPtr,SinNumeratorPtr,DenominatorPtr);
599   }
600
601
602 }
603 //=======================================================================
604 //function : BuildCosAndSin
605 //purpose  : 
606 //=======================================================================
607
608 void Convert_ConicToBSplineCurve::BuildCosAndSin(
609           const Convert_ParameterisationType     Parameterisation,
610           Handle(TColStd_HArray1OfReal)&         CosNumeratorPtr,
611           Handle(TColStd_HArray1OfReal)&         SinNumeratorPtr,
612           Handle(TColStd_HArray1OfReal)&         DenominatorPtr,
613           Standard_Integer&                      Degree,
614           Handle(TColStd_HArray1OfReal)&         KnotsPtr,
615           Handle(TColStd_HArray1OfInteger)&      MultsPtr)  const 
616 {
617   Standard_Real half_pi,
618   param,
619   first_param,
620   last_param,
621 //  direct,
622   inverse,
623   value1,
624   value2,
625   value3 ;
626
627   Standard_Integer 
628   ii,
629   jj,
630   index,
631   num_poles,
632   num_periodic_poles,
633   temp_degree,
634   pivot_index_problem,
635   num_flat_knots,
636   num_knots,
637   order ;
638
639   if (Parameterisation != Convert_TgtThetaOver2 &&
640       Parameterisation != Convert_RationalC1) {
641       Standard_ConstructionError::Raise() ;
642       }
643   Handle(TColStd_HArray1OfReal) temp_cos_ptr,
644   temp_sin_ptr,
645   temp_denominator_ptr,
646   temp_knots_ptr;
647   Handle(TColStd_HArray1OfInteger)  temp_mults_ptr;
648   if (Parameterisation == Convert_TgtThetaOver2) {
649     BuildCosAndSin(Convert_TgtThetaOver2_3,
650                    0.0e0,
651                    2 * M_PI,
652                    temp_cos_ptr,
653                    temp_sin_ptr,
654                    temp_denominator_ptr,
655                    Degree,
656                    KnotsPtr,
657                    MultsPtr) ;
658      CosNumeratorPtr =
659        new TColStd_HArray1OfReal(1,temp_cos_ptr->Length() -1) ;
660      SinNumeratorPtr =
661         new TColStd_HArray1OfReal(1,temp_cos_ptr->Length() -1) ;
662      DenominatorPtr =
663         new TColStd_HArray1OfReal(1,temp_cos_ptr->Length() -1) ; 
664      for (ii = temp_cos_ptr->Lower()  ; ii <= temp_cos_ptr->Upper()-1 ; ii++) {
665        CosNumeratorPtr->SetValue(ii,temp_cos_ptr->Value(ii)) ;
666        SinNumeratorPtr->SetValue(ii,temp_sin_ptr->Value(ii)) ;
667        DenominatorPtr->SetValue(ii,temp_denominator_ptr->Value(ii)) ;
668      }
669     for (ii = MultsPtr->Lower() ; ii <= MultsPtr->Upper() ; ii++) {
670       MultsPtr->SetValue(ii, Degree) ;
671     }
672   }
673   else if (Parameterisation == Convert_RationalC1) 
674     {
675      first_param = 0.0e0 ;
676      last_param  = M_PI ;
677      BuildCosAndSin(Convert_RationalC1,
678                    first_param,
679                    last_param,
680                    temp_cos_ptr,
681                    temp_sin_ptr,
682                    temp_denominator_ptr,
683                    temp_degree,
684                    temp_knots_ptr,
685                    temp_mults_ptr) ;
686
687
688      Degree = 4 ;
689      order = Degree + 1 ;
690      num_knots = 5 ;
691      num_flat_knots = (Degree -1) * num_knots + 2 * 2 ;
692      num_poles = num_flat_knots - order ;
693      num_periodic_poles = num_poles - 2 ;
694      TColStd_Array1OfReal  flat_knots(1,num_flat_knots) ;
695      CosNumeratorPtr = 
696       new TColStd_HArray1OfReal(1,num_periodic_poles) ;
697      SinNumeratorPtr = 
698       new TColStd_HArray1OfReal(1,num_periodic_poles) ;
699      DenominatorPtr = 
700       new TColStd_HArray1OfReal(1,num_periodic_poles) ;
701     
702      half_pi = M_PI * 0.5e0 ;
703      index = 1 ;
704      for (jj = 1 ; jj <= 2 ; jj++) {
705          flat_knots(index) = -  half_pi  ;
706          index += 1 ;
707        }
708      for (ii = 1 ; ii <= num_knots ; ii++) {
709        for (jj = 1 ; jj <= Degree -1 ; jj++) {
710          flat_knots(index) = (ii-1) * half_pi ;
711
712          index += 1 ;
713
714        }
715      }
716      for (jj = 1 ; jj <= 2 ; jj++) {
717          flat_knots(index) = 2 * M_PI +  half_pi  ;
718          index += 1 ;
719        }
720      KnotsPtr = 
721        new TColStd_HArray1OfReal(1,num_knots) ;
722      MultsPtr =
723        new TColStd_HArray1OfInteger(1,num_knots) ;
724      for ( ii = 1 ; ii <= num_knots  ; ii++) { 
725        KnotsPtr->SetValue(ii, (ii-1) * half_pi) ;
726        MultsPtr->SetValue(ii, Degree-1) ;
727      }
728      order = degree + 1 ;
729
730      TColStd_Array1OfReal      parameters(1,num_poles)  ;
731      TColgp_Array1OfPnt        poles_array(1,num_poles) ;
732      TColStd_Array1OfInteger   contact_order_array(1,num_poles) ;  
733      BSplCLib::BuildSchoenbergPoints(Degree,
734                                      flat_knots,
735                                      parameters) ;
736      inverse = 1.0e0 ;
737      for (ii = parameters.Lower() ; ii <= parameters.Upper() ; ii++) {
738        param = parameters(ii) ;
739        if (param > M_PI) {
740          inverse = -1.0e0 ;
741          param -= M_PI ;
742        }
743        BSplCLib::D0(param,
744                     0,
745                     temp_degree,
746                     Standard_False,
747                     temp_cos_ptr->Array1(),
748                     temp_denominator_ptr->Array1(),
749                     temp_knots_ptr->Array1(),
750                     temp_mults_ptr->Array1(),
751                     value1) ;
752
753        BSplCLib::D0(param,
754                     0,
755                     temp_degree,
756                     Standard_False,
757                     temp_sin_ptr->Array1(),
758                     temp_denominator_ptr->Array1(),
759                     temp_knots_ptr->Array1(),
760                     temp_mults_ptr->Array1(),
761                     value2) ;
762        BSplCLib::D0(param,
763                     0,
764                     temp_degree,
765                     Standard_False,
766                     temp_denominator_ptr->Array1(),
767                     BSplCLib::NoWeights(),
768                     temp_knots_ptr->Array1(),
769                     temp_mults_ptr->Array1(),
770                     value3) ;
771      contact_order_array(ii) = 0 ;
772      
773      poles_array(ii).SetCoord(1,
774                               value1 * value3 * inverse) ;
775      poles_array(ii).SetCoord(2,
776                               value2 * value3 * inverse) ;
777      poles_array(ii).SetCoord(3,
778                               value3) ;
779    }
780     BSplCLib::Interpolate(Degree,
781                           flat_knots,
782                           parameters,
783                           contact_order_array,
784                           poles_array,
785                           pivot_index_problem) ;
786     for (ii = 1 ; ii <= num_periodic_poles ; ii++) {
787       inverse = 1.0e0 / poles_array(ii).Coord(3) ;
788       CosNumeratorPtr->ChangeArray1()(ii) = poles_array(ii).Coord(1) * inverse ;
789       SinNumeratorPtr->ChangeArray1()(ii) = poles_array(ii).Coord(2) * inverse ;
790       DenominatorPtr->ChangeArray1()(ii)  = poles_array(ii).Coord(3) ;
791     }
792  }
793 }
794
795
796