106600b900e40269c01ca78920af6c49415e391e
[occt.git] / src / Geom2dAPI / Geom2dAPI_Interpolate.cxx
1 // Created on: 1994-08-18
2 // Created by: Laurent PAINNOT
3 // Copyright (c) 1994-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 // 8-Aug-95 : xab : interpolation uses BSplCLib::Interpolate
18
19 #include <BSplCLib.hxx>
20 #include <Geom2d_BSplineCurve.hxx>
21 #include <Geom2dAPI_Interpolate.hxx>
22 #include <gp_Vec2d.hxx>
23 #include <PLib.hxx>
24 #include <Standard_ConstructionError.hxx>
25 #include <StdFail_NotDone.hxx>
26 #include <TColgp_Array1OfPnt2d.hxx>
27 #include <TColStd_Array1OfBoolean.hxx>
28 #include <TColStd_Array1OfInteger.hxx>
29 #include <TColStd_HArray1OfBoolean.hxx>
30 #include <TColStd_HArray1OfReal.hxx>
31
32 //=======================================================================
33 //function : CheckPoints
34 //purpose  : 
35 //=======================================================================
36 static Standard_Boolean CheckPoints(const TColgp_Array1OfPnt2d& PointArray,
37                                     const Standard_Real    Tolerance) 
38 {
39   Standard_Integer ii ;
40   Standard_Real tolerance_squared = Tolerance * Tolerance,
41   distance_squared ;
42   Standard_Boolean result = Standard_True ;
43   for (ii = PointArray.Lower() ; result && ii < PointArray.Upper() ; ii++) {
44     distance_squared = 
45       PointArray.Value(ii).SquareDistance(PointArray.Value(ii+1)) ;
46     result = (distance_squared >= tolerance_squared) ;
47   }
48  return result ;
49
50 }
51 //=======================================================================
52 //function : CheckTangents
53 //purpose  : 
54 //=======================================================================
55 static Standard_Boolean CheckTangents(
56                                 const TColgp_Array1OfVec2d&      Tangents,
57                                 const TColStd_Array1OfBoolean& TangentFlags,
58                                 const Standard_Real            Tolerance) 
59 {
60   Standard_Integer ii,
61   index ;
62   Standard_Real tolerance_squared = Tolerance * Tolerance,
63   distance_squared ;
64   Standard_Boolean result = Standard_True ;
65   index = TangentFlags.Lower() ;
66   for (ii = Tangents.Lower(); result && ii <= Tangents.Upper() ; ii++) {
67     if(TangentFlags.Value(index)) {
68       distance_squared = 
69         Tangents.Value(ii).SquareMagnitude() ;
70       result = (distance_squared >= tolerance_squared) ;
71     }
72     index += 1 ;
73   }
74  return result ;
75
76 }
77 //=======================================================================
78 //function : CheckParameters
79 //purpose  : 
80 //=======================================================================
81 static Standard_Boolean CheckParameters(const 
82                                         TColStd_Array1OfReal&   Parameters) 
83 {
84   Standard_Integer ii ;
85   Standard_Real distance ;
86   Standard_Boolean result = Standard_True ;
87   for (ii = Parameters.Lower() ; result && ii < Parameters.Upper() ; ii++) {
88     distance = 
89       Parameters.Value(ii+1) - Parameters.Value(ii) ;
90     result = (distance >= RealSmall()) ;
91   }
92  return result ;
93 }       
94 //=======================================================================
95 //function : BuildParameters
96 //purpose  : 
97 //=======================================================================                     
98 static void  BuildParameters(const Standard_Boolean        PeriodicFlag,
99                              const TColgp_Array1OfPnt2d&     PointsArray,
100                              Handle(TColStd_HArray1OfReal)& ParametersPtr) 
101 {
102   Standard_Integer ii,
103   index ;
104   Standard_Real distance ;
105   Standard_Integer 
106     num_parameters = PointsArray.Length() ;
107   if (PeriodicFlag) {
108     num_parameters += 1 ;
109   }
110   ParametersPtr =
111     new TColStd_HArray1OfReal(1,
112                               num_parameters) ;
113   ParametersPtr->SetValue(1,0.0e0) ;
114   index = 2 ;
115   for (ii = PointsArray.Lower() ; ii < PointsArray.Upper() ; ii++) {
116     distance = 
117       PointsArray.Value(ii).Distance(PointsArray.Value(ii+1)) ;
118     ParametersPtr->SetValue(index,
119                             ParametersPtr->Value(ii) + distance) ;
120     index += 1 ;
121   }
122   if (PeriodicFlag) {
123     distance = 
124       PointsArray.Value(PointsArray.Upper()).
125         Distance(PointsArray.Value(PointsArray.Lower())) ;
126     ParametersPtr->SetValue(index,
127                             ParametersPtr->Value(ii) + distance) ;
128   }
129 }
130 //=======================================================================
131 //function : BuildPeriodicTangents
132 //purpose  : 
133 //=======================================================================
134                 
135 static void BuildPeriodicTangent(
136                       const TColgp_Array1OfPnt2d&      PointsArray,
137                       TColgp_Array1OfVec2d&            TangentsArray,
138                       TColStd_Array1OfBoolean&       TangentFlags,
139                       const TColStd_Array1OfReal&    ParametersArray)
140 {
141   Standard_Integer 
142     ii,
143     degree ;
144   Standard_Real *point_array,
145   *parameter_array,
146   eval_result[2][2] ;
147   
148   gp_Vec2d a_vector ;
149   
150   if (PointsArray.Length() < 3) {
151     throw Standard_ConstructionError();
152     }   
153  
154   if (!TangentFlags.Value(1)) {
155     degree = 3 ;
156     if (PointsArray.Length() == 3) {
157       degree = 2 ;
158     }
159     point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower()) ; 
160     parameter_array =
161       (Standard_Real *) &ParametersArray.Value(1) ;
162     TangentFlags.SetValue(1,Standard_True) ;
163     PLib::EvalLagrange(ParametersArray.Value(1),
164                        1,
165                        degree,
166                        2,
167                        point_array[0],
168                        parameter_array[0],
169                        eval_result[0][0]) ;
170     for (ii = 1 ; ii <= 2 ; ii++) {
171       a_vector.SetCoord(ii,eval_result[1][ii-1]) ;
172     }
173     TangentsArray.SetValue(1,a_vector) ;
174   }
175  } 
176 //=======================================================================
177 //function : BuildTangents
178 //purpose  : 
179 //=======================================================================
180                 
181 static void BuildTangents(const TColgp_Array1OfPnt2d&      PointsArray,
182                           TColgp_Array1OfVec2d&            TangentsArray,
183                           TColStd_Array1OfBoolean&       TangentFlags,
184                           const TColStd_Array1OfReal&    ParametersArray)
185 {
186  Standard_Integer ii,
187  degree ;
188  Standard_Real *point_array,
189  *parameter_array,
190  
191  eval_result[2][2] ;
192  gp_Vec2d a_vector ;
193  
194  degree = 3 ;
195
196  if ( PointsArray.Length() < 3) {
197    throw Standard_ConstructionError();
198    }   
199  if (PointsArray.Length() == 3) {
200    degree = 2 ;
201  }
202  if (!TangentFlags.Value(1)) {
203    point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower()) ; 
204    parameter_array =
205      (Standard_Real *) &ParametersArray.Value(1) ;
206    TangentFlags.SetValue(1,Standard_True) ;
207    PLib::EvalLagrange(ParametersArray.Value(1),
208                       1,
209                       degree,
210                       2,
211                       point_array[0],
212                       parameter_array[0],
213                       eval_result[0][0]) ;
214    for (ii = 1 ; ii <= 2 ; ii++) {
215      a_vector.SetCoord(ii,eval_result[1][ii-1]) ;
216    }
217    TangentsArray.SetValue(1,a_vector) ;
218  }
219  if (! TangentFlags.Value(TangentFlags.Upper())) {
220    point_array = 
221      (Standard_Real *) &PointsArray.Value(PointsArray.Upper() - degree) ;
222    TangentFlags.SetValue(TangentFlags.Upper(),Standard_True) ;
223    parameter_array =
224     (Standard_Real *)&ParametersArray.Value(ParametersArray.Upper() - degree) ;
225    PLib::EvalLagrange(ParametersArray.Value(ParametersArray.Upper()),
226                       1,
227                       degree,
228                       2,
229                       point_array[0],
230                       parameter_array[0],
231                       eval_result[0][0]) ;
232    for (ii = 1 ; ii <= 2 ; ii++) {
233      a_vector.SetCoord(ii,eval_result[1][ii-1]) ; 
234    }
235    TangentsArray.SetValue(TangentsArray.Upper(),a_vector) ;
236  }
237
238 //=======================================================================
239 //function : BuildTangents
240 //purpose  : scale the given tangent so that they have the length of
241 // the size of the derivative of the lagrange interpolation 
242 //
243 //=======================================================================    
244 static void ScaleTangents(const TColgp_Array1OfPnt2d&      PointsArray,
245                           TColgp_Array1OfVec2d&            TangentsArray,
246                           const TColStd_Array1OfBoolean& TangentFlags,
247                           const TColStd_Array1OfReal&     ParametersArray)
248 {
249  Standard_Integer ii,
250  jj,
251  degree=0,
252  index,
253  num_points ;
254
255  Standard_Real *point_array,
256  *parameter_array,
257  value[2],
258  ratio,
259  eval_result[2][2] ;
260
261  gp_Vec2d a_vector ;
262
263  num_points = PointsArray.Length() ; 
264  if (num_points == 2) {
265     degree = 1 ;
266   }
267  else if (num_points >= 3) {
268    degree = 2 ;
269  }
270  
271  index = PointsArray.Lower() ;
272  for (ii = TangentFlags.Lower()  ; ii <=  TangentFlags.Upper() ; ii++) {
273    if (TangentFlags.Value(ii)) {
274      point_array = 
275        (Standard_Real *) &PointsArray.Value(index) ; 
276      parameter_array =
277        (Standard_Real *) &ParametersArray.Value(index) ;
278      PLib::EvalLagrange(ParametersArray.Value(ii),
279                         1,
280                         degree,
281                         2,
282                         point_array[0],
283                         parameter_array[0],
284                         eval_result[0][0]) ;
285      value[0] = 
286        value[1] = 0.0e0 ;
287      for (jj = 1 ; jj <= 2 ; jj++) {
288        value[0] += Abs(TangentsArray.Value(ii).Coord(jj)) ;
289        value[1] += Abs(eval_result[1][jj-1]) ;
290      }
291      ratio = value[1] / value[0] ;
292      for (jj = 1 ; jj <= 2 ; jj++) {
293        a_vector.SetCoord(jj, ratio *
294                          TangentsArray.Value(ii).Coord(jj)) ;
295      }
296      TangentsArray.SetValue(ii, a_vector) ;
297      if (ii != TangentFlags.Lower()) {
298        index += 1 ;
299      }
300      if (index > PointsArray.Upper() - degree) {
301        index = PointsArray.Upper() - degree ;
302      } 
303      
304    }
305  }
306 }
307  
308 //=======================================================================
309 //function : Geom2dAPI_Interpolate
310 //purpose  : 
311 //=======================================================================
312
313 Geom2dAPI_Interpolate::Geom2dAPI_Interpolate
314    (const Handle(TColgp_HArray1OfPnt2d)& PointsPtr,
315     const Standard_Boolean            PeriodicFlag,
316     const Standard_Real               Tolerance) :
317 myTolerance(Tolerance),
318 myPoints(PointsPtr),
319 myIsDone(Standard_False),
320 myPeriodic(PeriodicFlag),
321 myTangentRequest(Standard_False) 
322
323 {
324  Standard_Integer ii ;
325  Standard_Boolean result = 
326    CheckPoints(PointsPtr->Array1(),
327                Tolerance) ;
328  myTangents = 
329      new TColgp_HArray1OfVec2d(myPoints->Lower(),
330                               myPoints->Upper()) ;
331  myTangentFlags =
332       new TColStd_HArray1OfBoolean(myPoints->Lower(),
333                                    myPoints->Upper()) ;
334
335  if (!result) {
336    throw Standard_ConstructionError();
337    }
338  BuildParameters(PeriodicFlag,
339                  PointsPtr->Array1(),
340                  myParameters) ;
341
342   for (ii = myPoints->Lower() ; ii <= myPoints->Upper() ; ii++) {
343     myTangentFlags->SetValue(ii,Standard_False) ;
344   }
345
346  
347                  
348 }
349
350 //=======================================================================
351 //function : Geom2dAPI_Interpolate
352 //purpose  : 
353 //=======================================================================
354
355 Geom2dAPI_Interpolate::Geom2dAPI_Interpolate
356    (const Handle(TColgp_HArray1OfPnt2d)&    PointsPtr,
357     const Handle(TColStd_HArray1OfReal)&  ParametersPtr,
358     const Standard_Boolean               PeriodicFlag,
359     const Standard_Real                  Tolerance) :
360 myTolerance(Tolerance),
361 myPoints(PointsPtr),
362 myIsDone(Standard_False),
363 myParameters(ParametersPtr),
364 myPeriodic(PeriodicFlag),
365 myTangentRequest(Standard_False) 
366 {
367  Standard_Integer ii ;
368     
369      
370  Standard_Boolean result = 
371    CheckPoints(PointsPtr->Array1(),
372                Tolerance) ;
373
374  if (PeriodicFlag) {
375    if ((PointsPtr->Length()) + 1 != ParametersPtr->Length()) {
376      throw Standard_ConstructionError();
377    }
378  }
379  myTangents = 
380      new TColgp_HArray1OfVec2d(myPoints->Lower(),
381                               myPoints->Upper()) ;
382  myTangentFlags =
383       new TColStd_HArray1OfBoolean(myPoints->Lower(),
384                                     myPoints->Upper()) ;
385  
386  if (!result) {
387    throw Standard_ConstructionError();
388    }
389                 
390  result =
391  CheckParameters(ParametersPtr->Array1()) ;
392  if (!result) {
393    throw Standard_ConstructionError();
394    }
395         
396  for (ii = myPoints->Lower() ; ii <= myPoints->Upper() ; ii++) {
397    myTangentFlags->SetValue(ii,Standard_False) ;
398  }
399          
400 }
401 //=======================================================================
402 //function : Load
403 //purpose  : 
404 //=======================================================================
405
406 void Geom2dAPI_Interpolate::Load( 
407    const TColgp_Array1OfVec2d&              Tangents,
408    const Handle(TColStd_HArray1OfBoolean)& TangentFlagsPtr,
409    const Standard_Boolean Scale)
410
411 {
412  Standard_Boolean result ;
413  Standard_Integer ii ;
414  myTangentRequest = Standard_True ;
415  myTangentFlags = TangentFlagsPtr ;
416  if (Tangents.Length() != myPoints->Length() ||
417      TangentFlagsPtr->Length() != myPoints->Length()) {
418    throw Standard_ConstructionError();
419    }
420  result  = 
421    CheckTangents(Tangents,
422                  TangentFlagsPtr->Array1(),
423                  myTolerance) ;
424  if (result) {
425     myTangents = 
426       new TColgp_HArray1OfVec2d(Tangents.Lower(),Tangents.Upper()) ;
427     for (ii = Tangents.Lower() ; ii <= Tangents.Upper() ; ii++ ) {
428       myTangents->SetValue(ii,Tangents.Value(ii)) ;
429     }
430     if (Scale) {
431       ScaleTangents(myPoints->Array1(),
432                     myTangents->ChangeArray1(),
433                     TangentFlagsPtr->Array1(),
434                     myParameters->Array1()) ;
435     } 
436   }
437  else {
438    throw Standard_ConstructionError();
439    }
440  
441   
442 }
443
444 //=======================================================================
445 //function : Load
446 //purpose  : 
447 //=======================================================================
448
449 void Geom2dAPI_Interpolate::Load(const gp_Vec2d& InitialTangent,
450                                const gp_Vec2d& FinalTangent,
451              const Standard_Boolean Scale)
452 {
453   Standard_Boolean result ;
454   myTangentRequest = Standard_True ;
455   myTangentFlags->SetValue(1,Standard_True) ;
456   myTangentFlags->SetValue(myPoints->Length(),Standard_True) ;
457   myTangents->SetValue(1,InitialTangent) ;
458   myTangents->SetValue(myPoints->Length(),FinalTangent);
459   result = 
460     CheckTangents(myTangents->Array1(),
461                   myTangentFlags->Array1(),
462                   myTolerance) ;
463   if (!result) {
464     throw Standard_ConstructionError();
465     }
466   if (Scale) {
467     ScaleTangents(myPoints->Array1(),
468                   myTangents->ChangeArray1(),
469                   myTangentFlags->Array1(),
470                   myParameters->Array1()) ;
471   }
472  
473 }
474 //=======================================================================
475 //function : Perform
476 //purpose  : 
477 //=======================================================================
478
479 void Geom2dAPI_Interpolate::Perform() 
480 {
481   if (myPeriodic) {
482     PerformPeriodic() ;
483   }
484   else {
485     PerformNonPeriodic() ;
486   }
487 }
488 //=======================================================================
489 //function : PerformPeriodic
490 //purpose  : 
491 //=======================================================================
492
493 void Geom2dAPI_Interpolate::PerformPeriodic()
494
495   Standard_Integer degree,
496   ii,
497   jj,
498   index,
499   index1,
500 //  index2,
501   mult_index,
502   half_order,
503   inversion_problem,
504   num_points,
505   num_distinct_knots,
506   num_poles  ;
507  
508   Standard_Real period ;
509
510   gp_Pnt2d a_point ;
511
512   num_points = myPoints->Length() ;
513   period = myParameters->Value(myParameters->Upper()) -
514    myParameters->Value(myParameters->Lower())  ;
515   num_poles = num_points + 1 ;
516   if (num_points == 2 && !myTangentRequest) {
517 //
518 // build a periodic curve of degree 1
519 //  
520
521     degree = 1 ;
522     TColStd_Array1OfInteger  deg1_mults(1,num_poles) ;    
523     for (ii = 1 ; ii <= num_poles ; ii++) {
524       deg1_mults.SetValue(ii,1) ;
525     }
526
527     myCurve =
528       new Geom2d_BSplineCurve(myPoints->Array1(),
529                             myParameters->Array1(),
530                             deg1_mults,
531                             degree,
532                             myPeriodic) ;
533     myIsDone = Standard_True ;
534
535   }
536   else {
537     num_distinct_knots = num_points + 1 ;
538     half_order = 2 ;  
539     degree = 3 ;
540     num_poles += 2 ;
541     if (myTangentRequest) 
542       for (ii = myTangentFlags->Lower() + 1 ; 
543            ii <= myTangentFlags->Upper() ; ii++) {
544         if (myTangentFlags->Value(ii)) {
545           num_poles += 1 ;
546         }
547       }
548   
549     TColStd_Array1OfReal     parameters(1,num_poles) ;  
550     TColStd_Array1OfReal     flatknots(1,num_poles + degree + 1) ;
551     TColStd_Array1OfInteger  mults(1,num_distinct_knots) ;
552     TColStd_Array1OfInteger  contact_order_array(1, num_poles) ;
553     TColgp_Array1OfPnt2d       poles(1,num_poles) ;
554
555     for (ii = 1 ; ii <= half_order ; ii++) {
556       flatknots.SetValue(ii,myParameters->Value(myParameters->Upper() -1) -
557                          period) ;
558       flatknots.SetValue(ii + half_order,myParameters->
559                          Value(myParameters->Lower())) ;
560       flatknots.SetValue(num_poles + ii,
561                          myParameters->Value(myParameters->Upper())) ;
562       flatknots.SetValue(num_poles + half_order + ii,
563                          myParameters->Value(half_order) + period)  ;
564     }
565     for (ii = 1 ; ii <= num_poles ; ii++) {
566       contact_order_array.SetValue(ii,0)  ;
567     }
568     for (ii = 2 ; ii < num_distinct_knots ; ii++) {
569       mults.SetValue(ii,1) ; 
570     }
571     mults.SetValue(1,half_order) ;
572     mults.SetValue(num_distinct_knots ,half_order) ;
573     if (num_points >= 3) {
574      
575 //
576 //   only enter here if there are more than 3 points otherwise
577 //   it means we have already the tangent
578 // 
579       BuildPeriodicTangent(myPoints->Array1(),
580                            myTangents->ChangeArray1(),
581                            myTangentFlags->ChangeArray1(),
582                            myParameters->Array1()) ;
583     }
584     contact_order_array.SetValue(2,1)  ;
585     parameters.SetValue(1,myParameters->Value(1)) ;
586     parameters.SetValue(2,myParameters->Value(1)) ;
587     poles.SetValue(1,myPoints->Value(1)) ;
588     for (jj = 1 ; jj <= 2 ; jj++) {
589       a_point.SetCoord(jj,myTangents->Value(1).Coord(jj)) ;
590     }
591     poles.SetValue(2,a_point) ;
592   
593     mult_index = 2 ;
594     index = 3 ;
595     index1 = degree + 2 ;
596     if (myTangentRequest) {
597       for (ii = myTangentFlags->Lower() + 1 ; 
598            ii <= myTangentFlags->Upper() ; ii++) {
599         parameters.SetValue(index,myParameters->Value(ii)) ;
600         flatknots.SetValue(index1,myParameters->Value(ii)) ;
601         poles.SetValue(index,myPoints->Value(ii)) ;
602         index += 1  ;
603         index1 += 1 ;
604         if (myTangentFlags->Value(ii)) {
605           mults.SetValue(mult_index,mults.Value(mult_index)  + 1)   ;
606           contact_order_array(index) = 1 ;
607
608           parameters.SetValue(index,
609                               myParameters->Value(ii)) ;
610           flatknots.SetValue(index1,myParameters->Value(ii)) ;
611           for (jj = 1 ; jj <= 2 ; jj++) {
612             a_point.SetCoord(jj,myTangents->Value(ii).Coord(jj)) ;
613           }
614           poles.SetValue(index,a_point) ;
615           index  += 1 ;
616           index1 += 1 ;
617         }
618        mult_index += 1 ;
619       }
620     }
621     else {
622       index = degree + 1 ;
623       index1 = 2 ;
624       for(ii = myParameters->Lower()  ; ii <= myParameters->Upper()  ; ii++) {
625         parameters.SetValue(index1, 
626                             myParameters->Value(ii)) ;
627         flatknots.SetValue(index,
628                            myParameters->Value(ii)) ;
629         index += 1 ;
630         index1 += 1 ;
631       }
632       index = 3 ;
633       for (ii = myPoints->Lower() + 1 ; ii <= myPoints->Upper() ; ii++) {
634 //
635 // copy all the given points since the last one will be initialized
636 // below by the first point in the array myPoints
637 //
638         poles.SetValue(index,
639                        myPoints->Value(ii)) ;
640         index += 1 ;
641       }
642    
643     }
644     contact_order_array.SetValue(num_poles - 1, 1)  ;
645     parameters.SetValue(num_poles-1,
646                         myParameters->Value(myParameters->Upper())) ;
647 //
648 // for the periodic curve ONLY  the tangent of the first point
649 // will be used since the curve should close itself at the first
650 // point See BuildPeriodicTangent
651 //
652     for (jj = 1 ; jj <= 2 ; jj++) {
653       a_point.SetCoord(jj,myTangents->Value(1).Coord(jj)) ;
654     }
655     poles.SetValue(num_poles-1,a_point) ;
656
657     parameters.SetValue(num_poles,
658                         myParameters->Value(myParameters->Upper())) ;
659
660     poles.SetValue(num_poles,
661                      myPoints->Value(1)) ;
662
663     
664     BSplCLib::Interpolate(degree,
665                           flatknots,
666                           parameters,
667                           contact_order_array,
668                           poles,
669                           inversion_problem) ;
670     if (!inversion_problem) {
671       TColgp_Array1OfPnt2d   newpoles(poles.Value(1),
672                                     1,
673                                     num_poles - 2) ;
674       myCurve =
675         new Geom2d_BSplineCurve(newpoles,
676                               myParameters->Array1(),
677                               mults,
678                               degree,
679                               myPeriodic) ;
680       myIsDone = Standard_True ;
681   }
682  }
683 }
684    
685
686 //=======================================================================
687 //function : PerformNonPeriodic
688 //purpose  : 
689 //=======================================================================
690
691 void Geom2dAPI_Interpolate::PerformNonPeriodic() 
692 {
693   Standard_Integer degree,
694   ii,
695   jj,
696   index,
697   index1,
698   index2,
699   index3,
700   mult_index,
701   inversion_problem,
702   num_points,
703   num_distinct_knots,
704   num_poles  ;
705   
706   gp_Pnt2d a_point ;
707
708   num_points =
709   num_distinct_knots =
710   num_poles = myPoints->Length() ;
711   if (num_poles == 2 &&   !myTangentRequest)  {
712     degree = 1 ;
713   } 
714   else if (num_poles == 3 && !myTangentRequest) {
715     degree = 2 ;
716     num_distinct_knots = 2 ;
717   }
718   else {
719     degree = 3 ;
720     num_poles += 2 ;
721     if (myTangentRequest) 
722       for (ii = myTangentFlags->Lower() + 1 ; 
723            ii < myTangentFlags->Upper() ; ii++) {
724         if (myTangentFlags->Value(ii)) {
725           num_poles += 1 ;
726         }
727       }
728     }
729   
730   
731   TColStd_Array1OfReal     parameters(1,num_poles) ;  
732   TColStd_Array1OfReal     flatknots(1,num_poles + degree + 1) ;
733   TColStd_Array1OfInteger  mults(1,num_distinct_knots) ;
734   TColStd_Array1OfReal     knots(1,num_distinct_knots) ;
735   TColStd_Array1OfInteger  contact_order_array(1, num_poles) ;
736   TColgp_Array1OfPnt2d       poles(1,num_poles) ;
737
738   for (ii = 1 ; ii <= degree + 1 ; ii++) {
739     flatknots.SetValue(ii,myParameters->Value(1)) ;
740     flatknots.SetValue(ii + num_poles,
741                        myParameters->Value(num_points)) ;
742   }
743   for (ii = 1 ; ii <= num_poles ; ii++) {
744     contact_order_array.SetValue(ii,0)  ;
745   }
746   for (ii = 2 ; ii < num_distinct_knots ; ii++) {
747     mults.SetValue(ii,1) ; 
748   }
749   mults.SetValue(1,degree + 1) ;
750   mults.SetValue(num_distinct_knots ,degree + 1) ;
751   
752   switch (degree) {
753   case 1:
754     for (ii = 1 ; ii <= num_poles ; ii++) {
755       poles.SetValue(ii ,myPoints->Value(ii)) ;
756     }
757     myCurve =
758       new Geom2d_BSplineCurve(poles,
759                             myParameters->Array1(),
760                             mults,
761                             degree) ;
762     myIsDone = Standard_True ;
763     break ;
764   case 2:
765     knots.SetValue(1,myParameters->Value(1)) ;
766     knots.SetValue(2,myParameters->Value(3)) ;
767     for (ii = 1 ; ii <= num_poles ; ii++) {
768       poles.SetValue(ii,myPoints->Value(ii)) ;
769       
770     }
771     BSplCLib::Interpolate(degree,
772                           flatknots,
773                           myParameters->Array1(),
774                           contact_order_array,
775                           poles,
776                           inversion_problem) ;
777     if (!inversion_problem) {
778       myCurve =
779         new Geom2d_BSplineCurve(poles,
780                               knots,
781                               mults,
782                               degree) ;
783       myIsDone = Standard_True ;
784     }
785     break ;
786   case 3:
787 //
788 // check if the boundary conditions are set
789 //
790     if (num_points >= 3) {
791 //
792 // cannot build the tangents with degree 3 with only 2 points
793 // if those where not given in advance 
794 //
795       BuildTangents(myPoints->Array1(),
796                     myTangents->ChangeArray1(),
797                     myTangentFlags->ChangeArray1(),
798                     myParameters->Array1()) ;
799     }
800     contact_order_array.SetValue(2,1)  ;
801     parameters.SetValue(1,myParameters->Value(1)) ;
802     parameters.SetValue(2,myParameters->Value(1)) ;
803     poles.SetValue(1,myPoints->Value(1)) ;
804     for (jj = 1 ; jj <= 2 ; jj++) {
805       a_point.SetCoord(jj,myTangents->Value(1).Coord(jj)) ;
806
807     }
808     poles.SetValue(2,a_point) ;
809     mult_index = 2 ;
810     index = 3 ;
811     index1 = 2 ;
812     index2 = myPoints->Lower() + 1 ;
813     index3 = degree + 2 ;
814     if (myTangentRequest) {
815       for (ii = myParameters->Lower() + 1 ; 
816            ii < myParameters->Upper() ; ii++) {
817         parameters.SetValue(index,myParameters->Value(ii)) ;
818         poles.SetValue(index,myPoints->Value(index2)) ;
819         flatknots.SetValue(index3,myParameters->Value(ii)) ; 
820         index += 1 ;
821         index3 += 1 ;
822         if (myTangentFlags->Value(index1)) {
823 //
824 // set the multiplicities, the order of the contact, the 
825 // the flatknots,
826 //
827           mults.SetValue(mult_index,mults.Value(mult_index) + 1) ;
828           contact_order_array(index) = 1 ;
829           flatknots.SetValue(index3, myParameters->Value(ii)) ;
830           parameters.SetValue(index,
831                               myParameters->Value(ii)) ;
832           for (jj = 1 ; jj <= 2 ; jj++) {
833             a_point.SetCoord(jj,myTangents->Value(ii).Coord(jj)) ;
834           }
835           poles.SetValue(index,a_point) ;
836           index += 1  ;
837           index3 += 1 ;
838         }
839         mult_index += 1 ;
840         index1 += 1 ;
841         index2 += 1 ;
842
843       }
844     }
845     else {
846       index1 = 2 ;
847       for(ii = myParameters->Lower()  ; ii <= myParameters->Upper()  ; ii++) {
848         parameters.SetValue(index1, 
849                             myParameters->Value(ii)) ;
850         index1 += 1 ;
851       }
852       index = 3 ;
853       for (ii = myPoints->Lower() + 1 ; ii <= myPoints->Upper() - 1 ; ii++) {
854         poles.SetValue(index,
855                        myPoints->Value(ii)) ;
856         index += 1 ;
857       }
858    
859    
860       index = degree + 1 ;
861       for(ii = myParameters->Lower()  ; ii <= myParameters->Upper()  ; ii++) {
862         flatknots.SetValue(index,
863                            myParameters->Value(ii)) ;
864         index += 1 ;
865       }
866     }
867     for (jj = 1 ; jj <= 2 ; jj++) {
868       a_point.SetCoord(jj,
869                        myTangents->Value(num_points).Coord(jj)) ;
870     }
871     poles.SetValue(num_poles-1 ,a_point) ;
872
873     contact_order_array.SetValue(num_poles - 1,1) ;
874     parameters.SetValue(num_poles,
875                         myParameters->Value(myParameters->Upper())) ;
876     parameters.SetValue(num_poles -1,
877                         myParameters->Value(myParameters->Upper())) ;
878
879     poles.SetValue(num_poles,
880                    myPoints->Value(num_points)) ;
881
882     BSplCLib::Interpolate(degree,
883                           flatknots,
884                           parameters,
885                           contact_order_array,
886                           poles,
887                           inversion_problem) ;
888     if (!inversion_problem) {
889       myCurve =
890         new Geom2d_BSplineCurve(poles,
891                               myParameters->Array1(),
892                               mults,
893                               degree) ;
894       myIsDone = Standard_True ;
895     }
896     break ;
897  
898   }
899 }
900 //=======================================================================
901 //function : Handle(Geom2d_BSplineCurve)&
902 //purpose  : 
903 //=======================================================================
904
905 const Handle(Geom2d_BSplineCurve)& Geom2dAPI_Interpolate::Curve() const 
906 {
907   if ( !myIsDone) 
908     throw StdFail_NotDone(" ");
909   return myCurve;
910 }
911
912
913
914 //=======================================================================
915 //function : Geom2d_BSplineCurve
916 //purpose  : 
917 //=======================================================================
918
919 Geom2dAPI_Interpolate::operator Handle(Geom2d_BSplineCurve)() const
920 {
921   return myCurve;
922 }
923
924
925 //=======================================================================
926 //function : IsDone
927 //purpose  : 
928 //=======================================================================
929
930 Standard_Boolean Geom2dAPI_Interpolate::IsDone() const
931 {
932   return myIsDone;
933 }