Commit | Line | Data |
---|---|---|
b311480e | 1 | // Copyright (c) 1995-1999 Matra Datavision |
973c2be1 | 2 | // Copyright (c) 1999-2014 OPEN CASCADE SAS |
b311480e | 3 | // |
973c2be1 | 4 | // This file is part of Open CASCADE Technology software library. |
b311480e | 5 | // |
d5f74e42 | 6 | // This library is free software; you can redistribute it and/or modify it under |
7 | // the terms of the GNU Lesser General Public License version 2.1 as published | |
973c2be1 | 8 | // by the Free Software Foundation, with special exception defined in the file |
9 | // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT | |
10 | // distribution for complete text of the license and disclaimer of any warranty. | |
b311480e | 11 | // |
973c2be1 | 12 | // Alternatively, this file may be used under the terms of Open CASCADE |
13 | // commercial license or contractual agreement. | |
b311480e | 14 | |
7fd59977 | 15 | // xab : modified 15-Mar-94 added EvalBSplineMatrix, FactorBandedMatrix, SolveBandedSystem |
16 | // Eval, BuildCache, cacheD0, cacheD1, cacheD2, cacheD3, LocalMatrix, | |
17 | // EvalPolynomial, RationalDerivatives. | |
18 | // xab : 22-Mar-94 : fixed rational problem in BuildCache | |
19 | // xab : 12-Mar-96 : added MovePointAndTangent | |
7fd59977 | 20 | #include <TColStd_Array1OfInteger.hxx> |
21 | #include <TColStd_Array1OfReal.hxx> | |
94f71cad | 22 | #include <TColStd_Array2OfReal.hxx> |
7fd59977 | 23 | #include <gp_Vec2d.hxx> |
24 | #include <Standard_ConstructionError.hxx> | |
25 | #include <PLib.hxx> | |
26 | #include <math_Matrix.hxx> | |
27 | ||
7fd59977 | 28 | //======================================================================= |
29 | //struct : BSplCLib_DataContainer | |
30 | //purpose: Auxiliary structure providing buffers for poles and knots used in | |
31 | // evaluation of bspline (allocated in the stack) | |
32 | //======================================================================= | |
33 | ||
34 | struct BSplCLib_DataContainer | |
35 | { | |
36 | BSplCLib_DataContainer(Standard_Integer Degree) | |
37 | { | |
105aae76 | 38 | (void)Degree; // avoid compiler warning |
41194117 K |
39 | Standard_OutOfRange_Raise_if (Degree > BSplCLib::MaxDegree() || |
40 | BSplCLib::MaxDegree() > 25, | |
41 | "BSplCLib: bspline degree is greater than maximum supported"); | |
7fd59977 | 42 | } |
43 | ||
44 | Standard_Real poles[(25+1)*(Dimension_gen+1)]; | |
45 | Standard_Real knots[2*25]; | |
46 | Standard_Real ders[Dimension_gen*4]; | |
47 | }; | |
48 | ||
49 | //======================================================================= | |
50 | //function : Reverse | |
51 | //purpose : | |
52 | //======================================================================= | |
53 | ||
54 | void BSplCLib::Reverse(Array1OfPoints& Poles, | |
55 | const Standard_Integer L) | |
56 | { | |
57 | Standard_Integer i, l = L; | |
58 | l = Poles.Lower() + (l-Poles.Lower()) % (Poles.Upper()-Poles.Lower()+1); | |
59 | ||
60 | Array1OfPoints temp(0,Poles.Length()-1); | |
61 | ||
62 | for (i = Poles.Lower(); i <= l; i++) | |
63 | temp(l-i) = Poles(i); | |
64 | ||
65 | for (i = l+1; i <= Poles.Upper(); i++) | |
66 | temp(l-Poles.Lower()+Poles.Upper()-i+1) = Poles(i); | |
67 | ||
68 | for (i = Poles.Lower(); i <= Poles.Upper(); i++) | |
69 | Poles(i) = temp(i-Poles.Lower()); | |
70 | } | |
71 | ||
72 | // | |
73 | // CURVES MODIFICATIONS | |
74 | // | |
75 | ||
76 | //======================================================================= | |
77 | //function : RemoveKnot | |
78 | //purpose : | |
79 | //======================================================================= | |
80 | ||
81 | Standard_Boolean BSplCLib::RemoveKnot | |
82 | (const Standard_Integer Index, | |
83 | const Standard_Integer Mult, | |
84 | const Standard_Integer Degree, | |
85 | const Standard_Boolean Periodic, | |
86 | const Array1OfPoints& Poles, | |
0e14656b | 87 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 88 | const TColStd_Array1OfReal& Knots, |
89 | const TColStd_Array1OfInteger& Mults, | |
90 | Array1OfPoints& NewPoles, | |
0e14656b | 91 | TColStd_Array1OfReal* NewWeights, |
7fd59977 | 92 | TColStd_Array1OfReal& NewKnots, |
93 | TColStd_Array1OfInteger& NewMults, | |
94 | const Standard_Real Tolerance) | |
95 | { | |
0e14656b | 96 | Standard_Boolean rational = Weights != NULL; |
7fd59977 | 97 | Standard_Integer dim; |
98 | dim = Dimension_gen; | |
99 | if (rational) dim++; | |
100 | ||
101 | TColStd_Array1OfReal poles(1,dim*Poles.Length()); | |
102 | TColStd_Array1OfReal newpoles(1,dim*NewPoles.Length()); | |
103 | ||
0e14656b | 104 | if (rational) PLib::SetPoles(Poles,*Weights,poles); |
7fd59977 | 105 | else PLib::SetPoles(Poles,poles); |
106 | ||
107 | if (!RemoveKnot(Index,Mult,Degree,Periodic,dim, | |
108 | poles,Knots,Mults,newpoles,NewKnots,NewMults,Tolerance)) | |
109 | return Standard_False; | |
110 | ||
0e14656b | 111 | if (rational) PLib::GetPoles(newpoles,NewPoles,*NewWeights); |
7fd59977 | 112 | else PLib::GetPoles(newpoles,NewPoles); |
113 | return Standard_True; | |
114 | } | |
115 | ||
116 | //======================================================================= | |
117 | //function : InsertKnots | |
118 | //purpose : insert an array of knots and multiplicities | |
119 | //======================================================================= | |
120 | ||
121 | void BSplCLib::InsertKnots | |
122 | (const Standard_Integer Degree, | |
123 | const Standard_Boolean Periodic, | |
124 | const Array1OfPoints& Poles, | |
0e14656b | 125 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 126 | const TColStd_Array1OfReal& Knots, |
127 | const TColStd_Array1OfInteger& Mults, | |
128 | const TColStd_Array1OfReal& AddKnots, | |
0e14656b | 129 | const TColStd_Array1OfInteger* AddMults, |
7fd59977 | 130 | Array1OfPoints& NewPoles, |
0e14656b | 131 | TColStd_Array1OfReal* NewWeights, |
7fd59977 | 132 | TColStd_Array1OfReal& NewKnots, |
133 | TColStd_Array1OfInteger& NewMults, | |
134 | const Standard_Real Epsilon, | |
135 | const Standard_Boolean Add) | |
136 | { | |
0e14656b | 137 | Standard_Boolean rational = Weights != NULL; |
7fd59977 | 138 | Standard_Integer dim; |
139 | dim = Dimension_gen; | |
140 | if (rational) dim++; | |
141 | ||
142 | TColStd_Array1OfReal poles(1,dim*Poles.Length()); | |
143 | TColStd_Array1OfReal newpoles(1,dim*NewPoles.Length()); | |
144 | ||
0e14656b | 145 | if (rational) PLib::SetPoles(Poles,*Weights,poles); |
7fd59977 | 146 | else PLib::SetPoles(Poles,poles); |
147 | ||
148 | InsertKnots(Degree,Periodic,dim,poles,Knots,Mults, | |
149 | AddKnots,AddMults,newpoles,NewKnots,NewMults,Epsilon,Add); | |
150 | ||
0e14656b | 151 | if (rational) PLib::GetPoles(newpoles,NewPoles,*NewWeights); |
7fd59977 | 152 | else PLib::GetPoles(newpoles,NewPoles); |
153 | } | |
154 | ||
155 | //======================================================================= | |
156 | //function : InsertKnot | |
157 | //purpose : | |
158 | //======================================================================= | |
159 | ||
160 | void BSplCLib::InsertKnot(const Standard_Integer , | |
161 | const Standard_Real U, | |
162 | const Standard_Integer UMult, | |
163 | const Standard_Integer Degree, | |
164 | const Standard_Boolean Periodic, | |
165 | const Array1OfPoints& Poles, | |
0e14656b | 166 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 167 | const TColStd_Array1OfReal& Knots, |
168 | const TColStd_Array1OfInteger& Mults, | |
169 | Array1OfPoints& NewPoles, | |
0e14656b | 170 | TColStd_Array1OfReal* NewWeights) |
7fd59977 | 171 | { |
172 | TColStd_Array1OfReal k(1,1); | |
173 | k(1) = U; | |
174 | TColStd_Array1OfInteger m(1,1); | |
175 | m(1) = UMult; | |
176 | TColStd_Array1OfReal nk(1,Knots.Length()+1); | |
177 | TColStd_Array1OfInteger nm(1,Knots.Length()+1); | |
178 | InsertKnots(Degree,Periodic,Poles,Weights,Knots,Mults, | |
0e14656b | 179 | k,&m,NewPoles,NewWeights,nk,nm,Epsilon(U)); |
7fd59977 | 180 | } |
181 | ||
182 | //======================================================================= | |
183 | //function : RaiseMultiplicity | |
184 | //purpose : | |
185 | //======================================================================= | |
186 | ||
187 | void BSplCLib::RaiseMultiplicity(const Standard_Integer KnotIndex, | |
188 | const Standard_Integer Mult, | |
189 | const Standard_Integer Degree, | |
190 | const Standard_Boolean Periodic, | |
191 | const Array1OfPoints& Poles, | |
0e14656b | 192 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 193 | const TColStd_Array1OfReal& Knots, |
194 | const TColStd_Array1OfInteger& Mults, | |
195 | Array1OfPoints& NewPoles, | |
0e14656b | 196 | TColStd_Array1OfReal* NewWeights) |
7fd59977 | 197 | { |
198 | TColStd_Array1OfReal k(1,1); | |
199 | k(1) = Knots(KnotIndex); | |
200 | TColStd_Array1OfInteger m(1,1); | |
201 | m(1) = Mult - Mults(KnotIndex); | |
202 | TColStd_Array1OfReal nk(1,Knots.Length()); | |
203 | TColStd_Array1OfInteger nm(1,Knots.Length()); | |
204 | InsertKnots(Degree,Periodic,Poles,Weights,Knots,Mults, | |
0e14656b | 205 | k,&m,NewPoles,NewWeights,nk,nm,Epsilon(k(1))); |
7fd59977 | 206 | } |
207 | ||
208 | //======================================================================= | |
209 | //function : IncreaseDegree | |
210 | //purpose : | |
211 | //======================================================================= | |
212 | ||
213 | void BSplCLib::IncreaseDegree | |
214 | (const Standard_Integer Degree, | |
215 | const Standard_Integer NewDegree, | |
216 | const Standard_Boolean Periodic, | |
217 | const Array1OfPoints& Poles, | |
0e14656b | 218 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 219 | const TColStd_Array1OfReal& Knots, |
220 | const TColStd_Array1OfInteger& Mults, | |
221 | Array1OfPoints& NewPoles, | |
0e14656b | 222 | TColStd_Array1OfReal* NewWeights, |
7fd59977 | 223 | TColStd_Array1OfReal& NewKnots, |
224 | TColStd_Array1OfInteger& NewMults) | |
225 | { | |
0e14656b | 226 | Standard_Boolean rational = Weights != NULL; |
7fd59977 | 227 | Standard_Integer dim; |
228 | dim = Dimension_gen; | |
229 | if (rational) dim++; | |
230 | ||
231 | TColStd_Array1OfReal poles(1,dim*Poles.Length()); | |
232 | TColStd_Array1OfReal newpoles(1,dim*NewPoles.Length()); | |
233 | ||
0e14656b | 234 | if (rational) PLib::SetPoles(Poles,*Weights,poles); |
7fd59977 | 235 | else PLib::SetPoles(Poles,poles); |
236 | ||
237 | IncreaseDegree(Degree,NewDegree,Periodic,dim,poles,Knots,Mults, | |
238 | newpoles,NewKnots,NewMults); | |
239 | ||
0e14656b | 240 | if (rational) PLib::GetPoles(newpoles,NewPoles,*NewWeights); |
7fd59977 | 241 | else PLib::GetPoles(newpoles,NewPoles); |
242 | } | |
243 | ||
244 | //======================================================================= | |
245 | //function : Unperiodize | |
246 | //purpose : | |
247 | //======================================================================= | |
248 | ||
249 | void BSplCLib::Unperiodize | |
250 | (const Standard_Integer Degree, | |
251 | const TColStd_Array1OfInteger& Mults, | |
252 | const TColStd_Array1OfReal& Knots, | |
253 | const Array1OfPoints& Poles, | |
0e14656b | 254 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 255 | TColStd_Array1OfInteger& NewMults, |
256 | TColStd_Array1OfReal& NewKnots, | |
257 | Array1OfPoints& NewPoles, | |
0e14656b | 258 | TColStd_Array1OfReal* NewWeights) |
7fd59977 | 259 | { |
0e14656b | 260 | Standard_Boolean rational = Weights != NULL; |
7fd59977 | 261 | Standard_Integer dim; |
262 | dim = Dimension_gen; | |
263 | if (rational) dim++; | |
264 | ||
265 | TColStd_Array1OfReal poles(1,dim*Poles.Length()); | |
266 | TColStd_Array1OfReal newpoles(1,dim*NewPoles.Length()); | |
267 | ||
0e14656b | 268 | if (rational) PLib::SetPoles(Poles,*Weights,poles); |
7fd59977 | 269 | else PLib::SetPoles(Poles,poles); |
270 | ||
271 | Unperiodize(Degree, dim, Mults, Knots, poles, | |
272 | NewMults,NewKnots, newpoles); | |
273 | ||
0e14656b | 274 | if (rational) PLib::GetPoles(newpoles,NewPoles,*NewWeights); |
7fd59977 | 275 | else PLib::GetPoles(newpoles,NewPoles); |
276 | } | |
277 | ||
278 | //======================================================================= | |
279 | //function : Trimming | |
280 | //purpose : | |
281 | //======================================================================= | |
282 | ||
283 | void BSplCLib::Trimming(const Standard_Integer Degree, | |
284 | const Standard_Boolean Periodic, | |
285 | const TColStd_Array1OfReal& Knots, | |
286 | const TColStd_Array1OfInteger& Mults, | |
287 | const Array1OfPoints& Poles, | |
0e14656b | 288 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 289 | const Standard_Real U1, |
290 | const Standard_Real U2, | |
291 | TColStd_Array1OfReal& NewKnots, | |
292 | TColStd_Array1OfInteger& NewMults, | |
293 | Array1OfPoints& NewPoles, | |
0e14656b | 294 | TColStd_Array1OfReal* NewWeights) |
7fd59977 | 295 | { |
0e14656b | 296 | Standard_Boolean rational = Weights != NULL; |
7fd59977 | 297 | Standard_Integer dim; |
298 | dim = Dimension_gen; | |
299 | if (rational) dim++; | |
300 | ||
301 | TColStd_Array1OfReal poles(1,dim*Poles.Length()); | |
302 | TColStd_Array1OfReal newpoles(1,dim*NewPoles.Length()); | |
303 | ||
0e14656b | 304 | if (rational) PLib::SetPoles(Poles,*Weights,poles); |
7fd59977 | 305 | else PLib::SetPoles(Poles,poles); |
306 | ||
307 | Trimming(Degree, Periodic, dim, Knots, Mults, poles, U1, U2, | |
308 | NewKnots, NewMults, newpoles); | |
309 | ||
0e14656b | 310 | if (rational) PLib::GetPoles(newpoles,NewPoles,*NewWeights); |
7fd59977 | 311 | else PLib::GetPoles(newpoles,NewPoles); |
312 | } | |
313 | ||
314 | //-------------------------------------------------------------------------- | |
315 | //ELEMENTARY COMPUTATIONS | |
316 | //-------------------------------------------------------------------------- | |
317 | // | |
318 | // All the following methods are methods of low level used in BSplCLib | |
319 | // but also in BSplSLib (but not in Geom) | |
320 | // At the creation time of this package there is no possibility | |
321 | // to declare private methods of package and to declare friend | |
322 | // methods of package. It could interesting to declare that BSplSLib | |
323 | // is a package friend to BSplCLib because it uses all the basis methods | |
324 | // of BSplCLib. | |
325 | //-------------------------------------------------------------------------- | |
326 | ||
327 | //======================================================================= | |
328 | //function : BuildEval | |
329 | //purpose : builds the local array for evaluation | |
330 | //======================================================================= | |
331 | ||
332 | void BSplCLib::BuildEval(const Standard_Integer Degree, | |
333 | const Standard_Integer Index, | |
334 | const Array1OfPoints& Poles, | |
0e14656b | 335 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 336 | Standard_Real& LP) |
337 | { | |
338 | Standard_Real w, *pole = &LP; | |
339 | Standard_Integer PLower = Poles.Lower(); | |
340 | Standard_Integer PUpper = Poles.Upper(); | |
341 | Standard_Integer i; | |
342 | Standard_Integer ip = PLower + Index - 1; | |
0e14656b | 343 | if (Weights == NULL) { |
7fd59977 | 344 | for (i = 0; i <= Degree; i++) { |
345 | ip++; | |
346 | if (ip > PUpper) ip = PLower; | |
347 | const Point& P = Poles(ip); | |
348 | PointToCoords (pole,P,+0); | |
349 | pole += Dimension_gen; | |
350 | } | |
351 | } | |
352 | else { | |
353 | for (i = 0; i <= Degree; i++) { | |
354 | ip++; | |
355 | if (ip > PUpper) ip = PLower; | |
356 | const Point& P = Poles(ip); | |
0e14656b | 357 | pole[Dimension_gen] = w = (*Weights)(ip); |
7fd59977 | 358 | PointToCoords (pole, P, * w); |
359 | pole += Dimension_gen + 1; | |
360 | } | |
361 | } | |
362 | } | |
363 | ||
364 | //======================================================================= | |
365 | //function : PrepareEval | |
366 | //purpose : stores data for Eval in the local arrays | |
367 | // dc.poles and dc.knots | |
368 | //======================================================================= | |
369 | ||
370 | static void PrepareEval | |
371 | (Standard_Real& u, | |
372 | Standard_Integer& index, | |
373 | Standard_Integer& dim, | |
374 | Standard_Boolean& rational, | |
375 | const Standard_Integer Degree, | |
376 | const Standard_Boolean Periodic, | |
377 | const Array1OfPoints& Poles, | |
0e14656b | 378 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 379 | const TColStd_Array1OfReal& Knots, |
0e14656b | 380 | const TColStd_Array1OfInteger* Mults, |
7fd59977 | 381 | BSplCLib_DataContainer& dc) |
382 | { | |
383 | // Set the Index | |
384 | BSplCLib::LocateParameter(Degree,Knots,Mults,u,Periodic,index,u); | |
385 | ||
386 | // make the knots | |
387 | BSplCLib::BuildKnots(Degree,index,Periodic,Knots,Mults,*dc.knots); | |
0e14656b | 388 | if (Mults == NULL) |
7fd59977 | 389 | index -= Knots.Lower() + Degree; |
390 | else | |
0e14656b | 391 | index = BSplCLib::PoleIndex(Degree,index,Periodic,*Mults); |
7fd59977 | 392 | |
393 | // check truly rational | |
0e14656b | 394 | rational = (Weights != NULL); |
7fd59977 | 395 | if (rational) { |
0e14656b | 396 | Standard_Integer WLower = Weights->Lower() + index; |
397 | rational = BSplCLib::IsRational(*Weights, WLower, WLower + Degree); | |
7fd59977 | 398 | } |
399 | ||
400 | // make the poles | |
401 | if (rational) { | |
402 | dim = Dimension_gen + 1; | |
403 | BSplCLib::BuildEval(Degree, index, Poles, Weights , *dc.poles); | |
404 | } | |
405 | else { | |
406 | dim = Dimension_gen; | |
407 | BSplCLib::BuildEval(Degree, index, Poles, BSplCLib::NoWeights(), *dc.poles); | |
408 | } | |
409 | } | |
410 | ||
411 | //======================================================================= | |
412 | //function : D0 | |
413 | //purpose : | |
414 | //======================================================================= | |
415 | ||
416 | void BSplCLib::D0 | |
417 | (const Standard_Real U, | |
418 | const Standard_Integer Index, | |
419 | const Standard_Integer Degree, | |
420 | const Standard_Boolean Periodic, | |
421 | const Array1OfPoints& Poles, | |
0e14656b | 422 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 423 | const TColStd_Array1OfReal& Knots, |
0e14656b | 424 | const TColStd_Array1OfInteger* Mults, |
7fd59977 | 425 | Point& P) |
426 | { | |
427 | // Standard_Integer k,dim,index = Index; | |
428 | Standard_Integer dim,index = Index; | |
429 | Standard_Real u = U; | |
430 | Standard_Boolean rational; | |
431 | BSplCLib_DataContainer dc(Degree); | |
432 | PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc); | |
433 | BSplCLib::Eval(u,Degree,*dc.knots,dim,*dc.poles); | |
434 | ||
435 | if (rational) { | |
436 | Standard_Real w = dc.poles[Dimension_gen]; | |
437 | CoordsToPoint (P, dc.poles, / w); | |
438 | } | |
439 | else | |
440 | CoordsToPoint (P, dc.poles, ); | |
441 | } | |
442 | ||
443 | //======================================================================= | |
444 | //function : D1 | |
445 | //purpose : | |
446 | //======================================================================= | |
447 | ||
448 | void BSplCLib::D1 | |
449 | (const Standard_Real U, | |
450 | const Standard_Integer Index, | |
451 | const Standard_Integer Degree, | |
452 | const Standard_Boolean Periodic, | |
453 | const Array1OfPoints& Poles, | |
0e14656b | 454 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 455 | const TColStd_Array1OfReal& Knots, |
0e14656b | 456 | const TColStd_Array1OfInteger* Mults, |
7fd59977 | 457 | Point& P, |
458 | Vector& V) | |
459 | { | |
460 | Standard_Integer dim,index = Index; | |
461 | Standard_Real u = U; | |
462 | Standard_Boolean rational; | |
463 | BSplCLib_DataContainer dc(Degree); | |
464 | PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc); | |
465 | BSplCLib::Bohm(u,Degree,1,*dc.knots,dim,*dc.poles); | |
466 | Standard_Real *result = dc.poles; | |
467 | if (rational) { | |
468 | PLib::RationalDerivative(Degree,1,Dimension_gen,*dc.poles,*dc.ders); | |
469 | result = dc.ders; | |
470 | } | |
471 | ||
472 | CoordsToPoint (P, result, ); | |
473 | CoordsToPoint (V, result + Dimension_gen, ); | |
474 | } | |
475 | ||
476 | //======================================================================= | |
477 | //function : D2 | |
478 | //purpose : | |
479 | //======================================================================= | |
480 | ||
481 | void BSplCLib::D2 | |
482 | (const Standard_Real U, | |
483 | const Standard_Integer Index, | |
484 | const Standard_Integer Degree, | |
485 | const Standard_Boolean Periodic, | |
486 | const Array1OfPoints& Poles, | |
0e14656b | 487 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 488 | const TColStd_Array1OfReal& Knots, |
0e14656b | 489 | const TColStd_Array1OfInteger* Mults, |
7fd59977 | 490 | Point& P, |
491 | Vector& V1, | |
492 | Vector& V2) | |
493 | { | |
494 | Standard_Integer dim,index = Index; | |
495 | Standard_Real u = U; | |
496 | Standard_Boolean rational; | |
497 | BSplCLib_DataContainer dc(Degree); | |
498 | PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc); | |
499 | BSplCLib::Bohm(u,Degree,2,*dc.knots,dim,*dc.poles); | |
500 | Standard_Real *result = dc.poles; | |
501 | if (rational) { | |
502 | PLib::RationalDerivative(Degree,2,Dimension_gen,*dc.poles,*dc.ders); | |
503 | result = dc.ders; | |
504 | } | |
505 | ||
506 | CoordsToPoint (P, result, ); | |
507 | CoordsToPoint (V1, result + Dimension_gen, ); | |
508 | if (!rational && (Degree < 2)) | |
509 | NullifyPoint (V2); | |
510 | else | |
511 | CoordsToPoint (V2, result + 2 * Dimension_gen, ); | |
512 | } | |
513 | ||
514 | //======================================================================= | |
515 | //function : D3 | |
516 | //purpose : | |
517 | //======================================================================= | |
518 | ||
519 | void BSplCLib::D3 | |
520 | (const Standard_Real U, | |
521 | const Standard_Integer Index, | |
522 | const Standard_Integer Degree, | |
523 | const Standard_Boolean Periodic, | |
524 | const Array1OfPoints& Poles, | |
0e14656b | 525 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 526 | const TColStd_Array1OfReal& Knots, |
0e14656b | 527 | const TColStd_Array1OfInteger* Mults, |
7fd59977 | 528 | Point& P, |
529 | Vector& V1, | |
530 | Vector& V2, | |
531 | Vector& V3) | |
532 | { | |
533 | Standard_Integer dim,index = Index; | |
534 | Standard_Real u = U; | |
535 | Standard_Boolean rational; | |
536 | BSplCLib_DataContainer dc(Degree); | |
537 | PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc); | |
538 | BSplCLib::Bohm(u,Degree,3,*dc.knots,dim,*dc.poles); | |
539 | Standard_Real *result = dc.poles; | |
540 | if (rational) { | |
541 | PLib::RationalDerivative(Degree,3,Dimension_gen,*dc.poles,*dc.ders); | |
542 | result = dc.ders; | |
543 | } | |
544 | ||
545 | CoordsToPoint (P, result, ); | |
546 | CoordsToPoint (V1, result + Dimension_gen, ); | |
547 | if (!rational && (Degree < 2)) | |
548 | NullifyPoint (V2); | |
549 | else | |
550 | CoordsToPoint (V2, result + 2 * Dimension_gen, ); | |
551 | if (!rational && (Degree < 3)) | |
552 | NullifyPoint (V3); | |
553 | else | |
554 | CoordsToPoint (V3, result + 3 * Dimension_gen, ); | |
555 | } | |
556 | ||
557 | //======================================================================= | |
558 | //function : DN | |
559 | //purpose : | |
560 | //======================================================================= | |
561 | ||
562 | void BSplCLib::DN | |
563 | (const Standard_Real U, | |
564 | const Standard_Integer N, | |
565 | const Standard_Integer Index, | |
566 | const Standard_Integer Degree, | |
567 | const Standard_Boolean Periodic, | |
568 | const Array1OfPoints& Poles, | |
0e14656b | 569 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 570 | const TColStd_Array1OfReal& Knots, |
0e14656b | 571 | const TColStd_Array1OfInteger* Mults, |
7fd59977 | 572 | Vector& VN) |
573 | { | |
574 | Standard_Integer dim,index = Index; | |
575 | Standard_Real u = U; | |
576 | Standard_Boolean rational; | |
577 | BSplCLib_DataContainer dc(Degree); | |
578 | PrepareEval(u,index,dim,rational,Degree,Periodic,Poles,Weights,Knots,Mults,dc); | |
579 | BSplCLib::Bohm(u,Degree,N,*dc.knots,dim,*dc.poles); | |
580 | ||
581 | if (rational) { | |
582 | Standard_Real v[Dimension_gen]; | |
583 | PLib::RationalDerivative(Degree,N,Dimension_gen,*dc.poles,v[0],Standard_False); | |
584 | CoordsToPoint (VN, v, ); | |
585 | } | |
586 | else { | |
587 | if (N > Degree) | |
588 | NullifyPoint (VN); | |
589 | else { | |
590 | Standard_Real *DN = dc.poles + N * Dimension_gen; | |
591 | CoordsToPoint (VN, DN, ); | |
592 | } | |
593 | } | |
594 | } | |
595 | ||
596 | //======================================================================= | |
597 | //function : Solves a LU factored Matrix | |
598 | //purpose : | |
599 | //======================================================================= | |
600 | ||
601 | Standard_Integer | |
602 | BSplCLib::SolveBandedSystem(const math_Matrix& Matrix, | |
603 | const Standard_Integer UpperBandWidth, | |
604 | const Standard_Integer LowerBandWidth, | |
605 | Array1OfPoints& PolesArray) | |
606 | { | |
607 | Standard_Real *PArray ; | |
608 | PArray = (Standard_Real *) &PolesArray(PolesArray.Lower()) ; | |
609 | ||
610 | return BSplCLib::SolveBandedSystem(Matrix, | |
611 | UpperBandWidth, | |
612 | LowerBandWidth, | |
613 | Dimension_gen, | |
614 | PArray[0]) ; | |
615 | } | |
616 | ||
617 | //======================================================================= | |
618 | //function : Solves a LU factored Matrix | |
619 | //purpose : if HomogeneousFlag is 1 then the input and the output | |
620 | // will be homogeneous that is no division or multiplication | |
621 | // by weigths will happen. On the contrary if HomogeneousFlag | |
622 | // is 0 then the poles will be multiplied first by the weights | |
623 | // and after interpolation they will be devided by the weights | |
624 | //======================================================================= | |
625 | ||
626 | Standard_Integer | |
627 | BSplCLib::SolveBandedSystem(const math_Matrix& Matrix, | |
628 | const Standard_Integer UpperBandWidth, | |
629 | const Standard_Integer LowerBandWidth, | |
630 | const Standard_Boolean HomogeneousFlag, | |
631 | Array1OfPoints& PolesArray, | |
632 | TColStd_Array1OfReal& WeightsArray) | |
633 | { | |
634 | Standard_Real *PArray, | |
635 | * WArray ; | |
636 | PArray = (Standard_Real *) &PolesArray(PolesArray.Lower()) ; | |
637 | WArray = (Standard_Real *) &WeightsArray(WeightsArray.Lower()) ; | |
638 | return BSplCLib::SolveBandedSystem(Matrix, | |
639 | UpperBandWidth, | |
640 | LowerBandWidth, | |
641 | HomogeneousFlag, | |
642 | Dimension_gen, | |
643 | PArray[0], | |
644 | WArray[0]) ; | |
645 | } | |
646 | ||
647 | //======================================================================= | |
648 | //function : Evaluates a Bspline function : uses the ExtrapMode | |
649 | //purpose : the function is extrapolated using the Taylor expansion | |
650 | // of degree ExtrapMode[0] to the left and the Taylor | |
651 | // expansion of degree ExtrapMode[1] to the right | |
652 | // if the HomogeneousFlag == True than the Poles are supposed | |
653 | // to be stored homogeneously and the result will also be homogeneous | |
654 | // Valid only if Weights | |
655 | //======================================================================= | |
656 | void BSplCLib::Eval(const Standard_Real Parameter, | |
657 | const Standard_Boolean PeriodicFlag, | |
658 | const Standard_Boolean HomogeneousFlag, | |
659 | Standard_Integer& ExtrapMode, | |
660 | const Standard_Integer Degree, | |
661 | const TColStd_Array1OfReal& FlatKnots, | |
662 | const Array1OfPoints& PolesArray, | |
663 | const TColStd_Array1OfReal& WeightsArray, | |
664 | Point& aPoint, | |
665 | Standard_Real& aWeight) | |
666 | { | |
667 | Standard_Real Inverse, | |
668 | P[Dimension_gen], | |
669 | *PArray, | |
670 | *WArray ; | |
671 | Standard_Integer kk ; | |
672 | PArray = (Standard_Real *) &PolesArray(PolesArray.Lower()) ; | |
673 | WArray = (Standard_Real *) &WeightsArray(WeightsArray.Lower()) ; | |
674 | if (HomogeneousFlag) { | |
675 | BSplCLib::Eval(Parameter, | |
676 | PeriodicFlag, | |
677 | 0, | |
678 | ExtrapMode, | |
679 | Degree, | |
680 | FlatKnots, | |
681 | Dimension_gen, | |
682 | PArray[0], | |
683 | P[0]) ; | |
684 | BSplCLib::Eval(Parameter, | |
685 | PeriodicFlag, | |
686 | 0, | |
687 | ExtrapMode, | |
688 | Degree, | |
689 | FlatKnots, | |
690 | 1, | |
691 | WArray[0], | |
692 | aWeight) ; | |
693 | } | |
694 | else { | |
695 | BSplCLib::Eval(Parameter, | |
696 | PeriodicFlag, | |
697 | 0, | |
698 | ExtrapMode, | |
699 | Degree, | |
700 | FlatKnots, | |
701 | Dimension_gen, | |
702 | PArray[0], | |
703 | WArray[0], | |
704 | P[0], | |
705 | aWeight) ; | |
706 | Inverse = 1.0e0 / aWeight ; | |
707 | ||
708 | for (kk = 0 ; kk < Dimension_gen ; kk++) { | |
709 | P[kk] *= Inverse ; | |
710 | } | |
711 | } | |
712 | ||
713 | for (kk = 0 ; kk < Dimension_gen ; kk++) | |
714 | aPoint.SetCoord(kk + 1, P[kk]) ; | |
715 | } | |
716 | ||
717 | //======================================================================= | |
718 | //function : CacheD0 | |
719 | //purpose : Evaluates the polynomial cache of the Bspline Curve | |
720 | // | |
721 | //======================================================================= | |
722 | void BSplCLib::CacheD0(const Standard_Real Parameter, | |
723 | const Standard_Integer Degree, | |
724 | const Standard_Real CacheParameter, | |
725 | const Standard_Real SpanLenght, | |
726 | const Array1OfPoints& PolesArray, | |
0e14656b | 727 | const TColStd_Array1OfReal* WeightsArray, |
7fd59977 | 728 | Point& aPoint) |
729 | { | |
730 | // | |
731 | // the CacheParameter is where the cache polynomial was evaluated in homogeneous | |
732 | // form | |
733 | // the SpanLenght is the normalizing factor so that the polynomial is between | |
734 | // 0 and 1 | |
735 | Standard_Real NewParameter, Inverse; | |
736 | Standard_Real * PArray = (Standard_Real *) &(PolesArray(PolesArray.Lower())) ; | |
737 | Standard_Real * myPoint = (Standard_Real *) &aPoint ; | |
738 | NewParameter = (Parameter - CacheParameter) / SpanLenght ; | |
739 | PLib::NoDerivativeEvalPolynomial(NewParameter, | |
740 | Degree, | |
741 | Dimension_gen, | |
742 | Degree * Dimension_gen, | |
743 | PArray[0], | |
744 | myPoint[0]) ; | |
0e14656b | 745 | if (WeightsArray != NULL) { |
746 | const TColStd_Array1OfReal& refWeights = *WeightsArray; | |
7fd59977 | 747 | Standard_Real * |
0e14656b | 748 | WArray = (Standard_Real *) &refWeights(refWeights.Lower()) ; |
7fd59977 | 749 | PLib::NoDerivativeEvalPolynomial(NewParameter, |
750 | Degree, | |
751 | 1, | |
752 | Degree, | |
753 | WArray[0], | |
754 | Inverse) ; | |
755 | ||
756 | Inverse = 1.0e0 / Inverse; | |
757 | ModifyCoords (myPoint, *= Inverse); | |
758 | } | |
759 | } | |
760 | ||
761 | //======================================================================= | |
762 | //function : CacheD1 | |
763 | //purpose : Evaluates the polynomial cache of the Bspline Curve | |
764 | // | |
765 | //======================================================================= | |
766 | void BSplCLib::CacheD1(const Standard_Real Parameter, | |
767 | const Standard_Integer Degree, | |
768 | const Standard_Real CacheParameter, | |
769 | const Standard_Real SpanLenght, | |
770 | const Array1OfPoints& PolesArray, | |
0e14656b | 771 | const TColStd_Array1OfReal* WeightsArray, |
7fd59977 | 772 | Point& aPoint, |
773 | Vector& aVector) | |
774 | { | |
775 | // | |
776 | // the CacheParameter is where the cache polynomial was evaluated in homogeneous | |
777 | // form | |
778 | // the SpanLenght is the normalizing factor so that the polynomial is between | |
779 | // 0 and 1 | |
780 | Standard_Real LocalPDerivatives[Dimension_gen << 1] ; | |
781 | // Standard_Real LocalWDerivatives[2], NewParameter, Inverse ; | |
782 | Standard_Real LocalWDerivatives[2], NewParameter ; | |
783 | ||
784 | Standard_Real * | |
785 | PArray = (Standard_Real *) &(PolesArray(PolesArray.Lower())) ; | |
786 | Standard_Real * | |
787 | myPoint = (Standard_Real *) &aPoint ; | |
788 | Standard_Real * | |
789 | myVector = (Standard_Real *) &aVector ; | |
790 | NewParameter = (Parameter - CacheParameter) / SpanLenght ; | |
791 | PLib::EvalPolynomial(NewParameter, | |
792 | 1, | |
793 | Degree, | |
794 | Dimension_gen, | |
795 | PArray[0], | |
796 | LocalPDerivatives[0]) ; | |
797 | // | |
798 | // unormalize derivatives since those are computed normalized | |
799 | // | |
800 | ||
801 | ModifyCoords (LocalPDerivatives + Dimension_gen, /= SpanLenght); | |
802 | ||
0e14656b | 803 | if (WeightsArray != NULL) { |
804 | const TColStd_Array1OfReal& refWeights = *WeightsArray; | |
7fd59977 | 805 | Standard_Real * |
0e14656b | 806 | WArray = (Standard_Real *) &refWeights(refWeights.Lower()) ; |
7fd59977 | 807 | PLib::EvalPolynomial(NewParameter, |
808 | 1, | |
809 | Degree, | |
810 | 1, | |
811 | WArray[0], | |
812 | LocalWDerivatives[0]) ; | |
813 | // | |
814 | // unormalize the result since the polynomial stored in the cache | |
815 | // is normalized between 0 and 1 | |
816 | // | |
817 | LocalWDerivatives[1] /= SpanLenght ; | |
818 | ||
819 | PLib::RationalDerivatives(1, | |
820 | Dimension_gen, | |
821 | LocalPDerivatives[0], | |
822 | LocalWDerivatives[0], | |
823 | LocalPDerivatives[0]) ; | |
824 | } | |
825 | ||
826 | CopyCoords (myPoint, LocalPDerivatives); | |
827 | CopyCoords (myVector, LocalPDerivatives + Dimension_gen); | |
828 | } | |
829 | ||
830 | //======================================================================= | |
831 | //function : CacheD2 | |
832 | //purpose : Evaluates the polynomial cache of the Bspline Curve | |
833 | // | |
834 | //======================================================================= | |
835 | void BSplCLib::CacheD2(const Standard_Real Parameter, | |
836 | const Standard_Integer Degree, | |
837 | const Standard_Real CacheParameter, | |
838 | const Standard_Real SpanLenght, | |
839 | const Array1OfPoints& PolesArray, | |
0e14656b | 840 | const TColStd_Array1OfReal* WeightsArray, |
7fd59977 | 841 | Point& aPoint, |
842 | Vector& aVector1, | |
843 | Vector& aVector2) | |
844 | { | |
845 | // | |
846 | // the CacheParameter is where the cache polynomial was evaluated in homogeneous | |
847 | // form | |
848 | // the SpanLenght is the normalizing factor so that the polynomial is between | |
849 | // 0 and 1 | |
850 | Standard_Integer ii,Index,EndIndex; | |
851 | Standard_Real LocalPDerivatives[(Dimension_gen << 1) + Dimension_gen] ; | |
852 | // Standard_Real LocalWDerivatives[3], NewParameter, Factor, Inverse ; | |
853 | Standard_Real LocalWDerivatives[3], NewParameter, Factor ; | |
854 | Standard_Real * PArray = (Standard_Real *) &(PolesArray(PolesArray.Lower())) ; | |
855 | Standard_Real * myPoint = (Standard_Real *) &aPoint ; | |
856 | Standard_Real * myVector1 = (Standard_Real *) &aVector1 ; | |
857 | Standard_Real * myVector2 = (Standard_Real *) &aVector2 ; | |
858 | NewParameter = (Parameter - CacheParameter) / SpanLenght ; | |
859 | PLib::EvalPolynomial(NewParameter, | |
860 | 2, | |
861 | Degree, | |
862 | Dimension_gen, | |
863 | PArray[0], | |
864 | LocalPDerivatives[0]) ; | |
865 | // | |
866 | // unormalize derivatives since those are computed normalized | |
867 | // | |
868 | Factor = 1.0e0/ SpanLenght ; | |
869 | Index = Dimension_gen ; | |
870 | EndIndex = Min (2, Degree) ; | |
871 | ||
872 | for (ii = 1 ; ii <= EndIndex ; ii++) { | |
873 | ModifyCoords (LocalPDerivatives + Index, *= Factor); | |
874 | Factor /= SpanLenght ; | |
875 | Index += Dimension_gen; | |
876 | } | |
877 | ||
878 | Index = (Degree + 1) * Dimension_gen; | |
879 | for (ii = Degree ; ii < 2 ; ii++) { | |
880 | NullifyCoords (LocalPDerivatives + Index); | |
881 | Index += Dimension_gen; | |
882 | } | |
883 | ||
0e14656b | 884 | if (WeightsArray != NULL) { |
885 | const TColStd_Array1OfReal& refWeights = *WeightsArray; | |
7fd59977 | 886 | Standard_Real * |
0e14656b | 887 | WArray = (Standard_Real *) &refWeights(refWeights.Lower()) ; |
7fd59977 | 888 | |
889 | PLib::EvalPolynomial(NewParameter, | |
890 | 2, | |
891 | Degree, | |
892 | 1, | |
893 | WArray[0], | |
894 | LocalWDerivatives[0]) ; | |
895 | ||
896 | for (ii = Degree + 1 ; ii <= 2 ; ii++) { | |
897 | LocalWDerivatives[ii] = 0.0e0 ; | |
898 | } | |
899 | // | |
900 | // unormalize the result since the polynomial stored in the cache | |
901 | // is normalized between 0 and 1 | |
902 | // | |
903 | Factor = 1.0e0 / SpanLenght ; | |
904 | ||
905 | for (ii = 1 ; ii <= EndIndex ; ii++) { | |
906 | LocalWDerivatives[ii] *= Factor ; | |
907 | Factor /= SpanLenght ; | |
908 | } | |
909 | PLib::RationalDerivatives(2, | |
910 | Dimension_gen, | |
911 | LocalPDerivatives[0], | |
912 | LocalWDerivatives[0], | |
913 | LocalPDerivatives[0]) ; | |
914 | } | |
915 | ||
916 | CopyCoords (myPoint, LocalPDerivatives); | |
917 | CopyCoords (myVector1, LocalPDerivatives + Dimension_gen); | |
918 | CopyCoords (myVector2, LocalPDerivatives + Dimension_gen * 2); | |
919 | } | |
920 | ||
921 | //======================================================================= | |
922 | //function : CacheD3 | |
923 | //purpose : Evaluates the polynomial cache of the Bspline Curve | |
924 | // | |
925 | //======================================================================= | |
926 | void BSplCLib::CacheD3(const Standard_Real Parameter, | |
927 | const Standard_Integer Degree, | |
928 | const Standard_Real CacheParameter, | |
929 | const Standard_Real SpanLenght, | |
930 | const Array1OfPoints& PolesArray, | |
0e14656b | 931 | const TColStd_Array1OfReal* WeightsArray, |
7fd59977 | 932 | Point& aPoint, |
933 | Vector& aVector1, | |
934 | Vector& aVector2, | |
935 | Vector& aVector3) | |
936 | { | |
937 | // | |
938 | // the CacheParameter is where the cache polynomial was evaluated in homogeneous | |
939 | // form | |
940 | // the SpanLenght is the normalizing factor so that the polynomial is between | |
941 | // 0 and 1 | |
942 | Standard_Integer ii, Index, EndIndex; | |
943 | Standard_Real LocalPDerivatives[Dimension_gen << 2] ; | |
944 | // Standard_Real LocalWDerivatives[4], Factor, NewParameter, Inverse ; | |
945 | Standard_Real LocalWDerivatives[4], Factor, NewParameter ; | |
946 | Standard_Real * PArray = (Standard_Real *) &(PolesArray(PolesArray.Lower())) ; | |
947 | Standard_Real * myPoint = (Standard_Real *) &aPoint ; | |
948 | Standard_Real * myVector1 = (Standard_Real *) &aVector1 ; | |
949 | Standard_Real * myVector2 = (Standard_Real *) &aVector2 ; | |
950 | Standard_Real * myVector3 = (Standard_Real *) &aVector3 ; | |
951 | NewParameter = (Parameter - CacheParameter) / SpanLenght ; | |
952 | PLib::EvalPolynomial(NewParameter, | |
953 | 3, | |
954 | Degree, | |
955 | Dimension_gen, | |
956 | PArray[0], | |
957 | LocalPDerivatives[0]) ; | |
958 | ||
959 | Index = (Degree + 1) * Dimension_gen; | |
960 | for (ii = Degree ; ii < 3 ; ii++) { | |
961 | NullifyCoords (LocalPDerivatives + Index); | |
962 | Index += Dimension_gen; | |
963 | } | |
964 | ||
965 | // | |
966 | // unormalize derivatives since those are computed normalized | |
967 | // | |
968 | Factor = 1.0e0 / SpanLenght ; | |
969 | Index = Dimension_gen ; | |
970 | EndIndex = Min(3,Degree) ; | |
971 | ||
972 | for (ii = 1 ; ii <= EndIndex ; ii++) { | |
973 | ModifyCoords (LocalPDerivatives + Index, *= Factor); | |
974 | Factor /= SpanLenght; | |
975 | Index += Dimension_gen; | |
976 | } | |
977 | ||
0e14656b | 978 | if (WeightsArray != NULL) { |
979 | const TColStd_Array1OfReal& refWeights = *WeightsArray; | |
7fd59977 | 980 | Standard_Real * |
0e14656b | 981 | WArray = (Standard_Real *) &refWeights(refWeights.Lower()) ; |
7fd59977 | 982 | |
983 | PLib::EvalPolynomial(NewParameter, | |
984 | 3, | |
985 | Degree, | |
986 | 1, | |
987 | WArray[0], | |
988 | LocalWDerivatives[0]) ; | |
989 | // | |
990 | // unormalize the result since the polynomial stored in the cache | |
991 | // is normalized between 0 and 1 | |
992 | // | |
993 | Factor = 1.0e0 / SpanLenght ; | |
994 | ||
995 | for (ii = 1 ; ii <= EndIndex ; ii++) { | |
996 | LocalWDerivatives[ii] *= Factor ; | |
997 | Factor /= SpanLenght ; | |
998 | } | |
999 | ||
1000 | for (ii = (Degree + 1) ; ii <= 3 ; ii++) { | |
1001 | LocalWDerivatives[ii] = 0.0e0 ; | |
1002 | } | |
1003 | PLib::RationalDerivatives(3, | |
1004 | Dimension_gen, | |
1005 | LocalPDerivatives[0], | |
1006 | LocalWDerivatives[0], | |
1007 | LocalPDerivatives[0]) ; | |
1008 | } | |
1009 | ||
1010 | CopyCoords (myPoint, LocalPDerivatives); | |
1011 | CopyCoords (myVector1, LocalPDerivatives + Dimension_gen); | |
1012 | CopyCoords (myVector2, LocalPDerivatives + Dimension_gen * 2); | |
1013 | CopyCoords (myVector3, LocalPDerivatives + Dimension_gen * 3); | |
1014 | } | |
1015 | ||
1016 | //======================================================================= | |
1017 | //function : BuildCache | |
1018 | //purpose : Stores theTaylor expansion normalized between 0,1 in the | |
1019 | // Cache : in case of a rational function the Poles are | |
1020 | // stored in homogeneous form | |
1021 | //======================================================================= | |
1022 | ||
1023 | void BSplCLib::BuildCache | |
1024 | (const Standard_Real U, | |
1025 | const Standard_Real SpanDomain, | |
1026 | const Standard_Boolean Periodic, | |
1027 | const Standard_Integer Degree, | |
1028 | const TColStd_Array1OfReal& FlatKnots, | |
1029 | const Array1OfPoints& Poles, | |
0e14656b | 1030 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 1031 | Array1OfPoints& CachePoles, |
0e14656b | 1032 | TColStd_Array1OfReal* CacheWeights) |
7fd59977 | 1033 | { |
1034 | Standard_Integer ii, | |
1035 | Dimension, | |
1036 | LocalIndex, | |
1037 | index = 0 ; | |
1038 | Standard_Real u = U, | |
1039 | LocalValue ; | |
1040 | Standard_Boolean rational; | |
1041 | ||
1042 | BSplCLib_DataContainer dc(Degree); | |
1043 | PrepareEval(u, | |
1044 | index, | |
1045 | Dimension, | |
1046 | rational, | |
1047 | Degree, | |
1048 | Periodic, | |
1049 | Poles, | |
1050 | Weights, | |
1051 | FlatKnots, | |
1052 | (BSplCLib::NoMults()), | |
1053 | dc); | |
1054 | // | |
1055 | // Watch out : PrepareEval checks if locally the Bspline is polynomial | |
1056 | // therefore rational flag could be set to False even if there are | |
1057 | // Weights. The Dimension is set accordingly. | |
1058 | // | |
1059 | ||
1060 | BSplCLib::Bohm(u,Degree,Degree,*dc.knots,Dimension,*dc.poles); | |
1061 | ||
1062 | LocalValue = 1.0e0 ; | |
1063 | LocalIndex = 0 ; | |
1064 | ||
1065 | if (rational) { | |
1066 | ||
1067 | for (ii = 1 ; ii <= Degree + 1 ; ii++) { | |
1068 | CoordsToPoint (CachePoles(ii), dc.poles + LocalIndex, * LocalValue); | |
1069 | LocalIndex += Dimension_gen + 1; | |
1070 | LocalValue *= SpanDomain / (Standard_Real) ii ; | |
1071 | } | |
1072 | ||
1073 | LocalIndex = Dimension_gen; | |
1074 | LocalValue = 1.0e0 ; | |
1075 | for (ii = 1 ; ii <= Degree + 1 ; ii++) { | |
0e14656b | 1076 | (*CacheWeights)(ii) = dc.poles[LocalIndex] * LocalValue ; |
7fd59977 | 1077 | LocalIndex += Dimension_gen + 1; |
1078 | LocalValue *= SpanDomain / (Standard_Real) ii ; | |
1079 | } | |
1080 | } | |
1081 | else { | |
1082 | ||
1083 | for (ii = 1 ; ii <= Degree + 1 ; ii++) { | |
1084 | CoordsToPoint (CachePoles(ii), dc.poles + LocalIndex, * LocalValue); | |
1085 | LocalIndex += Dimension_gen; | |
1086 | LocalValue *= SpanDomain / (Standard_Real) ii ; | |
1087 | } | |
1088 | ||
0e14656b | 1089 | if (Weights != NULL) { |
7fd59977 | 1090 | for (ii = 1 ; ii <= Degree + 1 ; ii++) |
0e14656b | 1091 | (*CacheWeights)(ii) = 0.0e0 ; |
1092 | (*CacheWeights)(1) = 1.0e0 ; | |
7fd59977 | 1093 | } |
1094 | } | |
1095 | } | |
1096 | ||
94f71cad | 1097 | void BSplCLib::BuildCache(const Standard_Real theParameter, |
1098 | const Standard_Real theSpanDomain, | |
1099 | const Standard_Boolean thePeriodicFlag, | |
1100 | const Standard_Integer theDegree, | |
a061150b | 1101 | const Standard_Integer theSpanIndex, |
94f71cad | 1102 | const TColStd_Array1OfReal& theFlatKnots, |
1103 | const Array1OfPoints& thePoles, | |
0e14656b | 1104 | const TColStd_Array1OfReal* theWeights, |
94f71cad | 1105 | TColStd_Array2OfReal& theCacheArray) |
1106 | { | |
1107 | Standard_Real aParam = theParameter; | |
a061150b | 1108 | Standard_Integer anIndex = theSpanIndex; |
94f71cad | 1109 | Standard_Integer aDimension; |
1110 | Standard_Boolean isRational; | |
1111 | ||
1112 | BSplCLib_DataContainer dc(theDegree); | |
1113 | PrepareEval(aParam, | |
1114 | anIndex, | |
1115 | aDimension, | |
1116 | isRational, | |
1117 | theDegree, | |
1118 | thePeriodicFlag, | |
1119 | thePoles, | |
1120 | theWeights, | |
1121 | theFlatKnots, | |
1122 | (BSplCLib::NoMults()), | |
1123 | dc); | |
1124 | // | |
1125 | // Watch out : PrepareEval checks if locally the Bspline is polynomial | |
1126 | // therefore rational flag could be set to False even if there are | |
1127 | // Weights. The Dimension is set accordingly. | |
1128 | // | |
1129 | ||
1130 | Standard_Integer aCacheShift = // helps to store weights when PrepareEval did not found that the curve is locally polynomial | |
0e14656b | 1131 | (theWeights != NULL && !isRational) ? aDimension + 1 : aDimension; |
94f71cad | 1132 | |
1133 | BSplCLib::Bohm(aParam, theDegree, theDegree, *dc.knots, aDimension, *dc.poles); | |
1134 | ||
1135 | Standard_Real aCoeff = 1.0; | |
1136 | Standard_Real* aCache = (Standard_Real *) &(theCacheArray(theCacheArray.LowerRow(), theCacheArray.LowerCol())); | |
1137 | Standard_Real* aPolyCoeffs = dc.poles; | |
1138 | ||
1139 | for (Standard_Integer i = 0; i <= theDegree; i++) | |
1140 | { | |
1141 | for (Standard_Integer j = 0; j< aDimension; j++) | |
1142 | aCache[j] = aPolyCoeffs[j] * aCoeff; | |
1143 | aCoeff *= theSpanDomain / (i+1); | |
1144 | aPolyCoeffs += aDimension; | |
1145 | aCache += aDimension; | |
1146 | if (aCacheShift > aDimension) | |
1147 | { | |
1148 | aCache[0] = 0.0; | |
1149 | aCache++; | |
1150 | } | |
1151 | } | |
1152 | ||
1153 | if (aCacheShift > aDimension) | |
1154 | theCacheArray.SetValue(theCacheArray.LowerRow(), theCacheArray.LowerCol() + aCacheShift - 1, 1.0); | |
1155 | } | |
1156 | ||
7fd59977 | 1157 | //======================================================================= |
1158 | //function : Interpolate | |
1159 | //purpose : | |
1160 | //======================================================================= | |
1161 | ||
1162 | void BSplCLib::Interpolate(const Standard_Integer Degree, | |
1163 | const TColStd_Array1OfReal& FlatKnots, | |
1164 | const TColStd_Array1OfReal& Parameters, | |
1165 | const TColStd_Array1OfInteger& ContactOrderArray, | |
1166 | Array1OfPoints& Poles, | |
1167 | Standard_Integer& InversionProblem) | |
1168 | ||
1169 | { | |
1170 | Standard_Real *PArray ; | |
1171 | ||
1172 | PArray = (Standard_Real *) &Poles(Poles.Lower()) ; | |
1173 | ||
1174 | BSplCLib::Interpolate(Degree, | |
1175 | FlatKnots, | |
1176 | Parameters, | |
1177 | ContactOrderArray, | |
1178 | Dimension_gen, | |
1179 | PArray[0], | |
1180 | InversionProblem) ; | |
1181 | } | |
1182 | ||
1183 | //======================================================================= | |
1184 | //function : Interpolate | |
1185 | //purpose : | |
1186 | //======================================================================= | |
1187 | ||
1188 | void BSplCLib::Interpolate(const Standard_Integer Degree, | |
1189 | const TColStd_Array1OfReal& FlatKnots, | |
1190 | const TColStd_Array1OfReal& Parameters, | |
1191 | const TColStd_Array1OfInteger& ContactOrderArray, | |
1192 | Array1OfPoints& Poles, | |
1193 | TColStd_Array1OfReal& Weights, | |
1194 | Standard_Integer& InversionProblem) | |
1195 | { | |
1196 | Standard_Real *PArray, | |
1197 | * WArray ; | |
1198 | PArray = (Standard_Real *) &Poles(Poles.Lower()) ; | |
1199 | WArray = (Standard_Real *) &Weights(Weights.Lower()) ; | |
1200 | BSplCLib::Interpolate(Degree, | |
1201 | FlatKnots, | |
1202 | Parameters, | |
1203 | ContactOrderArray, | |
1204 | Dimension_gen, | |
1205 | PArray[0], | |
1206 | WArray[0], | |
1207 | InversionProblem) ; | |
1208 | } | |
1209 | ||
1210 | //======================================================================= | |
1211 | //function : MovePoint | |
1212 | //purpose : Find the new poles which allows an old point (with a | |
1213 | // given u as parameter) to reach a new position | |
1214 | //======================================================================= | |
1215 | ||
1216 | void BSplCLib::MovePoint (const Standard_Real U, | |
1217 | const Vector& Displ, | |
1218 | const Standard_Integer Index1, | |
1219 | const Standard_Integer Index2, | |
1220 | const Standard_Integer Degree, | |
7fd59977 | 1221 | const Array1OfPoints& Poles, |
aff73fd5 | 1222 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 1223 | const TColStd_Array1OfReal& FlatKnots, |
1224 | Standard_Integer& FirstIndex, | |
1225 | Standard_Integer& LastIndex, | |
1226 | Array1OfPoints& NewPoles) | |
1227 | { | |
1228 | // calculate the BSplineBasis in the parameter U | |
1229 | Standard_Integer FirstNonZeroBsplineIndex; | |
1230 | math_Matrix BSplineBasis(1, 1, | |
1231 | 1, Degree+1); | |
1232 | Standard_Integer ErrorCode = | |
6143f12f | 1233 | BSplCLib::EvalBsplineBasis(0, |
7fd59977 | 1234 | Degree+1, |
1235 | FlatKnots, | |
1236 | U, | |
1237 | FirstNonZeroBsplineIndex, | |
1238 | BSplineBasis); | |
1239 | if (ErrorCode != 0) { | |
1240 | FirstIndex = 0; | |
1241 | LastIndex = 0; | |
1242 | ||
1243 | for (Standard_Integer i=Poles.Lower(); i<=Poles.Upper(); i++) { | |
1244 | NewPoles(i) = Poles(i); | |
1245 | } | |
1246 | return; | |
1247 | } | |
1248 | ||
1249 | // find the span which is predominant for parameter U | |
1250 | FirstIndex = FirstNonZeroBsplineIndex; | |
1251 | LastIndex = FirstNonZeroBsplineIndex + Degree ; | |
1252 | if (FirstIndex < Index1) FirstIndex = Index1; | |
1253 | if (LastIndex > Index2) LastIndex = Index2; | |
1254 | ||
1255 | Standard_Real maxValue = 0.0; | |
1256 | Standard_Integer i, kk1=0, kk2, ii; | |
1257 | ||
1258 | for (i = FirstIndex-FirstNonZeroBsplineIndex+1; | |
1259 | i <= LastIndex-FirstNonZeroBsplineIndex+1; i++) { | |
1260 | if (BSplineBasis(1,i) > maxValue) { | |
1261 | kk1 = i + FirstNonZeroBsplineIndex - 1; | |
1262 | maxValue = BSplineBasis(1,i); | |
1263 | } | |
1264 | } | |
1265 | ||
1266 | // find a kk2 if symetriy | |
1267 | kk2 = kk1; | |
1268 | i = kk1 - FirstNonZeroBsplineIndex + 2; | |
1269 | if ((kk1+1) <= LastIndex) { | |
1270 | if (Abs(BSplineBasis(1, kk1-FirstNonZeroBsplineIndex+2) - maxValue) < 1.e-10) { | |
1271 | kk2 = kk1+1; | |
1272 | } | |
1273 | } | |
1274 | ||
1275 | // compute the vector of displacement | |
1276 | Standard_Real D1 = 0.0; | |
1277 | Standard_Real D2 = 0.0; | |
1278 | Standard_Real hN, Coef, Dval; | |
1279 | ||
1280 | for (i = 1; i <= Degree+1; i++) { | |
1281 | ii = i + FirstNonZeroBsplineIndex - 1; | |
aff73fd5 | 1282 | if (Weights != NULL) { |
1283 | hN = Weights->Value(ii)*BSplineBasis(1, i); | |
7fd59977 | 1284 | D2 += hN; |
1285 | } | |
1286 | else { | |
1287 | hN = BSplineBasis(1, i); | |
1288 | } | |
1289 | if (ii >= FirstIndex && ii <= LastIndex) { | |
1290 | if (ii < kk1) { | |
1291 | Dval = kk1-ii; | |
1292 | } | |
1293 | else if (ii > kk2) { | |
1294 | Dval = ii - kk2; | |
1295 | } | |
1296 | else { | |
1297 | Dval = 0.0; | |
1298 | } | |
1299 | D1 += 1./(Dval + 1.) * hN; | |
1300 | } | |
1301 | } | |
1302 | ||
aff73fd5 | 1303 | if (Weights != NULL) { |
7fd59977 | 1304 | Coef = D2/D1; |
1305 | } | |
1306 | else { | |
1307 | Coef = 1./D1; | |
1308 | } | |
1309 | ||
1310 | // compute the new poles | |
1311 | ||
1312 | for (i=Poles.Lower(); i<=Poles.Upper(); i++) { | |
1313 | if (i >= FirstIndex && i <= LastIndex) { | |
1314 | if (i < kk1) { | |
1315 | Dval = kk1-i; | |
1316 | } | |
1317 | else if (i > kk2) { | |
1318 | Dval = i - kk2; | |
1319 | } | |
1320 | else { | |
1321 | Dval = 0.0; | |
1322 | } | |
1323 | NewPoles(i) = Poles(i).Translated((Coef/(Dval+1.))*Displ); | |
1324 | } | |
1325 | else { | |
1326 | NewPoles(i) = Poles(i); | |
1327 | } | |
1328 | } | |
1329 | } | |
1330 | ||
1331 | //======================================================================= | |
1332 | //function : MovePoint | |
1333 | //purpose : Find the new poles which allows an old point (with a | |
1334 | // given u as parameter) to reach a new position | |
1335 | //======================================================================= | |
1336 | ||
1337 | //======================================================================= | |
1338 | //function : MovePointAndTangent | |
1339 | //purpose : | |
1340 | //======================================================================= | |
1341 | ||
1342 | void BSplCLib::MovePointAndTangent (const Standard_Real U, | |
1343 | const Vector& Delta, | |
1344 | const Vector& DeltaDerivatives, | |
1345 | const Standard_Real Tolerance, | |
1346 | const Standard_Integer Degree, | |
7fd59977 | 1347 | const Standard_Integer StartingCondition, |
1348 | const Standard_Integer EndingCondition, | |
1349 | const Array1OfPoints& Poles, | |
aff73fd5 | 1350 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 1351 | const TColStd_Array1OfReal& FlatKnots, |
1352 | Array1OfPoints& NewPoles, | |
1353 | Standard_Integer & ErrorStatus) | |
1354 | { | |
1355 | Standard_Real *delta_array, | |
1356 | *delta_derivative_array, | |
1357 | *poles_array, | |
1358 | *new_poles_array ; | |
1359 | ||
1360 | Standard_Integer num_poles ; | |
1361 | num_poles = Poles.Length() ; | |
1362 | ||
1363 | if (NewPoles.Length() != num_poles) { | |
9775fa61 | 1364 | throw Standard_ConstructionError(); |
7fd59977 | 1365 | } |
1366 | delta_array = (Standard_Real *) &Delta ; | |
1367 | delta_derivative_array = (Standard_Real *)&DeltaDerivatives ; | |
1368 | poles_array = (Standard_Real *) &Poles(Poles.Lower()) ; | |
1369 | ||
1370 | new_poles_array = (Standard_Real *) &NewPoles(NewPoles.Lower()) ; | |
1371 | MovePointAndTangent(U, | |
1372 | Dimension_gen, | |
1373 | delta_array[0], | |
1374 | delta_derivative_array[0], | |
1375 | Tolerance, | |
1376 | Degree, | |
7fd59977 | 1377 | StartingCondition, |
1378 | EndingCondition, | |
1379 | poles_array[0], | |
1380 | Weights, | |
1381 | FlatKnots, | |
1382 | new_poles_array[0], | |
1383 | ErrorStatus) ; | |
1384 | } | |
1385 | ||
1386 | //======================================================================= | |
1387 | //function : Resolution | |
1388 | //purpose : | |
1389 | //======================================================================= | |
1390 | ||
1391 | void BSplCLib::Resolution(const Array1OfPoints& Poles, | |
0e14656b | 1392 | const TColStd_Array1OfReal* Weights, |
7fd59977 | 1393 | const Standard_Integer NumPoles, |
1394 | const TColStd_Array1OfReal& FlatKnots, | |
1395 | const Standard_Integer Degree, | |
1396 | const Standard_Real Tolerance3D, | |
1397 | Standard_Real& UTolerance) | |
1398 | { | |
1399 | Standard_Real *PolesArray ; | |
1400 | PolesArray = (Standard_Real *) &Poles(Poles.Lower()) ; | |
1401 | BSplCLib::Resolution(PolesArray[0], | |
1402 | Dimension_gen, | |
1403 | NumPoles, | |
1404 | Weights, | |
1405 | FlatKnots, | |
1406 | Degree, | |
1407 | Tolerance3D, | |
1408 | UTolerance) ; | |
1409 | } | |
1410 | ||
1411 | //======================================================================= | |
1412 | //function : FunctionMultiply | |
1413 | //purpose : | |
1414 | //======================================================================= | |
1415 | ||
1416 | void BSplCLib::FunctionMultiply | |
1417 | (const BSplCLib_EvaluatorFunction & FunctionPtr, | |
1418 | const Standard_Integer BSplineDegree, | |
1419 | const TColStd_Array1OfReal & BSplineFlatKnots, | |
1420 | const Array1OfPoints & Poles, | |
1421 | const TColStd_Array1OfReal & FlatKnots, | |
1422 | const Standard_Integer NewDegree, | |
1423 | Array1OfPoints & NewPoles, | |
9fd2d2c3 | 1424 | Standard_Integer & theStatus) |
7fd59977 | 1425 | { |
1426 | Standard_Integer num_bspline_poles = | |
1427 | BSplineFlatKnots.Length() - BSplineDegree - 1 ; | |
1428 | Standard_Integer num_new_poles = | |
1429 | FlatKnots.Length() - NewDegree - 1 ; | |
1430 | ||
1431 | if (Poles.Length() != num_bspline_poles || | |
1432 | NewPoles.Length() != num_new_poles) { | |
9775fa61 | 1433 | throw Standard_ConstructionError(); |
7fd59977 | 1434 | } |
1435 | Standard_Real * array_of_poles = | |
1436 | (Standard_Real *) &Poles(Poles.Lower()) ; | |
1437 | Standard_Real * array_of_new_poles = | |
1438 | (Standard_Real *) &NewPoles(NewPoles.Lower()) ; | |
1439 | BSplCLib::FunctionMultiply(FunctionPtr, | |
1440 | BSplineDegree, | |
1441 | BSplineFlatKnots, | |
1442 | Dimension_gen, | |
1443 | array_of_poles[0], | |
1444 | FlatKnots, | |
1445 | NewDegree, | |
1446 | array_of_new_poles[0], | |
9fd2d2c3 | 1447 | theStatus); |
7fd59977 | 1448 | } |
1449 | ||
1450 | //======================================================================= | |
1451 | //function : FunctionReparameterise | |
1452 | //purpose : | |
1453 | //======================================================================= | |
1454 | ||
1455 | void BSplCLib::FunctionReparameterise | |
1456 | (const BSplCLib_EvaluatorFunction & FunctionPtr, | |
1457 | const Standard_Integer BSplineDegree, | |
1458 | const TColStd_Array1OfReal & BSplineFlatKnots, | |
1459 | const Array1OfPoints & Poles, | |
1460 | const TColStd_Array1OfReal & FlatKnots, | |
1461 | const Standard_Integer NewDegree, | |
1462 | Array1OfPoints & NewPoles, | |
9fd2d2c3 | 1463 | Standard_Integer & theStatus) |
7fd59977 | 1464 | { |
1465 | Standard_Integer num_bspline_poles = | |
1466 | BSplineFlatKnots.Length() - BSplineDegree - 1 ; | |
1467 | Standard_Integer num_new_poles = | |
1468 | FlatKnots.Length() - NewDegree - 1 ; | |
1469 | ||
1470 | if (Poles.Length() != num_bspline_poles || | |
1471 | NewPoles.Length() != num_new_poles) { | |
9775fa61 | 1472 | throw Standard_ConstructionError(); |
7fd59977 | 1473 | } |
1474 | Standard_Real * array_of_poles = | |
1475 | (Standard_Real *) &Poles(Poles.Lower()) ; | |
1476 | Standard_Real * array_of_new_poles = | |
1477 | (Standard_Real *) &NewPoles(NewPoles.Lower()) ; | |
1478 | BSplCLib::FunctionReparameterise(FunctionPtr, | |
1479 | BSplineDegree, | |
1480 | BSplineFlatKnots, | |
1481 | Dimension_gen, | |
1482 | array_of_poles[0], | |
1483 | FlatKnots, | |
1484 | NewDegree, | |
1485 | array_of_new_poles[0], | |
9fd2d2c3 | 1486 | theStatus); |
7fd59977 | 1487 | } |