0031668: Visualization - WebGL sample doesn't work on Emscripten 1.39
[occt.git] / src / BSplSLib / BSplSLib.cxx
CommitLineData
b311480e 1// Created on: 1991-08-26
2// Created by: JCV
3// Copyright (c) 1991-1999 Matra Datavision
973c2be1 4// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 5//
973c2be1 6// This file is part of Open CASCADE Technology software library.
b311480e 7//
d5f74e42 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
973c2be1 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.
b311480e 13//
973c2be1 14// Alternatively, this file may be used under the terms of Open CASCADE
15// commercial license or contractual agreement.
7fd59977 16
17// Modifed RLE Aug 93 - Complete rewrite
18// xab 21-Mar-95 implemented cache mecanism
19// pmn 25-09-96 Interpolation
20// jct 25-09-96 : Correction de l'alloc de LocalArray dans RationalDerivative.
21// pmn 07-10-96 : Correction de DN dans le cas rationnal.
22// pmn 06-02-97 : Correction des poids dans RationalDerivative. (PRO700)
23
7fd59977 24#include <BSplCLib.hxx>
42cf5bc1 25#include <BSplSLib.hxx>
26#include <gp_Pnt.hxx>
27#include <gp_Vec.hxx>
28#include <math_Matrix.hxx>
29#include <NCollection_LocalArray.hxx>
30#include <PLib.hxx>
31#include <Standard_ConstructionError.hxx>
32#include <Standard_NotImplemented.hxx>
7fd59977 33#include <TColgp_Array1OfXYZ.hxx>
42cf5bc1 34#include <TColgp_Array2OfXYZ.hxx>
7fd59977 35#include <TColStd_HArray1OfInteger.hxx>
7fd59977 36
37// for null derivatives
41194117 38static Standard_Real BSplSLib_zero[3] = {0.0, 0.0, 0.0};
7fd59977 39
40//=======================================================================
41//struct : BSplCLib_DataContainer
42//purpose: Auxiliary structure providing buffers for poles and knots used in
43// evaluation of bspline (allocated in the stack)
44//=======================================================================
45
41194117 46struct BSplSLib_DataContainer
7fd59977 47{
41194117 48 BSplSLib_DataContainer (Standard_Integer UDegree, Standard_Integer VDegree)
7fd59977 49 {
105aae76 50 (void)UDegree; (void)VDegree; // just to avoid compiler warning in Release mode
41194117
K
51 Standard_OutOfRange_Raise_if (UDegree > BSplCLib::MaxDegree() ||
52 VDegree > BSplCLib::MaxDegree() || BSplCLib::MaxDegree() > 25,
53 "BSplSLib: bspline degree is greater than maximum supported");
7fd59977 54 }
7fd59977 55 Standard_Real poles[4*(25+1)*(25+1)];
56 Standard_Real knots1[2*25];
57 Standard_Real knots2[2*25];
58 Standard_Real ders[48];
59};
60
7fd59977 61//**************************************************************************
62// Evaluation methods
63//**************************************************************************
64
65//=======================================================================
66//function : RationalDerivative
67//purpose : computes the rational derivatives when whe have the
68// the derivatives of the homogeneous numerator and the
69// the derivatives of the denominator
70//=======================================================================
71
72void BSplSLib::RationalDerivative(const Standard_Integer UDeg,
73 const Standard_Integer VDeg,
74 const Standard_Integer N,
75 const Standard_Integer M,
76 Standard_Real& HDerivatives,
77 Standard_Real& RDerivatives,
78 const Standard_Boolean All)
79{
80 //
81 // if All is True all derivatives are computed. if Not only
82 // the requested N, M is computed
83 //
84 // Numerator(u,v)
85 // let f(u,v) be a rational function = ------------------
86 // Denominator(u,v)
87 //
88 //
89 // Let (N,M) the order of the derivatives we want : then since
90 // we have :
91 //
92 // Numerator = f * Denominator
93 //
94 // we derive :
95 //
96 // (N,M) 1 ( (N M) (p q) (N -p M-q) )
97 // f = ------------ ( Numerator - SUM SUM a * f * Denominator )
98 // (0,0) ( p<N q<M p q )
99 // Denominator
100 //
101 // with :
102 //
103 // ( N ) ( M )
104 // a = ( ) ( )
105 // p q ( p ) ( q )
106 //
107 //
108 // HDerivatives is an array where derivatives are stored in the following form
109 // Numerator is assumee to have 3 functions that is a vector of dimension
110 // 3
111 //
112 // (0,0) (0,0) (0, DegV) (0, DegV)
113 // Numerator Denominator ... Numerator Denominator
114 //
115 // (1,0) (1,0) (1, DegV) (1, DegV)
116 // Numerator Denominator ... Numerator Denominator
117 //
118 // ...........................................................
119 //
120 //
121 // (DegU,0) (DegU,0) (DegU, DegV) (DegU, DegV)
122 // Numerator Denominator ... Numerator Denominator
123 //
124 //
125 Standard_Integer ii,jj,pp,qq,index,index1,index2;
126 Standard_Integer M1,M3,M4,N1,iiM1,iiM3,jjM1,ppM1,ppM3;
127 Standard_Integer MinN,MinN1,MinM,MinM1;
128 Standard_Integer index_u,index_u1,index_v,index_v1,index_w;
129
130 M1 = M + 1;
131 N1 = N + 1;
132 ii = N1 * M1;
133 M3 = (M1 << 1) + M1;
134 M4 = (VDeg + 1) << 2;
135
f7b4312f 136 NCollection_LocalArray<Standard_Real> StoreDerivatives (All ? 0 : ii * 3);
7fd59977 137 Standard_Real *RArray = (All ? &RDerivatives : (Standard_Real*)StoreDerivatives);
f7b4312f 138 NCollection_LocalArray<Standard_Real> StoreW (ii);
7fd59977 139 Standard_Real *HomogeneousArray = &HDerivatives;
140 Standard_Real denominator,Pii,Pip,Pjq;
141
142 denominator = 1.0e0 / HomogeneousArray[3];
143 index_u = 0;
144 index_u1 = 0;
145 if (UDeg < N) MinN = UDeg;
146 else MinN = N;
147 if (VDeg < M) MinM = VDeg;
148 else MinM = M;
149 MinN1 = MinN + 1;
150 MinM1 = MinM + 1;
151 iiM1 = - M1;
152
153 for (ii = 0 ; ii < MinN1 ; ii++) {
154 iiM1 += M1;
155 index_v = index_u;
156 index_v1 = index_u1;
157 index_w = iiM1;
158
159 for (jj = 0 ; jj < MinM1 ; jj++) {
174178b9 160 RArray[index_v++] = HomogeneousArray[index_v1++];
161 RArray[index_v++] = HomogeneousArray[index_v1++];
162 RArray[index_v++] = HomogeneousArray[index_v1++];
163 StoreW[index_w++] = HomogeneousArray[index_v1++];
7fd59977 164 }
165
166 for (jj = MinM1 ; jj < M1 ; jj++) {
174178b9 167 RArray[index_v++] = 0.;
168 RArray[index_v++] = 0.;
169 RArray[index_v++] = 0.;
170 StoreW[index_w++] = 0.;
7fd59977 171 }
172 index_u1 += M4;
173 index_u += M3;
174 }
175 index_v = MinN1 * M3;
176 index_w = MinN1 * M1;
177
178 for (ii = MinN1 ; ii < N1 ; ii++) {
179
180 for (jj = 0 ; jj < M1 ; jj++) {
174178b9 181 RArray[index_v++] = 0.0e0;
182 RArray[index_v++] = 0.0e0;
183 RArray[index_v++] = 0.0e0;
184 StoreW[index_w++] = 0.0e0;
7fd59977 185 }
186 }
187
0d969553 188 // --------------- Calculation ----------------
7fd59977 189
190 iiM1 = - M1;
191 iiM3 = - M3;
7fd59977 192
193 for (ii = 0 ; ii <= N ; ii++) {
194 iiM1 += M1;
195 iiM3 += M3;
196 index1 = iiM3 - 3;
197 jjM1 = iiM1;
198
199 for (jj = 0 ; jj <= M ; jj++) {
200 jjM1 ++;
201 ppM1 = - M1;
202 ppM3 = - M3;
203 index1 += 3;
204
205 for (pp = 0 ; pp < ii ; pp++) {
206 ppM1 += M1;
207 ppM3 += M3;
208 index = ppM3;
209 index2 = jjM1 - ppM1;
210 Pip = PLib::Bin(ii,pp);
211
212 for (qq = 0 ; qq <= jj ; qq++) {
213 index2--;
214 Pjq = Pip * PLib::Bin(jj,qq) * StoreW[index2];
215 RArray[index1] -= Pjq * RArray[index]; index++; index1++;
216 RArray[index1] -= Pjq * RArray[index]; index++; index1++;
217 RArray[index1] -= Pjq * RArray[index]; index++;
218 index1 -= 2;
219 }
220 }
221 index = iiM3;
222 index2 = jj + 1;
223 Pii = PLib::Bin(ii,ii);
224
225 for (qq = 0 ; qq < jj ; qq++) {
226 index2--;
227 Pjq = Pii * PLib::Bin(jj,qq) * StoreW[index2];
228 RArray[index1] -= Pjq * RArray[index]; index++; index1++;
229 RArray[index1] -= Pjq * RArray[index]; index++; index1++;
230 RArray[index1] -= Pjq * RArray[index]; index++;
231 index1 -= 2;
232 }
233 RArray[index1] *= denominator; index1++;
234 RArray[index1] *= denominator; index1++;
235 RArray[index1] *= denominator;
236 index1 -= 2;
237 }
238 }
239 if (!All) {
240 RArray = &RDerivatives;
241 index = N * M1 + M;
242 index = (index << 1) + index;
243 RArray[0] = StoreDerivatives[index]; index++;
244 RArray[1] = StoreDerivatives[index]; index++;
245 RArray[2] = StoreDerivatives[index];
246 }
247}
248
249//=======================================================================
250//function : PrepareEval
251//purpose :
252//=======================================================================
253
254//
255// PrepareEval :
256//
0d969553 257// Prepare all data for computing points :
7fd59977 258// local arrays of knots
259// local array of poles (multiplied by the weights if rational)
260//
261// The first direction to compute (smaller degree) is returned
262// and the poles are stored according to this direction.
263
ab87e6fc 264static Standard_Boolean PrepareEval (const Standard_Real U,
265 const Standard_Real V,
266 const Standard_Integer Uindex,
267 const Standard_Integer Vindex,
268 const Standard_Integer UDegree,
269 const Standard_Integer VDegree,
270 const Standard_Boolean URat,
271 const Standard_Boolean VRat,
272 const Standard_Boolean UPer,
273 const Standard_Boolean VPer,
274 const TColgp_Array2OfPnt& Poles,
0e14656b 275 const TColStd_Array2OfReal* Weights,
ab87e6fc 276 const TColStd_Array1OfReal& UKnots,
277 const TColStd_Array1OfReal& VKnots,
0e14656b 278 const TColStd_Array1OfInteger* UMults,
279 const TColStd_Array1OfInteger* VMults,
ab87e6fc 280 Standard_Real& u1, // first parameter to use
281 Standard_Real& u2, // second parameter to use
282 Standard_Integer& d1, // first degree
283 Standard_Integer& d2, // second degree
284 Standard_Boolean& rational,
285 BSplSLib_DataContainer& dc)
286 {
7fd59977 287 rational = URat || VRat;
288 Standard_Integer uindex = Uindex;
289 Standard_Integer vindex = Vindex;
290 Standard_Integer UKLower = UKnots.Lower();
291 Standard_Integer UKUpper = UKnots.Upper();
292 Standard_Integer VKLower = VKnots.Lower();
293 Standard_Integer VKUpper = VKnots.Upper();
ab87e6fc 294
295 if (UDegree <= VDegree)
296 {
7fd59977 297 // compute the indices
298 if (uindex < UKLower || uindex > UKUpper)
299 BSplCLib::LocateParameter(UDegree,UKnots,UMults,U,UPer,uindex,u1);
ab87e6fc 300 else
301 u1 = U;
302
7fd59977 303 if (vindex < VKLower || vindex > VKUpper)
304 BSplCLib::LocateParameter(VDegree,VKnots,VMults,V,VPer,vindex,u2);
ab87e6fc 305 else
306 u2 = V;
307
7fd59977 308 // get the knots
309 d1 = UDegree;
310 d2 = VDegree;
311 BSplCLib::BuildKnots(UDegree,uindex,UPer,UKnots,UMults,*dc.knots1);
312 BSplCLib::BuildKnots(VDegree,vindex,VPer,VKnots,VMults,*dc.knots2);
ab87e6fc 313
0e14656b 314 if (UMults == NULL)
ab87e6fc 315 uindex -= UKLower + UDegree;
316 else
0e14656b 317 uindex = BSplCLib::PoleIndex(UDegree,uindex,UPer,*UMults);
ab87e6fc 318
0e14656b 319 if (VMults == NULL)
ab87e6fc 320 vindex -= VKLower + VDegree;
321 else
0e14656b 322 vindex = BSplCLib::PoleIndex(VDegree,vindex,VPer,*VMults);
ab87e6fc 323
7fd59977 324 // get the poles
7fd59977 325 Standard_Integer i,j,ip,jp;
326 Standard_Real w, *pole = dc.poles;
327 d1 = UDegree;
328 d2 = VDegree;
329 Standard_Integer PLowerRow = Poles.LowerRow();
330 Standard_Integer PUpperRow = Poles.UpperRow();
331 Standard_Integer PLowerCol = Poles.LowerCol();
332 Standard_Integer PUpperCol = Poles.UpperCol();
ab87e6fc 333
334 // verify if locally non rational
335 if (rational)
336 {
7fd59977 337 rational = Standard_False;
338 ip = PLowerRow + uindex;
339 jp = PLowerCol + vindex;
ab87e6fc 340
341 if(ip < PLowerRow) ip = PUpperRow;
342 if(jp < PLowerCol) jp = PUpperCol;
343
0e14656b 344 w = Weights->Value(ip,jp);
7fd59977 345 Standard_Real eps = Epsilon(w);
346 Standard_Real dw;
ab87e6fc 347
348 for (i = 0; i <= UDegree && !rational; i++)
349 {
350 jp = PLowerCol + vindex;
351
352 if(jp < PLowerCol)
353 jp = PUpperCol;
354
355 for (j = 0; j <= VDegree && !rational; j++)
356 {
0e14656b 357 dw = Weights->Value(ip,jp) - w;
ab87e6fc 358 if (dw < 0)
359 dw = - dw;
360
361 rational = (dw > eps);
362
363 jp++;
364
365 if (jp > PUpperCol)
366 jp = PLowerCol;
367 }
368
369 ip++;
370
371 if (ip > PUpperRow)
372 ip = PLowerRow;
373
374 }
7fd59977 375 }
ab87e6fc 376
7fd59977 377 // copy the poles
378 ip = PLowerRow + uindex;
ab87e6fc 379
380 if(ip < PLowerRow)
381 ip = PUpperRow;
382
383 if (rational)
384 {
385 for (i = 0; i <= d1; i++)
386 {
387 jp = PLowerCol + vindex;
388
389 if(jp < PLowerCol)
390 jp = PUpperCol;
391
392 for (j = 0; j <= d2; j++)
393 {
394 const gp_Pnt& P = Poles .Value(ip,jp);
0e14656b 395 pole[3] = w = Weights->Value(ip,jp);
ab87e6fc 396 pole[0] = P.X() * w;
397 pole[1] = P.Y() * w;
398 pole[2] = P.Z() * w;
399 pole += 4;
400 jp++;
401
402 if (jp > PUpperCol)
403 jp = PLowerCol;
404 }
405
406 ip++;
407
408 if (ip > PUpperRow)
409 ip = PLowerRow;
410
411 }
7fd59977 412 }
ab87e6fc 413 else
414 {
415 for (i = 0; i <= d1; i++)
416 {
417 jp = PLowerCol + vindex;
418
419 if(jp < PLowerCol)
420 jp = PUpperCol;
421
422 for (j = 0; j <= d2; j++)
423 {
424 const gp_Pnt& P = Poles.Value(ip,jp);
425 pole[0] = P.X();
426 pole[1] = P.Y();
427 pole[2] = P.Z();
428 pole += 3;
429 jp++;
430
431 if (jp > PUpperCol)
432 jp = PLowerCol;
433 }
434
435 ip++;
436
437 if (ip > PUpperRow)
438 ip = PLowerRow;
439 }
7fd59977 440 }
ab87e6fc 441
7fd59977 442 return Standard_True;
ab87e6fc 443 }
444 else
445 {
7fd59977 446 // compute the indices
447 if (uindex < UKLower || uindex > UKUpper)
448 BSplCLib::LocateParameter(UDegree,UKnots,UMults,U,UPer,uindex,u2);
ab87e6fc 449 else
450 u2 = U;
451
7fd59977 452 if (vindex < VKLower || vindex > VKUpper)
453 BSplCLib::LocateParameter(VDegree,VKnots,VMults,V,VPer,vindex,u1);
ab87e6fc 454 else
455 u1 = V;
456
7fd59977 457 // get the knots
ab87e6fc 458
7fd59977 459 d2 = UDegree;
460 d1 = VDegree;
ab87e6fc 461
7fd59977 462 BSplCLib::BuildKnots(UDegree,uindex,UPer,UKnots,UMults,*dc.knots2);
463 BSplCLib::BuildKnots(VDegree,vindex,VPer,VKnots,VMults,*dc.knots1);
ab87e6fc 464
0e14656b 465 if (UMults == NULL)
ab87e6fc 466 uindex -= UKLower + UDegree;
467 else
0e14656b 468 uindex = BSplCLib::PoleIndex(UDegree,uindex,UPer,*UMults);
ab87e6fc 469
0e14656b 470 if (VMults == NULL)
ab87e6fc 471 vindex -= VKLower + VDegree;
472 else
0e14656b 473 vindex = BSplCLib::PoleIndex(VDegree,vindex,VPer,*VMults);
ab87e6fc 474
7fd59977 475 // get the poles
7fd59977 476 Standard_Integer i,j,ip,jp;
477 Standard_Real w, *pole = dc.poles;
478 d1 = VDegree;
479 d2 = UDegree;
480 Standard_Integer PLowerRow = Poles.LowerRow();
481 Standard_Integer PUpperRow = Poles.UpperRow();
482 Standard_Integer PLowerCol = Poles.LowerCol();
483 Standard_Integer PUpperCol = Poles.UpperCol();
ab87e6fc 484
485 // verify if locally non rational
486 if (rational)
487 {
7fd59977 488 rational = Standard_False;
489 ip = PLowerRow + uindex;
490 jp = PLowerCol + vindex;
ab87e6fc 491
492 if(ip < PLowerRow)
493 ip = PUpperRow;
494
495 if(jp < PLowerCol)
496 jp = PUpperCol;
497
0e14656b 498 w = Weights->Value(ip,jp);
7fd59977 499 Standard_Real eps = Epsilon(w);
500 Standard_Real dw;
ab87e6fc 501
502 for (i = 0; i <= UDegree && !rational; i++)
503 {
504 jp = PLowerCol + vindex;
505
506 if(jp < PLowerCol)
507 jp = PUpperCol;
508
509 for (j = 0; j <= VDegree && !rational; j++)
510 {
0e14656b 511 dw = Weights->Value(ip,jp) - w;
ab87e6fc 512 if (dw < 0) dw = - dw;
513 rational = dw > eps;
514
515 jp++;
516
517 if (jp > PUpperCol)
518 jp = PLowerCol;
519 }
520
521 ip++;
522
523 if (ip > PUpperRow)
524 ip = PLowerRow;
525
526 }
7fd59977 527 }
ab87e6fc 528
7fd59977 529 // copy the poles
530 jp = PLowerCol + vindex;
ab87e6fc 531
532 if(jp < PLowerCol)
533 jp = PUpperCol;
534
535 if (rational)
536 {
537 for (i = 0; i <= d1; i++)
538 {
539 ip = PLowerRow + uindex;
540
541 if(ip < PLowerRow)
542 ip = PUpperRow;
543
544 for (j = 0; j <= d2; j++)
545 {
546 const gp_Pnt& P = Poles.Value(ip,jp);
0e14656b 547 pole[3] = w = Weights->Value(ip,jp);
ab87e6fc 548 pole[0] = P.X() * w;
549 pole[1] = P.Y() * w;
550 pole[2] = P.Z() * w;
551 pole += 4;
552 ip++;
553
554 if (ip > PUpperRow)
555 ip = PLowerRow;
556
557 }
558
559 jp++;
560
561 if (jp > PUpperCol)
562 jp = PLowerCol;
563
564 }
7fd59977 565 }
ab87e6fc 566 else
567 {
568 for (i = 0; i <= d1; i++)
569 {
570 ip = PLowerRow + uindex;
571
572 if(ip < PLowerRow)
573 ip = PUpperRow;
574
b26415fb 575 if(ip > PUpperRow)
576 ip = PLowerRow;
577
ab87e6fc 578 for (j = 0; j <= d2; j++)
579 {
580 const gp_Pnt& P = Poles.Value(ip,jp);
581 pole[0] = P.X();
582 pole[1] = P.Y();
583 pole[2] = P.Z();
584 pole += 3;
585 ip++;
586
587 if (ip > PUpperRow)
588 ip = PLowerRow;
589 }
590
591 jp++;
592
593 if (jp > PUpperCol)
594 jp = PLowerCol;
595
596 }
7fd59977 597 }
ab87e6fc 598
7fd59977 599 return Standard_False;
ab87e6fc 600 }
7fd59977 601 }
7fd59977 602
603//=======================================================================
604//function : D0
605//purpose :
606//=======================================================================
607
608void BSplSLib::D0
609(const Standard_Real U,
610 const Standard_Real V,
611 const Standard_Integer UIndex,
612 const Standard_Integer VIndex,
613 const TColgp_Array2OfPnt& Poles,
0e14656b 614 const TColStd_Array2OfReal* Weights,
7fd59977 615 const TColStd_Array1OfReal& UKnots,
616 const TColStd_Array1OfReal& VKnots,
0e14656b 617 const TColStd_Array1OfInteger* UMults,
618 const TColStd_Array1OfInteger* VMults,
7fd59977 619 const Standard_Integer UDegree,
620 const Standard_Integer VDegree,
621 const Standard_Boolean URat,
622 const Standard_Boolean VRat,
623 const Standard_Boolean UPer,
624 const Standard_Boolean VPer,
625 gp_Pnt& P)
626{
627// Standard_Integer k ;
628 Standard_Real W ;
629 HomogeneousD0(U,
630 V,
631 UIndex,
632 VIndex,
633 Poles,
634 Weights,
635 UKnots,
636 VKnots,
637 UMults,
638 VMults,
639 UDegree,
640 VDegree,
641 URat,
642 VRat,
643 UPer,
644 VPer,
645 W,
646 P) ;
647 P.SetX(P.X() / W);
648 P.SetY(P.Y() / W);
649 P.SetZ(P.Z() / W);
650}
651
652//=======================================================================
653//function : D0
654//purpose :
655//=======================================================================
656
657void BSplSLib::HomogeneousD0
658(const Standard_Real U,
659 const Standard_Real V,
660 const Standard_Integer UIndex,
661 const Standard_Integer VIndex,
662 const TColgp_Array2OfPnt& Poles,
0e14656b 663 const TColStd_Array2OfReal* Weights,
7fd59977 664 const TColStd_Array1OfReal& UKnots,
665 const TColStd_Array1OfReal& VKnots,
0e14656b 666 const TColStd_Array1OfInteger* UMults,
667 const TColStd_Array1OfInteger* VMults,
7fd59977 668 const Standard_Integer UDegree,
669 const Standard_Integer VDegree,
670 const Standard_Boolean URat,
671 const Standard_Boolean VRat,
672 const Standard_Boolean UPer,
673 const Standard_Boolean VPer,
674 Standard_Real & W,
675 gp_Pnt& P)
676{
677 Standard_Boolean rational;
678// Standard_Integer k,dim;
679 Standard_Integer dim;
680 Standard_Real u1,u2;
681 Standard_Integer d1,d2;
682 W = 1.0e0 ;
683
684 BSplSLib_DataContainer dc (UDegree, VDegree);
685 PrepareEval(U,V,UIndex,VIndex,UDegree,VDegree,URat,VRat,UPer,VPer,
686 Poles,Weights,UKnots,VKnots,UMults,VMults,
687 u1,u2,d1,d2,rational,dc);
688 if (rational) {
689 dim = 4;
690 BSplCLib::Eval(u1,d1,*dc.knots1,dim * (d2 + 1),*dc.poles);
691 BSplCLib::Eval(u2,d2,*dc.knots2,dim,*dc.poles);
692 W = dc.poles[3];
693 P.SetX(dc.poles[0]);
694 P.SetY(dc.poles[1]);
695 P.SetZ(dc.poles[2]);
696 }
697 else {
698 dim = 3;
699 BSplCLib::Eval(u1,d1,*dc.knots1,dim * (d2 + 1),*dc.poles);
700 BSplCLib::Eval(u2,d2,*dc.knots2,dim,*dc.poles);
701 P.SetX(dc.poles[0]);
702 P.SetY(dc.poles[1]);
703 P.SetZ(dc.poles[2]);
704 }
705}
706
707//=======================================================================
708//function : D1
709//purpose :
710//=======================================================================
711
712void BSplSLib::D1
713(const Standard_Real U,
714 const Standard_Real V,
715 const Standard_Integer UIndex,
716 const Standard_Integer VIndex,
717 const TColgp_Array2OfPnt& Poles,
0e14656b 718 const TColStd_Array2OfReal* Weights,
7fd59977 719 const TColStd_Array1OfReal& UKnots,
720 const TColStd_Array1OfReal& VKnots,
0e14656b 721 const TColStd_Array1OfInteger* UMults,
722 const TColStd_Array1OfInteger* VMults,
7fd59977 723 const Standard_Integer UDegree,
724 const Standard_Integer VDegree,
725 const Standard_Boolean URat,
726 const Standard_Boolean VRat,
727 const Standard_Boolean UPer,
728 const Standard_Boolean VPer,
729 gp_Pnt& P,
730 gp_Vec& Vu,
731 gp_Vec& Vv)
732{
733 Standard_Boolean rational;
734// Standard_Integer k,dim,dim2;
735 Standard_Integer dim,dim2;
736 Standard_Real u1,u2;
737 Standard_Integer d1,d2;
738 Standard_Real *result, *resVu, *resVv;
739 BSplSLib_DataContainer dc (UDegree, VDegree);
740 if (PrepareEval
741 (U,V,UIndex,VIndex,UDegree,VDegree,URat,VRat,UPer,VPer,
742 Poles,Weights,UKnots,VKnots,UMults,VMults,
743 u1,u2,d1,d2,rational,dc)) {
744 if (rational) {
745 dim = 4;
746 dim2 = (d2 + 1) << 2;
747 BSplCLib::Bohm(u1,d1,1,*dc.knots1,dim2,*dc.poles);
748 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*dc.poles);
749 BSplCLib::Eval(u2,d2, *dc.knots2,dim ,*(dc.poles + dim2));
750 BSplSLib::RationalDerivative(d1,d2,1,1,*dc.poles,*dc.ders);
751 result = dc.ders;
752 resVu = result + 6;
753 resVv = result + 3;
754 }
755 else {
756 dim = 3;
757 dim2 = d2 + 1;
758 dim2 = (dim2 << 1) + dim2;
759 BSplCLib::Bohm(u1,d1,1,*dc.knots1,dim2,*dc.poles);
760 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*dc.poles);
761 BSplCLib::Eval(u2,d2, *dc.knots2,dim ,*(dc.poles + dim2));
762 result = dc.poles;
763 resVu = result + dim2;
764 resVv = result + 3;
765 }
766 }
767 else {
768 if (rational) {
769 dim = 4;
770 dim2 = (d2 + 1) << 2;
771 BSplCLib::Bohm(u1,d1,1,*dc.knots1,dim2,*dc.poles);
772 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*dc.poles);
773 BSplCLib::Eval(u2,d2, *dc.knots2,dim ,*(dc.poles + dim2));
774 BSplSLib::RationalDerivative(d1,d2,1,1,*dc.poles,*dc.ders);
775 result = dc.ders;
776 resVu = result + 3;
777 resVv = result + 6;
778 }
779 else {
780 dim = 3;
781 dim2 = d2 + 1;
782 dim2 = (dim2 << 1) + dim2;
783 BSplCLib::Bohm(u1,d1,1,*dc.knots1,dim2,*dc.poles);
784 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*dc.poles);
785 BSplCLib::Eval(u2,d2 ,*dc.knots2,dim ,*(dc.poles + dim2));
786 result = dc.poles;
787 resVu = result + 3;
788 resVv = result + dim2;
789 }
790 }
791
792 P .SetX(result[0]);
793 Vu.SetX(resVu [0]);
794 Vv.SetX(resVv [0]);
795
796 P .SetY(result[1]);
797 Vu.SetY(resVu [1]);
798 Vv.SetY(resVv [1]);
799
800 P .SetZ(result[2]);
801 Vu.SetZ(resVu [2]);
802 Vv.SetZ(resVv [2]);
803}
804
805//=======================================================================
806//function : D1
807//purpose :
808//=======================================================================
809
810void BSplSLib::HomogeneousD1
811(const Standard_Real U,
812 const Standard_Real V,
813 const Standard_Integer UIndex,
814 const Standard_Integer VIndex,
815 const TColgp_Array2OfPnt& Poles,
0e14656b 816 const TColStd_Array2OfReal* Weights,
7fd59977 817 const TColStd_Array1OfReal& UKnots,
818 const TColStd_Array1OfReal& VKnots,
0e14656b 819 const TColStd_Array1OfInteger* UMults,
820 const TColStd_Array1OfInteger* VMults,
7fd59977 821 const Standard_Integer UDegree,
822 const Standard_Integer VDegree,
823 const Standard_Boolean URat,
824 const Standard_Boolean VRat,
825 const Standard_Boolean UPer,
826 const Standard_Boolean VPer,
827 gp_Pnt& N,
828 gp_Vec& Nu,
829 gp_Vec& Nv,
830 Standard_Real& D,
831 Standard_Real& Du,
832 Standard_Real& Dv)
833{
834 Standard_Boolean rational;
835// Standard_Integer k,dim;
836 Standard_Integer dim;
837 Standard_Real u1,u2;
838 Standard_Integer d1,d2;
839
840 D = 1.0e0 ;
841 Du = 0.0e0 ;
842 Dv = 0.0e0 ;
843 BSplSLib_DataContainer dc (UDegree, VDegree);
844 Standard_Boolean ufirst = PrepareEval
845 (U,V,UIndex,VIndex,UDegree,VDegree,URat,VRat,UPer,VPer,
846 Poles,Weights,UKnots,VKnots,UMults,VMults,
847 u1,u2,d1,d2,rational,dc);
848 dim = rational ? 4 : 3;
849
850 BSplCLib::Bohm(u1,d1,1,*dc.knots1,dim * (d2 + 1),*dc.poles);
851 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim,*dc.poles);
852 BSplCLib::Eval(u2,d2,*dc.knots2,dim,*(dc.poles+dim*(d2+1)));
853
854 Standard_Real *result, *resVu, *resVv;
855 result = dc.poles;
856 resVu = result + (ufirst ? dim*(d2+1) : dim);
857 resVv = result + (ufirst ? dim : dim*(d2+1));
858
859 N .SetX(result[0]);
860 Nu.SetX(resVu [0]);
861 Nv.SetX(resVv [0]);
862
863 N .SetY(result[1]);
864 Nu.SetY(resVu [1]);
865 Nv.SetY(resVv [1]);
866
867 N .SetZ(result[2]);
868 Nu.SetZ(resVu [2]);
869 Nv.SetZ(resVv [2]);
870
871 if (rational) {
872 D = result[3];
873 Du = resVu [3];
874 Dv = resVv [3];
875 }
876}
877
878//=======================================================================
879//function : D2
880//purpose :
881//=======================================================================
882
883void BSplSLib::D2
884(const Standard_Real U,
885 const Standard_Real V,
886 const Standard_Integer UIndex,
887 const Standard_Integer VIndex,
888 const TColgp_Array2OfPnt& Poles,
0e14656b 889 const TColStd_Array2OfReal* Weights,
7fd59977 890 const TColStd_Array1OfReal& UKnots,
891 const TColStd_Array1OfReal& VKnots,
0e14656b 892 const TColStd_Array1OfInteger* UMults,
893 const TColStd_Array1OfInteger* VMults,
7fd59977 894 const Standard_Integer UDegree,
895 const Standard_Integer VDegree,
896 const Standard_Boolean URat,
897 const Standard_Boolean VRat,
898 const Standard_Boolean UPer,
899 const Standard_Boolean VPer,
900 gp_Pnt& P,
901 gp_Vec& Vu,
902 gp_Vec& Vv,
903 gp_Vec& Vuu,
904 gp_Vec& Vvv,
905 gp_Vec& Vuv)
906{
907 Standard_Boolean rational;
908// Standard_Integer k,dim,dim2;
909 Standard_Integer dim,dim2;
910 Standard_Real u1,u2;
911 Standard_Integer d1,d2;
912 Standard_Real *result, *resVu, *resVv, *resVuu, *resVvv, *resVuv;
913 BSplSLib_DataContainer dc (UDegree, VDegree);
914 if (PrepareEval
915 (U,V,UIndex,VIndex,UDegree,VDegree,URat,VRat,UPer,VPer,
916 Poles,Weights,UKnots,VKnots,UMults,VMults,
917 u1,u2,d1,d2,rational,dc)) {
918 if (rational) {
919 dim = 4;
920 dim2 = (d2 + 1) << 2;
921 BSplCLib::Bohm(u1,d1,2,*dc.knots1,dim2,*dc.poles);
922 BSplCLib::Bohm(u2,d2,2,*dc.knots2,dim ,*dc.poles);
923 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*(dc.poles + dim2));
924 if (d1 > 1)
925 BSplCLib::Eval(u2,d2,*dc.knots2,dim ,*(dc.poles + (dim2 << 1)));
926 BSplSLib::RationalDerivative(d1,d2,2,2,*dc.poles,*dc.ders);
927 result = dc.ders;
928 resVu = result + 9;
929 resVv = result + 3;
930 resVuu = result + 18;
931 resVvv = result + 6;
932 resVuv = result + 12;
933 }
934 else {
935 dim = 3;
936 dim2 = d2 + 1;
937 dim2 = (dim2 << 1) + dim2;
938 BSplCLib::Bohm(u1,d1,2,*dc.knots1,dim2,*dc.poles);
939 BSplCLib::Bohm(u2,d2,2,*dc.knots2,dim ,*dc.poles);
940 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*(dc.poles + dim2));
941 if (d1 > 1)
942 BSplCLib::Eval(u2,d2,*dc.knots2,dim ,*(dc.poles + (dim2 << 1)));
943 result = dc.poles;
944 resVu = result + dim2;
945 resVv = result + 3;
946 if (UDegree <= 1) resVuu = BSplSLib_zero;
947 else resVuu = result + (dim2 << 1);
948 if (VDegree <= 1) resVvv = BSplSLib_zero;
949 else resVvv = result + 6;
950 resVuv = result + (d2 << 1) + d2 + 6;
951 }
952 }
953 else {
954 if (rational) {
955 dim = 4;
956 dim2 = (d2 + 1) << 2;
957 BSplCLib::Bohm(u1,d1,2,*dc.knots1,dim2,*dc.poles);
958 BSplCLib::Bohm(u2,d2,2,*dc.knots2,dim ,*dc.poles);
959 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*(dc.poles + dim2));
960 if (d1 > 1)
961 BSplCLib::Eval(u2,d2,*dc.knots2,dim ,*(dc.poles + (dim2 << 1)));
962 BSplSLib::RationalDerivative(d1,d2,2,2,*dc.poles,*dc.ders);
963 result = dc.ders;
964 resVu = result + 3;
965 resVv = result + 9;
966 resVuu = result + 6;
967 resVvv = result + 18;
968 resVuv = result + 12;
969 }
970 else {
971 dim = 3;
972 dim2 = d2 + 1;
973 dim2 = (dim2 << 1) + dim2;
974 BSplCLib::Bohm(u1,d1,2,*dc.knots1,dim2,*dc.poles);
975 BSplCLib::Bohm(u2,d2,2,*dc.knots2,dim ,*dc.poles);
976 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*(dc.poles + dim2));
977 if (d1 > 1)
978 BSplCLib::Eval(u2,d2,*dc.knots2,dim ,*(dc.poles + (dim2 << 1)));
979 result = dc.poles;
980 resVu = result + 3;
981 resVv = result + dim2;
982 if (UDegree <= 1) resVuu = BSplSLib_zero;
983 else resVuu = result + 6;
984 if (VDegree <= 1) resVvv = BSplSLib_zero;
985 else resVvv = result + (dim2 << 1);
986 resVuv = result + (d2 << 1) + d2 + 6;
987 }
988 }
989
990 P .SetX(result[0]);
991 Vu .SetX(resVu [0]);
992 Vv .SetX(resVv [0]);
993 Vuu.SetX(resVuu[0]);
994 Vvv.SetX(resVvv[0]);
995 Vuv.SetX(resVuv[0]);
996
997 P .SetY(result[1]);
998 Vu .SetY(resVu [1]);
999 Vv .SetY(resVv [1]);
1000 Vuu.SetY(resVuu[1]);
1001 Vvv.SetY(resVvv[1]);
1002 Vuv.SetY(resVuv[1]);
1003
1004 P .SetZ(result[2]);
1005 Vu .SetZ(resVu [2]);
1006 Vv .SetZ(resVv [2]);
1007 Vuu.SetZ(resVuu[2]);
1008 Vvv.SetZ(resVvv[2]);
1009 Vuv.SetZ(resVuv[2]);
1010}
1011
1012//=======================================================================
1013//function : D3
1014//purpose :
1015//=======================================================================
1016
1017void BSplSLib::D3
1018(const Standard_Real U,
1019 const Standard_Real V,
1020 const Standard_Integer UIndex,
1021 const Standard_Integer VIndex,
1022 const TColgp_Array2OfPnt& Poles,
0e14656b 1023 const TColStd_Array2OfReal* Weights,
7fd59977 1024 const TColStd_Array1OfReal& UKnots,
1025 const TColStd_Array1OfReal& VKnots,
0e14656b 1026 const TColStd_Array1OfInteger* UMults,
1027 const TColStd_Array1OfInteger* VMults,
7fd59977 1028 const Standard_Integer UDegree,
1029 const Standard_Integer VDegree,
1030 const Standard_Boolean URat,
1031 const Standard_Boolean VRat,
1032 const Standard_Boolean UPer,
1033 const Standard_Boolean VPer,
1034 gp_Pnt& P,
1035 gp_Vec& Vu,
1036 gp_Vec& Vv,
1037 gp_Vec& Vuu,
1038 gp_Vec& Vvv,
1039 gp_Vec& Vuv,
1040 gp_Vec& Vuuu,
1041 gp_Vec& Vvvv,
1042 gp_Vec& Vuuv,
1043 gp_Vec& Vuvv)
1044{
1045 Standard_Boolean rational;
1046// Standard_Integer k,dim,dim2;
1047 Standard_Integer dim,dim2;
1048 Standard_Real u1,u2;
1049 Standard_Integer d1,d2;
1050 Standard_Real *result, *resVu, *resVv, *resVuu, *resVvv, *resVuv,
1051 *resVuuu, *resVvvv, *resVuuv, *resVuvv;
1052 BSplSLib_DataContainer dc (UDegree, VDegree);
1053 if (PrepareEval
1054 (U,V,UIndex,VIndex,UDegree,VDegree,URat,VRat,UPer,VPer,
1055 Poles,Weights,UKnots,VKnots,UMults,VMults,
1056 u1,u2,d1,d2,rational,dc)) {
1057 if (rational) {
1058 dim = 4;
1059 dim2 = (d2 + 1) << 2;
1060 BSplCLib::Bohm (u1,d1,3,*dc.knots1,dim2,*dc.poles);
1061 BSplCLib::Bohm (u2,d2,3,*dc.knots2,dim ,*dc.poles);
1062 BSplCLib::Bohm (u2,d2,2,*dc.knots2,dim ,*(dc.poles + dim2));
1063 if (d1 > 1)
1064 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*(dc.poles + (dim2 << 1)));
1065 if (d1 > 2)
1066 BSplCLib::Eval(u2,d2 ,*dc.knots2,dim ,*(dc.poles + (dim2 << 1) + dim2));
1067 BSplSLib::RationalDerivative(d1,d2,3,3,*dc.poles,*dc.ders);
1068 result = dc.ders;
1069 resVu = result + 12;
1070 resVv = result + 3;
1071 resVuu = result + 24;
1072 resVvv = result + 6;
1073 resVuv = result + 15;
1074 resVuuu = result + 36;
1075 resVvvv = result + 9;
1076 resVuuv = result + 27;
1077 resVuvv = result + 18;
1078 }
1079 else {
1080 dim = 3;
1081 dim2 = (d2 + 1);
1082 dim2 = (dim2 << 1) + dim2;
1083 BSplCLib::Bohm (u1,d1,3,*dc.knots1,dim2,*dc.poles);
1084 BSplCLib::Bohm (u2,d2,3,*dc.knots2,dim ,*dc.poles);
1085 BSplCLib::Bohm (u2,d2,2,*dc.knots2,dim ,*(dc.poles + dim2));
1086 if (d1 > 1)
1087 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*(dc.poles + (dim2 << 1)));
1088 if (d1 > 2)
1089 BSplCLib::Eval(u2,d2 ,*dc.knots2,dim ,*(dc.poles + (dim2 << 1) + dim2));
1090 result = dc.poles;
1091 resVu = result + dim2;
1092 resVv = result + 3;
1093 if (UDegree <= 1) {
1094 resVuu = BSplSLib_zero;
1095 resVuuv = BSplSLib_zero;
1096 }
1097 else {
1098 resVuu = result + (dim2 << 1);
1099 resVuuv = result + (dim2 << 1) + 3;
1100 }
1101 if (VDegree <= 1) {
1102 resVvv = BSplSLib_zero;
1103 resVuvv = BSplSLib_zero;
1104 }
1105 else {
1106 resVvv = result + 6;
1107 resVuvv = result + dim2 + 6;
1108 }
1109 resVuv = result + (d2 << 1) + d2 + 6;
1110 if (UDegree <= 2) resVuuu = BSplSLib_zero;
1111 else resVuuu = result + (dim2 << 1) + dim2;
1112 if (VDegree <= 2) resVvvv = BSplSLib_zero;
1113 else resVvvv = result + 9;
1114 }
1115 }
1116 else {
1117 if (rational) {
1118 dim = 4;
1119 dim2 = (d2 + 1) << 2;
1120 BSplCLib::Bohm (u1,d1,3,*dc.knots1,dim2,*dc.poles);
1121 BSplCLib::Bohm (u2,d2,3,*dc.knots2,dim ,*dc.poles);
1122 BSplCLib::Bohm (u2,d2,2,*dc.knots2,dim ,*(dc.poles + dim2));
1123 if (d1 > 1)
1124 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*(dc.poles + (dim2 << 1)));
1125 if (d1 > 2)
1126 BSplCLib::Eval(u2,d2 ,*dc.knots2,dim ,*(dc.poles + (dim2 << 1) + dim2));
1127 BSplSLib::RationalDerivative(d1,d2,3,3,*dc.poles,*dc.ders);
1128 result = dc.ders;
1129 resVu = result + 3;
1130 resVv = result + 12;
1131 resVuu = result + 6;
1132 resVvv = result + 24;
1133 resVuv = result + 15;
1134 resVuuu = result + 9;
1135 resVvvv = result + 36;
1136 resVuuv = result + 18;
1137 resVuvv = result + 27;
1138 }
1139 else {
1140 dim = 3;
1141 dim2 = (d2 + 1);
1142 dim2 = (dim2 << 1) + dim2;
1143 BSplCLib::Bohm (u1,d1,3,*dc.knots1,dim2,*dc.poles);
1144 BSplCLib::Bohm (u2,d2,3,*dc.knots2,dim ,*dc.poles);
1145 BSplCLib::Bohm (u2,d2,2,*dc.knots2,dim ,*(dc.poles + dim2));
1146 if (d1 > 1)
1147 BSplCLib::Bohm(u2,d2,1,*dc.knots2,dim ,*(dc.poles + (dim2 << 1)));
1148 if (d1 > 2)
1149 BSplCLib::Eval(u2,d2 ,*dc.knots2,dim ,*(dc.poles + (dim2 << 1) + dim2));
1150 result = dc.poles;
1151 resVu = result + 3;
1152 resVv = result + dim2;
1153 if (UDegree <= 1) {
1154 resVuu = BSplSLib_zero;
1155 resVuuv = BSplSLib_zero;
1156 }
1157 else {
1158 resVuu = result + 6;
1159 resVuuv = result + dim2 + 6;
1160 }
1161 if (VDegree <= 1) {
1162 resVvv = BSplSLib_zero;
1163 resVuvv = BSplSLib_zero;
1164 }
1165 else {
1166 resVvv = result + (dim2 << 1);
1167 resVuvv = result + (dim2 << 1) + 3;
1168 }
1169 resVuv = result + (d2 << 1) + d2 + 6;
1170 if (UDegree <= 2) resVuuu = BSplSLib_zero;
1171 else resVuuu = result + 9;
1172 if (VDegree <= 2) resVvvv = BSplSLib_zero;
1173 else resVvvv = result + (dim2 << 1) + dim2;
1174 }
1175 }
1176
1177 P .SetX(result [0]);
1178 Vu .SetX(resVu [0]);
1179 Vv .SetX(resVv [0]);
1180 Vuu .SetX(resVuu [0]);
1181 Vvv .SetX(resVvv [0]);
1182 Vuv .SetX(resVuv [0]);
1183 Vuuu.SetX(resVuuu[0]);
1184 Vvvv.SetX(resVvvv[0]);
1185 Vuuv.SetX(resVuuv[0]);
1186 Vuvv.SetX(resVuvv[0]);
1187
1188 P .SetY(result [1]);
1189 Vu .SetY(resVu [1]);
1190 Vv .SetY(resVv [1]);
1191 Vuu .SetY(resVuu [1]);
1192 Vvv .SetY(resVvv [1]);
1193 Vuv .SetY(resVuv [1]);
1194 Vuuu.SetY(resVuuu[1]);
1195 Vvvv.SetY(resVvvv[1]);
1196 Vuuv.SetY(resVuuv[1]);
1197 Vuvv.SetY(resVuvv[1]);
1198
1199 P .SetZ(result [2]);
1200 Vu .SetZ(resVu [2]);
1201 Vv .SetZ(resVv [2]);
1202 Vuu .SetZ(resVuu [2]);
1203 Vvv .SetZ(resVvv [2]);
1204 Vuv .SetZ(resVuv [2]);
1205 Vuuu.SetZ(resVuuu[2]);
1206 Vvvv.SetZ(resVvvv[2]);
1207 Vuuv.SetZ(resVuuv[2]);
1208 Vuvv.SetZ(resVuvv[2]);
1209}
1210
1211//=======================================================================
1212//function : DN
1213//purpose :
1214//=======================================================================
1215
1216void BSplSLib::DN
1217(const Standard_Real U,
1218 const Standard_Real V,
1219 const Standard_Integer Nu,
1220 const Standard_Integer Nv,
1221 const Standard_Integer UIndex,
1222 const Standard_Integer VIndex,
1223 const TColgp_Array2OfPnt& Poles,
0e14656b 1224 const TColStd_Array2OfReal* Weights,
7fd59977 1225 const TColStd_Array1OfReal& UKnots,
1226 const TColStd_Array1OfReal& VKnots,
0e14656b 1227 const TColStd_Array1OfInteger* UMults,
1228 const TColStd_Array1OfInteger* VMults,
7fd59977 1229 const Standard_Integer UDegree,
1230 const Standard_Integer VDegree,
1231 const Standard_Boolean URat,
1232 const Standard_Boolean VRat,
1233 const Standard_Boolean UPer,
1234 const Standard_Boolean VPer,
1235 gp_Vec& Vn)
1236{
1237 Standard_Boolean rational;
1238 Standard_Integer k,dim;
1239 Standard_Real u1,u2;
1240 Standard_Integer d1,d2;
1241
1242 BSplSLib_DataContainer dc (UDegree, VDegree);
1243 Standard_Boolean ufirst = PrepareEval
1244 (U,V,UIndex,VIndex,UDegree,VDegree,URat,VRat,UPer,VPer,
1245 Poles,Weights,UKnots,VKnots,UMults,VMults,
1246 u1,u2,d1,d2,rational,dc);
1247 dim = rational ? 4 : 3;
1248
1249 if (!rational) {
1250 if ((Nu > UDegree) || (Nv > VDegree)) {
1251 Vn.SetX(0.);
1252 Vn.SetY(0.);
1253 Vn.SetZ(0.);
1254 return;
1255 }
1256 }
1257
1258 Standard_Integer n1 = ufirst ? Nu : Nv;
1259 Standard_Integer n2 = ufirst ? Nv : Nu;
1260
1261 BSplCLib::Bohm(u1,d1,n1,*dc.knots1,dim * (d2 + 1),*dc.poles);
1262
1263 for (k = 0; k <= Min(n1,d1); k++)
1264 BSplCLib::Bohm(u2,d2,n2,*dc.knots2,dim,*(dc.poles+k*dim*(d2+1)));
1265
1266 Standard_Real *result;
1267 if (rational) {
1268 BSplSLib::RationalDerivative(d1,d2,n1,n2,*dc.poles,*dc.ders,Standard_False);
1269 result = dc.ders; // because Standard_False ci-dessus.
1270
1271 }
1272 else {
1273 result = dc.poles + (n1 * (d2+1) + n2) * dim ;
1274 }
1275
1276 Vn.SetX(result[0]);
1277 Vn.SetY(result[1]);
1278 Vn.SetZ(result[2]);
1279}
1280
1281//
1282// Surface modifications
1283//
1284// a surface is processed as a curve of curves.
1285// i.e : a curve of parameter u where the current point is the set of poles
1286// of the iso.
1287//
1288
1289//=======================================================================
1290//function : Iso
1291//purpose :
1292//=======================================================================
1293
1294void BSplSLib::Iso(const Standard_Real Param,
1295 const Standard_Boolean IsU,
1296 const TColgp_Array2OfPnt& Poles,
0e14656b 1297 const TColStd_Array2OfReal* Weights,
7fd59977 1298 const TColStd_Array1OfReal& Knots,
0e14656b 1299 const TColStd_Array1OfInteger* Mults,
7fd59977 1300 const Standard_Integer Degree,
1301 const Standard_Boolean Periodic,
1302 TColgp_Array1OfPnt& CPoles,
0e14656b 1303 TColStd_Array1OfReal* CWeights)
7fd59977 1304{
1305 Standard_Integer index = 0;
1306 Standard_Real u = Param;
0e14656b 1307 Standard_Boolean rational = Weights != NULL;
7fd59977 1308 Standard_Integer dim = rational ? 4 : 3;
1309
1310 // compute local knots
1311
f7b4312f 1312 NCollection_LocalArray<Standard_Real> locknots1 (2*Degree);
7fd59977 1313 BSplCLib::LocateParameter(Degree,Knots,Mults,u,Periodic,index,u);
1314 BSplCLib::BuildKnots(Degree,index,Periodic,Knots,Mults,*locknots1);
0e14656b 1315 if (Mults == NULL)
7fd59977 1316 index -= Knots.Lower() + Degree;
1317 else
0e14656b 1318 index = BSplCLib::PoleIndex(Degree,index,Periodic,*Mults);
7fd59977 1319
1320
1321 // copy the local poles
1322
1323// Standard_Integer f1,l1,f2,l2,i,j,k;
1324 Standard_Integer f1,l1,f2,l2,i,j;
1325
1326 if (IsU) {
1327 f1 = Poles.LowerRow();
1328 l1 = Poles.UpperRow();
1329 f2 = Poles.LowerCol();
1330 l2 = Poles.UpperCol();
1331 }
1332 else {
1333 f1 = Poles.LowerCol();
1334 l1 = Poles.UpperCol();
1335 f2 = Poles.LowerRow();
1336 l2 = Poles.UpperRow();
1337 }
1338
f7b4312f 1339 NCollection_LocalArray<Standard_Real> locpoles ((Degree+1) * (l2-f2+1) * dim);
7fd59977 1340
1341 Standard_Real w, *pole = locpoles;
1342 index += f1;
1343
1344 for (i = 0; i <= Degree; i++) {
1345
1346 for (j = f2; j <= l2; j++) {
1347
1348 const gp_Pnt& P = IsU ? Poles(index,j) : Poles(j,index);
1349 if (rational) {
0e14656b 1350 pole[3] = w = IsU ? (*Weights)(index,j) : (*Weights)(j,index);
7fd59977 1351 pole[0] = P.X() * w;
1352 pole[1] = P.Y() * w;
1353 pole[2] = P.Z() * w;
1354 }
1355 else {
1356 pole[0] = P.X();
1357 pole[1] = P.Y();
1358 pole[2] = P.Z();
1359 }
1360 pole += dim;
1361 }
1362 index++;
1363 if (index > l1) index = f1;
1364 }
1365
1366 // compute the iso
1367 BSplCLib::Eval(u,Degree,*locknots1,(l2-f2+1)*dim,*locpoles);
1368
1369 // get the result
1370 pole = locpoles;
1371
1372 for (i = CPoles.Lower(); i <= CPoles.Upper(); i++) {
1373 gp_Pnt& P = CPoles(i);
1374 if (rational) {
0e14656b 1375 (*CWeights)(i) = w = pole[3];
7fd59977 1376 P.SetX( pole[0] / w);
1377 P.SetY( pole[1] / w);
1378 P.SetZ( pole[2] / w);
1379 }
1380 else {
1381 P.SetX( pole[0]);
1382 P.SetY( pole[1]);
1383 P.SetZ( pole[2]);
1384 }
1385 pole += dim;
1386 }
1387
1388 // if the input is not rational but weights are wanted
0e14656b 1389 if (!rational && (CWeights != NULL)) {
7fd59977 1390
0e14656b 1391 for (i = CWeights->Lower(); i <= CWeights->Upper(); i++)
1392 (*CWeights)(i) = 1.;
7fd59977 1393 }
1394}
1395
1396//=======================================================================
1397//function : Reverse
1398//purpose :
1399//=======================================================================
1400
1401void BSplSLib::Reverse( TColgp_Array2OfPnt& Poles,
1402 const Standard_Integer Last,
1403 const Standard_Boolean UDirection)
1404{
1405 Standard_Integer i,j, l = Last;
1406 if ( UDirection) {
1407 l = Poles.LowerRow() +
1408 (l - Poles.LowerRow())%(Poles.ColLength());
1409 TColgp_Array2OfPnt temp(0, Poles.ColLength()-1,
1410 Poles.LowerCol(), Poles.UpperCol());
1411
1412 for (i = Poles.LowerRow(); i <= l; i++) {
1413
1414 for (j = Poles.LowerCol(); j <= Poles.UpperCol(); j++) {
1415 temp(l-i,j) = Poles(i,j);
1416 }
1417 }
1418
1419 for (i = l+1; i <= Poles.UpperRow(); i++) {
1420
1421 for (j = Poles.LowerCol(); j <= Poles.UpperCol(); j++) {
1422 temp(l+Poles.ColLength()-i,j) = Poles(i,j);
1423 }
1424 }
1425
1426 for (i = Poles.LowerRow(); i <= Poles.UpperRow(); i++) {
1427
1428 for (j = Poles.LowerCol(); j <= Poles.UpperCol(); j++) {
1429 Poles(i,j) = temp (i-Poles.LowerRow(),j);
1430 }
1431 }
1432 }
1433 else {
1434 l = Poles.LowerCol() +
1435 (l - Poles.LowerCol())%(Poles.RowLength());
1436 TColgp_Array2OfPnt temp(Poles.LowerRow(), Poles.UpperRow(),
1437 0, Poles.RowLength()-1);
1438
1439 for (j = Poles.LowerCol(); j <= l; j++) {
1440
1441 for (i = Poles.LowerRow(); i <= Poles.UpperRow(); i++) {
1442 temp(i,l-j) = Poles(i,j);
1443 }
1444 }
1445
1446 for (j = l+1; j <= Poles.UpperCol(); j++) {
1447
1448 for (i = Poles.LowerRow(); i <= Poles.UpperRow(); i++) {
1449 temp(i,l+Poles.RowLength()-j) = Poles(i,j);
1450 }
1451 }
1452
1453 for (i = Poles.LowerRow(); i <= Poles.UpperRow(); i++) {
1454
1455 for (j = Poles.LowerCol(); j <= Poles.UpperCol(); j++) {
1456 Poles(i,j) = temp (i,j-Poles.LowerCol());
1457 }
1458 }
1459 }
1460}
1461
1462//=======================================================================
1463//function : Reverse
1464//purpose :
1465//=======================================================================
1466
1467void BSplSLib::Reverse( TColStd_Array2OfReal& Weights,
1468 const Standard_Integer Last,
1469 const Standard_Boolean UDirection)
1470{
1471 Standard_Integer i,j, l = Last;
1472 if ( UDirection) {
1473 l = Weights.LowerRow() +
1474 (l - Weights.LowerRow())%(Weights.ColLength());
1475 TColStd_Array2OfReal temp(0, Weights.ColLength()-1,
1476 Weights.LowerCol(), Weights.UpperCol());
1477
1478 for (i = Weights.LowerRow(); i <= l; i++) {
1479
1480 for (j = Weights.LowerCol(); j <= Weights.UpperCol(); j++) {
1481 temp(l-i,j) = Weights(i,j);
1482 }
1483 }
1484
1485 for (i = l+1; i <= Weights.UpperRow(); i++) {
1486
1487 for (j = Weights.LowerCol(); j <= Weights.UpperCol(); j++) {
1488 temp(l+Weights.ColLength()-i,j) = Weights(i,j);
1489 }
1490 }
1491
1492 for (i = Weights.LowerRow(); i <= Weights.UpperRow(); i++) {
1493
1494 for (j = Weights.LowerCol(); j <= Weights.UpperCol(); j++) {
1495 Weights(i,j) = temp (i-Weights.LowerRow(),j);
1496 }
1497 }
1498 }
1499 else {
1500 l = Weights.LowerCol() +
1501 (l - Weights.LowerCol())%(Weights.RowLength());
1502 TColStd_Array2OfReal temp(Weights.LowerRow(), Weights.UpperRow(),
1503 0, Weights.RowLength()-1);
1504
1505 for (j = Weights.LowerCol(); j <= l; j++) {
1506
1507 for (i = Weights.LowerRow(); i <= Weights.UpperRow(); i++) {
1508 temp(i,l-j) = Weights(i,j);
1509 }
1510 }
1511
1512 for (j = l+1; j <= Weights.UpperCol(); j++) {
1513
1514 for (i = Weights.LowerRow(); i <= Weights.UpperRow(); i++) {
1515 temp(i,l+Weights.RowLength()-j) = Weights(i,j);
1516 }
1517 }
1518
1519 for (i = Weights.LowerRow(); i <= Weights.UpperRow(); i++) {
1520
1521 for (j = Weights.LowerCol(); j <= Weights.UpperCol(); j++) {
1522 Weights(i,j) = temp (i,j-Weights.LowerCol());
1523 }
1524 }
1525 }
1526}
1527
1528//=======================================================================
1529//function : IsRational
1530//purpose :
1531//=======================================================================
1532
1533Standard_Boolean BSplSLib::IsRational
1534(const TColStd_Array2OfReal& Weights,
1535 const Standard_Integer I1,
1536 const Standard_Integer I2,
1537 const Standard_Integer J1,
1538 const Standard_Integer J2,
1539 const Standard_Real Epsi)
1540{
1541 Standard_Real eps = (Epsi > 0.0) ? Epsi : Epsilon(Weights(I1,I2));
1542 Standard_Integer i,j;
1543 Standard_Integer fi = Weights.LowerRow(), li = Weights.ColLength();
1544 Standard_Integer fj = Weights.LowerCol(), lj = Weights.RowLength();
1545
1546 for (i = I1 - fi; i < I2 - fi; i++) {
1547
1548 for (j = J1 - fj; j < J2 - fj; j++) {
1549 if (Abs(Weights(fi+i%li,fj+j%lj)-Weights(fi+(i+1)%li,fj+j%lj))>eps)
1550 return Standard_True;
1551 }
1552 }
1553 return Standard_False;
1554}
1555
1556//=======================================================================
1557//function : SetPoles
1558//purpose :
1559//=======================================================================
1560
1561void BSplSLib::SetPoles(const TColgp_Array2OfPnt& Poles,
1562 TColStd_Array1OfReal& FP,
1563 const Standard_Boolean UDirection)
1564{
1565 Standard_Integer i, j, l = FP.Lower();
1566 Standard_Integer PLowerRow = Poles.LowerRow();
1567 Standard_Integer PUpperRow = Poles.UpperRow();
1568 Standard_Integer PLowerCol = Poles.LowerCol();
1569 Standard_Integer PUpperCol = Poles.UpperCol();
1570 if (UDirection) {
1571
1572 for ( i = PLowerRow; i <= PUpperRow; i++) {
1573
1574 for ( j = PLowerCol; j <= PUpperCol; j++) {
1575 const gp_Pnt& P = Poles.Value(i,j);
1576 FP(l) = P.X(); l++;
1577 FP(l) = P.Y(); l++;
1578 FP(l) = P.Z(); l++;
1579 }
1580 }
1581 }
1582 else {
1583
1584 for ( j = PLowerCol; j <= PUpperCol; j++) {
1585
1586 for ( i = PLowerRow; i <= PUpperRow; i++) {
1587 const gp_Pnt& P = Poles.Value(i,j);
1588 FP(l) = P.X(); l++;
1589 FP(l) = P.Y(); l++;
1590 FP(l) = P.Z(); l++;
1591 }
1592 }
1593 }
1594}
1595
1596//=======================================================================
1597//function : SetPoles
1598//purpose :
1599//=======================================================================
1600
1601void BSplSLib::SetPoles(const TColgp_Array2OfPnt& Poles,
1602 const TColStd_Array2OfReal& Weights,
1603 TColStd_Array1OfReal& FP,
1604 const Standard_Boolean UDirection)
1605{
1606 Standard_Integer i, j, l = FP.Lower();
1607 Standard_Integer PLowerRow = Poles.LowerRow();
1608 Standard_Integer PUpperRow = Poles.UpperRow();
1609 Standard_Integer PLowerCol = Poles.LowerCol();
1610 Standard_Integer PUpperCol = Poles.UpperCol();
1611 if (UDirection) {
1612
1613 for ( i = PLowerRow; i <= PUpperRow; i++) {
1614
1615 for ( j = PLowerCol; j <= PUpperCol; j++) {
1616 const gp_Pnt& P = Poles .Value(i,j);
1617 Standard_Real w = Weights.Value(i,j);
1618 FP(l) = P.X() * w; l++;
1619 FP(l) = P.Y() * w; l++;
1620 FP(l) = P.Z() * w; l++;
1621 FP(l) = w; l++;
1622 }
1623 }
1624 }
1625 else {
1626
1627 for ( j = PLowerCol; j <= PUpperCol; j++) {
1628
1629 for ( i = PLowerRow; i <= PUpperRow; i++) {
1630 const gp_Pnt& P = Poles .Value(i,j);
1631 Standard_Real w = Weights.Value(i,j);
1632 FP(l) = P.X() * w; l++;
1633 FP(l) = P.Y() * w; l++;
1634 FP(l) = P.Z() * w; l++;
1635 FP(l) = w; l++;
1636 }
1637 }
1638 }
1639}
1640
1641//=======================================================================
1642//function : GetPoles
1643//purpose :
1644//=======================================================================
1645
1646void BSplSLib::GetPoles(const TColStd_Array1OfReal& FP,
1647 TColgp_Array2OfPnt& Poles,
1648 const Standard_Boolean UDirection)
1649{
1650 Standard_Integer i, j, l = FP.Lower();
1651 Standard_Integer PLowerRow = Poles.LowerRow();
1652 Standard_Integer PUpperRow = Poles.UpperRow();
1653 Standard_Integer PLowerCol = Poles.LowerCol();
1654 Standard_Integer PUpperCol = Poles.UpperCol();
1655 if (UDirection) {
1656
1657 for ( i = PLowerRow; i <= PUpperRow; i++) {
1658
1659 for ( j = PLowerCol; j <= PUpperCol; j++) {
1660 gp_Pnt& P = Poles.ChangeValue(i,j);
1661 P.SetX(FP(l)); l++;
1662 P.SetY(FP(l)); l++;
1663 P.SetZ(FP(l)); l++;
1664 }
1665 }
1666 }
1667 else {
1668
1669 for ( j = PLowerCol; j <= PUpperCol; j++) {
1670
1671 for ( i = PLowerRow; i <= PUpperRow; i++) {
1672 gp_Pnt& P = Poles.ChangeValue(i,j);
1673 P.SetX(FP(l)); l++;
1674 P.SetY(FP(l)); l++;
1675 P.SetZ(FP(l)); l++;
1676 }
1677 }
1678 }
1679}
1680
1681//=======================================================================
1682//function : GetPoles
1683//purpose :
1684//=======================================================================
1685
1686void BSplSLib::GetPoles(const TColStd_Array1OfReal& FP,
1687 TColgp_Array2OfPnt& Poles,
1688 TColStd_Array2OfReal& Weights,
1689 const Standard_Boolean UDirection)
1690{
1691 Standard_Integer i, j, l = FP.Lower();
1692 Standard_Integer PLowerRow = Poles.LowerRow();
1693 Standard_Integer PUpperRow = Poles.UpperRow();
1694 Standard_Integer PLowerCol = Poles.LowerCol();
1695 Standard_Integer PUpperCol = Poles.UpperCol();
1696 if (UDirection) {
1697
1698 for ( i = PLowerRow; i <= PUpperRow; i++) {
1699
1700 for ( j = PLowerCol; j <= PUpperCol; j++) {
1701 Standard_Real w = FP( l + 3);
1702 Weights(i,j) = w;
1703 gp_Pnt& P = Poles.ChangeValue(i,j);
1704 P.SetX(FP(l) / w); l++;
1705 P.SetY(FP(l) / w); l++;
1706 P.SetZ(FP(l) / w); l++;
1707 l++;
1708 }
1709 }
1710 }
1711 else {
1712
1713 for ( j = PLowerCol; j <= PUpperCol; j++) {
1714
1715 for ( i = PLowerRow; i <= PUpperRow; i++) {
1716 Standard_Real w = FP( l + 3);
1717 Weights(i,j) = w;
1718 gp_Pnt& P = Poles.ChangeValue(i,j);
1719 P.SetX(FP(l) / w); l++;
1720 P.SetY(FP(l) / w); l++;
1721 P.SetZ(FP(l) / w); l++;
1722 l++;
1723 }
1724 }
1725 }
1726}
1727
1728//=======================================================================
1729//function : InsertKnots
1730//purpose :
1731//=======================================================================
1732
1733void BSplSLib::InsertKnots(const Standard_Boolean UDirection,
1734 const Standard_Integer Degree,
1735 const Standard_Boolean Periodic,
1736 const TColgp_Array2OfPnt& Poles,
0e14656b 1737 const TColStd_Array2OfReal* Weights,
7fd59977 1738 const TColStd_Array1OfReal& Knots,
1739 const TColStd_Array1OfInteger& Mults,
1740 const TColStd_Array1OfReal& AddKnots,
0e14656b 1741 const TColStd_Array1OfInteger* AddMults,
7fd59977 1742 TColgp_Array2OfPnt& NewPoles,
0e14656b 1743 TColStd_Array2OfReal* NewWeights,
7fd59977 1744 TColStd_Array1OfReal& NewKnots,
1745 TColStd_Array1OfInteger& NewMults,
1746 const Standard_Real Epsilon,
1747 const Standard_Boolean Add )
1748{
0e14656b 1749 Standard_Boolean rational = Weights != NULL;
7fd59977 1750 Standard_Integer dim = 3;
1751 if (rational) dim++;
1752
1753 TColStd_Array1OfReal poles( 1, dim*Poles.RowLength()*Poles.ColLength());
1754 TColStd_Array1OfReal
1755 newpoles( 1, dim*NewPoles.RowLength()*NewPoles.ColLength());
1756
0e14656b 1757 if (rational) SetPoles(Poles,*Weights,poles,UDirection);
7fd59977 1758 else SetPoles(Poles,poles,UDirection);
1759
1760 if (UDirection) {
1761 dim *= Poles.RowLength();
1762 }
1763 else {
1764 dim *= Poles.ColLength();
1765 }
1766 BSplCLib::InsertKnots(Degree,Periodic,dim,poles,Knots,Mults,
1767 AddKnots,AddMults,newpoles,NewKnots,NewMults,
1768 Epsilon,Add);
1769
0e14656b 1770 if (rational) GetPoles(newpoles,NewPoles,*NewWeights,UDirection);
7fd59977 1771 else GetPoles(newpoles,NewPoles,UDirection);
1772}
1773
1774//=======================================================================
1775//function : RemoveKnot
1776//purpose :
1777//=======================================================================
1778
1779Standard_Boolean BSplSLib::RemoveKnot
1780(const Standard_Boolean UDirection,
1781 const Standard_Integer Index,
1782 const Standard_Integer Mult,
1783 const Standard_Integer Degree,
1784 const Standard_Boolean Periodic,
1785 const TColgp_Array2OfPnt& Poles,
0e14656b 1786 const TColStd_Array2OfReal* Weights,
7fd59977 1787 const TColStd_Array1OfReal& Knots,
1788 const TColStd_Array1OfInteger& Mults,
1789 TColgp_Array2OfPnt& NewPoles,
0e14656b 1790 TColStd_Array2OfReal* NewWeights,
7fd59977 1791 TColStd_Array1OfReal& NewKnots,
1792 TColStd_Array1OfInteger& NewMults,
1793 const Standard_Real Tolerance)
1794{
0e14656b 1795 Standard_Boolean rational = Weights != NULL;
7fd59977 1796 Standard_Integer dim = 3;
1797 if (rational) dim++;
1798
1799 TColStd_Array1OfReal poles( 1, dim*Poles.RowLength()*Poles.ColLength());
1800 TColStd_Array1OfReal
1801 newpoles( 1, dim*NewPoles.RowLength()*NewPoles.ColLength());
1802
0e14656b 1803 if (rational) SetPoles(Poles,*Weights,poles,UDirection);
7fd59977 1804 else SetPoles(Poles,poles,UDirection);
1805
1806 if (UDirection) {
1807 dim *= Poles.RowLength();
1808 }
1809 else {
1810 dim *= Poles.ColLength();
1811 }
1812
1813 if ( !BSplCLib::RemoveKnot(Index,Mult,Degree,Periodic,dim,
1814 poles,Knots,Mults,newpoles,NewKnots,NewMults,
1815 Tolerance))
1816 return Standard_False;
1817
0e14656b 1818 if (rational) GetPoles(newpoles,NewPoles,*NewWeights,UDirection);
7fd59977 1819 else GetPoles(newpoles,NewPoles,UDirection);
1820 return Standard_True;
1821}
1822
1823//=======================================================================
1824//function : IncreaseDegree
1825//purpose :
1826//=======================================================================
1827
1828void BSplSLib::IncreaseDegree
1829(const Standard_Boolean UDirection,
1830 const Standard_Integer Degree,
1831 const Standard_Integer NewDegree,
1832 const Standard_Boolean Periodic,
1833 const TColgp_Array2OfPnt& Poles,
0e14656b 1834 const TColStd_Array2OfReal* Weights,
7fd59977 1835 const TColStd_Array1OfReal& Knots,
1836 const TColStd_Array1OfInteger& Mults,
1837 TColgp_Array2OfPnt& NewPoles,
0e14656b 1838 TColStd_Array2OfReal* NewWeights,
7fd59977 1839 TColStd_Array1OfReal& NewKnots,
1840 TColStd_Array1OfInteger& NewMults)
1841{
0e14656b 1842 Standard_Boolean rational = Weights != NULL;
7fd59977 1843 Standard_Integer dim = 3;
1844 if (rational) dim++;
1845
1846 TColStd_Array1OfReal poles( 1, dim*Poles.RowLength()*Poles.ColLength());
1847 TColStd_Array1OfReal
1848 newpoles( 1, dim*NewPoles.RowLength()*NewPoles.ColLength());
1849
0e14656b 1850 if (rational) SetPoles(Poles,*Weights,poles,UDirection);
7fd59977 1851 else SetPoles(Poles,poles,UDirection);
1852
1853 if (UDirection) {
1854 dim *= Poles.RowLength();
1855 }
1856 else {
1857 dim *= Poles.ColLength();
1858 }
1859
1860 BSplCLib::IncreaseDegree(Degree,NewDegree,Periodic,dim,poles,Knots,Mults,
1861 newpoles,NewKnots,NewMults);
1862
0e14656b 1863 if (rational) GetPoles(newpoles,NewPoles,*NewWeights,UDirection);
7fd59977 1864 else GetPoles(newpoles,NewPoles,UDirection);
1865}
1866
1867//=======================================================================
1868//function : Unperiodize
1869//purpose :
1870//=======================================================================
1871
1872void BSplSLib::Unperiodize
1873(const Standard_Boolean UDirection,
1874 const Standard_Integer Degree,
1875 const TColStd_Array1OfInteger& Mults,
1876 const TColStd_Array1OfReal& Knots,
1877 const TColgp_Array2OfPnt& Poles,
0e14656b 1878 const TColStd_Array2OfReal* Weights,
7fd59977 1879 TColStd_Array1OfInteger& NewMults,
1880 TColStd_Array1OfReal& NewKnots,
1881 TColgp_Array2OfPnt& NewPoles,
0e14656b 1882 TColStd_Array2OfReal* NewWeights)
7fd59977 1883{
0e14656b 1884 Standard_Boolean rational = Weights != NULL;
7fd59977 1885 Standard_Integer dim = 3;
1886 if (rational) dim++;
1887
1888 TColStd_Array1OfReal poles( 1, dim*Poles.RowLength()*Poles.ColLength());
1889 TColStd_Array1OfReal
1890 newpoles( 1, dim*NewPoles.RowLength()*NewPoles.ColLength());
1891
0e14656b 1892 if (rational) SetPoles(Poles,*Weights,poles,UDirection);
7fd59977 1893 else SetPoles(Poles,poles,UDirection);
1894
1895 if (UDirection) {
1896 dim *= Poles.RowLength();
1897 }
1898 else {
1899 dim *= Poles.ColLength();
1900 }
1901
1902 BSplCLib::Unperiodize(Degree, dim, Mults, Knots, poles,
1903 NewMults, NewKnots, newpoles);
1904
0e14656b 1905 if (rational) GetPoles(newpoles,NewPoles,*NewWeights,UDirection);
7fd59977 1906 else GetPoles(newpoles,NewPoles,UDirection);
1907}
1908
1909//=======================================================================
1910//function : BuildCache
1911//purpose : Stores theTaylor expansion normalized between 0,1 in the
1912// Cache : in case of a rational function the Poles are
1913// stored in homogeneous form
1914//=======================================================================
1915
1916void BSplSLib::BuildCache
1917(const Standard_Real U,
1918 const Standard_Real V,
1919 const Standard_Real USpanDomain,
1920 const Standard_Real VSpanDomain,
1921 const Standard_Boolean UPeriodic,
1922 const Standard_Boolean VPeriodic,
1923 const Standard_Integer UDegree,
1924 const Standard_Integer VDegree,
1925 const Standard_Integer UIndex,
1926 const Standard_Integer VIndex,
1927 const TColStd_Array1OfReal& UFlatKnots,
1928 const TColStd_Array1OfReal& VFlatKnots,
1929 const TColgp_Array2OfPnt& Poles,
0e14656b 1930 const TColStd_Array2OfReal* Weights,
7fd59977 1931 TColgp_Array2OfPnt& CachePoles,
0e14656b 1932 TColStd_Array2OfReal* CacheWeights)
7fd59977 1933{
1934 Standard_Boolean rational,rational_u,rational_v,flag_u_or_v;
1935 Standard_Integer kk,d1,d1p1,d2,d2p1,ii,jj,iii,jjj,Index;
1936 Standard_Real u1,min_degree_domain,max_degree_domain,f,factor[2],u2;
0e14656b 1937 if (Weights != NULL)
7fd59977 1938 rational_u = rational_v = Standard_True;
1939 else
1940 rational_u = rational_v = Standard_False;
1941 BSplSLib_DataContainer dc (UDegree, VDegree);
1942 flag_u_or_v =
1943 PrepareEval (U,
1944 V,
1945 UIndex,
1946 VIndex,
1947 UDegree,
1948 VDegree,
1949 rational_u,
1950 rational_v,
1951 UPeriodic,
1952 VPeriodic,
1953 Poles,
1954 Weights,
1955 UFlatKnots,
1956 VFlatKnots,
1957 (BSplCLib::NoMults()),
1958 (BSplCLib::NoMults()),
1959 u1,
1960 u2,
1961 d1,
1962 d2,
1963 rational,
1964 dc);
1965 d1p1 = d1 + 1;
1966 d2p1 = d2 + 1;
1967 if (rational) {
1968 BSplCLib::Bohm(u1,d1,d1,*dc.knots1,4 * d2p1,*dc.poles);
1969
1970 for (kk = 0; kk <= d1 ; kk++)
1971 BSplCLib::Bohm(u2,d2,d2,*dc.knots2,4,*(dc.poles + kk * 4 * d2p1));
1972 if (flag_u_or_v) {
1973 min_degree_domain = USpanDomain ;
1974 max_degree_domain = VSpanDomain ;
1975 }
1976 else {
1977 min_degree_domain = VSpanDomain ;
1978 max_degree_domain = USpanDomain ;
1979 }
1980 factor[0] = 1.0e0 ;
1981
1982 for (ii = 0 ; ii <= d2 ; ii++) {
1983 iii = ii + 1;
1984 factor[1] = 1.0e0 ;
1985
1986 for (jj = 0 ; jj <= d1 ; jj++) {
1987 jjj = jj + 1;
1988 Index = jj * d2p1 + ii ;
1989 Index = Index << 2;
1990 gp_Pnt& P = CachePoles(iii,jjj);
1991 f = factor[0] * factor[1];
1992 P.SetX( f * dc.poles[Index]); Index++;
1993 P.SetY( f * dc.poles[Index]); Index++;
1994 P.SetZ( f * dc.poles[Index]); Index++;
0e14656b 1995 (*CacheWeights)(iii ,jjj) = f * dc.poles[Index] ;
7fd59977 1996 factor[1] *= min_degree_domain / (Standard_Real) (jjj) ;
1997 }
1998 factor[0] *= max_degree_domain / (Standard_Real) (iii) ;
1999 }
2000 }
2001 else {
2002 BSplCLib::Bohm(u1,d1,d1,*dc.knots1,3 * d2p1,*dc.poles);
2003
2004 for (kk = 0; kk <= d1 ; kk++)
2005 BSplCLib::Bohm(u2,d2,d2,*dc.knots2,3,*(dc.poles + kk * 3 * d2p1));
2006 if (flag_u_or_v) {
2007 min_degree_domain = USpanDomain ;
2008 max_degree_domain = VSpanDomain ;
2009 }
2010 else {
2011 min_degree_domain = VSpanDomain ;
2012 max_degree_domain = USpanDomain ;
2013 }
2014 factor[0] = 1.0e0 ;
2015
2016 for (ii = 0 ; ii <= d2 ; ii++) {
2017 iii = ii + 1;
2018 factor[1] = 1.0e0 ;
2019
2020 for (jj = 0 ; jj <= d1 ; jj++) {
2021 jjj = jj + 1;
2022 Index = jj * d2p1 + ii ;
2023 Index = (Index << 1) + Index;
2024 gp_Pnt& P = CachePoles(iii,jjj);
2025 f = factor[0] * factor[1];
2026 P.SetX( f * dc.poles[Index]); Index++;
2027 P.SetY( f * dc.poles[Index]); Index++;
2028 P.SetZ( f * dc.poles[Index]);
2029 factor[1] *= min_degree_domain / (Standard_Real) (jjj) ;
2030 }
2031 factor[0] *= max_degree_domain / (Standard_Real) (iii) ;
2032 }
0e14656b 2033 if (Weights != NULL) {
7fd59977 2034 //
2035 // means that PrepareEval did found out that the surface was
2036 // locally polynomial but since the surface is constructed
2037 // with some weights we need to set the weight polynomial to constant
2038 //
2039
2040 for (ii = 1 ; ii <= d2p1 ; ii++) {
2041
2042 for (jj = 1 ; jj <= d1p1 ; jj++) {
0e14656b 2043 (*CacheWeights)(ii,jj) = 0.0e0 ;
7fd59977 2044 }
2045 }
0e14656b 2046 (*CacheWeights)(1,1) = 1.0e0 ;
7fd59977 2047 }
2048 }
2049}
2050
94f71cad 2051void BSplSLib::BuildCache(const Standard_Real theU,
2052 const Standard_Real theV,
2053 const Standard_Real theUSpanDomain,
2054 const Standard_Real theVSpanDomain,
2055 const Standard_Boolean theUPeriodicFlag,
2056 const Standard_Boolean theVPeriodicFlag,
2057 const Standard_Integer theUDegree,
2058 const Standard_Integer theVDegree,
2059 const Standard_Integer theUIndex,
2060 const Standard_Integer theVIndex,
2061 const TColStd_Array1OfReal& theUFlatKnots,
2062 const TColStd_Array1OfReal& theVFlatKnots,
2063 const TColgp_Array2OfPnt& thePoles,
0e14656b 2064 const TColStd_Array2OfReal* theWeights,
94f71cad 2065 TColStd_Array2OfReal& theCacheArray)
2066{
2067 Standard_Boolean flag_u_or_v;
2068 Standard_Integer d1, d2;
2069 Standard_Real u1, u2;
0e14656b 2070 Standard_Boolean isRationalOnParam = (theWeights != NULL);
94f71cad 2071 Standard_Boolean isRational;
2072
2073 BSplSLib_DataContainer dc(theUDegree, theVDegree);
2074 flag_u_or_v =
2075 PrepareEval(theU, theV, theUIndex, theVIndex, theUDegree, theVDegree,
2076 isRationalOnParam, isRationalOnParam,
2077 theUPeriodicFlag, theVPeriodicFlag,
2078 thePoles, theWeights,
2079 theUFlatKnots, theVFlatKnots,
2080 (BSplCLib::NoMults()), (BSplCLib::NoMults()),
2081 u1, u2, d1, d2, isRational, dc);
2082
2083 Standard_Integer d2p1 = d2 + 1;
2084 Standard_Integer aDimension = isRational ? 4 : 3;
2085 Standard_Integer aCacheShift = // helps to store weights when PrepareEval did not found that the surface is locally polynomial
2086 (isRationalOnParam && !isRational) ? aDimension + 1 : aDimension;
2087
2088 Standard_Real aDomains[2];
2089 // aDomains[0] corresponds to variable with minimal degree
2090 // aDomains[1] corresponds to variable with maximal degree
2091 if (flag_u_or_v)
2092 {
2093 aDomains[0] = theUSpanDomain;
2094 aDomains[1] = theVSpanDomain;
2095 }
2096 else
2097 {
2098 aDomains[0] = theVSpanDomain;
2099 aDomains[1] = theUSpanDomain;
2100 }
2101
2102 BSplCLib::Bohm(u1, d1, d1, *dc.knots1, aDimension * d2p1, *dc.poles);
2103 for (Standard_Integer kk = 0; kk <= d1 ; kk++)
2104 BSplCLib::Bohm(u2, d2, d2, *dc.knots2, aDimension, *(dc.poles + kk * aDimension * d2p1));
2105
2106 Standard_Real* aCache = (Standard_Real *) &(theCacheArray(theCacheArray.LowerRow(), theCacheArray.LowerCol()));
2107 Standard_Real* aPolyCoeffs = dc.poles;
2108
2109 Standard_Real aFactors[2];
2110 // aFactors[0] corresponds to variable with minimal degree
2111 // aFactors[1] corresponds to variable with maximal degree
2112 aFactors[1] = 1.0;
2113 Standard_Integer aRow, aCol, i;
2114 Standard_Real aCoeff;
2115 for (aRow = 0; aRow <= d2; aRow++)
2116 {
2117 aFactors[0] = 1.0;
2118 for (aCol = 0; aCol <= d1; aCol++)
2119 {
2120 aPolyCoeffs = dc.poles + (aCol * d2p1 + aRow) * aDimension;
2121 aCoeff = aFactors[0] * aFactors[1];
2122 for (i = 0; i < aDimension; i++)
2123 aCache[i] = aPolyCoeffs[i] * aCoeff;
2124 aCache += aCacheShift;
2125 aFactors[0] *= aDomains[0] / (aCol + 1);
2126 }
2127 aFactors[1] *= aDomains[1] / (aRow + 1);
2128 }
2129
2130 // Fill the weights for the surface which is not locally polynomial
2131 if (aCacheShift > aDimension)
2132 {
2133 aCache = (Standard_Real *) &(theCacheArray(theCacheArray.LowerRow(), theCacheArray.LowerCol()));
2134 aCache += aCacheShift - 1;
2135 for (aRow = 0; aRow <= d2; aRow++)
2136 for (aCol = 0; aCol <= d1; aCol++)
2137 {
2138 *aCache = 0.0;
2139 aCache += aCacheShift;
2140 }
2141 theCacheArray.SetValue(theCacheArray.LowerRow(), theCacheArray.LowerCol() + aCacheShift - 1, 1.0);
2142 }
2143}
2144
2145
7fd59977 2146//=======================================================================
2147//function : CacheD0
2148//purpose : Evaluates the polynomial cache of the Bspline Curve
2149//
2150//=======================================================================
2151
2152void BSplSLib::CacheD0(const Standard_Real UParameter,
2153 const Standard_Real VParameter,
2154 const Standard_Integer UDegree,
2155 const Standard_Integer VDegree,
2156 const Standard_Real UCacheParameter,
2157 const Standard_Real VCacheParameter,
2158 const Standard_Real USpanLenght,
2159 const Standard_Real VSpanLenght,
2160 const TColgp_Array2OfPnt& PolesArray,
0e14656b 2161 const TColStd_Array2OfReal* WeightsArray,
7fd59977 2162 gp_Pnt& aPoint)
2163{
2164 //
2165 // the CacheParameter is where the cache polynomial was evaluated in homogeneous
2166 // form
2167 // the SpanLenght is the normalizing factor so that the polynomial is between
2168 // 0 and 1
2169 Standard_Integer
2170// ii,
2171 dimension,
2172 min_degree,
2173 max_degree ;
2174
2175 Standard_Real
2176 new_parameter[2] ,
2177 inverse ;
2178
2179 Standard_Real *
2180 PArray = (Standard_Real *)
2181 &(PolesArray(PolesArray.LowerCol(), PolesArray.LowerRow())) ;
2182 Standard_Real *
2183 myPoint = (Standard_Real *) &aPoint ;
2184 if (UDegree <= VDegree) {
2185 min_degree = UDegree ;
2186 max_degree = VDegree ;
2187 new_parameter[1] = (UParameter - UCacheParameter) / USpanLenght ;
2188 new_parameter[0] = (VParameter - VCacheParameter) / VSpanLenght ;
2189 dimension = 3 * (UDegree + 1) ;
2190 }
2191 else {
2192 min_degree = VDegree ;
2193 max_degree = UDegree ;
2194 new_parameter[0] = (UParameter - UCacheParameter) / USpanLenght ;
2195 new_parameter[1] = (VParameter - VCacheParameter) / VSpanLenght ;
2196 dimension = 3 * (VDegree + 1) ;
2197 }
f7b4312f 2198 NCollection_LocalArray<Standard_Real> locpoles(dimension);
7fd59977 2199
2200 PLib::NoDerivativeEvalPolynomial(new_parameter[0],
2201 max_degree,
2202 dimension,
2203 max_degree*dimension,
2204 PArray[0],
2205 locpoles[0]) ;
2206
2207 PLib::NoDerivativeEvalPolynomial(new_parameter[1],
2208 min_degree,
2209 3,
2210 (min_degree << 1) + min_degree,
2211 locpoles[0],
2212 myPoint[0]) ;
0e14656b 2213 if (WeightsArray != NULL) {
7fd59977 2214 dimension = min_degree + 1 ;
0e14656b 2215 const TColStd_Array2OfReal& refWeights = *WeightsArray;
7fd59977 2216 Standard_Real *
2217 WArray = (Standard_Real *)
0e14656b 2218 &refWeights(refWeights.LowerCol(),refWeights.LowerRow()) ;
7fd59977 2219 PLib::NoDerivativeEvalPolynomial(new_parameter[0],
2220 max_degree,
2221 dimension,
2222 max_degree*dimension,
2223 WArray[0],
2224 locpoles[0]) ;
2225
2226 PLib::NoDerivativeEvalPolynomial(new_parameter[1],
2227 min_degree,
2228 1,
2229 min_degree,
2230 locpoles[0],
2231 inverse) ;
2232 inverse = 1.0e0/ inverse ;
2233
2234 myPoint[0] *= inverse ;
2235 myPoint[1] *= inverse ;
2236 myPoint[2] *= inverse ;
2237 }
2238}
2239
2240//=======================================================================
2241//function : CacheD1
2242//purpose : Evaluates the polynomial cache of the Bspline Curve
2243//
2244//=======================================================================
2245
2246void BSplSLib::CacheD1(const Standard_Real UParameter,
2247 const Standard_Real VParameter,
2248 const Standard_Integer UDegree,
2249 const Standard_Integer VDegree,
2250 const Standard_Real UCacheParameter,
2251 const Standard_Real VCacheParameter,
2252 const Standard_Real USpanLenght,
2253 const Standard_Real VSpanLenght,
2254 const TColgp_Array2OfPnt& PolesArray,
0e14656b 2255 const TColStd_Array2OfReal* WeightsArray,
7fd59977 2256 gp_Pnt& aPoint,
2257 gp_Vec& aVecU,
2258 gp_Vec& aVecV)
2259{
2260 //
2261 // the CacheParameter is where the cache polynomial was evaluated in homogeneous
2262 // form
2263 // the SpanLenght is the normalizing factor so that the polynomial is between
2264 // 0 and 1
2265 Standard_Integer
2266// ii,
2267// jj,
2268// kk,
2269 dimension,
2270 min_degree,
2271 max_degree ;
2272
2273 Standard_Real
2274 inverse_min,
2275 inverse_max,
2276 new_parameter[2] ;
2277
2278 Standard_Real *
2279 PArray = (Standard_Real *)
2280 &(PolesArray(PolesArray.LowerCol(), PolesArray.LowerRow())) ;
2281 Standard_Real local_poles_array[2][2][3],
2282 local_poles_and_weights_array[2][2][4],
2283 local_weights_array[2][2] ;
2284 Standard_Real * my_vec_min,
2285 * my_vec_max,
2286 * my_point ;
2287 my_point = (Standard_Real *) &aPoint ;
2288 //
2289 // initialize in case of rational evaluation
2290 // because RationalDerivative will use all
2291 // the coefficients
2292 //
2293 //
0e14656b 2294 if (WeightsArray != NULL) {
7fd59977 2295
2296 local_poles_array [0][0][0] = 0.0e0 ;
2297 local_poles_array [0][0][1] = 0.0e0 ;
2298 local_poles_array [0][0][2] = 0.0e0 ;
2299 local_weights_array [0][0] = 0.0e0 ;
2300 local_poles_and_weights_array[0][0][0] = 0.0e0 ;
2301 local_poles_and_weights_array[0][0][1] = 0.0e0 ;
2302 local_poles_and_weights_array[0][0][2] = 0.0e0 ;
2303 local_poles_and_weights_array[0][0][3] = 0.0e0 ;
2304
2305 local_poles_array [0][1][0] = 0.0e0 ;
2306 local_poles_array [0][1][1] = 0.0e0 ;
2307 local_poles_array [0][1][2] = 0.0e0 ;
2308 local_weights_array [0][1] = 0.0e0 ;
2309 local_poles_and_weights_array[0][1][0] = 0.0e0 ;
2310 local_poles_and_weights_array[0][1][1] = 0.0e0 ;
2311 local_poles_and_weights_array[0][1][2] = 0.0e0 ;
2312 local_poles_and_weights_array[0][1][3] = 0.0e0 ;
2313
2314 local_poles_array [1][0][0] = 0.0e0 ;
2315 local_poles_array [1][0][1] = 0.0e0 ;
2316 local_poles_array [1][0][2] = 0.0e0 ;
2317 local_weights_array [1][0] = 0.0e0 ;
2318 local_poles_and_weights_array[1][0][0] = 0.0e0 ;
2319 local_poles_and_weights_array[1][0][1] = 0.0e0 ;
2320 local_poles_and_weights_array[1][0][2] = 0.0e0 ;
2321 local_poles_and_weights_array[1][0][3] = 0.0e0 ;
2322
2323 local_poles_array [1][1][0] = 0.0e0 ;
2324 local_poles_array [1][1][1] = 0.0e0 ;
2325 local_poles_array [1][1][2] = 0.0e0 ;
2326 local_weights_array [1][1] = 0.0e0 ;
2327 local_poles_and_weights_array[1][1][0] = 0.0e0 ;
2328 local_poles_and_weights_array[1][1][1] = 0.0e0 ;
2329 local_poles_and_weights_array[1][1][2] = 0.0e0 ;
2330 local_poles_and_weights_array[1][1][3] = 0.0e0 ;
2331 }
2332
2333 if (UDegree <= VDegree) {
2334 min_degree = UDegree ;
2335 max_degree = VDegree ;
2336 inverse_min = 1.0e0/USpanLenght ;
2337 inverse_max = 1.0e0/VSpanLenght ;
2338 new_parameter[0] = (VParameter - VCacheParameter) * inverse_max ;
2339 new_parameter[1] = (UParameter - UCacheParameter) * inverse_min ;
2340
2341 dimension = 3 * (UDegree + 1) ;
2342 my_vec_min = (Standard_Real *) &aVecU ;
2343 my_vec_max = (Standard_Real *) &aVecV ;
2344 }
2345 else {
2346 min_degree = VDegree ;
2347 max_degree = UDegree ;
2348 inverse_min = 1.0e0/VSpanLenght ;
2349 inverse_max = 1.0e0/USpanLenght ;
2350 new_parameter[0] = (UParameter - UCacheParameter) * inverse_max ;
2351 new_parameter[1] = (VParameter - VCacheParameter) * inverse_min ;
2352 dimension = 3 * (VDegree + 1) ;
2353 my_vec_min = (Standard_Real *) &aVecV ;
2354 my_vec_max = (Standard_Real *) &aVecU ;
2355 }
2356
f7b4312f 2357 NCollection_LocalArray<Standard_Real> locpoles (2 * dimension);
7fd59977 2358
2359 PLib::EvalPolynomial(new_parameter[0],
2360 1,
2361 max_degree,
2362 dimension,
2363 PArray[0],
2364 locpoles[0]) ;
2365
2366 PLib::EvalPolynomial(new_parameter[1],
2367 1,
2368 min_degree,
2369 3,
2370 locpoles[0],
2371 local_poles_array[0][0][0]) ;
2372 PLib::NoDerivativeEvalPolynomial(new_parameter[1],
2373 min_degree,
2374 3,
2375 (min_degree << 1) + min_degree,
2376 locpoles[dimension],
2377 local_poles_array[1][0][0]) ;
2378
0e14656b 2379 if (WeightsArray != NULL) {
7fd59977 2380 dimension = min_degree + 1 ;
0e14656b 2381 const TColStd_Array2OfReal& refWeights = *WeightsArray;
7fd59977 2382 Standard_Real *
2383 WArray = (Standard_Real *)
0e14656b 2384 &refWeights(refWeights.LowerCol(),refWeights.LowerRow()) ;
7fd59977 2385 PLib::EvalPolynomial(new_parameter[0],
2386 1,
2387 max_degree,
2388 dimension,
2389 WArray[0],
2390 locpoles[0]) ;
2391
2392 PLib::EvalPolynomial(new_parameter[1],
2393 1,
2394 min_degree,
2395 1,
2396 locpoles[0],
2397 local_weights_array[0][0]) ;
2398 PLib::NoDerivativeEvalPolynomial(new_parameter[1],
2399 min_degree,
2400 1,
2401 min_degree,
2402 locpoles[dimension],
2403 local_weights_array[1][0]) ;
2404
2405 local_poles_and_weights_array[0][0][0] = local_poles_array [0][0][0] ;
2406 local_poles_and_weights_array[0][0][1] = local_poles_array [0][0][1] ;
2407 local_poles_and_weights_array[0][0][2] = local_poles_array [0][0][2] ;
2408 local_poles_and_weights_array[0][0][3] = local_weights_array[0][0] ;
2409
2410 local_poles_and_weights_array[0][1][0] = local_poles_array [0][1][0] ;
2411 local_poles_and_weights_array[0][1][1] = local_poles_array [0][1][1] ;
2412 local_poles_and_weights_array[0][1][2] = local_poles_array [0][1][2] ;
2413 local_poles_and_weights_array[0][1][3] = local_weights_array[0][1] ;
2414
2415 local_poles_and_weights_array[1][0][0] = local_poles_array [1][0][0] ;
2416 local_poles_and_weights_array[1][0][1] = local_poles_array [1][0][1] ;
2417 local_poles_and_weights_array[1][0][2] = local_poles_array [1][0][2] ;
2418 local_poles_and_weights_array[1][0][3] = local_weights_array[1][0] ;
2419
2420 local_poles_and_weights_array[1][1][0] = local_poles_array [1][1][0] ;
2421 local_poles_and_weights_array[1][1][1] = local_poles_array [1][1][1] ;
2422 local_poles_and_weights_array[1][1][2] = local_poles_array [1][1][2] ;
2423 local_poles_and_weights_array[1][1][3] = local_weights_array[1][1] ;
2424
2425 BSplSLib::RationalDerivative(1,
2426 1,
2427 1,
2428 1,
2429 local_poles_and_weights_array[0][0][0],
2430 local_poles_array[0][0][0]) ;
2431 }
2432
2433 my_point [0] = local_poles_array [0][0][0] ;
2434 my_vec_min[0] = inverse_min * local_poles_array[0][1][0] ;
2435 my_vec_max[0] = inverse_max * local_poles_array[1][0][0] ;
2436
2437 my_point [1] = local_poles_array [0][0][1] ;
2438 my_vec_min[1] = inverse_min * local_poles_array[0][1][1] ;
2439 my_vec_max[1] = inverse_max * local_poles_array[1][0][1] ;
2440
2441 my_point [2] = local_poles_array [0][0][2] ;
2442 my_vec_min[2] = inverse_min * local_poles_array[0][1][2] ;
2443 my_vec_max[2] = inverse_max * local_poles_array[1][0][2] ;
2444}
2445
2446//=======================================================================
2447//function : CacheD2
2448//purpose : Evaluates the polynomial cache of the Bspline Curve
2449//
2450//=======================================================================
2451
2452void BSplSLib::CacheD2(const Standard_Real UParameter,
2453 const Standard_Real VParameter,
2454 const Standard_Integer UDegree,
2455 const Standard_Integer VDegree,
2456 const Standard_Real UCacheParameter,
2457 const Standard_Real VCacheParameter,
2458 const Standard_Real USpanLenght,
2459 const Standard_Real VSpanLenght,
2460 const TColgp_Array2OfPnt& PolesArray,
0e14656b 2461 const TColStd_Array2OfReal* WeightsArray,
7fd59977 2462 gp_Pnt& aPoint,
2463 gp_Vec& aVecU,
2464 gp_Vec& aVecV,
2465 gp_Vec& aVecUU,
2466 gp_Vec& aVecUV,
2467 gp_Vec& aVecVV)
2468{
2469 //
2470 // the CacheParameter is where the cache polynomial was evaluated in homogeneous
2471 // form
2472 // the SpanLenght is the normalizing factor so that the polynomial is between
2473 // 0 and 1
2474 Standard_Integer
2475 ii,
2476// jj,
2477 kk,
2478 index,
2479 dimension,
2480 min_degree,
2481 max_degree ;
2482
2483 Standard_Real
2484 inverse_min,
2485 inverse_max,
2486 new_parameter[2] ;
2487
2488 Standard_Real *
2489 PArray = (Standard_Real *)
2490 &(PolesArray(PolesArray.LowerCol(), PolesArray.LowerRow())) ;
2491 Standard_Real local_poles_array[3][3][3],
2492 local_poles_and_weights_array[3][3][4],
2493 local_weights_array[3][3] ;
2494 Standard_Real * my_vec_min,
2495 * my_vec_max,
2496 * my_vec_min_min,
2497 * my_vec_max_max,
2498 * my_vec_min_max,
2499 * my_point ;
2500 my_point = (Standard_Real *) &aPoint ;
2501
2502 //
2503 // initialize in case the min and max degree are less than 2
2504 //
2505 local_poles_array[0][0][0] = 0.0e0 ;
2506 local_poles_array[0][0][1] = 0.0e0 ;
2507 local_poles_array[0][0][2] = 0.0e0 ;
2508 local_poles_array[0][1][0] = 0.0e0 ;
2509 local_poles_array[0][1][1] = 0.0e0 ;
2510 local_poles_array[0][1][2] = 0.0e0 ;
2511 local_poles_array[0][2][0] = 0.0e0 ;
2512 local_poles_array[0][2][1] = 0.0e0 ;
2513 local_poles_array[0][2][2] = 0.0e0 ;
2514
2515 local_poles_array[1][0][0] = 0.0e0 ;
2516 local_poles_array[1][0][1] = 0.0e0 ;
2517 local_poles_array[1][0][2] = 0.0e0 ;
2518 local_poles_array[1][1][0] = 0.0e0 ;
2519 local_poles_array[1][1][1] = 0.0e0 ;
2520 local_poles_array[1][1][2] = 0.0e0 ;
2521 local_poles_array[1][2][0] = 0.0e0 ;
2522 local_poles_array[1][2][1] = 0.0e0 ;
2523 local_poles_array[1][2][2] = 0.0e0 ;
2524
2525 local_poles_array[2][0][0] = 0.0e0 ;
2526 local_poles_array[2][0][1] = 0.0e0 ;
2527 local_poles_array[2][0][2] = 0.0e0 ;
2528 local_poles_array[2][1][0] = 0.0e0 ;
2529 local_poles_array[2][1][1] = 0.0e0 ;
2530 local_poles_array[2][1][2] = 0.0e0 ;
2531 local_poles_array[2][2][0] = 0.0e0 ;
2532 local_poles_array[2][2][1] = 0.0e0 ;
2533 local_poles_array[2][2][2] = 0.0e0 ;
2534 //
2535 // initialize in case of rational evaluation
2536 // because RationalDerivative will use all
2537 // the coefficients
2538 //
2539 //
0e14656b 2540 if (WeightsArray != NULL) {
7fd59977 2541
2542 local_poles_and_weights_array[0][0][0] = 0.0e0 ;
2543 local_poles_and_weights_array[0][0][1] = 0.0e0 ;
2544 local_poles_and_weights_array[0][0][2] = 0.0e0 ;
2545 local_poles_and_weights_array[0][1][0] = 0.0e0 ;
2546 local_poles_and_weights_array[0][1][1] = 0.0e0 ;
2547 local_poles_and_weights_array[0][1][2] = 0.0e0 ;
2548 local_poles_and_weights_array[0][2][0] = 0.0e0 ;
2549 local_poles_and_weights_array[0][2][1] = 0.0e0 ;
2550 local_poles_and_weights_array[0][2][2] = 0.0e0 ;
2551
2552 local_poles_and_weights_array[1][0][0] = 0.0e0 ;
2553 local_poles_and_weights_array[1][0][1] = 0.0e0 ;
2554 local_poles_and_weights_array[1][0][2] = 0.0e0 ;
2555 local_poles_and_weights_array[1][1][0] = 0.0e0 ;
2556 local_poles_and_weights_array[1][1][1] = 0.0e0 ;
2557 local_poles_and_weights_array[1][1][2] = 0.0e0 ;
2558 local_poles_and_weights_array[1][2][0] = 0.0e0 ;
2559 local_poles_and_weights_array[1][2][1] = 0.0e0 ;
2560 local_poles_and_weights_array[1][2][2] = 0.0e0 ;
2561
2562 local_poles_and_weights_array[2][0][0] = 0.0e0 ;
2563 local_poles_and_weights_array[2][0][1] = 0.0e0 ;
2564 local_poles_and_weights_array[2][0][2] = 0.0e0 ;
2565 local_poles_and_weights_array[2][1][0] = 0.0e0 ;
2566 local_poles_and_weights_array[2][1][1] = 0.0e0 ;
2567 local_poles_and_weights_array[2][1][2] = 0.0e0 ;
2568 local_poles_and_weights_array[2][2][0] = 0.0e0 ;
2569 local_poles_and_weights_array[2][2][1] = 0.0e0 ;
2570 local_poles_and_weights_array[2][2][2] = 0.0e0 ;
2571
2572 local_poles_and_weights_array[0][0][3] =
2573 local_weights_array[0][0] = 0.0e0 ;
2574 local_poles_and_weights_array[0][1][3] =
2575 local_weights_array[0][1] = 0.0e0 ;
2576 local_poles_and_weights_array[0][2][3] =
2577 local_weights_array[0][2] = 0.0e0 ;
2578 local_poles_and_weights_array[1][0][3] =
2579 local_weights_array[1][0] = 0.0e0 ;
2580 local_poles_and_weights_array[1][1][3] =
2581 local_weights_array[1][1] = 0.0e0 ;
2582 local_poles_and_weights_array[1][2][3] =
2583 local_weights_array[1][2] = 0.0e0 ;
2584 local_poles_and_weights_array[2][0][3] =
2585 local_weights_array[2][0] = 0.0e0 ;
2586 local_poles_and_weights_array[2][1][3] =
2587 local_weights_array[2][1] = 0.0e0 ;
2588 local_poles_and_weights_array[2][2][3] =
2589 local_weights_array[2][2] = 0.0e0 ;
2590 }
2591
2592 if (UDegree <= VDegree) {
2593 min_degree = UDegree ;
2594 max_degree = VDegree ;
2595 inverse_min = 1.0e0/USpanLenght ;
2596 inverse_max = 1.0e0/VSpanLenght ;
2597 new_parameter[0] = (VParameter - VCacheParameter) * inverse_max ;
2598 new_parameter[1] = (UParameter - UCacheParameter) * inverse_min ;
2599
2600 dimension = 3 * (UDegree + 1) ;
2601 my_vec_min = (Standard_Real *) &aVecU ;
2602 my_vec_max = (Standard_Real *) &aVecV ;
2603 my_vec_min_min = (Standard_Real *) &aVecUU ;
2604 my_vec_min_max = (Standard_Real *) &aVecUV ;
2605 my_vec_max_max = (Standard_Real *) &aVecVV ;
2606 }
2607 else {
2608 min_degree = VDegree ;
2609 max_degree = UDegree ;
2610 inverse_min = 1.0e0/VSpanLenght ;
2611 inverse_max = 1.0e0/USpanLenght ;
2612 new_parameter[0] = (UParameter - UCacheParameter) * inverse_max ;
2613 new_parameter[1] = (VParameter - VCacheParameter) * inverse_min ;
2614 dimension = 3 * (VDegree + 1) ;
2615 my_vec_min = (Standard_Real *) &aVecV ;
2616 my_vec_max = (Standard_Real *) &aVecU ;
2617 my_vec_min_min = (Standard_Real *) &aVecVV ;
2618 my_vec_min_max = (Standard_Real *) &aVecUV ;
2619 my_vec_max_max = (Standard_Real *) &aVecUU ;
2620 }
2621
f7b4312f 2622 NCollection_LocalArray<Standard_Real> locpoles (3 * dimension);
7fd59977 2623
2624 //
2625 // initialize in case min or max degree are less than 2
2626 //
2627 Standard_Integer MinIndMax = 2;
2628 if ( max_degree < 2) MinIndMax = max_degree;
2629 Standard_Integer MinIndMin = 2;
2630 if ( min_degree < 2) MinIndMin = min_degree;
2631
2632 index = MinIndMax * dimension ;
2633
2634 for (ii = MinIndMax ; ii < 3 ; ii++) {
2635
2636 for (kk = 0 ; kk < dimension ; kk++) {
2637 locpoles[index] = 0.0e0 ;
2638 index += 1 ;
2639 }
2640 }
2641
2642 PLib::EvalPolynomial(new_parameter[0],
2643 MinIndMax,
2644 max_degree,
2645 dimension,
2646 PArray[0],
2647 locpoles[0]) ;
2648
2649 PLib::EvalPolynomial(new_parameter[1],
2650 MinIndMin,
2651 min_degree,
2652 3,
2653 locpoles[0],
2654 local_poles_array[0][0][0]) ;
2655 PLib::EvalPolynomial(new_parameter[1],
2656 1,
2657 min_degree,
2658 3,
2659 locpoles[dimension],
2660 local_poles_array[1][0][0]) ;
2661
2662 PLib::NoDerivativeEvalPolynomial(new_parameter[1],
2663 min_degree,
2664 3,
2665 (min_degree << 1) + min_degree,
2666 locpoles[dimension + dimension],
2667 local_poles_array[2][0][0]) ;
2668
0e14656b 2669 if (WeightsArray != NULL) {
7fd59977 2670 dimension = min_degree + 1 ;
0e14656b 2671 const TColStd_Array2OfReal& refWeights = *WeightsArray;
7fd59977 2672 Standard_Real *
2673 WArray = (Standard_Real *)
0e14656b 2674 &refWeights(refWeights.LowerCol(),refWeights.LowerRow()) ;
7fd59977 2675 PLib::EvalPolynomial(new_parameter[0],
2676 MinIndMax,
2677 max_degree,
2678 dimension,
2679 WArray[0],
2680 locpoles[0]) ;
2681
2682 PLib::EvalPolynomial(new_parameter[1],
2683 MinIndMin,
2684 min_degree,
2685 1,
2686 locpoles[0],
2687 local_weights_array[0][0]) ;
2688 PLib::EvalPolynomial(new_parameter[1],
2689 1,
2690 min_degree,
2691 1,
2692 locpoles[dimension],
2693 local_weights_array[1][0]) ;
2694 PLib::NoDerivativeEvalPolynomial(new_parameter[1],
2695 min_degree,
2696 1,
2697 min_degree,
2698 locpoles[dimension + dimension],
2699 local_weights_array[2][0]) ;
2700
2701
2702 local_poles_and_weights_array[0][0][0] = local_poles_array[0][0][0];
2703 local_poles_and_weights_array[0][0][1] = local_poles_array[0][0][1];
2704 local_poles_and_weights_array[0][0][2] = local_poles_array[0][0][2];
2705 local_poles_and_weights_array[0][1][0] = local_poles_array[0][1][0];
2706 local_poles_and_weights_array[0][1][1] = local_poles_array[0][1][1];
2707 local_poles_and_weights_array[0][1][2] = local_poles_array[0][1][2];
2708 local_poles_and_weights_array[0][2][0] = local_poles_array[0][2][0];
2709 local_poles_and_weights_array[0][2][1] = local_poles_array[0][2][1];
2710 local_poles_and_weights_array[0][2][2] = local_poles_array[0][2][2];
2711
2712 local_poles_and_weights_array[1][0][0] = local_poles_array[1][0][0];
2713 local_poles_and_weights_array[1][0][1] = local_poles_array[1][0][1];
2714 local_poles_and_weights_array[1][0][2] = local_poles_array[1][0][2];
2715 local_poles_and_weights_array[1][1][0] = local_poles_array[1][1][0];
2716 local_poles_and_weights_array[1][1][1] = local_poles_array[1][1][1];
2717 local_poles_and_weights_array[1][1][2] = local_poles_array[1][1][2];
2718 local_poles_and_weights_array[1][2][0] = local_poles_array[1][2][0];
2719 local_poles_and_weights_array[1][2][1] = local_poles_array[1][2][1];
2720 local_poles_and_weights_array[1][2][2] = local_poles_array[1][2][2];
2721
2722 local_poles_and_weights_array[2][0][0] = local_poles_array[2][0][0];
2723 local_poles_and_weights_array[2][0][1] = local_poles_array[2][0][1];
2724 local_poles_and_weights_array[2][0][2] = local_poles_array[2][0][2];
2725 local_poles_and_weights_array[2][1][0] = local_poles_array[2][1][0];
2726 local_poles_and_weights_array[2][1][1] = local_poles_array[2][1][1];
2727 local_poles_and_weights_array[2][1][2] = local_poles_array[2][1][2];
2728 local_poles_and_weights_array[2][2][0] = local_poles_array[2][2][0];
2729 local_poles_and_weights_array[2][2][1] = local_poles_array[2][2][1];
2730 local_poles_and_weights_array[2][2][2] = local_poles_array[2][2][2];
2731
2732
2733 local_poles_and_weights_array[0][0][3] = local_weights_array[0][0];
2734 local_poles_and_weights_array[0][1][3] = local_weights_array[0][1];
2735 local_poles_and_weights_array[0][2][3] = local_weights_array[0][2];
2736 local_poles_and_weights_array[1][0][3] = local_weights_array[1][0];
2737 local_poles_and_weights_array[1][1][3] = local_weights_array[1][1];
2738 local_poles_and_weights_array[1][2][3] = local_weights_array[1][2];
2739 local_poles_and_weights_array[2][0][3] = local_weights_array[2][0];
2740 local_poles_and_weights_array[2][1][3] = local_weights_array[2][1];
2741 local_poles_and_weights_array[2][2][3] = local_weights_array[2][2];
2742
2743 BSplSLib::RationalDerivative(2,
2744 2,
2745 2,
2746 2,
2747 local_poles_and_weights_array[0][0][0],
2748 local_poles_array[0][0][0]) ;
2749 }
2750
2751
2752 Standard_Real minmin = inverse_min * inverse_min;
2753 Standard_Real minmax = inverse_min * inverse_max;
2754 Standard_Real maxmax = inverse_max * inverse_max;
2755
2756 my_point [0] = local_poles_array [0][0][0] ;
2757 my_vec_min [0] = inverse_min * local_poles_array[0][1][0] ;
2758 my_vec_max [0] = inverse_max * local_poles_array[1][0][0] ;
2759 my_vec_min_min[0] = minmin * local_poles_array [0][2][0] ;
2760 my_vec_min_max[0] = minmax * local_poles_array [1][1][0] ;
2761 my_vec_max_max[0] = maxmax * local_poles_array [2][0][0] ;
2762
2763 my_point [1] = local_poles_array [0][0][1] ;
2764 my_vec_min [1] = inverse_min * local_poles_array[0][1][1] ;
2765 my_vec_max [1] = inverse_max * local_poles_array[1][0][1] ;
2766 my_vec_min_min[1] = minmin * local_poles_array [0][2][1] ;
2767 my_vec_min_max[1] = minmax * local_poles_array [1][1][1] ;
2768 my_vec_max_max[1] = maxmax * local_poles_array [2][0][1] ;
2769
2770 my_point [2] = local_poles_array [0][0][2] ;
2771 my_vec_min [2] = inverse_min * local_poles_array[0][1][2] ;
2772 my_vec_max [2] = inverse_max * local_poles_array[1][0][2] ;
2773 my_vec_min_min[2] = minmin * local_poles_array [0][2][2] ;
2774 my_vec_min_max[2] = minmax * local_poles_array [1][1][2] ;
2775 my_vec_max_max[2] = maxmax * local_poles_array [2][0][2] ;
2776}
2777
2778//=======================================================================
2779//function : MovePoint
2780//purpose : Find the new poles which allows an old point (with a
2781// given u and v as parameters) to reach a new position
2782//=======================================================================
2783
2784void BSplSLib::MovePoint (const Standard_Real U,
2785 const Standard_Real V,
2786 const gp_Vec& Displ,
2787 const Standard_Integer UIndex1,
2788 const Standard_Integer UIndex2,
2789 const Standard_Integer VIndex1,
2790 const Standard_Integer VIndex2,
2791 const Standard_Integer UDegree,
2792 const Standard_Integer VDegree,
2793 const Standard_Boolean Rational,
2794 const TColgp_Array2OfPnt& Poles,
2795 const TColStd_Array2OfReal& Weights,
2796 const TColStd_Array1OfReal& UFlatKnots,
2797 const TColStd_Array1OfReal& VFlatKnots,
2798 Standard_Integer& UFirstIndex,
2799 Standard_Integer& ULastIndex,
2800 Standard_Integer& VFirstIndex,
2801 Standard_Integer& VLastIndex,
2802 TColgp_Array2OfPnt& NewPoles)
2803{
2804 // calculate the UBSplineBasis in the parameter U
2805 Standard_Integer UFirstNonZeroBsplineIndex;
2806 math_Matrix UBSplineBasis(1, 1,
2807 1, UDegree+1);
6143f12f 2808 Standard_Integer ErrorCod1 = BSplCLib::EvalBsplineBasis(0,
7fd59977 2809 UDegree+1,
2810 UFlatKnots,
2811 U,
2812 UFirstNonZeroBsplineIndex,
2813 UBSplineBasis);
2814 // calculate the VBSplineBasis in the parameter V
2815 Standard_Integer VFirstNonZeroBsplineIndex;
2816 math_Matrix VBSplineBasis(1, 1,
2817 1, VDegree+1);
6143f12f 2818 Standard_Integer ErrorCod2 = BSplCLib::EvalBsplineBasis(0,
7fd59977 2819 VDegree+1,
2820 VFlatKnots,
2821 V,
2822 VFirstNonZeroBsplineIndex,
2823 VBSplineBasis);
2824 if (ErrorCod1 || ErrorCod2) {
2825 UFirstIndex = 0;
2826 ULastIndex = 0;
2827 VFirstIndex = 0;
2828 VLastIndex = 0;
2829 return;
2830 }
2831
2832 // find the span which is predominant for parameter U
2833 UFirstIndex = UFirstNonZeroBsplineIndex;
2834 ULastIndex = UFirstNonZeroBsplineIndex + UDegree ;
2835 if (UFirstIndex < UIndex1) UFirstIndex = UIndex1;
2836 if (ULastIndex > UIndex2) ULastIndex = UIndex2;
2837
2838 Standard_Real maxValue = 0.0;
2839 Standard_Integer i, ukk1=0, ukk2;
2840
2841 for (i = UFirstIndex-UFirstNonZeroBsplineIndex+1; i <= ULastIndex-UFirstNonZeroBsplineIndex+1; i++) {
2842 if (UBSplineBasis(1,i) > maxValue) {
2843 ukk1 = i + UFirstNonZeroBsplineIndex - 1;
2844 maxValue = UBSplineBasis(1,i);
2845 }
2846 }
2847
2848 // find a ukk2 if symetriy
2849 ukk2 = ukk1;
2850 i = ukk1 - UFirstNonZeroBsplineIndex + 2;
2851 if ((ukk1+1) <= ULastIndex) {
2852 if (Abs(UBSplineBasis(1, ukk1-UFirstNonZeroBsplineIndex+2) - maxValue) < 1.e-10) {
2853 ukk2 = ukk1+1;
2854 }
2855 }
2856
2857 // find the span which is predominant for parameter V
2858 VFirstIndex = VFirstNonZeroBsplineIndex;
2859 VLastIndex = VFirstNonZeroBsplineIndex + VDegree ;
2860
2861 if (VFirstIndex < VIndex1) VFirstIndex = VIndex1;
2862 if (VLastIndex > VIndex2) VLastIndex = VIndex2;
2863
2864 maxValue = 0.0;
2865 Standard_Integer j, vkk1=0, vkk2;
2866
2867 for (j = VFirstIndex-VFirstNonZeroBsplineIndex+1; j <= VLastIndex-VFirstNonZeroBsplineIndex+1; j++) {
2868 if (VBSplineBasis(1,j) > maxValue) {
2869 vkk1 = j + VFirstNonZeroBsplineIndex - 1;
2870 maxValue = VBSplineBasis(1,j);
2871 }
2872 }
2873
2874 // find a vkk2 if symetriy
2875 vkk2 = vkk1;
2876 j = vkk1 - VFirstNonZeroBsplineIndex + 2;
2877 if ((vkk1+1) <= VLastIndex) {
2878 if (Abs(VBSplineBasis(1, vkk1-VFirstNonZeroBsplineIndex+2) - maxValue) < 1.e-10) {
2879 vkk2 = vkk1+1;
2880 }
2881 }
2882
2883 // compute the vector of displacement
2884 Standard_Real D1 = 0.0;
2885 Standard_Real D2 = 0.0;
2886 Standard_Real hN, Coef, DvalU, DvalV;
2887
2888 Standard_Integer ii, jj;
2889
2890 for (i = 1; i <= UDegree+1; i++) {
2891 ii = i + UFirstNonZeroBsplineIndex - 1;
2892 if (ii < ukk1) {
2893 DvalU = ukk1-ii;
2894 }
2895 else if (ii > ukk2) {
2896 DvalU = ii - ukk2;
2897 }
2898 else {
2899 DvalU = 0.0;
2900 }
2901
2902 for (j = 1; j <= VDegree+1; j++) {
2903 jj = j + VFirstNonZeroBsplineIndex - 1;
2904 if (Rational) {
2905 hN = Weights(ii, jj)*UBSplineBasis(1, i)*VBSplineBasis(1,j);
2906 D2 += hN;
2907 }
2908 else {
2909 hN = UBSplineBasis(1, i)*VBSplineBasis(1,j);
2910 }
2911 if (ii >= UFirstIndex && ii <= ULastIndex && jj >= VFirstIndex && jj <= VLastIndex) {
2912 if (jj < vkk1) {
2913 DvalV = vkk1-jj;
2914 }
2915 else if (jj > vkk2) {
2916 DvalV = jj - vkk2;
2917 }
2918 else {
2919 DvalV = 0.0;
2920 }
2921 D1 += 1./(DvalU + DvalV + 1.) * hN;
2922 }
2923 }
2924 }
2925
2926 if (Rational) {
2927 Coef = D2/D1;
2928 }
2929 else {
2930 Coef = 1./D1;
2931 }
2932
2933 // compute the new poles
2934
2935 for (i=Poles.LowerRow(); i<=Poles.UpperRow(); i++) {
2936 if (i < ukk1) {
2937 DvalU = ukk1-i;
2938 }
2939 else if (i > ukk2) {
2940 DvalU = i - ukk2;
2941 }
2942 else {
2943 DvalU = 0.0;
2944 }
2945
2946 for (j=Poles.LowerCol(); j<=Poles.UpperCol(); j++) {
2947 if (i >= UFirstIndex && i <= ULastIndex && j >= VFirstIndex && j <= VLastIndex) {
2948 if (j < vkk1) {
2949 DvalV = vkk1-j;
2950 }
2951 else if (j > vkk2) {
2952 DvalV = j - vkk2;
2953 }
2954 else {
2955 DvalV = 0.0;
2956 }
2957 NewPoles(i,j) = Poles(i,j).Translated((Coef/(DvalU + DvalV + 1.))*Displ);
2958 }
2959 else {
2960 NewPoles(i,j) = Poles(i,j);
2961 }
2962 }
2963 }
2964}
2965
2966//=======================================================================
0d969553
Y
2967// function : Resolution
2968// purpose : this computes an estimate for the maximum of the
7fd59977 2969// partial derivatives both in U and in V
2970//
2971//
0d969553
Y
2972// The calculation resembles at the calculation of curves with
2973// additional index for the control point. Let Si,j be the
2974// control points for ls surface and Di,j the weights.
2975// The checking of upper bounds for the partial derivatives
2976// will be omitted and Su is the next upper bound in the polynomial case :
7fd59977 2977//
2978//
2979//
2980// | Si,j - Si-1,j |
2981// d * Max | ------------- |
2982// i = 2,n | ti+d - ti |
2983// i=1.m
2984//
2985//
0d969553 2986// and in the rational case :
7fd59977 2987//
2988//
2989//
2990// Di,j * (Si,j - Sk,j) - Di-1,j * (Si-1,j - Sk,j)
2991// Max Max d * -----------------------------------------------
2992// k=1,n i dans Rj ti+d - ti
2993// j=1,m
2994// ----------------------------------------------------------------------
2995//
2996// Min Di,j
2997// i=1,n
2998// j=1,m
2999//
3000//
3001//
0d969553 3002// with Rj = {j-d, ...., j+d+d+1}.
7fd59977 3003//
3004//
3005//=======================================================================
3006
3007void BSplSLib::Resolution(const TColgp_Array2OfPnt& Poles,
0e14656b 3008 const TColStd_Array2OfReal* Weights,
7fd59977 3009 const TColStd_Array1OfReal& UKnots,
3010 const TColStd_Array1OfReal& VKnots,
3011 const TColStd_Array1OfInteger& UMults,
3012 const TColStd_Array1OfInteger& VMults,
3013 const Standard_Integer UDegree,
3014 const Standard_Integer VDegree,
3015 const Standard_Boolean URational,
3016 const Standard_Boolean VRational,
3017 const Standard_Boolean UPeriodic,
3018 const Standard_Boolean VPeriodic,
3019 const Standard_Real Tolerance3D,
3020 Standard_Real& UTolerance,
3021 Standard_Real& VTolerance)
3022{
3023 Standard_Real Wij,Wmj,Wji,Wjm;
3024 Standard_Real Xij,Xmj,Xji,Xjm,Xpq,Xqp;
3025 Standard_Real Yij,Ymj,Yji,Yjm,Ypq,Yqp;
3026 Standard_Real Zij,Zmj,Zji,Zjm,Zpq,Zqp;
3027 Standard_Real factor,value,min,min_weights=0,inverse,max_derivative[2];
3028
3029 max_derivative[0] = max_derivative[1] = 0.0e0 ;
3030
3031 Standard_Integer PRowLength, PColLength;
3032 Standard_Integer ii,jj,pp,qq,ii_index,jj_index,pp_index,qq_index;
3033 Standard_Integer ii_minus,upper[2],lower[2],poles_length[2];
3034 Standard_Integer num_poles[2],num_flat_knots[2];
3035
3036 num_flat_knots[0] =
3037 BSplCLib::KnotSequenceLength(UMults,
3038 UDegree,
3039 UPeriodic) ;
3040 num_flat_knots[1] =
3041 BSplCLib::KnotSequenceLength(VMults,
3042 VDegree,
3043 VPeriodic) ;
3044 TColStd_Array1OfReal flat_knots_in_u(1,num_flat_knots[0]) ;
3045 TColStd_Array1OfReal flat_knots_in_v(1,num_flat_knots[1]) ;
3046 BSplCLib::KnotSequence(UKnots,
3047 UMults,
3048 UDegree,
3049 UPeriodic,
3050 flat_knots_in_u) ;
3051 BSplCLib::KnotSequence(VKnots,
3052 VMults,
3053 VDegree,
3054 VPeriodic,
3055 flat_knots_in_v) ;
3056 PRowLength = Poles.RowLength();
3057 PColLength = Poles.ColLength();
3058 if (URational || VRational) {
3059 Standard_Integer Wsize = PRowLength * PColLength;
0e14656b 3060 const TColStd_Array2OfReal& refWights = *Weights;
3061 const Standard_Real * WG = &refWights(refWights.LowerRow(), refWights.LowerCol());
7fd59977 3062 min_weights = WG[0];
3063
3064 for (ii = 1 ; ii < Wsize ; ii++) {
3065 min = WG[ii];
3066 if (min_weights > min) min_weights = min;
3067 }
3068 }
3069 Standard_Integer UD1 = UDegree + 1;
3070 Standard_Integer VD1 = VDegree + 1;
3071 num_poles[0] = num_flat_knots[0] - UD1;
3072 num_poles[1] = num_flat_knots[1] - VD1;
3073 poles_length[0] = PColLength;
3074 poles_length[1] = PRowLength;
3075 if(URational) {
3076 Standard_Integer UD2 = UDegree << 1;
3077 Standard_Integer VD2 = VDegree << 1;
3078
3079 for (ii = 2 ; ii <= num_poles[0] ; ii++) {
3080 ii_index = (ii - 1) % poles_length[0] + 1 ;
3081 ii_minus = (ii - 2) % poles_length[0] + 1 ;
3082 inverse = flat_knots_in_u(ii + UDegree) - flat_knots_in_u(ii) ;
3083 inverse = 1.0e0 / inverse ;
3084 lower[0] = ii - UD1;
3085 if (lower[0] < 1) lower[0] = 1;
3086 upper[0] = ii + UD2 + 1;
3087 if (upper[0] > num_poles[0]) upper[0] = num_poles[0];
3088
3089 for ( jj = 1 ; jj <= num_poles[1] ; jj++) {
3090 jj_index = (jj - 1) % poles_length[1] + 1 ;
3091 lower[1] = jj - VD1;
3092 if (lower[1] < 1) lower[1] = 1;
3093 upper[1] = jj + VD2 + 1;
3094 if (upper[1] > num_poles[1]) upper[1] = num_poles[1];
3095 const gp_Pnt& Pij = Poles .Value(ii_index,jj_index);
0e14656b 3096 Wij = Weights->Value(ii_index,jj_index);
7fd59977 3097 const gp_Pnt& Pmj = Poles .Value(ii_minus,jj_index);
0e14656b 3098 Wmj = Weights->Value(ii_minus,jj_index);
7fd59977 3099 Xij = Pij.X();
3100 Yij = Pij.Y();
3101 Zij = Pij.Z();
3102 Xmj = Pmj.X();
3103 Ymj = Pmj.Y();
3104 Zmj = Pmj.Z();
3105
3106 for (pp = lower[0] ; pp <= upper[0] ; pp++) {
3107 pp_index = (pp - 1) % poles_length[0] + 1 ;
3108
3109 for (qq = lower[1] ; qq <= upper[1] ; qq++) {
3110 value = 0.0e0 ;
3111 qq_index = (qq - 1) % poles_length[1] + 1 ;
3112 const gp_Pnt& Ppq = Poles.Value(pp_index,qq_index);
3113 Xpq = Ppq.X();
3114 Ypq = Ppq.Y();
3115 Zpq = Ppq.Z();
3116 factor = (Xpq - Xij) * Wij;
3117 factor -= (Xpq - Xmj) * Wmj;
3118 if (factor < 0) factor = - factor;
3119 value += factor ;
3120 factor = (Ypq - Yij) * Wij;
3121 factor -= (Ypq - Ymj) * Wmj;
3122 if (factor < 0) factor = - factor;
3123 value += factor ;
3124 factor = (Zpq - Zij) * Wij;
3125 factor -= (Zpq - Zmj) * Wmj;
3126 if (factor < 0) factor = - factor;
3127 value += factor ;
3128 value *= inverse ;
3129 if (max_derivative[0] < value) max_derivative[0] = value ;
3130 }
3131 }
3132 }
3133 }
3134 max_derivative[0] /= min_weights ;
3135 }
3136 else {
3137
3138 for (ii = 2 ; ii <= num_poles[0] ; ii++) {
3139 ii_index = (ii - 1) % poles_length[0] + 1 ;
3140 ii_minus = (ii - 2) % poles_length[0] + 1 ;
3141 inverse = flat_knots_in_u(ii + UDegree) - flat_knots_in_u(ii) ;
3142 inverse = 1.0e0 / inverse ;
3143
3144 for ( jj = 1 ; jj <= num_poles[1] ; jj++) {
3145 jj_index = (jj - 1) % poles_length[1] + 1 ;
3146 value = 0.0e0 ;
3147 const gp_Pnt& Pij = Poles.Value(ii_index,jj_index);
3148 const gp_Pnt& Pmj = Poles.Value(ii_minus,jj_index);
3149 factor = Pij.X() - Pmj.X();
3150 if (factor < 0) factor = - factor;
3151 value += factor;
3152 factor = Pij.Y() - Pmj.Y();
3153 if (factor < 0) factor = - factor;
3154 value += factor;
3155 factor = Pij.Z() - Pmj.Z();
3156 if (factor < 0) factor = - factor;
3157 value += factor;
3158 value *= inverse ;
3159 if (max_derivative[0] < value) max_derivative[0] = value ;
3160 }
3161 }
3162 }
3163 max_derivative[0] *= UDegree ;
3164 if(VRational) {
3165 Standard_Integer UD2 = UDegree << 1;
3166 Standard_Integer VD2 = VDegree << 1;
3167
3168 for (ii = 2 ; ii <= num_poles[1] ; ii++) {
3169 ii_index = (ii - 1) % poles_length[1] + 1 ;
3170 ii_minus = (ii - 2) % poles_length[1] + 1 ;
3171 inverse = flat_knots_in_v(ii + VDegree) - flat_knots_in_v(ii) ;
3172 inverse = 1.0e0 / inverse ;
3173 lower[0] = ii - VD1;
3174 if (lower[0] < 1) lower[0] = 1;
3175 upper[0] = ii + VD2 + 1;
3176 if (upper[0] > num_poles[1]) upper[0] = num_poles[1];
3177
3178 for ( jj = 1 ; jj <= num_poles[0] ; jj++) {
3179 jj_index = (jj - 1) % poles_length[0] + 1 ;
3180 lower[1] = jj - UD1;
3181 if (lower[1] < 1) lower[1] = 1;
3182 upper[1] = jj + UD2 + 1;
3183 if (upper[1] > num_poles[0]) upper[1] = num_poles[0];
3184 const gp_Pnt& Pji = Poles .Value(jj_index,ii_index);
0e14656b 3185 Wji = Weights->Value(jj_index,ii_index);
7fd59977 3186 const gp_Pnt& Pjm = Poles .Value(jj_index,ii_minus);
0e14656b 3187 Wjm = Weights->Value(jj_index,ii_minus);
7fd59977 3188 Xji = Pji.X();
3189 Yji = Pji.Y();
3190 Zji = Pji.Z();
3191 Xjm = Pjm.X();
3192 Yjm = Pjm.Y();
3193 Zjm = Pjm.Z();
3194
3195 for (pp = lower[1] ; pp <= upper[1] ; pp++) {
3196 pp_index = (pp - 1) % poles_length[1] + 1 ;
3197
3198 for (qq = lower[0] ; qq <= upper[0] ; qq++) {
3199 value = 0.0e0 ;
3200 qq_index = (qq - 1) % poles_length[0] + 1 ;
3201 const gp_Pnt& Pqp = Poles.Value(qq_index,pp_index);
3202 Xqp = Pqp.X();
3203 Yqp = Pqp.Y();
3204 Zqp = Pqp.Z();
3205 factor = (Xqp - Xji) * Wji;
3206 factor -= (Xqp - Xjm) * Wjm;
3207 if (factor < 0) factor = - factor;
3208 value += factor ;
3209 factor = (Yqp - Yji) * Wji;
3210 factor -= (Yqp - Yjm) * Wjm;
3211 if (factor < 0) factor = - factor;
3212 value += factor ;
3213 factor = (Zqp - Zji) * Wji;
3214 factor -= (Zqp - Zjm) * Wjm;
3215 if (factor < 0) factor = - factor;
3216 value += factor ;
3217 value *= inverse ;
3218 if (max_derivative[1] < value) max_derivative[1] = value ;
3219 }
3220 }
3221 }
3222 }
3223 max_derivative[1] /= min_weights ;
3224 }
3225 else {
3226
3227 for (ii = 2 ; ii <= num_poles[1] ; ii++) {
3228 ii_index = (ii - 1) % poles_length[1] + 1 ;
3229 ii_minus = (ii - 2) % poles_length[1] + 1 ;
3230 inverse = flat_knots_in_v(ii + VDegree) - flat_knots_in_v(ii) ;
3231 inverse = 1.0e0 / inverse ;
3232
3233 for ( jj = 1 ; jj <= num_poles[0] ; jj++) {
3234 jj_index = (jj - 1) % poles_length[0] + 1 ;
3235 value = 0.0e0 ;
3236 const gp_Pnt& Pji = Poles.Value(jj_index,ii_index);
3237 const gp_Pnt& Pjm = Poles.Value(jj_index,ii_minus);
3238 factor = Pji.X() - Pjm.X() ;
3239 if (factor < 0) factor = - factor;
3240 value += factor;
3241 factor = Pji.Y() - Pjm.Y() ;
3242 if (factor < 0) factor = - factor;
3243 value += factor;
3244 factor = Pji.Z() - Pjm.Z() ;
3245 if (factor < 0) factor = - factor;
3246 value += factor;
3247 value *= inverse ;
3248 if (max_derivative[1] < value) max_derivative[1] = value ;
3249 }
3250 }
3251 }
3252 max_derivative[1] *= VDegree ;
3253 max_derivative[0] *= M_SQRT2 ;
3254 max_derivative[1] *= M_SQRT2 ;
3255 if(max_derivative[0] && max_derivative[1]) {
3256 UTolerance = Tolerance3D / max_derivative[0] ;
3257 VTolerance = Tolerance3D / max_derivative[1] ;
3258 }
3259 else {
3260 UTolerance=VTolerance=0.0;
0797d9d3 3261#ifdef OCCT_DEBUG
04232180 3262 std::cout<<"ElSLib.cxx : maxderivative = 0.0 "<<std::endl;
7fd59977 3263#endif
3264 }
3265}
3266
3267//=======================================================================
3268//function : Interpolate
3269//purpose :
3270//=======================================================================
3271
3272void BSplSLib::Interpolate(const Standard_Integer UDegree,
3273 const Standard_Integer VDegree,
3274 const TColStd_Array1OfReal& UFlatKnots,
3275 const TColStd_Array1OfReal& VFlatKnots,
3276 const TColStd_Array1OfReal& UParameters,
3277 const TColStd_Array1OfReal& VParameters,
3278 TColgp_Array2OfPnt& Poles,
3279 TColStd_Array2OfReal& Weights,
3280 Standard_Integer& InversionProblem)
3281{
3282 Standard_Integer ii, jj, ll, kk, dimension;
3283 Standard_Integer ULength = UParameters.Length();
3284 Standard_Integer VLength = VParameters.Length();
3285 Standard_Real * poles_array;
3286
0d969553 3287 // extraction of iso u
7fd59977 3288 dimension = 4*ULength;
3289 TColStd_Array2OfReal Points(1, VLength,
3290 1, dimension);
3291
3292 Handle(TColStd_HArray1OfInteger) ContactOrder =
3293 new (TColStd_HArray1OfInteger)(1, VLength);
3294 ContactOrder->Init(0);
3295
3296 for (ii=1; ii <= VLength; ii++) {
3297
3298 for (jj=1, ll=1; jj<= ULength; jj++, ll+=4) {
3299 Points(ii,ll) = Poles(jj, ii).X();
3300 Points(ii,ll+1) = Poles(jj, ii).Y();
3301 Points(ii,ll+2) = Poles(jj, ii).Z();
3302 Points(ii,ll+3) = Weights(jj,ii) ;
3303 }
3304 }
3305
0d969553 3306 // interpolation of iso u
7fd59977 3307 poles_array = (Standard_Real *) &Points.ChangeValue(1,1) ;
3308 BSplCLib::Interpolate(VDegree,
3309 VFlatKnots,
3310 VParameters,
3311 ContactOrder->Array1(),
3312 dimension,
3313 poles_array[0],
3314 InversionProblem) ;
3315 if (InversionProblem != 0) return;
3316
0d969553 3317 // extraction of iso v
7fd59977 3318
3319 dimension = VLength*4;
3320 TColStd_Array2OfReal IsoPoles(1, ULength,
3321 1, dimension);
3322
3323 ContactOrder = new (TColStd_HArray1OfInteger)(1, ULength);
3324 ContactOrder->Init(0);
3325 poles_array = (Standard_Real *) &IsoPoles.ChangeValue(1,1) ;
3326
3327 for (ii=1, kk=1; ii <= ULength; ii++, kk+=4) {
3328
3329 for (jj=1, ll=1; jj<= VLength; jj++, ll+=4) {
3330 IsoPoles (ii,ll) = Points(jj, kk);
3331 IsoPoles (ii,ll+1) = Points(jj, kk+1);
3332 IsoPoles (ii,ll+2) = Points(jj, kk+2);
3333 IsoPoles (ii,ll+3) = Points(jj, kk+3);
3334 }
3335 }
0d969553 3336 // interpolation of iso v
7fd59977 3337 BSplCLib::Interpolate(UDegree,
3338 UFlatKnots,
3339 UParameters,
3340 ContactOrder->Array1(),
3341 dimension,
3342 poles_array[0],
3343 InversionProblem);
3344
0d969553 3345 // return results
7fd59977 3346
3347 for (ii=1; ii <= ULength; ii++) {
3348
3349 for (jj=1, ll=1; jj<= VLength; jj++, ll+=4) {
3350 gp_Pnt Pnt(IsoPoles(ii,ll), IsoPoles(ii,ll+1), IsoPoles(ii,ll+2));
3351 Poles.SetValue(ii, jj, Pnt);
3352 Weights.SetValue(ii,jj,IsoPoles(ii,ll+3)) ;
3353 }
3354 }
3355}
3356
3357//=======================================================================
3358//function : Interpolate
3359//purpose :
3360//=======================================================================
3361
3362void BSplSLib::Interpolate(const Standard_Integer UDegree,
3363 const Standard_Integer VDegree,
3364 const TColStd_Array1OfReal& UFlatKnots,
3365 const TColStd_Array1OfReal& VFlatKnots,
3366 const TColStd_Array1OfReal& UParameters,
3367 const TColStd_Array1OfReal& VParameters,
3368 TColgp_Array2OfPnt& Poles,
3369 Standard_Integer& InversionProblem)
3370{
3371 Standard_Integer ii, jj, ll, kk, dimension;
3372 Standard_Integer ULength = UParameters.Length();
3373 Standard_Integer VLength = VParameters.Length();
3374 Standard_Real * poles_array;
3375
0d969553 3376 // extraction of iso u
7fd59977 3377 dimension = 3*ULength;
3378 TColStd_Array2OfReal Points(1, VLength,
3379 1, dimension);
3380
3381 Handle(TColStd_HArray1OfInteger) ContactOrder =
3382 new (TColStd_HArray1OfInteger)(1, VLength);
3383 ContactOrder->Init(0);
3384
3385 for (ii=1; ii <= VLength; ii++) {
3386
3387 for (jj=1, ll=1; jj<= ULength; jj++, ll+=3) {
3388 Points(ii,ll) = Poles(jj, ii).X();
3389 Points(ii,ll+1) = Poles(jj, ii).Y();
3390 Points(ii,ll+2) = Poles(jj, ii).Z();
3391 }
3392 }
3393
0d969553 3394 // interpolation of iso u
7fd59977 3395 poles_array = (Standard_Real *) &Points.ChangeValue(1,1) ;
3396 BSplCLib::Interpolate(VDegree,
3397 VFlatKnots,
3398 VParameters,
3399 ContactOrder->Array1(),
3400 dimension,
3401 poles_array[0],
3402 InversionProblem) ;
3403 if (InversionProblem != 0) return;
3404
0d969553 3405 // extraction of iso v
7fd59977 3406
3407 dimension = VLength*3;
3408 TColStd_Array2OfReal IsoPoles(1, ULength,
3409 1, dimension);
3410
3411 ContactOrder = new (TColStd_HArray1OfInteger)(1, ULength);
3412 ContactOrder->Init(0);
3413 poles_array = (Standard_Real *) &IsoPoles.ChangeValue(1,1) ;
3414
3415 for (ii=1, kk=1; ii <= ULength; ii++, kk+=3) {
3416
3417 for (jj=1, ll=1; jj<= VLength; jj++, ll+=3) {
3418 IsoPoles (ii,ll) = Points(jj, kk);
3419 IsoPoles (ii,ll+1) = Points(jj, kk+1);
3420 IsoPoles (ii,ll+2) = Points(jj, kk+2);
3421 }
3422 }
0d969553 3423 // interpolation of iso v
7fd59977 3424 BSplCLib::Interpolate(UDegree,
3425 UFlatKnots,
3426 UParameters,
3427 ContactOrder->Array1(),
3428 dimension,
3429 poles_array[0],
3430 InversionProblem);
3431
0d969553 3432 // return results
7fd59977 3433
3434 for (ii=1; ii <= ULength; ii++) {
3435
3436 for (jj=1, ll=1; jj<= VLength; jj++, ll+=3) {
3437 gp_Pnt Pnt(IsoPoles(ii,ll), IsoPoles(ii,ll+1), IsoPoles(ii,ll+2));
3438 Poles.SetValue(ii, jj, Pnt);
3439 }
3440 }
3441}
3442
3443//=======================================================================
3444//function : FunctionMultiply
3445//purpose :
3446//=======================================================================
3447
3448void BSplSLib::FunctionMultiply
3449(const BSplSLib_EvaluatorFunction& Function,
3450 const Standard_Integer UBSplineDegree,
3451 const Standard_Integer VBSplineDegree,
3452 const TColStd_Array1OfReal& UBSplineKnots,
3453 const TColStd_Array1OfReal& VBSplineKnots,
0e14656b 3454 const TColStd_Array1OfInteger * UMults,
3455 const TColStd_Array1OfInteger * VMults,
7fd59977 3456 const TColgp_Array2OfPnt& Poles,
0e14656b 3457 const TColStd_Array2OfReal* Weights,
7fd59977 3458 const TColStd_Array1OfReal& UFlatKnots,
3459 const TColStd_Array1OfReal& VFlatKnots,
3460 const Standard_Integer UNewDegree,
3461 const Standard_Integer VNewDegree,
3462 TColgp_Array2OfPnt& NewNumerator,
3463 TColStd_Array2OfReal& NewDenominator,
9fd2d2c3 3464 Standard_Integer& theStatus)
7fd59977 3465{
3466 Standard_Integer num_uparameters,
3467// ii,jj,kk,
3468 ii,jj,
3469 error_code,
3470 num_vparameters ;
3471 Standard_Real result ;
3472
3473 num_uparameters = UFlatKnots.Length() - UNewDegree - 1 ;
3474 num_vparameters = VFlatKnots.Length() - VNewDegree - 1 ;
3475 TColStd_Array1OfReal UParameters(1,num_uparameters) ;
3476 TColStd_Array1OfReal VParameters(1,num_vparameters) ;
3477
3478 if ((NewNumerator.ColLength() == num_uparameters) &&
3479 (NewNumerator.RowLength() == num_vparameters) &&
3480 (NewDenominator.ColLength() == num_uparameters) &&
3481 (NewDenominator.RowLength() == num_vparameters)) {
3482
3483
3484 BSplCLib::BuildSchoenbergPoints(UNewDegree,
3485 UFlatKnots,
3486 UParameters) ;
3487
3488 BSplCLib::BuildSchoenbergPoints(VNewDegree,
3489 VFlatKnots,
3490 VParameters) ;
3491
3492 for (ii = 1 ; ii <= num_uparameters ; ii++) {
3493
3494 for (jj = 1 ; jj <= num_vparameters ; jj++) {
3495 HomogeneousD0(UParameters(ii),
3496 VParameters(jj),
3497 0,
3498 0,
3499 Poles,
3500 Weights,
3501 UBSplineKnots,
3502 VBSplineKnots,
3503 UMults,
3504 VMults,
3505 UBSplineDegree,
3506 VBSplineDegree,
3507 Standard_True,
3508 Standard_True,
3509 Standard_False,
3510 Standard_False,
3511 NewDenominator(ii,jj),
3512 NewNumerator(ii,jj)) ;
3513
41194117 3514 Function.Evaluate (0,
7fd59977 3515 UParameters(ii),
3516 VParameters(jj),
3517 result,
3518 error_code) ;
3519 if (error_code) {
9775fa61 3520 throw Standard_ConstructionError();
7fd59977 3521 }
3522 gp_Pnt& P = NewNumerator(ii,jj);
3523 P.SetX(P.X() * result);
3524 P.SetY(P.Y() * result);
3525 P.SetZ(P.Z() * result);
3526 NewDenominator(ii,jj) *= result ;
3527 }
3528 }
3529 Interpolate(UNewDegree,
3530 VNewDegree,
3531 UFlatKnots,
3532 VFlatKnots,
3533 UParameters,
3534 VParameters,
3535 NewNumerator,
3536 NewDenominator,
9fd2d2c3 3537 theStatus);
7fd59977 3538 }
3539 else {
9775fa61 3540 throw Standard_ConstructionError();
7fd59977 3541 }
3542}