1 // Created on: 1995-01-17
2 // Created by: Laurent BOURESCHE
3 // Copyright (c) 1995-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 // xab : modified 15-Mar-95 : added BuildBSplMatrix,FactorBandedMatrix,
18 // EvalBsplineBasis,
19 // EvalPolynomial : Horners method
20
21 #include <Standard_Stream.hxx>
22
23 #include <BSplCLib.hxx>
24 #include <gp_Mat2d.hxx>
25 #include <PLib.hxx>
26 #include <Standard_ConstructionError.hxx>
27 #include <TColStd_Array1OfReal.hxx>
28 #include <TColStd_Array1OfInteger.hxx>
29 #include <TColStd_HArray1OfReal.hxx>
30 #include <TColStd_HArray1OfInteger.hxx>
31
32 #include <math_Matrix.hxx>
33
34 //=======================================================================
35 //struct : BSplCLib_DataContainer
36 //purpose: Auxiliary structure providing buffers for poles and knots used in
37 // evaluation of bspline (allocated in the stack)
38 //=======================================================================
39
40 struct BSplCLib_DataContainer
41 {
BSplCLib_DataContainerBSplCLib_DataContainer42 BSplCLib_DataContainer(Standard_Integer Degree)
43 {
44 (void)Degree; // avoid compiler warning
45 Standard_OutOfRange_Raise_if (Degree > BSplCLib::MaxDegree(),
46 "BSplCLib: bspline degree is greater than maximum supported");
47 }
48
49 Standard_Real poles[2*(25+1)];
50 Standard_Real knots[2*25];
51 Standard_Real ders[4];
52 };
53
54 // methods for 1 dimensional BSplines
55
56 //=======================================================================
57 //function : BuildEval
58 //purpose : builds the local array for evaluation
59 //=======================================================================
60
BuildEval(const Standard_Integer Degree,const Standard_Integer Index,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal * Weights,Standard_Real & LP)61 void BSplCLib::BuildEval(const Standard_Integer Degree,
62 const Standard_Integer Index,
63 const TColStd_Array1OfReal& Poles,
64 const TColStd_Array1OfReal* Weights,
65 Standard_Real& LP)
66 {
67 Standard_Integer PLower = Poles.Lower();
68 Standard_Integer PUpper = Poles.Upper();
69 Standard_Integer i;
70 Standard_Integer ip = PLower + Index - 1;
71 Standard_Real w, *pole = &LP;
72 if (Weights == NULL) {
73
74 for (i = 0; i <= Degree; i++) {
75 ip++;
76 if (ip > PUpper) ip = PLower;
77 pole[0] = Poles(ip);
78 pole += 1;
79 }
80 }
81 else {
82
83 for (i = 0; i <= Degree; i++) {
84 ip++;
85 if (ip > PUpper) ip = PLower;
86 pole[1] = w = (*Weights)(ip);
87 pole[0] = Poles(ip) * w;
88 pole += 2;
89 }
90 }
91 }
92
93 //=======================================================================
94 //function : PrepareEval
95 //purpose : stores data for Eval in the local arrays
96 // dc.poles and dc.knots
97 //=======================================================================
98
PrepareEval(Standard_Real & u,Standard_Integer & index,Standard_Integer & dim,Standard_Boolean & rational,const Standard_Integer Degree,const Standard_Boolean Periodic,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal * Weights,const TColStd_Array1OfReal & Knots,const TColStd_Array1OfInteger * Mults,BSplCLib_DataContainer & dc)99 static void PrepareEval
100 (Standard_Real& u,
101 Standard_Integer& index,
102 Standard_Integer& dim,
103 Standard_Boolean& rational,
104 const Standard_Integer Degree,
105 const Standard_Boolean Periodic,
106 const TColStd_Array1OfReal& Poles,
107 const TColStd_Array1OfReal* Weights,
108 const TColStd_Array1OfReal& Knots,
109 const TColStd_Array1OfInteger* Mults,
110 BSplCLib_DataContainer& dc)
111 {
112 // Set the Index
113 BSplCLib::LocateParameter(Degree,Knots,Mults,u,Periodic,index,u);
114
115 // make the knots
116 BSplCLib::BuildKnots(Degree,index,Periodic,Knots,Mults,*dc.knots);
117 if (Mults == NULL)
118 index -= Knots.Lower() + Degree;
119 else
120 index = BSplCLib::PoleIndex(Degree,index,Periodic,*Mults);
121
122 // check truly rational
123 rational = (Weights != NULL);
124 if (rational) {
125 Standard_Integer WLower = Weights->Lower() + index;
126 rational = BSplCLib::IsRational(*Weights, WLower, WLower + Degree);
127 }
128
129 // make the poles
130 if(rational) {
131 dim = 2;
132 BSplCLib::BuildEval(Degree, index, Poles, Weights , *dc.poles);
133 }
134 else {
135 dim = 1;
136 BSplCLib::BuildEval(Degree, index, Poles, BSplCLib::NoWeights(), *dc.poles);
137 }
138 }
139
140 //=======================================================================
141 //function : D0
142 //purpose :
143 //=======================================================================
144
D0(const Standard_Real U,const Standard_Integer Index,const Standard_Integer Degree,const Standard_Boolean Periodic,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal * Weights,const TColStd_Array1OfReal & Knots,const TColStd_Array1OfInteger * Mults,Standard_Real & P)145 void BSplCLib::D0
146 (const Standard_Real U,
147 const Standard_Integer Index,
148 const Standard_Integer Degree,
149 const Standard_Boolean Periodic,
150 const TColStd_Array1OfReal& Poles,
151 const TColStd_Array1OfReal* Weights,
152 const TColStd_Array1OfReal& Knots,
153 const TColStd_Array1OfInteger* Mults,
154 Standard_Real& P)
155 {
156 Standard_Integer dim,index = Index;
157 Standard_Real u = U;
158 Standard_Boolean rational;
159 BSplCLib_DataContainer dc(Degree);
160 PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc);
161 BSplCLib::Eval(u,Degree,*dc.knots,dim,*dc.poles);
162 if (rational) P = dc.poles[0] / dc.poles[1];
163 else P = dc.poles[0];
164 }
165
166 //=======================================================================
167 //function : D1
168 //purpose :
169 //=======================================================================
170
D1(const Standard_Real U,const Standard_Integer Index,const Standard_Integer Degree,const Standard_Boolean Periodic,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal * Weights,const TColStd_Array1OfReal & Knots,const TColStd_Array1OfInteger * Mults,Standard_Real & P,Standard_Real & V)171 void BSplCLib::D1
172 (const Standard_Real U,
173 const Standard_Integer Index,
174 const Standard_Integer Degree,
175 const Standard_Boolean Periodic,
176 const TColStd_Array1OfReal& Poles,
177 const TColStd_Array1OfReal* Weights,
178 const TColStd_Array1OfReal& Knots,
179 const TColStd_Array1OfInteger* Mults,
180 Standard_Real& P,
181 Standard_Real& V)
182 {
183 Standard_Integer dim,index = Index;
184 Standard_Real u = U;
185 Standard_Boolean rational;
186 BSplCLib_DataContainer dc(Degree);
187 PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc);
188 BSplCLib::Bohm(u,Degree,1,*dc.knots,dim,*dc.poles);
189 Standard_Real *result = dc.poles;
190 if (rational) {
191 PLib::RationalDerivative(Degree,1,1,*dc.poles,*dc.ders);
192 result = dc.ders;
193 }
194 P = result[0];
195 V = result[1];
196 }
197
198 //=======================================================================
199 //function : D2
200 //purpose :
201 //=======================================================================
202
D2(const Standard_Real U,const Standard_Integer Index,const Standard_Integer Degree,const Standard_Boolean Periodic,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal * Weights,const TColStd_Array1OfReal & Knots,const TColStd_Array1OfInteger * Mults,Standard_Real & P,Standard_Real & V1,Standard_Real & V2)203 void BSplCLib::D2
204 (const Standard_Real U,
205 const Standard_Integer Index,
206 const Standard_Integer Degree,
207 const Standard_Boolean Periodic,
208 const TColStd_Array1OfReal& Poles,
209 const TColStd_Array1OfReal* Weights,
210 const TColStd_Array1OfReal& Knots,
211 const TColStd_Array1OfInteger* Mults,
212 Standard_Real& P,
213 Standard_Real& V1,
214 Standard_Real& V2)
215 {
216 Standard_Integer dim,index = Index;
217 Standard_Real u = U;
218 Standard_Boolean rational;
219 BSplCLib_DataContainer dc(Degree);
220 PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc);
221 BSplCLib::Bohm(u,Degree,2,*dc.knots,dim,*dc.poles);
222 Standard_Real *result = dc.poles;
223 if (rational) {
224 PLib::RationalDerivative(Degree,2,1,*dc.poles,*dc.ders);
225 result = dc.ders;
226 }
227 P = result[0];
228 V1 = result[1];
229 if (!rational && (Degree < 2)) V2 = 0.;
230 else V2 = result[2];
231 }
232
233 //=======================================================================
234 //function : D3
235 //purpose :
236 //=======================================================================
237
D3(const Standard_Real U,const Standard_Integer Index,const Standard_Integer Degree,const Standard_Boolean Periodic,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal * Weights,const TColStd_Array1OfReal & Knots,const TColStd_Array1OfInteger * Mults,Standard_Real & P,Standard_Real & V1,Standard_Real & V2,Standard_Real & V3)238 void BSplCLib::D3
239 (const Standard_Real U,
240 const Standard_Integer Index,
241 const Standard_Integer Degree,
242 const Standard_Boolean Periodic,
243 const TColStd_Array1OfReal& Poles,
244 const TColStd_Array1OfReal* Weights,
245 const TColStd_Array1OfReal& Knots,
246 const TColStd_Array1OfInteger* Mults,
247 Standard_Real& P,
248 Standard_Real& V1,
249 Standard_Real& V2,
250 Standard_Real& V3)
251 {
252 Standard_Integer dim,index = Index;
253 Standard_Real u = U;
254 Standard_Boolean rational;
255 BSplCLib_DataContainer dc(Degree);
256 PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc);
257 BSplCLib::Bohm(u,Degree,3,*dc.knots,dim,*dc.poles);
258 Standard_Real *result = dc.poles;
259 if (rational) {
260 PLib::RationalDerivative(Degree,3,1,*dc.poles,*dc.ders);
261 result = dc.ders;
262 }
263 P = result[0];
264 V1 = result[1];
265 if (!rational && (Degree < 2)) V2 = 0.;
266 else V2 = result[2];
267 if (!rational && (Degree < 3)) V3 = 0.;
268 else V3 = result[3];
269 }
270
271 //=======================================================================
272 //function : DN
273 //purpose :
274 //=======================================================================
275
DN(const Standard_Real U,const Standard_Integer N,const Standard_Integer Index,const Standard_Integer Degree,const Standard_Boolean Periodic,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal * Weights,const TColStd_Array1OfReal & Knots,const TColStd_Array1OfInteger * Mults,Standard_Real & VN)276 void BSplCLib::DN
277 (const Standard_Real U,
278 const Standard_Integer N,
279 const Standard_Integer Index,
280 const Standard_Integer Degree,
281 const Standard_Boolean Periodic,
282 const TColStd_Array1OfReal& Poles,
283 const TColStd_Array1OfReal* Weights,
284 const TColStd_Array1OfReal& Knots,
285 const TColStd_Array1OfInteger* Mults,
286 Standard_Real& VN)
287 {
288 Standard_Integer dim,index = Index;
289 Standard_Real u = U;
290 Standard_Boolean rational;
291 BSplCLib_DataContainer dc(Degree);
292 PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc);
293 BSplCLib::Bohm(u,Degree,N,*dc.knots,dim,*dc.poles);
294 if (rational) {
295 Standard_Real v;
296 PLib::RationalDerivative(Degree,N,1,*dc.poles,v,Standard_False);
297 VN = v;
298 }
299 else {
300 if (N > Degree) VN = 0.;
301 else VN = dc.poles[N];
302 }
303 }
304
305 //=======================================================================
306 //function : Build BSpline Matrix
307 //purpose : Builds the Bspline Matrix
308 //=======================================================================
309
310 Standard_Integer
BuildBSpMatrix(const TColStd_Array1OfReal & Parameters,const TColStd_Array1OfInteger & ContactOrderArray,const TColStd_Array1OfReal & FlatKnots,const Standard_Integer Degree,math_Matrix & Matrix,Standard_Integer & UpperBandWidth,Standard_Integer & LowerBandWidth)311 BSplCLib::BuildBSpMatrix(const TColStd_Array1OfReal& Parameters,
312 const TColStd_Array1OfInteger& ContactOrderArray,
313 const TColStd_Array1OfReal& FlatKnots,
314 const Standard_Integer Degree,
315 math_Matrix& Matrix,
316 Standard_Integer& UpperBandWidth,
317 Standard_Integer& LowerBandWidth)
318 {
319 Standard_Integer ii,
320 jj,
321 Index,
322 ErrorCode,
323 ReturnCode = 0,
324 FirstNonZeroBsplineIndex,
325 BandWidth,
326 MaxOrder = BSplCLib::MaxDegree() + 1,
327 Order ;
328
329 math_Matrix BSplineBasis(1, MaxOrder,
330 1, MaxOrder) ;
331
332 Order = Degree + 1 ;
333 UpperBandWidth = Degree ;
334 LowerBandWidth = Degree ;
335 BandWidth = UpperBandWidth + LowerBandWidth + 1 ;
336 if (Matrix.LowerRow() != Parameters.Lower() ||
337 Matrix.UpperRow() != Parameters.Upper() ||
338 Matrix.LowerCol() != 1 ||
339 Matrix.UpperCol() != BandWidth) {
340 ReturnCode = 1;
341 goto FINISH ;
342 }
343
344 for (ii = Parameters.Lower() ; ii <= Parameters.Upper() ; ii++) {
345 ErrorCode =
346 BSplCLib::EvalBsplineBasis(ContactOrderArray(ii),
347 Order,
348 FlatKnots,
349 Parameters(ii),
350
351 FirstNonZeroBsplineIndex,
352 BSplineBasis) ;
353 if (ErrorCode != 0) {
354 ReturnCode = 2 ;
355 goto FINISH ;
356 }
357 Index = LowerBandWidth + 1 + FirstNonZeroBsplineIndex - ii ;
358
359 for (jj = 1 ; jj < Index ; jj++) {
360 Matrix.Value(ii,jj) = 0.0e0 ;
361 }
362
363 for (jj = 1 ; jj <= Order ; jj++) {
364 Matrix.Value(ii,Index) = BSplineBasis(ContactOrderArray(ii) + 1, jj) ;
365 Index += 1 ;
366 }
367
368 for (jj = Index ; jj <= BandWidth ; jj++) {
369 Matrix.Value(ii,jj) = 0.0e0 ;
370 }
371 }
372 FINISH : ;
373 return (ReturnCode) ;
374 }
375
376 //=======================================================================
377 //function : Makes LU decompositiomn without Pivoting
378 //purpose : Builds the Bspline Matrix
379 //=======================================================================
380
381 Standard_Integer
FactorBandedMatrix(math_Matrix & Matrix,const Standard_Integer UpperBandWidth,const Standard_Integer LowerBandWidth,Standard_Integer & PivotIndexProblem)382 BSplCLib::FactorBandedMatrix(math_Matrix& Matrix,
383 const Standard_Integer UpperBandWidth,
384 const Standard_Integer LowerBandWidth,
385 Standard_Integer& PivotIndexProblem)
386 {
387 Standard_Integer ii,
388 jj,
389 kk,
390 Index,
391 MinIndex,
392 MaxIndex,
393 ReturnCode = 0,
394 BandWidth = UpperBandWidth + LowerBandWidth + 1 ;
395
396 Standard_Real Inverse ;
397 PivotIndexProblem = 0 ;
398
399 for (ii = Matrix.LowerRow() + 1 ; ii <= Matrix.UpperRow() ; ii++) {
400 MinIndex = ( LowerBandWidth - ii + 2 >= 1 ? LowerBandWidth - ii + 2 : 1) ;
401
402 for (jj = MinIndex ; jj <= LowerBandWidth ; jj++) {
403 Index = ii - LowerBandWidth + jj - 1 ;
404 Inverse = Matrix(Index ,LowerBandWidth + 1) ;
405 if (Abs(Inverse) > RealSmall()) {
406 Inverse = -1.0e0/Inverse ;
407 }
408 else {
409 ReturnCode = 1 ;
410 PivotIndexProblem = Index ;
411 goto FINISH ;
412 }
413 Matrix(ii,jj) = Matrix(ii,jj) * Inverse ;
414 MaxIndex = BandWidth + Index - ii ;
415
416 for (kk = jj + 1 ; kk <= MaxIndex ; kk++) {
417 Matrix(ii,kk) += Matrix(ii,jj) * Matrix(Index, kk + ii - Index) ;
418 }
419 }
420 }
421 FINISH :
422 return (ReturnCode) ;
423 }
424
425 //=======================================================================
426 //function : Build BSpline Matrix
427 //purpose : Builds the Bspline Matrix
428 //=======================================================================
429
430 Standard_Integer
EvalBsplineBasis(const Standard_Integer DerivativeRequest,const Standard_Integer Order,const TColStd_Array1OfReal & FlatKnots,const Standard_Real Parameter,Standard_Integer & FirstNonZeroBsplineIndex,math_Matrix & BsplineBasis,Standard_Boolean isPeriodic)431 BSplCLib::EvalBsplineBasis
432 (const Standard_Integer DerivativeRequest,
433 const Standard_Integer Order,
434 const TColStd_Array1OfReal& FlatKnots,
435 const Standard_Real Parameter,
436 Standard_Integer& FirstNonZeroBsplineIndex,
437 math_Matrix& BsplineBasis,
438 Standard_Boolean isPeriodic)
439 {
440 // the matrix must have at least DerivativeRequest + 1
441 // row and Order columns
442 // the result are stored in the following way in
443 // the Bspline matrix
444 // Let i be the FirstNonZeroBsplineIndex and
445 // t be the parameter value, k the order of the
446 // knot vector, r the DerivativeRequest :
447 //
448 // B (t) B (t) B (t)
449 // i i+1 i+k-1
450 //
451 // (1) (1) (1)
452 // B (t) B (t) B (t)
453 // i i+1 i+k-1
454 //
455 //
456 //
457 //
458 // (r) (r) (r)
459 // B (t) B (t) B (t)
460 // i i+1 i+k-1
461 //
462 Standard_Integer
463 ReturnCode,
464 ii,
465 pp,
466 qq,
467 ss,
468 NumPoles,
469 LocalRequest ;
470 // ,Index ;
471
472 Standard_Real NewParameter,
473 Inverse,
474 Factor,
475 LocalInverse,
476 Saved ;
477 // , *FlatKnotsArray ;
478
479 ReturnCode = 0 ;
480 FirstNonZeroBsplineIndex = 0 ;
481 LocalRequest = DerivativeRequest ;
482 if (DerivativeRequest >= Order) {
483 LocalRequest = Order - 1 ;
484 }
485
486 if (BsplineBasis.LowerCol() != 1 ||
487 BsplineBasis.UpperCol() < Order ||
488 BsplineBasis.LowerRow() != 1 ||
489 BsplineBasis.UpperRow() <= LocalRequest) {
490 ReturnCode = 1;
491 goto FINISH ;
492 }
493 NumPoles = FlatKnots.Upper() - FlatKnots.Lower() + 1 - Order ;
494 BSplCLib::LocateParameter(Order - 1,
495 FlatKnots,
496 Parameter,
497 isPeriodic,
498 Order,
499 NumPoles+1,
500 ii,
501 NewParameter) ;
502
503 FirstNonZeroBsplineIndex = ii - Order + 1 ;
504
505 BsplineBasis(1,1) = 1.0e0 ;
506 LocalRequest = DerivativeRequest ;
507 if (DerivativeRequest >= Order) {
508 LocalRequest = Order - 1 ;
509 }
510
511 for (qq = 2 ; qq <= Order - LocalRequest ; qq++) {
512 BsplineBasis(1,qq) = 0.0e0 ;
513
514 for (pp = 1 ; pp <= qq - 1 ; pp++) {
515 //
516 // this should be always invertible if ii is correctly computed
517 //
518 Factor = (Parameter - FlatKnots(ii - qq + pp + 1))
519 / (FlatKnots(ii + pp) - FlatKnots(ii - qq + pp + 1)) ;
520 Saved = Factor * BsplineBasis(1,pp) ;
521 BsplineBasis(1,pp) *= (1.0e0 - Factor) ;
522 BsplineBasis(1,pp) += BsplineBasis(1,qq) ;
523 BsplineBasis(1,qq) = Saved ;
524 }
525 }
526
527 for (qq = Order - LocalRequest + 1 ; qq <= Order ; qq++) {
528
529 for (pp = 1 ; pp <= qq - 1 ; pp++) {
530 BsplineBasis(Order - qq + 2,pp) = BsplineBasis(1,pp) ;
531 }
532 BsplineBasis(1,qq) = 0.0e0 ;
533
534 for (ss = Order - LocalRequest + 1 ; ss <= qq ; ss++) {
535 BsplineBasis(Order - ss + 2,qq) = 0.0e0 ;
536 }
537
538 for (pp = 1 ; pp <= qq - 1 ; pp++) {
539 Inverse = 1.0e0 / (FlatKnots(ii + pp) - FlatKnots(ii - qq + pp + 1)) ;
540 Factor = (Parameter - FlatKnots(ii - qq + pp + 1)) * Inverse ;
541 Saved = Factor * BsplineBasis(1,pp) ;
542 BsplineBasis(1,pp) *= (1.0e0 - Factor) ;
543 BsplineBasis(1,pp) += BsplineBasis(1,qq) ;
544 BsplineBasis(1,qq) = Saved ;
545 LocalInverse = (Standard_Real) (qq - 1) * Inverse ;
546
547 for (ss = Order - LocalRequest + 1 ; ss <= qq ; ss++) {
548 Saved = LocalInverse * BsplineBasis(Order - ss + 2, pp) ;
549 BsplineBasis(Order - ss + 2, pp) *= - LocalInverse ;
550 BsplineBasis(Order - ss + 2, pp) += BsplineBasis(Order - ss + 2,qq) ;
551 BsplineBasis(Order - ss + 2,qq) = Saved ;
552 }
553 }
554 }
555 FINISH :
556 return (ReturnCode) ;
557 }
558
559 //=======================================================================
560 //function : MovePointAndTangent
561 //purpose :
562 //=======================================================================
563
MovePointAndTangent(const Standard_Real U,const Standard_Integer ArrayDimension,Standard_Real & Delta,Standard_Real & DeltaDerivatives,const Standard_Real Tolerance,const Standard_Integer Degree,const Standard_Integer StartingCondition,const Standard_Integer EndingCondition,Standard_Real & Poles,const TColStd_Array1OfReal * Weights,const TColStd_Array1OfReal & FlatKnots,Standard_Real & NewPoles,Standard_Integer & ErrorStatus)564 void BSplCLib::MovePointAndTangent(const Standard_Real U,
565 const Standard_Integer ArrayDimension,
566 Standard_Real &Delta,
567 Standard_Real &DeltaDerivatives,
568 const Standard_Real Tolerance,
569 const Standard_Integer Degree,
570 const Standard_Integer StartingCondition,
571 const Standard_Integer EndingCondition,
572 Standard_Real& Poles,
573 const TColStd_Array1OfReal* Weights,
574 const TColStd_Array1OfReal& FlatKnots,
575 Standard_Real& NewPoles,
576 Standard_Integer& ErrorStatus)
577 {
578 Standard_Integer num_poles,
579 num_knots,
580 ii,
581 jj,
582 conditions,
583 start_num_poles,
584 end_num_poles,
585 index,
586 start_index,
587 end_index,
588 other_index,
589 type,
590 order ;
591
592 Standard_Real new_parameter,
593 value,
594 divide,
595 end_value,
596 start_value,
597 *poles_array,
598 *new_poles_array,
599 *delta_array,
600 *derivatives_array,
601 *weights_array ;
602
603 ErrorStatus = 0 ;
604 weights_array = NULL ;
605 if (Weights != NULL) {
606 weights_array = const_cast<Standard_Real*>(&Weights->First());
607 }
608
609 poles_array = &Poles ;
610 new_poles_array = &NewPoles ;
611 delta_array = &Delta ;
612 derivatives_array = &DeltaDerivatives ;
613 order = Degree + 1 ;
614 num_knots = FlatKnots.Length() ;
615 num_poles = num_knots - order ;
616 conditions = StartingCondition + EndingCondition + 4 ;
617 //
618 // check validity of input data
619 //
620 if (StartingCondition >= -1 &&
621 StartingCondition <= Degree &&
622 EndingCondition >= -1 &&
623 EndingCondition <= Degree &&
624 conditions <= num_poles) {
625 //
626 // check the parameter is within bounds
627 //
628 start_index = FlatKnots.Lower() + Degree ;
629 end_index = FlatKnots.Upper() - Degree ;
630 //
631 // check if there is enough room to move the poles
632 //
633 conditions = 1 ;
634 if (StartingCondition == -1) {
635 conditions = conditions && (FlatKnots(start_index) <= U) ;
636 }
637 else {
638 conditions = conditions && (FlatKnots(start_index) + Tolerance < U) ;
639 }
640 if (EndingCondition == -1) {
641 conditions = conditions && (FlatKnots(end_index) >= U) ;
642 }
643 else {
644 conditions = conditions && (FlatKnots(end_index) - Tolerance > U) ;
645 }
646
647 if (conditions) {
648 //
649 // build 2 auxiliary functions
650 //
651 TColStd_Array1OfReal schoenberg_points(1,num_poles) ;
652 TColStd_Array1OfReal first_function (1,num_poles) ;
653 TColStd_Array1OfReal second_function (1,num_poles) ;
654
655 BuildSchoenbergPoints(Degree,
656 FlatKnots,
657 schoenberg_points) ;
658 start_index = StartingCondition + 2 ;
659 end_index = num_poles - EndingCondition - 1 ;
660 LocateParameter(schoenberg_points,
661 U,
662 Standard_False,
663 start_index,
664 end_index,
665 index,
666 new_parameter,
667 0, 1) ;
668
669 if (index == start_index) {
670 other_index = index + 1 ;
671 }
672 else if (index == end_index) {
673 other_index = index -1 ;
674 }
675 else if (U - FlatKnots(index) < FlatKnots(index + 1) - U ) {
676 other_index = index - 1 ;
677 }
678 else {
679 other_index = index + 1 ;
680 }
681 type = 3 ;
682
683 start_num_poles = StartingCondition + 2 ;
684 end_num_poles = num_poles - EndingCondition - 1 ;
685 if (start_num_poles == 1) {
686 start_value = schoenberg_points(num_poles) - schoenberg_points(1) ;
687 start_value = schoenberg_points(1) - start_value ;
688 }
689 else {
690 start_value = schoenberg_points(start_num_poles - 1) ;
691 }
692 if (end_num_poles == num_poles) {
693 end_value = schoenberg_points(num_poles) - schoenberg_points(1) ;
694 end_value = schoenberg_points(num_poles) + end_value ;
695 }
696 else {
697 end_value = schoenberg_points(end_num_poles + 1) ;
698 }
699
700 for (ii = 1 ; ii < start_num_poles ; ii++) {
701 first_function(ii) = 0.e0 ;
702 second_function(ii) = 0.0e0 ;
703 }
704
705 for (ii = end_num_poles + 1 ; ii <= num_poles ; ii++) {
706 first_function(ii) = 0.e0 ;
707 second_function(ii) = 0.0e0 ;
708 }
709 divide = 1.0e0 / (schoenberg_points(index) - start_value) ;
710
711 for (ii = start_num_poles ; ii <= index ; ii++) {
712 value = schoenberg_points(ii) - start_value ;
713 value *= divide ;
714 first_function(ii) = 1.0e0 ;
715
716 for (jj = 0 ; jj < type ; jj++) {
717 first_function(ii) *= value ;
718 }
719 }
720 divide = 1.0e0 /(end_value - schoenberg_points(index)) ;
721
722 for (ii = index ; ii <= end_num_poles ; ii++) {
723 value = end_value - schoenberg_points(ii) ;
724 value *= divide ;
725 first_function(ii) = 1.0e0 ;
726
727 for (jj = 0 ; jj < type ; jj++) {
728 first_function(ii) *= value ;
729 }
730 }
731
732 divide = 1.0e0 / (schoenberg_points(other_index) - start_value) ;
733
734 for (ii = start_num_poles ; ii <= other_index ; ii++) {
735 value = schoenberg_points(ii) - start_value ;
736 value *= divide ;
737 second_function(ii) = 1.0e0 ;
738
739 for (jj = 0 ; jj < type ; jj++) {
740 second_function(ii) *= value ;
741 }
742 }
743 divide = 1.0e0/( end_value - schoenberg_points(other_index)) ;
744
745 for (ii = other_index ; ii <= end_num_poles ; ii++) {
746 value = end_value - schoenberg_points(ii) ;
747 value *= divide ;
748 second_function(ii) = 1.0e0 ;
749
750 for (jj = 0 ; jj < type ; jj++) {
751 second_function(ii) *= value ;
752 }
753 }
754
755 //
756 // compute the point and derivatives of both functions
757 //
758 Standard_Real results[2][2],
759 weights_results[2][2];
760 Standard_Integer extrap_mode[2],
761 derivative_request = 1,
762 dimension = 1 ;
763 Standard_Boolean periodic_flag = Standard_False ;
764
765 extrap_mode[0] = Degree ;
766 extrap_mode[1] = Degree ;
767 if (Weights != NULL) {
768 //
769 // evaluate in homogenised form
770 //
771 Eval(U,
772 periodic_flag,
773 derivative_request,
774 extrap_mode[0],
775 Degree,
776 FlatKnots,
777 dimension,
778 first_function(1),
779 weights_array[0],
780 results[0][0],
781 weights_results[0][0]) ;
782
783 Eval(U,
784 periodic_flag,
785 derivative_request,
786 extrap_mode[0],
787 Degree,
788 FlatKnots,
789 dimension,
790 second_function(1),
791 weights_array[0],
792 results[1][0],
793 weights_results[1][0]) ;
794 //
795 // compute the rational derivatives values
796 //
797
798 for (ii = 0 ; ii < 2 ; ii++) {
799 PLib::RationalDerivatives(1,
800 1,
801 results[ii][0],
802 weights_results[ii][0],
803 results[ii][0]) ;
804 }
805 }
806 else {
807 Eval(U,
808 Standard_False,
809 1,
810 extrap_mode[0],
811 Degree,
812 FlatKnots,
813 1,
814 first_function(1),
815 results[0][0]) ;
816
817 Eval(U,
818 Standard_False,
819 1,
820 extrap_mode[0],
821 Degree,
822 FlatKnots,
823 1,
824 second_function(1),
825 results[1][0]) ;
826 }
827 gp_Mat2d a_matrix ;
828
829 for (ii = 0 ; ii < 2 ; ii++) {
830
831 for (jj = 0 ; jj < 2 ; jj++) {
832 a_matrix.SetValue(ii+1,jj+1,results[ii][jj]) ;
833 }
834 }
835 a_matrix.Invert() ;
836 TColStd_Array1OfReal the_a_vector(0,ArrayDimension-1) ;
837 TColStd_Array1OfReal the_b_vector(0,ArrayDimension-1) ;
838
839 for( ii = 0 ; ii < ArrayDimension ; ii++) {
840 the_a_vector(ii) =
841 a_matrix.Value(1,1) * delta_array[ii] +
842 a_matrix.Value(2,1) * derivatives_array[ii] ;
843 the_b_vector(ii) =
844 a_matrix.Value(1,2) * delta_array[ii] +
845 a_matrix.Value(2,2) * derivatives_array[ii] ;
846 }
847 index = 0 ;
848
849 for (ii = 0 ; ii < num_poles ; ii++) {
850
851 for (jj = 0 ; jj < ArrayDimension ; jj++) {
852 new_poles_array[index] = poles_array[index] ;
853 new_poles_array[index] +=
854 first_function(ii+1) * the_a_vector(jj) ;
855 new_poles_array[index] +=
856 second_function(ii+1) * the_b_vector(jj) ;
857 index += 1 ;
858 }
859 }
860 }
861 else {
862 ErrorStatus = 1 ;
863 }
864 }
865 else {
866 ErrorStatus = 2 ;
867 }
868 }
869
870 //=======================================================================
871 //function : FunctionMultiply
872 //purpose :
873 //=======================================================================
874
FunctionMultiply(const BSplCLib_EvaluatorFunction & FunctionPtr,const Standard_Integer BSplineDegree,const TColStd_Array1OfReal & BSplineFlatKnots,const Standard_Integer PolesDimension,Standard_Real & Poles,const TColStd_Array1OfReal & FlatKnots,const Standard_Integer NewDegree,Standard_Real & NewPoles,Standard_Integer & theStatus)875 void BSplCLib::FunctionMultiply
876 (const BSplCLib_EvaluatorFunction & FunctionPtr,
877 const Standard_Integer BSplineDegree,
878 const TColStd_Array1OfReal & BSplineFlatKnots,
879 const Standard_Integer PolesDimension,
880 Standard_Real & Poles,
881 const TColStd_Array1OfReal & FlatKnots,
882 const Standard_Integer NewDegree,
883 Standard_Real & NewPoles,
884 Standard_Integer & theStatus)
885 {
886 Standard_Integer ii,
887 jj,
888 index ;
889 Standard_Integer extrap_mode[2],
890 error_code,
891 num_new_poles,
892 derivative_request = 0 ;
893 Standard_Boolean periodic_flag = Standard_False ;
894 Standard_Real result,
895 start_end[2],
896 *array_of_poles,
897 *array_of_new_poles ;
898
899 array_of_poles = (Standard_Real *) &NewPoles ;
900 extrap_mode[0] =
901 extrap_mode[1] = BSplineDegree ;
902 num_new_poles =
903 FlatKnots.Length() - NewDegree - 1 ;
904 start_end[0] = FlatKnots(NewDegree+1) ;
905 start_end[1] = FlatKnots(num_new_poles+1) ;
906 TColStd_Array1OfReal parameters(1,num_new_poles) ;
907 TColStd_Array1OfInteger contact_order_array(1,num_new_poles) ;
908 TColStd_Array1OfReal new_poles_array(1, num_new_poles * PolesDimension) ;
909
910 array_of_new_poles =
911 (Standard_Real *) &new_poles_array(1) ;
912 BuildSchoenbergPoints(NewDegree,
913 FlatKnots,
914 parameters) ;
915 //
916 // on recadre sur les bornes
917 //
918 if (parameters(1) < start_end[0]) {
919 parameters(1) = start_end[0] ;
920 }
921 if (parameters(num_new_poles) > start_end[1]) {
922 parameters(num_new_poles) = start_end[1] ;
923 }
924 index = 0 ;
925
926 for (ii = 1 ; ii <= num_new_poles ; ii++) {
927 contact_order_array(ii) = 0 ;
928 FunctionPtr.Evaluate (contact_order_array(ii),
929 start_end,
930 parameters(ii),
931 result,
932 error_code);
933 if (error_code) {
934 theStatus = 1;
935 goto FINISH ;
936 }
937
938 Eval(parameters(ii),
939 periodic_flag,
940 derivative_request,
941 extrap_mode[0],
942 BSplineDegree,
943 BSplineFlatKnots,
944 PolesDimension,
945 Poles,
946 array_of_new_poles[index]) ;
947
948 for (jj = 0 ; jj < PolesDimension ; jj++) {
949 array_of_new_poles[index] *= result ;
950 index += 1 ;
951 }
952 }
953 Interpolate(NewDegree,
954 FlatKnots,
955 parameters,
956 contact_order_array,
957 PolesDimension,
958 array_of_new_poles[0],
959 theStatus);
960
961 for (ii = 0 ; ii < num_new_poles * PolesDimension ; ii++) {
962 array_of_poles[ii] = array_of_new_poles[ii] ;
963
964 }
965 FINISH :
966 ;
967 }
968
969 //=======================================================================
970 // function : FunctionMultiply
971 //purpose :
972 //=======================================================================
973
FunctionReparameterise(const BSplCLib_EvaluatorFunction & FunctionPtr,const Standard_Integer BSplineDegree,const TColStd_Array1OfReal & BSplineFlatKnots,const Standard_Integer PolesDimension,Standard_Real & Poles,const TColStd_Array1OfReal & FlatKnots,const Standard_Integer NewDegree,Standard_Real & NewPoles,Standard_Integer & theStatus)974 void BSplCLib::FunctionReparameterise
975 (const BSplCLib_EvaluatorFunction & FunctionPtr,
976 const Standard_Integer BSplineDegree,
977 const TColStd_Array1OfReal & BSplineFlatKnots,
978 const Standard_Integer PolesDimension,
979 Standard_Real & Poles,
980 const TColStd_Array1OfReal & FlatKnots,
981 const Standard_Integer NewDegree,
982 Standard_Real & NewPoles,
983 Standard_Integer & theStatus)
984 {
985 Standard_Integer ii,
986 // jj,
987 index ;
988 Standard_Integer extrap_mode[2],
989 error_code,
990 num_new_poles,
991 derivative_request = 0 ;
992 Standard_Boolean periodic_flag = Standard_False ;
993 Standard_Real result,
994 start_end[2],
995 *array_of_poles,
996 *array_of_new_poles ;
997
998 array_of_poles = (Standard_Real *) &NewPoles ;
999 extrap_mode[0] =
1000 extrap_mode[1] = BSplineDegree ;
1001 num_new_poles =
1002 FlatKnots.Length() - NewDegree - 1 ;
1003 start_end[0] = FlatKnots(NewDegree+1) ;
1004 start_end[1] = FlatKnots(num_new_poles+1) ;
1005 TColStd_Array1OfReal parameters(1,num_new_poles) ;
1006 TColStd_Array1OfInteger contact_order_array(1,num_new_poles) ;
1007 TColStd_Array1OfReal new_poles_array(1, num_new_poles * PolesDimension) ;
1008
1009 array_of_new_poles =
1010 (Standard_Real *) &new_poles_array(1) ;
1011 BuildSchoenbergPoints(NewDegree,
1012 FlatKnots,
1013 parameters) ;
1014 index = 0 ;
1015
1016 for (ii = 1 ; ii <= num_new_poles ; ii++) {
1017 contact_order_array(ii) = 0 ;
1018 FunctionPtr.Evaluate (contact_order_array(ii),
1019 start_end,
1020 parameters(ii),
1021 result,
1022 error_code);
1023 if (error_code) {
1024 theStatus = 1;
1025 goto FINISH ;
1026 }
1027
1028 Eval(result,
1029 periodic_flag,
1030 derivative_request,
1031 extrap_mode[0],
1032 BSplineDegree,
1033 BSplineFlatKnots,
1034 PolesDimension,
1035 Poles,
1036 array_of_new_poles[index]) ;
1037 index += PolesDimension ;
1038 }
1039 Interpolate(NewDegree,
1040 FlatKnots,
1041 parameters,
1042 contact_order_array,
1043 PolesDimension,
1044 array_of_new_poles[0],
1045 theStatus);
1046
1047 for (ii = 0 ; ii < num_new_poles * PolesDimension ; ii++) {
1048 array_of_poles[ii] = array_of_new_poles[ii] ;
1049
1050 }
1051 FINISH :
1052 ;
1053 }
1054
1055 //=======================================================================
1056 //function : FunctionMultiply
1057 //purpose :
1058 //=======================================================================
1059
FunctionMultiply(const BSplCLib_EvaluatorFunction & FunctionPtr,const Standard_Integer BSplineDegree,const TColStd_Array1OfReal & BSplineFlatKnots,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal & FlatKnots,const Standard_Integer NewDegree,TColStd_Array1OfReal & NewPoles,Standard_Integer & theStatus)1060 void BSplCLib::FunctionMultiply
1061 (const BSplCLib_EvaluatorFunction & FunctionPtr,
1062 const Standard_Integer BSplineDegree,
1063 const TColStd_Array1OfReal & BSplineFlatKnots,
1064 const TColStd_Array1OfReal & Poles,
1065 const TColStd_Array1OfReal & FlatKnots,
1066 const Standard_Integer NewDegree,
1067 TColStd_Array1OfReal & NewPoles,
1068 Standard_Integer & theStatus)
1069 {
1070 Standard_Integer num_bspline_poles =
1071 BSplineFlatKnots.Length() - BSplineDegree - 1 ;
1072 Standard_Integer num_new_poles =
1073 FlatKnots.Length() - NewDegree - 1 ;
1074
1075 if (Poles.Length() != num_bspline_poles ||
1076 NewPoles.Length() != num_new_poles) {
1077 throw Standard_ConstructionError();
1078 }
1079 Standard_Real * array_of_poles =
1080 (Standard_Real *) &Poles(Poles.Lower()) ;
1081 Standard_Real * array_of_new_poles =
1082 (Standard_Real *) &NewPoles(NewPoles.Lower()) ;
1083 BSplCLib::FunctionMultiply(FunctionPtr,
1084 BSplineDegree,
1085 BSplineFlatKnots,
1086 1,
1087 array_of_poles[0],
1088 FlatKnots,
1089 NewDegree,
1090 array_of_new_poles[0],
1091 theStatus);
1092 }
1093
1094 //=======================================================================
1095 //function : FunctionReparameterise
1096 //purpose :
1097 //=======================================================================
1098
FunctionReparameterise(const BSplCLib_EvaluatorFunction & FunctionPtr,const Standard_Integer BSplineDegree,const TColStd_Array1OfReal & BSplineFlatKnots,const TColStd_Array1OfReal & Poles,const TColStd_Array1OfReal & FlatKnots,const Standard_Integer NewDegree,TColStd_Array1OfReal & NewPoles,Standard_Integer & theStatus)1099 void BSplCLib::FunctionReparameterise
1100 (const BSplCLib_EvaluatorFunction & FunctionPtr,
1101 const Standard_Integer BSplineDegree,
1102 const TColStd_Array1OfReal & BSplineFlatKnots,
1103 const TColStd_Array1OfReal & Poles,
1104 const TColStd_Array1OfReal & FlatKnots,
1105 const Standard_Integer NewDegree,
1106 TColStd_Array1OfReal & NewPoles,
1107 Standard_Integer & theStatus)
1108 {
1109 Standard_Integer num_bspline_poles =
1110 BSplineFlatKnots.Length() - BSplineDegree - 1 ;
1111 Standard_Integer num_new_poles =
1112 FlatKnots.Length() - NewDegree - 1 ;
1113
1114 if (Poles.Length() != num_bspline_poles ||
1115 NewPoles.Length() != num_new_poles) {
1116 throw Standard_ConstructionError();
1117 }
1118 Standard_Real * array_of_poles =
1119 (Standard_Real *) &Poles(Poles.Lower()) ;
1120 Standard_Real * array_of_new_poles =
1121 (Standard_Real *) &NewPoles(NewPoles.Lower()) ;
1122 BSplCLib::FunctionReparameterise(
1123 FunctionPtr,
1124 BSplineDegree,
1125 BSplineFlatKnots,
1126 1,
1127 array_of_poles[0],
1128 FlatKnots,
1129 NewDegree,
1130 array_of_new_poles[0],
1131 theStatus);
1132 }
1133
1134 //=======================================================================
1135 //function : MergeBSplineKnots
1136 //purpose :
1137 //=======================================================================
MergeBSplineKnots(const Standard_Real Tolerance,const Standard_Real StartValue,const Standard_Real EndValue,const Standard_Integer Degree1,const TColStd_Array1OfReal & Knots1,const TColStd_Array1OfInteger & Mults1,const Standard_Integer Degree2,const TColStd_Array1OfReal & Knots2,const TColStd_Array1OfInteger & Mults2,Standard_Integer & NumPoles,Handle (TColStd_HArray1OfReal)& NewKnots,Handle (TColStd_HArray1OfInteger)& NewMults)1138 void BSplCLib::MergeBSplineKnots
1139 (const Standard_Real Tolerance,
1140 const Standard_Real StartValue,
1141 const Standard_Real EndValue,
1142 const Standard_Integer Degree1,
1143 const TColStd_Array1OfReal & Knots1,
1144 const TColStd_Array1OfInteger & Mults1,
1145 const Standard_Integer Degree2,
1146 const TColStd_Array1OfReal & Knots2,
1147 const TColStd_Array1OfInteger & Mults2,
1148 Standard_Integer & NumPoles,
1149 Handle(TColStd_HArray1OfReal) & NewKnots,
1150 Handle(TColStd_HArray1OfInteger) & NewMults)
1151 {
1152 Standard_Integer ii,
1153 jj,
1154 continuity,
1155 set_mults_flag,
1156 degree,
1157 index,
1158 num_knots ;
1159 if (StartValue < EndValue - Tolerance) {
1160 TColStd_Array1OfReal knots1(1,Knots1.Length()) ;
1161 TColStd_Array1OfReal knots2(1,Knots2.Length()) ;
1162 degree = Degree1 + Degree2 ;
1163 index = 1 ;
1164
1165 for (ii = Knots1.Lower() ; ii <= Knots1.Upper() ; ii++) {
1166 knots1(index) = Knots1(ii) ;
1167 index += 1 ;
1168 }
1169 index = 1 ;
1170
1171 for (ii = Knots2.Lower() ; ii <= Knots2.Upper() ; ii++) {
1172 knots2(index) = Knots2(ii) ;
1173 index += 1 ;
1174 }
1175 BSplCLib::Reparametrize(StartValue,
1176 EndValue,
1177 knots1) ;
1178
1179 BSplCLib::Reparametrize(StartValue,
1180 EndValue,
1181 knots2) ;
1182 num_knots = 0 ;
1183 jj = 1 ;
1184
1185 for (ii = 1 ; ii <= knots1.Length() ; ii++) {
1186
1187 while (jj <= knots2.Length() && knots2(jj) <= knots1(ii) - Tolerance) {
1188 jj += 1 ;
1189 num_knots += 1 ;
1190 }
1191
1192 while (jj <= knots2.Length() && knots2(jj) <= knots1(ii) + Tolerance) {
1193 jj += 1 ;
1194 }
1195 num_knots += 1 ;
1196 }
1197 NewKnots =
1198 new TColStd_HArray1OfReal(1,num_knots) ;
1199 NewMults =
1200 new TColStd_HArray1OfInteger(1,num_knots) ;
1201 num_knots = 1 ;
1202 jj = 1 ;
1203
1204 for (ii = 1 ; ii <= knots1.Length() ; ii++) {
1205
1206 while (jj <= knots2.Length() && knots2(jj) <= knots1(ii) - Tolerance) {
1207 NewKnots->ChangeArray1()(num_knots) = knots2(jj) ;
1208 NewMults->ChangeArray1()(num_knots) = Mults2(jj) + Degree1 ;
1209 jj += 1 ;
1210 num_knots += 1 ;
1211 }
1212 set_mults_flag = 0 ;
1213
1214 while (jj <= knots2.Length() && knots2(jj) <= knots1(ii) + Tolerance) {
1215 continuity = Min(Degree1 - Mults1(ii), Degree2 - Mults2(jj)) ;
1216 set_mults_flag = 1 ;
1217 NewMults->ChangeArray1()(num_knots) = degree - continuity ;
1218 jj += 1 ;
1219 }
1220
1221 NewKnots->ChangeArray1()(num_knots) = knots1(ii) ;
1222 if (! set_mults_flag) {
1223 NewMults->ChangeArray1()(num_knots) = Mults1(ii) + Degree2 ;
1224 }
1225 num_knots += 1 ;
1226 }
1227 num_knots -= 1 ;
1228 NewMults->ChangeArray1()(1) = degree + 1 ;
1229 NewMults->ChangeArray1()(num_knots) = degree + 1 ;
1230 index = 0 ;
1231
1232 for (ii = 1 ; ii <= num_knots ; ii++) {
1233 index += NewMults->Value(ii) ;
1234 }
1235 NumPoles = index - degree - 1 ;
1236 }
1237 }
1238
1239