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