00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 #ifndef __itkOptLinearInterpolateImageFunction_h
00018 #define __itkOptLinearInterpolateImageFunction_h
00019
00020 #include "itkInterpolateImageFunction.h"
00021
00022 namespace itk
00023 {
00024
00042 template <class TInputImage, class TCoordRep = float>
00043 class ITK_EXPORT LinearInterpolateImageFunction :
00044 public InterpolateImageFunction<TInputImage,TCoordRep>
00045 {
00046 public:
00048 typedef LinearInterpolateImageFunction Self;
00049 typedef InterpolateImageFunction<TInputImage,TCoordRep> Superclass;
00050 typedef SmartPointer<Self> Pointer;
00051 typedef SmartPointer<const Self> ConstPointer;
00052
00054 itkTypeMacro(LinearInterpolateImageFunction, InterpolateImageFunction);
00055
00057 itkNewMacro(Self);
00058
00060 typedef typename Superclass::OutputType OutputType;
00061
00063 typedef typename Superclass::InputImageType InputImageType;
00064
00066 typedef typename Superclass::InputPixelType InputPixelType;
00067
00069 typedef typename Superclass::RealType RealType;
00070
00072 itkStaticConstMacro(ImageDimension, unsigned int,Superclass::ImageDimension);
00073
00075 typedef typename Superclass::IndexType IndexType;
00076
00078 typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
00079
00088 virtual inline OutputType EvaluateAtContinuousIndex( const
00089 ContinuousIndexType &
00090 index ) const
00091 {
00092 return this->EvaluateOptimized( Dispatch< ImageDimension >(), index );
00093 }
00094
00095 protected:
00096 LinearInterpolateImageFunction();
00097 ~LinearInterpolateImageFunction();
00098 void PrintSelf(std::ostream& os, Indent indent) const;
00099
00100 private:
00101 LinearInterpolateImageFunction( const Self& );
00102 void operator=( const Self& );
00103
00105 static const unsigned long m_Neighbors;
00106
00107 struct DispatchBase {};
00108 template< unsigned int > struct Dispatch : DispatchBase {};
00109
00110 inline OutputType EvaluateOptimized( const Dispatch<0> &,
00111 const ContinuousIndexType & index) const
00112 {
00113 return 0;
00114 }
00115
00116 inline OutputType EvaluateOptimized( const Dispatch<1>&,
00117 const ContinuousIndexType & index) const
00118 {
00119 IndexType basei;
00120
00121 double i = index[0];
00122 basei[0] = (long)i;
00123 if( i < 0.0 && double(basei[0]) != i)
00124 {
00125 basei[0]--;
00126 }
00127
00128 double distance = i - double(basei[0]);
00129
00130 double val0 = this->GetInputImage()->GetPixel( basei );
00131 double val1 = val0;
00132 ++basei[0];
00133 val1 = this->GetInputImage()->GetPixel( basei );
00134
00135 return( static_cast<OutputType>( val0 + distance * ( val1 - val0 ) ) );
00136 }
00137
00138 inline OutputType EvaluateOptimized( const Dispatch<2>&,
00139 const ContinuousIndexType & index) const
00140 {
00141 IndexType basei;
00142
00143 double i = index[0];
00144 basei[0] = (long)i;
00145 if( i < 0.0 && double(basei[0]) != i)
00146 {
00147 basei[0]--;
00148 }
00149 double distance0 = i - double(basei[0]);
00150
00151 i = index[1];
00152 basei[1] = (long)i;
00153 if( i < 0.0 && double(basei[1]) != i)
00154 {
00155 basei[1]--;
00156 }
00157 double distance1 = i - double(basei[1]);
00158
00159
00160 double val00 = this->GetInputImage()->GetPixel( basei );
00161 if(distance0+distance1 == 0)
00162 {
00163 return( static_cast<OutputType>( val00 ) );
00164 }
00165 else if(distance1 == 0)
00166 {
00167 ++basei[0];
00168 if(basei[0]>this->m_EndIndex[0])
00169 {
00170 return( static_cast<OutputType>( val00 ) );
00171 }
00172 double val10 = this->GetInputImage()->GetPixel( basei );
00173 return( static_cast<OutputType>(val00 + distance0 * (val10 - val00)) );
00174 }
00175 else if(distance0 == 0)
00176 {
00177 ++basei[1];
00178 if(basei[1]>this->m_EndIndex[1])
00179 {
00180 return( static_cast<OutputType>( val00 ) );
00181 }
00182 double val01 = this->GetInputImage()->GetPixel( basei );
00183 return( static_cast<OutputType>(val00 + distance1 * (val01 - val00)) );
00184 }
00185 else
00186 {
00187 ++basei[0];
00188 double val10 = this->GetInputImage()->GetPixel( basei );
00189 ++basei[1];
00190 double val11 = this->GetInputImage()->GetPixel( basei );
00191 --basei[0];
00192 double val01 = this->GetInputImage()->GetPixel( basei );
00193
00194 double val0 = val00 + distance1 * ( val01 - val00 );
00195 double val1 = val10 + distance1 * ( val11 - val10 );
00196
00197 return( static_cast<OutputType>( val0 + distance0 * (val1-val0) ) );
00198 }
00199 }
00200
00201 inline OutputType EvaluateOptimized( const Dispatch<3>&,
00202 const ContinuousIndexType & index) const
00203 {
00204 IndexType basei;
00205 double val[8];
00206
00207 unsigned long min0;
00208 unsigned long max0;
00209 unsigned long min1;
00210 unsigned long max1;
00211 unsigned long min2;
00212 unsigned long max2;
00213
00214 double i = index[0];
00215 basei[0] = (long)i;
00216 if( i < 0.0 && double(basei[0]) != i)
00217 {
00218 basei[0]--;
00219 }
00220 double distance0 = i - double(basei[0]);
00221
00222 i = index[1];
00223 basei[1] = (long)i;
00224 if( i < 0.0 && double(basei[1]) != i)
00225 {
00226 basei[1]--;
00227 }
00228 double distance1 = i - double(basei[1]);
00229
00230 i = index[2];
00231 basei[2] = (long)i;
00232 if( i < 0.0 && double(basei[2]) != i)
00233 {
00234 basei[2]--;
00235 }
00236 double distance2 = i - double(basei[2]);
00237
00238 val[0] = this->GetInputImage()->GetPixel( basei );
00239 if(distance0+distance1+distance2 == 0)
00240 {
00241 return( static_cast<OutputType>( val[0] ) );
00242 }
00243 if(distance0 > 0.0)
00244 {
00245 min0 = basei[0];
00246 max0 = basei[0]+1;
00247 if(max0>this->m_EndIndex[0])
00248 {
00249 max0 = this->m_EndIndex[0];
00250 }
00251 }
00252 if(distance1 > 0.0)
00253 {
00254 min1 = basei[1];
00255 max1 = basei[1]+1;
00256 if(max1>this->m_EndIndex[1])
00257 {
00258 max1 = this->m_EndIndex[1];
00259 }
00260 }
00261 if(distance2 > 0.0)
00262 {
00263 min2 = basei[2];
00264 max2 = basei[2]+1;
00265 if(max2>this->m_EndIndex[2])
00266 {
00267 max2 = this->m_EndIndex[2];
00268 }
00269 }
00270 if(distance2 == 0)
00271 {
00272 if(distance1 == 0)
00273 {
00274 basei[0] = max0;
00275 val[1] = this->GetInputImage()->GetPixel( basei );
00276
00277 double val0 = val[0] + distance0 * (val[1]-val[0]);
00278
00279 return( static_cast<OutputType>( val0 ) );
00280 }
00281 else if(distance0 == 0)
00282 {
00283 basei[1] = max1;
00284 val[2] = this->GetInputImage()->GetPixel( basei );
00285
00286 double val0 = val[0] + distance1 * (val[2]-val[0]);
00287
00288 return( static_cast<OutputType>( val0 ) );
00289 }
00290 else
00291 {
00292 basei[0] = max0;
00293 val[1] = this->GetInputImage()->GetPixel( basei );
00294 basei[1] = max1;
00295 val[3] = this->GetInputImage()->GetPixel( basei );
00296 basei[0] = min0;
00297 val[2] = this->GetInputImage()->GetPixel( basei );
00298
00299 double val0 = val[0] + distance0 * (val[1]-val[0]);
00300 double val1 = val[2] + distance0 * (val[3]-val[2]);
00301
00302 return( static_cast<OutputType>( val0 + distance1 * (val1 - val0) ) );
00303 }
00304 }
00305 else
00306 {
00307 basei[2] = max2;
00308 val[4] = this->GetInputImage()->GetPixel( basei );
00309 if(distance1 == 0)
00310 {
00311 if(distance0 == 0)
00312 {
00313 return( static_cast<OutputType>( val[0] + distance2
00314 * (val[4] - val[0]) ) );
00315 }
00316 else
00317 {
00318 basei[0] = max0;
00319 val[5] = this->GetInputImage()->GetPixel( basei );
00320 basei[2] = min2;
00321 val[1] = this->GetInputImage()->GetPixel( basei );
00322
00323 double val0 = val[0] + distance0 * (val[1]-val[0]);
00324 double val1 = val[4] + distance0 * (val[5]-val[4]);
00325
00326 return( static_cast<OutputType>( val0 + distance2 * (val1 - val0) ) );
00327 }
00328 }
00329 else if(distance0 == 0)
00330 {
00331 basei[1] = max1;
00332 val[6] = this->GetInputImage()->GetPixel( basei );
00333 basei[2] = min2;
00334 val[2] = this->GetInputImage()->GetPixel( basei );
00335
00336 double val0 = val[0] + distance1 * (val[2]-val[0]);
00337 double val1 = val[4] + distance1 * (val[6]-val[4]);
00338
00339 return( static_cast<OutputType>( val0 + distance2 * (val1 - val0) ) );
00340 }
00341 else
00342 {
00343 basei[0] = max0;
00344 val[5] = this->GetInputImage()->GetPixel( basei );
00345 basei[1] = max1;
00346 val[7] = this->GetInputImage()->GetPixel( basei );
00347 basei[0] = min0;
00348 val[6] = this->GetInputImage()->GetPixel( basei );
00349 basei[2] = min2;
00350 val[2] = this->GetInputImage()->GetPixel( basei );
00351 basei[0] = max0;
00352 val[3] = this->GetInputImage()->GetPixel( basei );
00353 basei[1] = min1;
00354 val[1] = this->GetInputImage()->GetPixel( basei );
00355
00356 double val00 = val[0] + distance0 * (val[1]-val[0]);
00357 double val01 = val[2] + distance0 * (val[3]-val[2]);
00358 double val0 = val00 + distance1 * (val01-val00);
00359
00360 double val10 = val[4] + distance0 * (val[5]-val[4]);
00361 double val11 = val[6] + distance0 * (val[7]-val[6]);
00362 double val1 = val10 + distance1 * (val11-val10);
00363
00364 return( static_cast<OutputType>( val0 + distance2 * (val1-val0) ) );
00365 }
00366 }
00367 }
00368
00369 inline OutputType EvaluateOptimized( const DispatchBase &,
00370 const ContinuousIndexType & index) const
00371 {
00372 return this->EvaluateUnoptimized( index );
00373 }
00374
00375 virtual inline OutputType EvaluateUnoptimized(
00376 const ContinuousIndexType & index) const;
00377 };
00378
00379 }
00380
00381
00382 #define ITK_TEMPLATE_LinearInterpolateImageFunction(_, EXPORT, x, y) namespace itk { \
00383 _(2(class EXPORT LinearInterpolateImageFunction< ITK_TEMPLATE_2 x >)) \
00384 namespace Templates { typedef LinearInterpolateImageFunction< ITK_TEMPLATE_2 x > \
00385 LinearInterpolateImageFunction##y; } \
00386 }
00387
00388 #if ITK_TEMPLATE_EXPLICIT
00389 # include "Templates/itkLinearInterpolateImageFunction+-.h"
00390 #endif
00391
00392 #if ITK_TEMPLATE_TXX
00393 # include "itkOptLinearInterpolateImageFunction.txx"
00394 #endif
00395
00396 #endif
00397