Library of Assembled Shared Sources
spline_bezier_path.inl
Go to the documentation of this file.
1/** @file
2 * @author Bram de Greve (bram@cocamware.com)
3 * @author Tom De Muer (tom@cocamware.com)
4 *
5 * Distributed under the terms of the GPL (GNU Public License)
6 *
7 * The LASS License:
8 *
9 * Copyright 2004-2006 Bram de Greve and Tom De Muer
10 *
11 * LASS is free software; you can redistribute it and/or modify
12 * it under the terms of the GNU General Public License as published by
13 * the Free Software Foundation; either version 2 of the License, or
14 * (at your option) any later version.
15 *
16 * This program is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 * GNU General Public License for more details.
20 *
21 * You should have received a copy of the GNU General Public License
22 * along with this program; if not, write to the Free Software
23 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
24 */
25
26
27
28#ifndef LASS_GUARDIAN_OF_INCLUSION_NUM_SPLINE_BEZIER_PATH_INL
29#define LASS_GUARDIAN_OF_INCLUSION_NUM_SPLINE_BEZIER_PATH_INL
30
31#include "num_common.h"
32#include "spline_bezier_path.h"
34
35namespace lass
36{
37namespace num
38{
39
40// --- public --------------------------------------------------------------------------------------
41
42template <typename S, typename D, typename T>
43SplineBezierPath<S, D, T>::SplineBezierPath()
44{
45}
46
47
48
49template <typename S, typename D, typename T>
50template <typename PairInputIterator>
51SplineBezierPath<S, D, T>::SplineBezierPath(
52 PairInputIterator first, PairInputIterator last)
53{
54 init(first, last, first->second);
55}
56
57
58
59template <typename S, typename D, typename T>
60template <typename ScalarInputIterator, typename DataInputIterator>
61SplineBezierPath<S, D, T>::SplineBezierPath(ScalarInputIterator firstControl,
62 ScalarInputIterator lastControl,
63 DataInputIterator firstData)
64{
65 init(firstControl, lastControl, firstData, *firstData);
66}
67
68
69
70template <typename S, typename D, typename T>
71const typename SplineBezierPath<S, D, T>::TData
73{
74 LASS_ASSERT(!isEmpty());
75
76 const typename TNodes::const_iterator first = findNode(x);
77 LASS_ASSERT(first != nodes_.end());
78 const typename TNodes::const_iterator second = stde::next(first);
79 LASS_ASSERT(second != nodes_.end());
80
81 const TData& a = first->knot();
82 const TData& b = first->right();
83 const TData& c = second->left();
84 const TData& d = second->knot();
85
86 const TScalar t = (x - first->x) / (second->x - first->x);
87 const TScalar s = 1 - t;
88
89 TData y(a);
90 TDataTraits::scale(y, s * s * s);
91 TDataTraits::multiplyAccumulate(y, b, 3 * s * s * t);
92 TDataTraits::multiplyAccumulate(y, c, 3 * s * t * t);
93 TDataTraits::multiplyAccumulate(y, d, t * t * t);
94 return y;
95
96}
97
98
99
100template <typename S, typename D, typename T>
101const typename SplineBezierPath<S, D, T>::TData
102SplineBezierPath<S, D, T>::derivative(TScalar x) const
103{
104 LASS_ASSERT(!isEmpty());
105
106 const typename TNodes::const_iterator first = findNode(x);
107 LASS_ASSERT(first != nodes_.end());
108 const typename TNodes::const_iterator second = stde::next(first);
109 LASS_ASSERT(second != nodes_.end());
110
111 const TData& a = first->knot();
112 const TData& b = first->right();
113 const TData& c = second->left();
114 const TData& d = second->knot();
115 const TScalar t = (x - first->x) / (second->x - first->x);
116 const TScalar s = 1 - t;
117
118 TData dy = a;
119 TDataTraits::scale(dy, -3 * s * s);
120 TDataTraits::multiplyAccumulate(dy, b, 3 * s * (s - 2 * t));
121 TDataTraits::multiplyAccumulate(dy, c, 3 * t * (2 * s - t));
122 TDataTraits::multiplyAccumulate(dy, d, 3 * t * t);
123 TDataTraits::scale(dy, num::inv(second->x - first->x));
124 return dy;
125}
126
127
128
129template <typename S, typename D, typename T>
130const typename SplineBezierPath<S, D, T>::TData
131SplineBezierPath<S, D, T>::derivative2(TScalar x) const
132{
133 LASS_ASSERT(!isEmpty());
134
135 const typename TNodes::const_iterator first = findNode(x);
136 LASS_ASSERT(first != nodes_.end());
137 const typename TNodes::const_iterator second = stde::next(first);
138 LASS_ASSERT(second != nodes_.end());
139
140 const TData& a = first->knot();
141 const TData& b = first->right();
142 const TData& c = second->left();
143 const TData& d = second->knot();
144 const TScalar dx = second->x - first->x;
145 const TScalar t = (x - first->x) / dx;
146
147 // (6-6*t)*A+(18*t-12)*B+(-18*t+6)*C+6*t*D
148
149 TData dy = a;
150 TDataTraits::scale(dy, 6 - 6 * t);
151 TDataTraits::multiplyAccumulate(dy, b, 18 * t - 12);
152 TDataTraits::multiplyAccumulate(dy, c, -18 * t + 6);
153 TDataTraits::multiplyAccumulate(dy, d, 6 * t);
154 TDataTraits::scale(dy, num::inv(num::sqr(dx)));
155 return dy;
156}
157
158
159
160template <typename S, typename D, typename T>
161const typename SplineBezierPath<S, D, T>::TData
162SplineBezierPath<S, D, T>::integral(TScalar begin, TScalar end) const
163{
164 // -1/4*(1-t)^4*A+3*B*(1/4*t^4-2/3*t^3+1/2*t^2)+3*C*(-1/4*t^4+1/3*t^3)+1/4*t^4*D
165
166 LASS_ASSERT(!isEmpty());
167
168 typename TNodes::const_iterator first = findNode(begin);
169 LASS_ASSERT(first != nodes_.end());
170 typename TNodes::const_iterator last = findNode(end);
171 LASS_ASSERT(last != nodes_.end());
172
173 if (first == last)
174 {
175 const typename TNodes::const_iterator second = stde::next(first);
176 LASS_ASSERT(second != nodes_.end());
177
178 const TData& a = first->knot();
179 const TData& b = first->right();
180 const TData& c = second->left();
181 const TData& d = second->knot();
182
183 const TScalar dx = second->x - first->x;
184 const TScalar invDx = num::inv(dx);
185 const TScalar s = invDx * (begin - first->x);
186 const TScalar t = invDx * (end - first->x);
187
188 const TScalar st = t - s;
189 const TScalar st2 = num::sqr(t) - num::sqr(s);
190 const TScalar st3 = num::cubic(t) - num::cubic(s);
191 const TScalar st4 = num::sqr(num::sqr(t)) - num::sqr(num::sqr(s));
192 //const TScalar s4 = num::sqr(num::sqr(1 - end)) - num::sqr(num::sqr(1 - begin));
193
194 TData inty(a);
195 TDataTraits::scale(inty, dx * (-.25f * st4 + st3 - 1.5f * st2 + st));
196 TDataTraits::multiplyAccumulate(inty, b, dx * (.75f * st4 - 2 * st3 + 1.5f * st2));
197 TDataTraits::multiplyAccumulate(inty, c, dx * (-.75f * st4 + st3));
198 TDataTraits::multiplyAccumulate(inty, d, dx * .25f * st4);
199 TDataTraits::scale(inty, second->x - first->x);
200 return inty;
201 }
202 else
203 {
204 TScalar multiplier = 1;
205 if (end < begin)
206 {
207 std::swap(begin, end);
208 std::swap(first, last);
209 multiplier = -1;
210 }
211
212 typename TNodes::const_iterator second = stde::next(first);
213 LASS_ASSERT(second != nodes_.end());
214
215 TData inty(first->knot());
216 {
217 const TScalar dx = second->x - first->x;
218 const TScalar s = (begin - first->x) / dx;
219 const TScalar s2 = num::sqr(s);
220 const TScalar s3 = num::cubic(s);
221 const TScalar s4 = num::sqr(s2);
222
223 TDataTraits::scale(inty, dx * (.25f + .25f * s4 - s3 + 1.5f * s2 - s));
224 TDataTraits::multiplyAccumulate(inty, first->right(), dx * (.25f - .75f * s4 + 2 * s3 - 1.5f * s2));
225 TDataTraits::multiplyAccumulate(inty, second->left(), dx * (.25f + .75f * s4 - s3));
226 TDataTraits::multiplyAccumulate(inty, second->knot(), dx * .25f * (1 - s4));
227 }
228
229 first = second;
230 second += 1;
231 LASS_ASSERT(second != nodes_.end());
232
233 while (first != last)
234 {
235 const TScalar dx = second->x - first->x;
236 TDataTraits::multiplyAccumulate(inty, first->knot(), dx * .25f);
237 TDataTraits::multiplyAccumulate(inty, first->right(), dx * .25f);
238 TDataTraits::multiplyAccumulate(inty, second->left(), dx * .25f);
239 TDataTraits::multiplyAccumulate(inty, second->knot(), dx * .25f);
240 first = second;
241 second += 1;
242 LASS_ASSERT(second != nodes_.end());
243 }
244
245 {
246 const TScalar dx = second->x - first->x;
247 const TScalar t = (end - first->x) / dx;
248 const TScalar t2 = num::sqr(t);
249 const TScalar t3 = num::cubic(t);
250 const TScalar t4 = num::sqr(t2);
251
252 TDataTraits::multiplyAccumulate(inty, first->knot(),
253 dx * (-.25f * t4 + t3 - 1.5f * t2 + t));
254 TDataTraits::multiplyAccumulate(inty, first->right(),
255 dx * (.75f * t4 - 2 * t3 + 1.5f * t2));
256 TDataTraits::multiplyAccumulate(inty, second->left(),
257 dx * (-.75f * t4 + t3));
258 TDataTraits::multiplyAccumulate(inty, second->knot(),
259 dx * .25f * t4);
260 }
261
262 TDataTraits::scale(inty, multiplier);
263 return inty;
264 }
265}
266
267
268
269/** return true if the spline contains any nodes at all.
270 * @par complexity:
271 * O(1)
272 */
273template <typename S, typename D, typename T>
275{
276 return nodes_.empty();
277}
278
279
280
281/** return the range of control values for which the spline can interpolate.
282 * @par complexity:
283 * O(1)
284 */
285template <typename S, typename D, typename T>
286const typename SplineBezierPath<S, D, T>::TControlRange
288{
289 return TControlRange(nodes_.front().x, nodes_.back().x);
290}
291
292
293
294// --- private -------------------------------------------------------------------------------------
295
296template <typename S, typename D, typename T>
297template <typename PairInputIterator>
298void SplineBezierPath<S, D, T>::init(
299 PairInputIterator first, PairInputIterator last, const DataTriplet& /* dummy */)
300{
301 while (first != last)
302 {
303 Node node;
304 node.x = first->first;
305 node.triplet = first->second;
306 nodes_.push_back(node);
307 ++first;
308 }
309 finalInit();
310}
311
312template <typename S, typename D, typename T>
313template <typename PairInputIterator>
314void SplineBezierPath<S, D, T>::init(
315 PairInputIterator first, PairInputIterator last, const TData& /* dummy */)
316{
317 TSimpleNodes simpleNodes;
318 std::copy(first, last, std::back_inserter(simpleNodes));
319 makeFullNodes(simpleNodes);
320 finalInit();
321}
322
323template <typename S, typename D, typename T>
324template <typename ScalarInputIterator, typename DataInputIterator>
325void SplineBezierPath<S, D, T>::init(
326 ScalarInputIterator firstControl, ScalarInputIterator lastControl,
327 DataInputIterator firstData, const DataTriplet& /* dummy */)
328{
329 while (firstControl != lastControl)
330 {
331 Node node;
332 node.x = *firstControl++;
333 node.triplet = *firstData++;
334 nodes_.push_back(node);
335 }
336 finalInit();
337}
338
339template <typename S, typename D, typename T>
340template <typename ScalarInputIterator, typename DataInputIterator>
341void SplineBezierPath<S, D, T>::init(
342 ScalarInputIterator firstControl, ScalarInputIterator lastControl,
343 DataInputIterator firstData, const TData& /* dummy */)
344{
345 TSimpleNodes simpleNodes;
346 while (firstControl != lastControl)
347 {
348 simpleNodes.push_back(std::make_pair(*firstControl++, *firstData++));
349 }
350 makeFullNodes(simpleNodes);
351 finalInit();
352}
353
354template <typename S, typename D, typename T>
355void SplineBezierPath<S, D, T>::makeFullNodes(const TSimpleNodes& simpleNodes)
356{
357 const TScalar dt = TScalar(1) / 3;
358 const typename TNodes::size_type size = simpleNodes.size();
359
360 TNodes nodes;
361 switch (size)
362 {
363 case 0:
364 break;
365
366 case 1:
367 {
368 const TData& knot = simpleNodes[0].second;
369 TData null;
370 TDataTraits::zero(null, TDataTraits::dimension(knot));
371 nodes.push_back(Node(DataTriplet(null, knot, null), simpleNodes[0].first));
372 }
373 break;
374
375 default:
376 {
377 const TData& knot = simpleNodes[0].second;
378 TData dy = simpleNodes[1].second;
379 TDataTraits::multiplyAccumulate(dy, knot, -1);
380 TData right = knot;
381 TDataTraits::multiplyAccumulate(right, dy, dt);
382 TData left = knot;
383 TDataTraits::multiplyAccumulate(left, dy, -dt);
384 nodes.push_back(Node(DataTriplet(left, knot, right), simpleNodes[0].first));
385 }
386
387 for (size_t i = 1; i < size - 1; ++i)
388 {
389 const TData& knot = simpleNodes[i].second;
390 TData dy = simpleNodes[i + 1].second;
391 TDataTraits::multiplyAccumulate(dy, simpleNodes[i - 1].second, -1);
392 TDataTraits::scale(dy, .5);
393 TData right = knot;
394 TDataTraits::multiplyAccumulate(right, dy, dt);
395 TData left = knot;
396 TDataTraits::multiplyAccumulate(left, dy, -dt);
397 nodes.push_back(Node(DataTriplet(left, knot, right), simpleNodes[i].first));
398 }
399
400 {
401 const TData& knot = simpleNodes[size - 1].second;
402 TData dy = knot;
403 TDataTraits::multiplyAccumulate(dy, simpleNodes[size - 2].second, -1);
404 TData right = knot;
405 TDataTraits::multiplyAccumulate(right, dy, dt);
406 TData left = knot;
407 TDataTraits::multiplyAccumulate(left, dy, -dt);
408 nodes.push_back(Node(DataTriplet(left, knot, right), simpleNodes[size - 1].first));
409 }
410 }
411
412 nodes_.swap(nodes);
413}
414
415template <typename S, typename D, typename T>
416void SplineBezierPath<S, D, T>::finalInit()
417{
418 // are there any nodes at all? we need at least one!
419 //
420 const typename TNodes::size_type size = nodes_.size();
421 if (size < 1)
422 {
423 LASS_THROW("A bezier path interpolator needs at least one node!");
424 }
425
426 typename TNodes::iterator prev = nodes_.begin();
427 dataDimension_ = TDataTraits::dimension(prev->knot());
428 const typename TNodes::iterator end = nodes_.end();
429 for (typename TNodes::iterator i = stde::next(nodes_.begin()); i != end; prev = i++)
430 {
431 // check for increasing control components
432 //
433 if (prev->x >= i->x)
434 {
435 LASS_THROW("Nodes in bezier path interpolator must have absolutely increasing control "
436 "values.");
437 }
438
439 // check dimension of data
440 //
441 if (TDataTraits::dimension(i->left()) != dataDimension_ ||
442 TDataTraits::dimension(i->knot()) != dataDimension_ ||
443 TDataTraits::dimension(i->right()) != dataDimension_)
444 {
445 LASS_THROW("All data in linear interpolator must be of same dimension.");
446 }
447 }
448}
449
450
451
452/** binary search to find node that belongs to x
453 *
454 * @return
455 * @arg the iterator to node such that @a x is in the interval [@c node->x, @c next(node)->x)
456 * @arg nodes_.begin() @a x is smaller than @c nodes_.front().x</tt>
457 * @arg prev(nodes_.end(), 2) if @a x is greater than @c nodes_.back().x
458 *
459 * complexity: O(ln N)
460 */
461template <typename S, typename D, typename T>
462const typename SplineBezierPath<S, D, T>::TNodeConstIterator
463SplineBezierPath<S, D, T>::findNode(TScalar x) const
464{
465 LASS_ASSERT(nodes_.size() >= 2);
466
467 // early outs
468 //
469 if (x < nodes_.front().x)
470 {
471 return nodes_.begin();
472 }
473 if (x >= nodes_.back().x)
474 {
475 return stde::prev(nodes_.end(), 2);
476 }
477
478 // binary search
479 //
480 TNodeConstIterator first = nodes_.begin();
481 TNodeConstIterator last = nodes_.end();
482 while (stde::next(first) != last)
483 {
484 TNodeConstIterator middle = stde::next(first, std::distance(first, last) / 2);
485 LASS_ASSERT(middle != first && middle != last);
486
487 if (middle->x <= x)
488 {
489 first = middle;
490 }
491 else
492 {
493 last = middle;
494 }
495 LASS_ASSERT(first != last);
496 }
497
498 LASS_ASSERT(first->x <= x && last->x > x);
499 return first;
500}
501
502
503}
504
505}
506
507#endif
508
509// EOF
connects the data nodes with cubic bezier splines.
bool isEmpty() const override
return true if the spline contains any nodes at all.
const TControlRange controlRange() const override
return the range of control values for which the spline can interpolate.
T sqr(const T &x)
return x * x
Definition basic_ops.h:162
T inv(const T &x)
return x ^ -1
Definition basic_ops.h:178
T cubic(const T &x)
return x * x * x
Definition basic_ops.h:170
numeric types and traits.
Definition basic_ops.h:70
Library for Assembled Shared Sources.
Definition config.h:53