Orfeo Toolbox  3.16
itkOptLinearInterpolateImageFunction.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Insight Segmentation & Registration Toolkit
4  Module: $RCSfile: itkOptLinearInterpolateImageFunction.h,v $
5  Language: C++
6  Date: $Date: 2009-10-29 11:19:19 $
7  Version: $Revision: 1.11 $
8 
9  Copyright (c) Insight Software Consortium. All rights reserved.
10  See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
11 
12  This software is distributed WITHOUT ANY WARRANTY; without even
13  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
14  PURPOSE. See the above copyright notices for more information.
15 
16 =========================================================================*/
17 #ifndef __itkOptLinearInterpolateImageFunction_h
18 #define __itkOptLinearInterpolateImageFunction_h
19 
20 #include "itkInterpolateImageFunction.h"
21 
22 namespace itk
23 {
24 
42 template <class TInputImage, class TCoordRep = double>
43 class ITK_EXPORT LinearInterpolateImageFunction :
44  public InterpolateImageFunction<TInputImage,TCoordRep>
45 {
46 public:
48  typedef LinearInterpolateImageFunction Self;
49  typedef InterpolateImageFunction<TInputImage,TCoordRep> Superclass;
50  typedef SmartPointer<Self> Pointer;
51  typedef SmartPointer<const Self> ConstPointer;
52 
54  itkTypeMacro(LinearInterpolateImageFunction, InterpolateImageFunction);
55 
57  itkNewMacro(Self);
58 
60  typedef typename Superclass::OutputType OutputType;
61 
63  typedef typename Superclass::InputImageType InputImageType;
64 
66  typedef typename Superclass::InputPixelType InputPixelType;
67 
69  typedef typename Superclass::RealType RealType;
70 
72  itkStaticConstMacro(ImageDimension, unsigned int,Superclass::ImageDimension);
73 
75  typedef typename Superclass::IndexType IndexType;
76  typedef typename Superclass::IndexValueType IndexValueType;
77 
79  typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
80 
89  virtual inline OutputType EvaluateAtContinuousIndex( const
90  ContinuousIndexType &
91  index ) const
92  {
93  return this->EvaluateOptimized( Dispatch< ImageDimension >(), index );
94  }
95 
96 protected:
97  LinearInterpolateImageFunction();
98  ~LinearInterpolateImageFunction();
99  void PrintSelf(std::ostream& os, Indent indent) const;
100 
101 private:
102  LinearInterpolateImageFunction( const Self& ); //purposely not implemented
103  void operator=( const Self& ); //purposely not implemented
104 
106  static const unsigned long m_Neighbors;
107 
108  struct DispatchBase {};
109  template< unsigned int > struct Dispatch : DispatchBase {};
110 
111  inline OutputType EvaluateOptimized( const Dispatch<0> &,
112  const ContinuousIndexType & index) const
113  {
114  return 0;
115  }
116 
117  inline OutputType EvaluateOptimized( const Dispatch<1>&,
118  const ContinuousIndexType & index) const
119  {
120  IndexType basei;
121 
122  basei[0] = Math::Floor<IndexValueType>(index[0]);
123  if( basei[0] < this->m_StartIndex[0] )
124  {
125  basei[0] = this->m_StartIndex[0];
126  }
127 
128  const double distance = index[0] - static_cast<double>(basei[0]);
129 
130  const RealType val0 = this->GetInputImage()->GetPixel( basei );
131 
132  if(distance <= 0.)
133  {
134  return( static_cast<OutputType>( val0 ) );
135  }
136 
137  ++basei[0];
138  if(basei[0]>this->m_EndIndex[0])
139  {
140  return( static_cast<OutputType>( val0 ) );
141  }
142  const RealType val1 = this->GetInputImage()->GetPixel( basei );
143 
144  return( static_cast<OutputType>( val0 + ( val1 - val0 ) * distance ) );
145  }
146 
147  inline OutputType EvaluateOptimized( const Dispatch<2>&,
148  const ContinuousIndexType & index) const
149  {
150  IndexType basei;
151 
152  basei[0] = Math::Floor<IndexValueType>(index[0]);
153  if( basei[0] < this->m_StartIndex[0] )
154  {
155  basei[0] = this->m_StartIndex[0];
156  }
157  const double distance0 = index[0] - static_cast<double>(basei[0]);
158 
159  basei[1] = Math::Floor<IndexValueType>(index[1]);
160  if( basei[1] < this->m_StartIndex[1] )
161  {
162  basei[1] = this->m_StartIndex[1];
163  }
164  const double distance1 = index[1] - static_cast<double>(basei[1]);
165 
166 
167  const RealType val00 = this->GetInputImage()->GetPixel( basei );
168  if(distance0 <= 0. && distance1 <= 0.)
169  {
170  return( static_cast<OutputType>( val00 ) );
171  }
172  else if(distance1 <= 0.) // if they have the same "y"
173  {
174  ++basei[0]; // then interpolate across "x"
175  if(basei[0]>this->m_EndIndex[0])
176  {
177  return( static_cast<OutputType>( val00 ) );
178  }
179  const RealType val10 = this->GetInputImage()->GetPixel( basei );
180  return( static_cast<OutputType>(val00 + (val10 - val00) * distance0) );
181  }
182  else if(distance0 <= 0.) // if they have the same "x"
183  {
184  ++basei[1]; // then interpolate across "y"
185  if(basei[1]>this->m_EndIndex[1])
186  {
187  return( static_cast<OutputType>( val00 ) );
188  }
189  const RealType val01 = this->GetInputImage()->GetPixel( basei );
190  return( static_cast<OutputType>(val00 + (val01 - val00) * distance1) );
191  }
192  else // interpolate across "xy"
193  {
194  ++basei[0];
195  if(basei[0]>this->m_EndIndex[0]) // interpolate across "y"
196  {
197  --basei[0];
198  ++basei[1];
199  if(basei[1]>this->m_EndIndex[1])
200  {
201  return( static_cast<OutputType>( val00 ) );
202  }
203  const RealType val01 = this->GetInputImage()->GetPixel( basei );
204  return( static_cast<OutputType>(val00 + (val01 - val00) * distance1) );
205  }
206  const RealType val10 = this->GetInputImage()->GetPixel( basei );
207 
208  const RealType valx0 = val00 + (val10 - val00) * distance0;
209 
210  ++basei[1];
211  if(basei[1]>this->m_EndIndex[1]) // interpolate across "x"
212  {
213  return( static_cast<OutputType>( valx0 ) );
214  }
215  const RealType val11 = this->GetInputImage()->GetPixel( basei );
216  --basei[0];
217  const RealType val01 = this->GetInputImage()->GetPixel( basei );
218 
219  const RealType valx1 = val01 + (val11 - val01) * distance0;
220 
221  return( static_cast<OutputType>( valx0 + (valx1-valx0) * distance1 ) );
222  }
223  }
224 
225  inline OutputType EvaluateOptimized( const Dispatch<3>&,
226  const ContinuousIndexType & index) const
227  {
228  IndexType basei;
229 
230  basei[0] = Math::Floor<IndexValueType>(index[0]);
231  if( basei[0] < this->m_StartIndex[0] )
232  {
233  basei[0] = this->m_StartIndex[0];
234  }
235  const double distance0 = index[0] - static_cast<double>(basei[0]);
236 
237  basei[1] = Math::Floor<IndexValueType>(index[1]);
238  if( basei[1] < this->m_StartIndex[1] )
239  {
240  basei[1] = this->m_StartIndex[1];
241  }
242  const double distance1 = index[1] - static_cast<double>(basei[1]);
243 
244  basei[2] = Math::Floor<IndexValueType>(index[2]);
245  if( basei[2] < this->m_StartIndex[2] )
246  {
247  basei[2] = this->m_StartIndex[2];
248  }
249  const double distance2 = index[2] - static_cast<double>(basei[2]);
250 
251  if(distance0<=0. && distance1<=0. && distance2<=0.)
252  {
253  return( static_cast<OutputType>( this->GetInputImage()->GetPixel( basei ) ) );
254  }
255 
256 
257  const RealType val000 = this->GetInputImage()->GetPixel( basei );
258 
259  if(distance2 <= 0.)
260  {
261  if(distance1 <= 0.) // interpolate across "x"
262  {
263  ++basei[0];
264  if(basei[0]>this->m_EndIndex[0])
265  {
266  return( static_cast<OutputType>( val000 ) );
267  }
268  const RealType val100 = this->GetInputImage()->GetPixel( basei );
269 
270  return static_cast<OutputType>( val000 + (val100-val000) * distance0 );
271  }
272  else if(distance0 <= 0.) // interpolate across "y"
273  {
274  ++basei[1];
275  if(basei[1]>this->m_EndIndex[1])
276  {
277  return( static_cast<OutputType>( val000 ) );
278  }
279  const RealType val010 = this->GetInputImage()->GetPixel( basei );
280 
281  return static_cast<OutputType>( val000 + (val010-val000) * distance1 );
282  }
283  else // interpolate across "xy"
284  {
285  ++basei[0];
286  if(basei[0]>this->m_EndIndex[0]) // interpolate across "y"
287  {
288  --basei[0];
289  ++basei[1];
290  if(basei[1]>this->m_EndIndex[1])
291  {
292  return( static_cast<OutputType>( val000 ) );
293  }
294  const RealType val010 = this->GetInputImage()->GetPixel( basei );
295 
296  return static_cast<OutputType>( val000 + (val010-val000) * distance1 );
297  }
298  const RealType val100 = this->GetInputImage()->GetPixel( basei );
299 
300  const RealType valx00 = val000 + (val100-val000) * distance0;
301 
302  ++basei[1];
303  if(basei[1]>this->m_EndIndex[1]) // interpolate across "x"
304  {
305  return( static_cast<OutputType>( valx00 ) );
306  }
307  const RealType val110 = this->GetInputImage()->GetPixel( basei );
308 
309  --basei[0];
310  const RealType val010 = this->GetInputImage()->GetPixel( basei );
311 
312  const RealType valx10 = val010 + (val110-val010) * distance0;
313 
314  return static_cast<OutputType>( valx00 + (valx10-valx00) * distance1 );
315  }
316  }
317  else
318  {
319  if(distance1 <= 0.)
320  {
321  if(distance0 <= 0.) // interpolate across "z"
322  {
323  ++basei[2];
324  if(basei[2]>this->m_EndIndex[2])
325  {
326  return( static_cast<OutputType>( val000 ) );
327  }
328  const RealType val001 = this->GetInputImage()->GetPixel( basei );
329 
330  return static_cast<OutputType>( val000 + (val001-val000) * distance2 );
331  }
332  else // interpolate across "xz"
333  {
334  ++basei[0];
335  if(basei[0]>this->m_EndIndex[0]) // interpolate across "z"
336  {
337  --basei[0];
338  ++basei[2];
339  if(basei[2]>this->m_EndIndex[2])
340  {
341  return( static_cast<OutputType>( val000 ) );
342  }
343  const RealType val001 = this->GetInputImage()->GetPixel( basei );
344 
345  return static_cast<OutputType>( val000 + (val001-val000) * distance2 );
346  }
347  const RealType val100 = this->GetInputImage()->GetPixel( basei );
348 
349  const RealType valx00 = val000 + (val100-val000) * distance0;
350 
351  ++basei[2];
352  if(basei[2]>this->m_EndIndex[2]) // interpolate across "x"
353  {
354  return( static_cast<OutputType>( valx00 ) );
355  }
356  const RealType val101 = this->GetInputImage()->GetPixel( basei );
357 
358  --basei[0];
359  const RealType val001 = this->GetInputImage()->GetPixel( basei );
360 
361  const RealType valx01 = val001 + (val101-val001) * distance0;
362 
363  return static_cast<OutputType>( valx00 + (valx01-valx00) * distance2 );
364  }
365  }
366  else if(distance0 <= 0.) // interpolate across "yz"
367  {
368  ++basei[1];
369  if(basei[1]>this->m_EndIndex[1]) // interpolate across "z"
370  {
371  --basei[1];
372  ++basei[2];
373  if(basei[2]>this->m_EndIndex[2])
374  {
375  return( static_cast<OutputType>( val000 ) );
376  }
377  const RealType val001 = this->GetInputImage()->GetPixel( basei );
378 
379  return static_cast<OutputType>( val000 + (val001-val000) * distance2 );
380  }
381  const RealType val010 = this->GetInputImage()->GetPixel( basei );
382 
383  const RealType val0x0 = val000 + (val010-val000) * distance1;
384 
385  ++basei[2];
386  if(basei[2]>this->m_EndIndex[2]) // interpolate across "y"
387  {
388  return( static_cast<OutputType>( val0x0 ) );
389  }
390  const RealType val011 = this->GetInputImage()->GetPixel( basei );
391 
392  --basei[1];
393  const RealType val001 = this->GetInputImage()->GetPixel( basei );
394 
395  const RealType val0x1 = val001 + (val011-val001) * distance1;
396 
397  return static_cast<OutputType>( val0x0 + (val0x1-val0x0) * distance2 );
398  }
399  else // interpolate across "xyz"
400  {
401  ++basei[0];
402  if(basei[0]>this->m_EndIndex[0]) // interpolate across "yz"
403  {
404  --basei[0];
405  ++basei[1];
406  if(basei[1]>this->m_EndIndex[1]) // interpolate across "z"
407  {
408  --basei[1];
409  ++basei[2];
410  if(basei[2]>this->m_EndIndex[2])
411  {
412  return( static_cast<OutputType>( val000 ) );
413  }
414  const RealType val001 = this->GetInputImage()->GetPixel( basei );
415 
416  return static_cast<OutputType>( val000 + (val001-val000) * distance2 );
417  }
418  const RealType val010 = this->GetInputImage()->GetPixel( basei );
419 
420  const RealType val0x0 = val000 + (val010-val000) * distance1;
421 
422  ++basei[2];
423  if(basei[2]>this->m_EndIndex[2]) // interpolate across "y"
424  {
425  return( static_cast<OutputType>( val0x0 ) );
426  }
427  const RealType val011 = this->GetInputImage()->GetPixel( basei );
428 
429  --basei[1];
430  const RealType val001 = this->GetInputImage()->GetPixel( basei );
431 
432  const RealType val0x1 = val001 + (val011-val001) * distance1;
433 
434  return static_cast<OutputType>( val0x0 + (val0x1-val0x0) * distance2 );
435  }
436  const RealType val100 = this->GetInputImage()->GetPixel( basei );
437 
438  const RealType valx00 = val000 + (val100-val000) * distance0;
439 
440  ++basei[1];
441  if(basei[1]>this->m_EndIndex[1]) // interpolate across "xz"
442  {
443  --basei[1];
444  ++basei[2];
445  if(basei[2]>this->m_EndIndex[2]) // interpolate across "x"
446  {
447  return( static_cast<OutputType>( valx00 ) );
448  }
449  const RealType val101 = this->GetInputImage()->GetPixel( basei );
450 
451  --basei[0];
452  const RealType val001 = this->GetInputImage()->GetPixel( basei );
453 
454  const RealType valx01 = val001 + (val101-val001) * distance0;
455 
456  return static_cast<OutputType>( valx00 + (valx01-valx00) * distance2 );
457  }
458  const RealType val110 = this->GetInputImage()->GetPixel( basei );
459 
460  --basei[0];
461  const RealType val010 = this->GetInputImage()->GetPixel( basei );
462 
463  const RealType valx10 = val010 + (val110-val010) * distance0;
464 
465  const RealType valxx0 = valx00 + (valx10-valx00) * distance1;
466 
467 
468  ++basei[2];
469  if(basei[2]>this->m_EndIndex[2]) // interpolate across "xy"
470  {
471  return( static_cast<OutputType>( valxx0 ) );
472  }
473  const RealType val011 = this->GetInputImage()->GetPixel( basei );
474 
475  ++basei[0];
476  const RealType val111 = this->GetInputImage()->GetPixel( basei );
477 
478  --basei[1];
479  const RealType val101 = this->GetInputImage()->GetPixel( basei );
480 
481  --basei[0];
482  const RealType val001 = this->GetInputImage()->GetPixel( basei );
483 
484  const RealType valx01 = val001 + (val101-val001) * distance0;
485  const RealType valx11 = val011 + (val111-val011) * distance0;
486  const RealType valxx1 = valx01 + (valx11-valx01) * distance1;
487 
488  return( static_cast<OutputType>( valxx0 + (valxx1-valxx0) * distance2 ) );
489  }
490  }
491  }
492 
493  inline OutputType EvaluateOptimized( const DispatchBase &,
494  const ContinuousIndexType & index) const
495  {
496  return this->EvaluateUnoptimized( index );
497  }
498 
499  virtual inline OutputType EvaluateUnoptimized(
500  const ContinuousIndexType & index) const;
501 };
502 
503 } // end namespace itk
504 
505 // Define instantiation macro for this template.
506 #define ITK_TEMPLATE_LinearInterpolateImageFunction(_, EXPORT, x, y) namespace itk { \
507  _(2(class EXPORT LinearInterpolateImageFunction< ITK_TEMPLATE_2 x >)) \
508  namespace Templates { typedef LinearInterpolateImageFunction< ITK_TEMPLATE_2 x > \
509  LinearInterpolateImageFunction##y; } \
510  }
511 
512 #if ITK_TEMPLATE_EXPLICIT
513 # include "Templates/itkLinearInterpolateImageFunction+-.h"
514 #endif
515 
516 #if ITK_TEMPLATE_TXX
517 # include "itkOptLinearInterpolateImageFunction.txx"
518 #endif
519 
520 #endif
Generated at Sat Feb 2 2013 23:57:17 for Orfeo Toolbox with doxygen 1.8.1.1