Point Cloud Library (PCL)  1.7.1
bspline_data.hpp
1 /*
2 Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
3 All rights reserved.
4 
5 Redistribution and use in source and binary forms, with or without modification,
6 are permitted provided that the following conditions are met:
7 
8 Redistributions of source code must retain the above copyright notice, this list of
9 conditions and the following disclaimer. Redistributions in binary form must reproduce
10 the above copyright notice, this list of conditions and the following disclaimer
11 in the documentation and/or other materials provided with the distribution.
12 
13 Neither the name of the Johns Hopkins University nor the names of its contributors
14 may be used to endorse or promote products derived from this software without specific
15 prior written permission.
16 
17 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
18 EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
19 OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
20 SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
21 INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
22 TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
23 BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
24 CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
25 ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
26 DAMAGE.
27 */
28 
29 
30 namespace pcl
31 {
32  namespace poisson
33  {
34 
35  /////////////////
36  // BSplineData //
37  /////////////////
38  // Support[i]:
39  // Odd: i +/- 0.5 * ( 1 + Degree )
40  // i - 0.5 * ( 1 + Degree ) < 0
41  // <=> i < 0.5 * ( 1 + Degree )
42  // i + 0.5 * ( 1 + Degree ) > 0
43  // <=> i > - 0.5 * ( 1 + Degree )
44  // i + 0.5 * ( 1 + Degree ) > r
45  // <=> i > r - 0.5 * ( 1 + Degree )
46  // i - 0.5 * ( 1 + Degree ) < r
47  // <=> i < r + 0.5 * ( 1 + Degree )
48  // Even: i + 0.5 +/- 0.5 * ( 1 + Degree )
49  // i - 0.5 * Degree < 0
50  // <=> i < 0.5 * Degree
51  // i + 1 + 0.5 * Degree > 0
52  // <=> i > -1 - 0.5 * Degree
53  // i + 1 + 0.5 * Degree > r
54  // <=> i > r - 1 - 0.5 * Degree
55  // i - 0.5 * Degree < r
56  // <=> i < r + 0.5 * Degree
57  template< int Degree > inline bool LeftOverlap( unsigned int depth , int offset )
58  {
59  offset <<= 1;
60  if( Degree & 1 ) return (offset < 1+Degree) && (offset > -1-Degree );
61  else return (offset < Degree) && (offset > -2-Degree );
62  }
63  template< int Degree > inline bool RightOverlap( unsigned int depth , int offset )
64  {
65  offset <<= 1;
66  int r = 1<<(depth+1);
67  if( Degree & 1 ) return (offset > 2-1-Degree) && (offset < 2+1+Degree );
68  else return (offset > 2-2-Degree) && (offset < 2+ Degree );
69  }
70  template< int Degree > inline int ReflectLeft( unsigned int depth , int offset )
71  {
72  if( Degree&1 ) return -offset;
73  else return -1-offset;
74  }
75  template< int Degree > inline int ReflectRight( unsigned int depth , int offset )
76  {
77  int r = 1<<(depth+1);
78  if( Degree&1 ) return r -offset;
79  else return r-1-offset;
80  }
81 
82  template< int Degree , class Real >
84  {
85  vvDotTable = dvDotTable = ddDotTable = NULL;
86  valueTables = dValueTables = NULL;
87  functionCount = sampleCount = 0;
88  }
89 
90  template< int Degree , class Real >
92  {
93  if( functionCount )
94  {
95  if( vvDotTable ) delete[] vvDotTable;
96  if( dvDotTable ) delete[] dvDotTable;
97  if( ddDotTable ) delete[] ddDotTable;
98 
99  if( valueTables ) delete[] valueTables;
100  if( dValueTables ) delete[] dValueTables;
101  }
102  vvDotTable = dvDotTable = ddDotTable = NULL;
103  valueTables = dValueTables=NULL;
104  functionCount = 0;
105  }
106 
107  template<int Degree,class Real>
108  void BSplineData<Degree,Real>::set( int maxDepth , bool useDotRatios , bool reflectBoundary )
109  {
110  this->useDotRatios = useDotRatios;
111  this->reflectBoundary = reflectBoundary;
112 
113  depth = maxDepth;
114  // [Warning] This assumes that the functions spacing is dual
115  functionCount = BinaryNode< double >::CumulativeCenterCount( depth );
117  baseFunctions = new PPolynomial<Degree>[functionCount];
118  baseBSplines = new BSplineComponents[functionCount];
119 
120  baseFunction = PPolynomial< Degree >::BSpline();
121  for( int i=0 ; i<=Degree ; i++ ) baseBSpline[i] = Polynomial< Degree >::BSplineComponent( i ).shift( double(-(Degree+1)/2) + i - 0.5 );
122  dBaseFunction = baseFunction.derivative();
123  StartingPolynomial< Degree > sPolys[Degree+3];
124 
125  for( int i=0 ; i<Degree+3 ; i++ )
126  {
127  sPolys[i].start = double(-(Degree+1)/2) + i - 1.5;
128  sPolys[i].p *= 0;
129  if( i<=Degree ) sPolys[i].p += baseBSpline[i ].shift( -1 );
130  if( i>=1 && i<=Degree+1 ) sPolys[i].p += baseBSpline[i-1];
131  for( int j=0 ; j<i ; j++ ) sPolys[i].p -= sPolys[j].p;
132  }
133  leftBaseFunction.set( sPolys , Degree+3 );
134  for( int i=0 ; i<Degree+3 ; i++ )
135  {
136  sPolys[i].start = double(-(Degree+1)/2) + i - 0.5;
137  sPolys[i].p *= 0;
138  if( i<=Degree ) sPolys[i].p += baseBSpline[i ];
139  if( i>=1 && i<=Degree+1 ) sPolys[i].p += baseBSpline[i-1].shift( 1 );
140  for( int j=0 ; j<i ; j++ ) sPolys[i].p -= sPolys[j].p;
141  }
142  rightBaseFunction.set( sPolys , Degree+3 );
143  dLeftBaseFunction = leftBaseFunction.derivative();
144  dRightBaseFunction = rightBaseFunction.derivative();
145  leftBSpline = rightBSpline = baseBSpline;
146  leftBSpline [1] += leftBSpline[2].shift( -1 ) , leftBSpline[0] *= 0;
147  rightBSpline[1] += rightBSpline[0].shift( 1 ) , rightBSpline[2] *= 0;
148  double c , w;
149  for( int i=0 ; i<functionCount ; i++ )
150  {
152  baseFunctions[i] = baseFunction.scale(w).shift(c);
153  baseBSplines[i] = baseBSpline.scale(w).shift(c);
154  if( reflectBoundary )
155  {
156  int d , off , r;
158  r = 1<<d;
159  if ( off==0 ) baseFunctions[i] = leftBaseFunction.scale(w).shift(c);
160  else if( off==r-1 ) baseFunctions[i] = rightBaseFunction.scale(w).shift(c);
161  if ( off==0 ) baseBSplines[i] = leftBSpline.scale(w).shift(c);
162  else if( off==r-1 ) baseBSplines[i] = rightBSpline.scale(w).shift(c);
163  }
164  }
165  }
166  template<int Degree,class Real>
168  {
169  clearDotTables( flags );
170  int size = ( functionCount*functionCount + functionCount )>>1;
171  int fullSize = functionCount*functionCount;
172  if( flags & VV_DOT_FLAG )
173  {
174  vvDotTable = new Real[size];
175  memset( vvDotTable , 0 , sizeof(Real)*size );
176  }
177  if( flags & DV_DOT_FLAG )
178  {
179  dvDotTable = new Real[fullSize];
180  memset( dvDotTable , 0 , sizeof(Real)*fullSize );
181  }
182  if( flags & DD_DOT_FLAG )
183  {
184  ddDotTable = new Real[size];
185  memset( ddDotTable , 0 , sizeof(Real)*size );
186  }
187  double vvIntegrals[Degree+1][Degree+1];
188  double vdIntegrals[Degree+1][Degree ];
189  double dvIntegrals[Degree ][Degree+1];
190  double ddIntegrals[Degree ][Degree ];
191  int vvSums[Degree+1][Degree+1];
192  int vdSums[Degree+1][Degree ];
193  int dvSums[Degree ][Degree+1];
194  int ddSums[Degree ][Degree ];
195  SetBSplineElementIntegrals< Degree , Degree >( vvIntegrals );
196  SetBSplineElementIntegrals< Degree , Degree-1 >( vdIntegrals );
197  SetBSplineElementIntegrals< Degree-1 , Degree >( dvIntegrals );
198  SetBSplineElementIntegrals< Degree-1 , Degree-1 >( ddIntegrals );
199 
200  for( int d1=0 ; d1<=depth ; d1++ )
201  for( int off1=0 ; off1<(1<<d1) ; off1++ )
202  {
203  int ii = BinaryNode< Real >::CenterIndex( d1 , off1 );
205  BSplineElements< Degree-1 > db1;
206  b1.differentiate( db1 );
207 
208  int start1 , end1;
209 
210  start1 = -1;
211  for( int i=0 ; i<int(b1.size()) ; i++ ) for( int j=0 ; j<=Degree ; j++ )
212  {
213  if( b1[i][j] && start1==-1 ) start1 = i;
214  if( b1[i][j] ) end1 = i+1;
215  }
216  for( int d2=d1 ; d2<=depth ; d2++ )
217  {
218  for( int off2=0 ; off2<(1<<d2) ; off2++ )
219  {
220  int start2 = off2-Degree;
221  int end2 = off2+Degree+1;
222  if( start2>=end1 || start1>=end2 ) continue;
223  start2 = std::max< int >( start1 , start2 );
224  end2 = std::min< int >( end1 , end2 );
225  if( d1==d2 && off2<off1 ) continue;
226  int jj = BinaryNode< Real >::CenterIndex( d2 , off2 );
227  BSplineElements< Degree > b2( 1<<d2 , off2 , reflectBoundary ? BSplineElements<Degree>::NEUMANN : BSplineElements< Degree>::NONE );
228  BSplineElements< Degree-1 > db2;
229  b2.differentiate( db2 );
230 
231  int idx = SymmetricIndex( ii , jj );
232  int idx1 = Index( ii , jj ) , idx2 = Index( jj , ii );
233 
234  memset( vvSums , 0 , sizeof( int ) * ( Degree+1 ) * ( Degree+1 ) );
235  memset( vdSums , 0 , sizeof( int ) * ( Degree+1 ) * ( Degree ) );
236  memset( dvSums , 0 , sizeof( int ) * ( Degree ) * ( Degree+1 ) );
237  memset( ddSums , 0 , sizeof( int ) * ( Degree ) * ( Degree ) );
238  for( int i=start2 ; i<end2 ; i++ )
239  {
240  for( int j=0 ; j<=Degree ; j++ ) for( int k=0 ; k<=Degree ; k++ ) vvSums[j][k] += b1[i][j] * b2[i][k];
241  for( int j=0 ; j<=Degree ; j++ ) for( int k=0 ; k< Degree ; k++ ) vdSums[j][k] += b1[i][j] * db2[i][k];
242  for( int j=0 ; j< Degree ; j++ ) for( int k=0 ; k<=Degree ; k++ ) dvSums[j][k] += db1[i][j] * b2[i][k];
243  for( int j=0 ; j< Degree ; j++ ) for( int k=0 ; k< Degree ; k++ ) ddSums[j][k] += db1[i][j] * db2[i][k];
244  }
245  double vvDot = 0 , dvDot = 0 , vdDot = 0 , ddDot = 0;
246  for( int j=0 ; j<=Degree ; j++ ) for( int k=0 ; k<=Degree ; k++ ) vvDot += vvIntegrals[j][k] * vvSums[j][k];
247  for( int j=0 ; j<=Degree ; j++ ) for( int k=0 ; k< Degree ; k++ ) vdDot += vdIntegrals[j][k] * vdSums[j][k];
248  for( int j=0 ; j< Degree ; j++ ) for( int k=0 ; k<=Degree ; k++ ) dvDot += dvIntegrals[j][k] * dvSums[j][k];
249  for( int j=0 ; j< Degree ; j++ ) for( int k=0 ; k< Degree ; k++ ) ddDot += ddIntegrals[j][k] * ddSums[j][k];
250  vvDot /= (1<<d2);
251  ddDot *= (1<<d2);
252  vvDot /= ( b1.denominator * b2.denominator );
253  dvDot /= ( b1.denominator * b2.denominator );
254  vdDot /= ( b1.denominator * b2.denominator );
255  ddDot /= ( b1.denominator * b2.denominator );
256  if( fabs(vvDot)<1e-15 ) continue;
257  if( flags & VV_DOT_FLAG ) vvDotTable [idx] = Real( vvDot );
258  if( useDotRatios )
259  {
260  if( flags & DV_DOT_FLAG ) dvDotTable[idx1] = Real( dvDot / vvDot );
261  if( flags & DV_DOT_FLAG ) dvDotTable[idx2] = Real( vdDot / vvDot );
262  if( flags & DD_DOT_FLAG ) ddDotTable[idx ] = Real( ddDot / vvDot );
263  }
264  else
265  {
266  if( flags & DV_DOT_FLAG ) dvDotTable[idx1] = Real( dvDot );
267  if( flags & DV_DOT_FLAG ) dvDotTable[idx2] = Real( dvDot );
268  if( flags & DD_DOT_FLAG ) ddDotTable[idx ] = Real( ddDot );
269  }
270  }
272  b = b1;
273  b.upSample( b1 );
274  b1.differentiate( db1 );
275  start1 = -1;
276  for( int i=0 ; i<int(b1.size()) ; i++ ) for( int j=0 ; j<=Degree ; j++ )
277  {
278  if( b1[i][j] && start1==-1 ) start1 = i;
279  if( b1[i][j] ) end1 = i+1;
280  }
281  }
282  }
283  }
284  template<int Degree,class Real>
286  {
287  if( (flags & VV_DOT_FLAG) && vvDotTable ) delete[] vvDotTable , vvDotTable = NULL;
288  if( (flags & DV_DOT_FLAG) && dvDotTable ) delete[] dvDotTable , dvDotTable = NULL;
289  if( (flags & DD_DOT_FLAG) && ddDotTable ) delete[] ddDotTable , ddDotTable = NULL;
290  }
291  template< int Degree , class Real >
292  void BSplineData< Degree , Real >::setSampleSpan( int idx , int& start , int& end , double smooth ) const
293  {
294  int d , off , res;
295  BinaryNode< double >::DepthAndOffset( idx , d , off );
296  res = 1<<d;
297  double _start = ( off + 0.5 - 0.5*(Degree+1) ) / res - smooth;
298  double _end = ( off + 0.5 + 0.5*(Degree+1) ) / res + smooth;
299  // (start)/(sampleCount-1) >_start && (start-1)/(sampleCount-1)<=_start
300  // => start > _start * (sampleCount-1 ) && start <= _start*(sampleCount-1) + 1
301  // => _start * (sampleCount-1) + 1 >= start > _start * (sampleCount-1)
302  start = int( floor( _start * (sampleCount-1) + 1 ) );
303  if( start<0 ) start = 0;
304  // (end)/(sampleCount-1)<_end && (end+1)/(sampleCount-1)>=_end
305  // => end < _end * (sampleCount-1 ) && end >= _end*(sampleCount-1) - 1
306  // => _end * (sampleCount-1) > end >= _end * (sampleCount-1) - 1
307  end = int( ceil( _end * (sampleCount-1) - 1 ) );
308  if( end>=sampleCount ) end = sampleCount-1;
309  }
310  template<int Degree,class Real>
311  void BSplineData<Degree,Real>::setValueTables( int flags , double smooth )
312  {
313  clearValueTables();
314  if( flags & VALUE_FLAG ) valueTables = new Real[functionCount*sampleCount];
315  if( flags & D_VALUE_FLAG ) dValueTables = new Real[functionCount*sampleCount];
316  PPolynomial<Degree+1> function;
317  PPolynomial<Degree> dFunction;
318  for( int i=0 ; i<functionCount ; i++ )
319  {
320  if(smooth>0)
321  {
322  function = baseFunctions[i].MovingAverage(smooth);
323  dFunction = baseFunctions[i].derivative().MovingAverage(smooth);
324  }
325  else
326  {
327  function = baseFunctions[i];
328  dFunction = baseFunctions[i].derivative();
329  }
330  for( int j=0 ; j<sampleCount ; j++ )
331  {
332  double x=double(j)/(sampleCount-1);
333  if(flags & VALUE_FLAG){ valueTables[j*functionCount+i]=Real( function(x));}
334  if(flags & D_VALUE_FLAG){dValueTables[j*functionCount+i]=Real(dFunction(x));}
335  }
336  }
337  }
338  template<int Degree,class Real>
339  void BSplineData<Degree,Real>::setValueTables(int flags,double valueSmooth,double normalSmooth){
340  clearValueTables();
341  if(flags & VALUE_FLAG){ valueTables=new Real[functionCount*sampleCount];}
342  if(flags & D_VALUE_FLAG){dValueTables=new Real[functionCount*sampleCount];}
343  PPolynomial<Degree+1> function;
344  PPolynomial<Degree> dFunction;
345  for(int i=0;i<functionCount;i++){
346  if(valueSmooth>0) { function=baseFunctions[i].MovingAverage(valueSmooth);}
347  else { function=baseFunctions[i];}
348  if(normalSmooth>0) {dFunction=baseFunctions[i].derivative().MovingAverage(normalSmooth);}
349  else {dFunction=baseFunctions[i].derivative();}
350 
351  for(int j=0;j<sampleCount;j++){
352  double x=double(j)/(sampleCount-1);
353  if(flags & VALUE_FLAG){ valueTables[j*functionCount+i]=Real( function(x));}
354  if(flags & D_VALUE_FLAG){dValueTables[j*functionCount+i]=Real(dFunction(x));}
355  }
356  }
357  }
358 
359 
360  template<int Degree,class Real>
362  if( valueTables){delete[] valueTables;}
363  if(dValueTables){delete[] dValueTables;}
364  valueTables=dValueTables=NULL;
365  }
366 
367  template<int Degree,class Real>
368  inline int BSplineData<Degree,Real>::Index( int i1 , int i2 ) const { return i1*functionCount+i2; }
369  template<int Degree,class Real>
370  inline int BSplineData<Degree,Real>::SymmetricIndex( int i1 , int i2 )
371  {
372  if( i1>i2 ) return ((i1*i1+i1)>>1)+i2;
373  else return ((i2*i2+i2)>>1)+i1;
374  }
375  template<int Degree,class Real>
376  inline int BSplineData<Degree,Real>::SymmetricIndex( int i1 , int i2 , int& index )
377  {
378  if( i1<i2 )
379  {
380  index = ((i2*i2+i2)>>1)+i1;
381  return 1;
382  }
383  else
384  {
385  index = ((i1*i1+i1)>>1)+i2;
386  return 0;
387  }
388  }
389 
390 
391  ////////////////////////
392  // BSplineElementData //
393  ////////////////////////
394  template< int Degree >
395  BSplineElements< Degree >::BSplineElements( int res , int offset , int boundary )
396  {
397  denominator = 1;
398  this->resize( res , BSplineElementCoefficients<Degree>() );
399 
400  for( int i=0 ; i<=Degree ; i++ )
401  {
402  int idx = -_off + offset + i;
403  if( idx>=0 && idx<res ) (*this)[idx][i] = 1;
404  }
405  if( boundary!=0 )
406  {
407  _addLeft( offset-2*res , boundary ) , _addRight( offset+2*res , boundary );
408  if( Degree&1 ) _addLeft( offset-res , boundary ) , _addRight( offset+res , boundary );
409  else _addLeft( -offset-1 , boundary ) , _addRight( -offset-1+2*res , boundary );
410  }
411  }
412  template< int Degree >
413  void BSplineElements< Degree >::_addLeft( int offset , int boundary )
414  {
415  int res = int( this->size() );
416  bool set = false;
417  for( int i=0 ; i<=Degree ; i++ )
418  {
419  int idx = -_off + offset + i;
420  if( idx>=0 && idx<res ) (*this)[idx][i] += boundary , set = true;
421  }
422  if( set ) _addLeft( offset-2*res , boundary );
423  }
424  template< int Degree >
425  void BSplineElements< Degree >::_addRight( int offset , int boundary )
426  {
427  int res = int( this->size() );
428  bool set = false;
429  for( int i=0 ; i<=Degree ; i++ )
430  {
431  int idx = -_off + offset + i;
432  if( idx>=0 && idx<res ) (*this)[idx][i] += boundary , set = true;
433  }
434  if( set ) _addRight( offset+2*res , boundary );
435  }
436  template< int Degree >
438  {
439  fprintf( stderr , "[ERROR] B-spline up-sampling not supported for degree %d\n" , Degree );
440  exit( 0 );
441  }
442  template<>
444  {
445  high.resize( size()*2 );
446  high.assign( high.size() , BSplineElementCoefficients<1>() );
447  for( int i=0 ; i<int(size()) ; i++ )
448  {
449  high[2*i+0][0] += 1 * (*this)[i][0];
450  high[2*i+0][1] += 0 * (*this)[i][0];
451  high[2*i+1][0] += 2 * (*this)[i][0];
452  high[2*i+1][1] += 1 * (*this)[i][0];
453 
454  high[2*i+0][0] += 1 * (*this)[i][1];
455  high[2*i+0][1] += 2 * (*this)[i][1];
456  high[2*i+1][0] += 0 * (*this)[i][1];
457  high[2*i+1][1] += 1 * (*this)[i][1];
458  }
459  high.denominator = denominator * 2;
460  }
461  template<>
463  {
464  // /----\
465  // / \
466  // / \ = 1 /--\ +3 /--\ +3 /--\ +1 /--\
467  // / \ / \ / \ / \ / \
468  // |----------| |----------| |----------| |----------| |----------|
469 
470  high.resize( size()*2 );
471  high.assign( high.size() , BSplineElementCoefficients<2>() );
472  for( int i=0 ; i<int(size()) ; i++ )
473  {
474  high[2*i+0][0] += 1 * (*this)[i][0];
475  high[2*i+0][1] += 0 * (*this)[i][0];
476  high[2*i+0][2] += 0 * (*this)[i][0];
477  high[2*i+1][0] += 3 * (*this)[i][0];
478  high[2*i+1][1] += 1 * (*this)[i][0];
479  high[2*i+1][2] += 0 * (*this)[i][0];
480 
481  high[2*i+0][0] += 3 * (*this)[i][1];
482  high[2*i+0][1] += 3 * (*this)[i][1];
483  high[2*i+0][2] += 1 * (*this)[i][1];
484  high[2*i+1][0] += 1 * (*this)[i][1];
485  high[2*i+1][1] += 3 * (*this)[i][1];
486  high[2*i+1][2] += 3 * (*this)[i][1];
487 
488  high[2*i+0][0] += 0 * (*this)[i][2];
489  high[2*i+0][1] += 1 * (*this)[i][2];
490  high[2*i+0][2] += 3 * (*this)[i][2];
491  high[2*i+1][0] += 0 * (*this)[i][2];
492  high[2*i+1][1] += 0 * (*this)[i][2];
493  high[2*i+1][2] += 1 * (*this)[i][2];
494  }
495  high.denominator = denominator * 4;
496  }
497 
498  template< int Degree >
500  {
501  d.resize( this->size() );
502  d.assign( d.size() , BSplineElementCoefficients< Degree-1 >() );
503  for( int i=0 ; i<int(this->size()) ; i++ ) for( int j=0 ; j<=Degree ; j++ )
504  {
505  if( j-1>=0 ) d[i][j-1] -= (*this)[i][j];
506  if( j<Degree ) d[i][j ] += (*this)[i][j];
507  }
508  d.denominator = denominator;
509  }
510  // If we were really good, we would implement this integral table to store
511  // rational values to improve precision...
512  template< int Degree1 , int Degree2 >
513  void SetBSplineElementIntegrals( double integrals[Degree1+1][Degree2+1] )
514  {
515  for( int i=0 ; i<=Degree1 ; i++ )
516  {
518  for( int j=0 ; j<=Degree2 ; j++ )
519  {
521  integrals[i][j] = ( p1 * p2 ).integral( 0 , 1 );
522  }
523  }
524  }
525 
526 
527  }
528 }
virtual void clearDotTables(int flags)
StartingPolynomial shift(double t) const
Definition: ppolynomial.hpp:58
Polynomial< Degree > p
Definition: ppolynomial.h:46
bool RightOverlap(unsigned int depth, int offset)
static int CenterIndex(int depth, int offSet)
Definition: binary_node.h:47
static int SymmetricIndex(int i1, int i2)
virtual void clearValueTables(void)
void _addLeft(int offset, int boundary)
static int CenterCount(int depth)
Definition: binary_node.h:42
static int CumulativeCenterCount(int maxDepth)
Definition: binary_node.h:44
PPolynomial< Degree-1 > derivative(void) const
bool LeftOverlap(unsigned int depth, int offset)
virtual void setDotTables(int flags)
void upSample(BSplineElements &high) const
void differentiate(BSplineElements< Degree-1 > &d) const
void _addRight(int offset, int boundary)
int Index(int i1, int i2) const
static void DepthAndOffset(int idx, int &depth, int &offset)
Definition: binary_node.h:64
void set(int maxDepth, bool useDotRatios=true, bool reflectBoundary=false)
static Polynomial BSplineComponent(int i)
Definition: polynomial.hpp:307
void setSampleSpan(int idx, int &start, int &end, double smooth=0) const
void SetBSplineElementIntegrals(double integrals[Degree1+1][Degree2+1])
virtual void setValueTables(int flags, double smooth=0)
int ReflectRight(unsigned int depth, int offset)
int ReflectLeft(unsigned int depth, int offset)
static PPolynomial BSpline(double radius=0.5)
static int CornerCount(int depth)
Definition: binary_node.h:43
PPolynomial< Degree+1 > MovingAverage(double radius)
static void CenterAndWidth(int depth, int offset, Real &center, Real &width)
Definition: binary_node.h:53