ITK  6.0.0
Insight Toolkit
itkLinearInterpolateImageFunction.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright NumFOCUS
4 *
5 * Licensed under the Apache License, Version 2.0 (the "License");
6 * you may not use this file except in compliance with the License.
7 * You may obtain a copy of the License at
8 *
9 * https://www.apache.org/licenses/LICENSE-2.0.txt
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 *
17 *=========================================================================*/
18#ifndef itkLinearInterpolateImageFunction_h
19#define itkLinearInterpolateImageFunction_h
20
23#include <algorithm> // For max.
24
25namespace itk
26{
50template <typename TInputImage, typename TCoordRep = double>
51class ITK_TEMPLATE_EXPORT LinearInterpolateImageFunction : public InterpolateImageFunction<TInputImage, TCoordRep>
52{
53public:
54 ITK_DISALLOW_COPY_AND_MOVE(LinearInterpolateImageFunction);
55
61
63 itkOverrideGetNameOfClassMacro(LinearInterpolateImageFunction);
64
66 itkNewMacro(Self);
67
69 using typename Superclass::OutputType;
70
72 using typename Superclass::InputImageType;
73
75 using typename Superclass::InputPixelType;
76
78 using typename Superclass::RealType;
79
81 static constexpr unsigned int ImageDimension = Superclass::ImageDimension;
82
84 using typename Superclass::IndexType;
85
87 using typename Superclass::SizeType;
88
90 using typename Superclass::ContinuousIndexType;
92
102 EvaluateAtContinuousIndex(const ContinuousIndexType & index) const override
103 {
104 return this->EvaluateOptimized(Dispatch<ImageDimension>(), index);
105 }
106
108 GetRadius() const override
109 {
110 return SizeType::Filled(1);
111 }
112
113protected:
116 void
117 PrintSelf(std::ostream & os, Indent indent) const override;
118
119private:
121 {};
122 template <unsigned int>
123 struct Dispatch : public DispatchBase
124 {};
125
126 inline OutputType
128 {
129 return 0;
130 }
131
132 inline OutputType
134 {
135 IndexType basei;
136 basei[0] = std::max(Math::Floor<IndexValueType>(index[0]), this->m_StartIndex[0]);
137
138 const InternalComputationType distance = index[0] - static_cast<InternalComputationType>(basei[0]);
139
140 const TInputImage * const inputImagePtr = this->GetInputImage();
141 const RealType & val0 = inputImagePtr->GetPixel(basei);
142 if (distance <= 0.)
143 {
144 return (static_cast<OutputType>(val0));
145 }
146
147 ++basei[0];
148 if (basei[0] > this->m_EndIndex[0])
149 {
150 return (static_cast<OutputType>(val0));
151 }
152 const RealType & val1 = inputImagePtr->GetPixel(basei);
153
154 return (static_cast<OutputType>(val0 + (val1 - val0) * distance));
155 }
156
157 inline OutputType
159 {
160 IndexType basei;
161
162 basei[0] = std::max(Math::Floor<IndexValueType>(index[0]), this->m_StartIndex[0]);
163 const InternalComputationType distance0 = index[0] - static_cast<InternalComputationType>(basei[0]);
164
165 basei[1] = std::max(Math::Floor<IndexValueType>(index[1]), this->m_StartIndex[1]);
166 const InternalComputationType distance1 = index[1] - static_cast<InternalComputationType>(basei[1]);
167
168 const TInputImage * const inputImagePtr = this->GetInputImage();
169 const RealType & val00 = inputImagePtr->GetPixel(basei);
170 if (distance0 <= 0. && distance1 <= 0.)
171 {
172 return (static_cast<OutputType>(val00));
173 }
174 else if (distance1 <= 0.) // if they have the same "y"
175 {
176 ++basei[0]; // then interpolate across "x"
177 if (basei[0] > this->m_EndIndex[0])
178 {
179 return (static_cast<OutputType>(val00));
180 }
181 const RealType & val10 = inputImagePtr->GetPixel(basei);
182 return (static_cast<OutputType>(val00 + (val10 - val00) * distance0));
183 }
184 else if (distance0 <= 0.) // if they have the same "x"
185 {
186 ++basei[1]; // then interpolate across "y"
187 if (basei[1] > this->m_EndIndex[1])
188 {
189 return (static_cast<OutputType>(val00));
190 }
191 const RealType & val01 = inputImagePtr->GetPixel(basei);
192 return (static_cast<OutputType>(val00 + (val01 - val00) * distance1));
193 }
194 // fall-through case:
195 // interpolate across "xy"
196 ++basei[0];
197 if (basei[0] > this->m_EndIndex[0]) // interpolate across "y"
198 {
199 --basei[0];
200 ++basei[1];
201 if (basei[1] > this->m_EndIndex[1])
202 {
203 return (static_cast<OutputType>(val00));
204 }
205 const RealType & val01 = inputImagePtr->GetPixel(basei);
206 return (static_cast<OutputType>(val00 + (val01 - val00) * distance1));
207 }
208 const RealType & val10 = inputImagePtr->GetPixel(basei);
209
210 const RealType & valx0 = val00 + (val10 - val00) * distance0;
211
212 ++basei[1];
213 if (basei[1] > this->m_EndIndex[1]) // interpolate across "x"
214 {
215 return (static_cast<OutputType>(valx0));
216 }
217 const RealType & val11 = inputImagePtr->GetPixel(basei);
218 --basei[0];
219 const RealType & val01 = inputImagePtr->GetPixel(basei);
220
221 const RealType & valx1 = val01 + (val11 - val01) * distance0;
222
223 return (static_cast<OutputType>(valx0 + (valx1 - valx0) * distance1));
224 }
225
226 inline OutputType
228 {
229 IndexType basei;
230 basei[0] = std::max(Math::Floor<IndexValueType>(index[0]), this->m_StartIndex[0]);
231 const InternalComputationType distance0 = index[0] - static_cast<InternalComputationType>(basei[0]);
232
233 basei[1] = std::max(Math::Floor<IndexValueType>(index[1]), this->m_StartIndex[1]);
234 const InternalComputationType distance1 = index[1] - static_cast<InternalComputationType>(basei[1]);
235
236 basei[2] = std::max(Math::Floor<IndexValueType>(index[2]), this->m_StartIndex[2]);
237 const InternalComputationType distance2 = index[2] - static_cast<InternalComputationType>(basei[2]);
238
239 const TInputImage * const inputImagePtr = this->GetInputImage();
240 const RealType & val000 = inputImagePtr->GetPixel(basei);
241 if (distance0 <= 0. && distance1 <= 0. && distance2 <= 0.)
242 {
243 return (static_cast<OutputType>(val000));
244 }
245
246 if (distance2 <= 0.)
247 {
248 if (distance1 <= 0.) // interpolate across "x"
249 {
250 ++basei[0];
251 if (basei[0] > this->m_EndIndex[0])
252 {
253 return (static_cast<OutputType>(val000));
254 }
255 const RealType & val100 = inputImagePtr->GetPixel(basei);
256
257 return static_cast<OutputType>(val000 + (val100 - val000) * distance0);
258 }
259 else if (distance0 <= 0.) // interpolate across "y"
260 {
261 ++basei[1];
262 if (basei[1] > this->m_EndIndex[1])
263 {
264 return (static_cast<OutputType>(val000));
265 }
266 const RealType & val010 = inputImagePtr->GetPixel(basei);
267
268 return static_cast<OutputType>(val000 + (val010 - val000) * distance1);
269 }
270 else // interpolate across "xy"
271 {
272 ++basei[0];
273 if (basei[0] > this->m_EndIndex[0]) // interpolate across "y"
274 {
275 --basei[0];
276 ++basei[1];
277 if (basei[1] > this->m_EndIndex[1])
278 {
279 return (static_cast<OutputType>(val000));
280 }
281 const RealType & val010 = inputImagePtr->GetPixel(basei);
282 return static_cast<OutputType>(val000 + (val010 - val000) * distance1);
283 }
284 const RealType & val100 = inputImagePtr->GetPixel(basei);
285 const RealType & valx00 = val000 + (val100 - val000) * distance0;
286
287 ++basei[1];
288 if (basei[1] > this->m_EndIndex[1]) // interpolate across "x"
289 {
290 return (static_cast<OutputType>(valx00));
291 }
292 const RealType & val110 = inputImagePtr->GetPixel(basei);
293
294 --basei[0];
295 const RealType & val010 = inputImagePtr->GetPixel(basei);
296 const RealType & valx10 = val010 + (val110 - val010) * distance0;
297
298 return static_cast<OutputType>(valx00 + (valx10 - valx00) * distance1);
299 }
300 }
301 else
302 {
303 if (distance1 <= 0.)
304 {
305 if (distance0 <= 0.) // interpolate across "z"
306 {
307 ++basei[2];
308 if (basei[2] > this->m_EndIndex[2])
309 {
310 return (static_cast<OutputType>(val000));
311 }
312 const RealType & val001 = inputImagePtr->GetPixel(basei);
313
314 return static_cast<OutputType>(val000 + (val001 - val000) * distance2);
315 }
316 else // interpolate across "xz"
317 {
318 ++basei[0];
319 if (basei[0] > this->m_EndIndex[0]) // interpolate across "z"
320 {
321 --basei[0];
322 ++basei[2];
323 if (basei[2] > this->m_EndIndex[2])
324 {
325 return (static_cast<OutputType>(val000));
326 }
327 const RealType & val001 = inputImagePtr->GetPixel(basei);
328
329 return static_cast<OutputType>(val000 + (val001 - val000) * distance2);
330 }
331 const RealType & val100 = inputImagePtr->GetPixel(basei);
332
333 const RealType & valx00 = val000 + (val100 - val000) * distance0;
334
335 ++basei[2];
336 if (basei[2] > this->m_EndIndex[2]) // interpolate across "x"
337 {
338 return (static_cast<OutputType>(valx00));
339 }
340 const RealType & val101 = inputImagePtr->GetPixel(basei);
341
342 --basei[0];
343 const RealType & val001 = inputImagePtr->GetPixel(basei);
344
345 const RealType & valx01 = val001 + (val101 - val001) * distance0;
346
347 return static_cast<OutputType>(valx00 + (valx01 - valx00) * distance2);
348 }
349 }
350 else if (distance0 <= 0.) // interpolate across "yz"
351 {
352 ++basei[1];
353 if (basei[1] > this->m_EndIndex[1]) // interpolate across "z"
354 {
355 --basei[1];
356 ++basei[2];
357 if (basei[2] > this->m_EndIndex[2])
358 {
359 return (static_cast<OutputType>(val000));
360 }
361 const RealType & val001 = inputImagePtr->GetPixel(basei);
362
363 return static_cast<OutputType>(val000 + (val001 - val000) * distance2);
364 }
365 const RealType & val010 = inputImagePtr->GetPixel(basei);
366
367 const RealType & val0x0 = val000 + (val010 - val000) * distance1;
368
369 ++basei[2];
370 if (basei[2] > this->m_EndIndex[2]) // interpolate across "y"
371 {
372 return (static_cast<OutputType>(val0x0));
373 }
374 const RealType & val011 = inputImagePtr->GetPixel(basei);
375
376 --basei[1];
377 const RealType & val001 = inputImagePtr->GetPixel(basei);
378
379 const RealType & val0x1 = val001 + (val011 - val001) * distance1;
380
381 return static_cast<OutputType>(val0x0 + (val0x1 - val0x0) * distance2);
382 }
383 else // interpolate across "xyz"
384 {
385 ++basei[0];
386 if (basei[0] > this->m_EndIndex[0]) // interpolate across "yz"
387 {
388 --basei[0];
389 ++basei[1];
390 if (basei[1] > this->m_EndIndex[1]) // interpolate across "z"
391 {
392 --basei[1];
393 ++basei[2];
394 if (basei[2] > this->m_EndIndex[2])
395 {
396 return (static_cast<OutputType>(val000));
397 }
398 const RealType & val001 = inputImagePtr->GetPixel(basei);
399
400 return static_cast<OutputType>(val000 + (val001 - val000) * distance2);
401 }
402 const RealType & val010 = inputImagePtr->GetPixel(basei);
403 const RealType & val0x0 = val000 + (val010 - val000) * distance1;
404
405 ++basei[2];
406 if (basei[2] > this->m_EndIndex[2]) // interpolate across "y"
407 {
408 return (static_cast<OutputType>(val0x0));
409 }
410 const RealType & val011 = inputImagePtr->GetPixel(basei);
411
412 --basei[1];
413 const RealType & val001 = inputImagePtr->GetPixel(basei);
414
415 const RealType & val0x1 = val001 + (val011 - val001) * distance1;
416
417 return static_cast<OutputType>(val0x0 + (val0x1 - val0x0) * distance2);
418 }
419 const RealType & val100 = inputImagePtr->GetPixel(basei);
420
421 const RealType & valx00 = val000 + (val100 - val000) * distance0;
422
423 ++basei[1];
424 if (basei[1] > this->m_EndIndex[1]) // interpolate across "xz"
425 {
426 --basei[1];
427 ++basei[2];
428 if (basei[2] > this->m_EndIndex[2]) // interpolate across "x"
429 {
430 return (static_cast<OutputType>(valx00));
431 }
432 const RealType & val101 = inputImagePtr->GetPixel(basei);
433
434 --basei[0];
435 const RealType & val001 = inputImagePtr->GetPixel(basei);
436
437 const RealType & valx01 = val001 + (val101 - val001) * distance0;
438
439 return static_cast<OutputType>(valx00 + (valx01 - valx00) * distance2);
440 }
441 const RealType & val110 = inputImagePtr->GetPixel(basei);
442
443 --basei[0];
444 const RealType & val010 = inputImagePtr->GetPixel(basei);
445
446 const RealType & valx10 = val010 + (val110 - val010) * distance0;
447
448 const RealType & valxx0 = valx00 + (valx10 - valx00) * distance1;
449
450 ++basei[2];
451 if (basei[2] > this->m_EndIndex[2]) // interpolate across "xy"
452 {
453 return (static_cast<OutputType>(valxx0));
454 }
455 const RealType & val011 = inputImagePtr->GetPixel(basei);
456
457 ++basei[0];
458 const RealType & val111 = inputImagePtr->GetPixel(basei);
459
460 --basei[1];
461 const RealType & val101 = inputImagePtr->GetPixel(basei);
462
463 --basei[0];
464 const RealType & val001 = inputImagePtr->GetPixel(basei);
465
466 const RealType & valx01 = val001 + (val101 - val001) * distance0;
467 const RealType & valx11 = val011 + (val111 - val011) * distance0;
468 const RealType & valxx1 = valx01 + (valx11 - valx01) * distance1;
469
470 return (static_cast<OutputType>(valxx0 + (valxx1 - valxx0) * distance2));
471 }
472 }
473 }
474
475 inline OutputType
477 {
478 return this->EvaluateUnoptimized(index);
479 }
480
482 virtual inline OutputType
484
487 template <typename RealTypeScalarRealType>
488 void
489 MakeZeroInitializer(const TInputImage * const inputImagePtr,
491 {
492 // Variable length vector version to get the size of the pixel correct.
493 constexpr typename TInputImage::IndexType idx = { { 0 } };
494 const typename TInputImage::PixelType & tempPixel = inputImagePtr->GetPixel(idx);
495 const unsigned int sizeOfVarLengthVector = tempPixel.GetSize();
496 tempZeros.SetSize(sizeOfVarLengthVector);
497 tempZeros.Fill(RealTypeScalarRealType{});
498 }
499
500 template <typename RealTypeScalarRealType>
501 void
502 MakeZeroInitializer(const TInputImage * const itkNotUsed(inputImagePtr), RealTypeScalarRealType & tempZeros) const
503 {
504 // All other cases
505 tempZeros = RealTypeScalarRealType{};
506 }
507};
508} // end namespace itk
509
510#ifndef ITK_MANUAL_INSTANTIATION
511# include "itkLinearInterpolateImageFunction.hxx"
512#endif
513
514#endif
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Base class for all image interpolators.
typename NumericTraits< typename TInputImage::PixelType >::RealType RealType
Light weight base class for most itk classes.
Linearly interpolate an image at specified positions.
void PrintSelf(std::ostream &os, Indent indent) const override
~LinearInterpolateImageFunction() override=default
OutputType EvaluateOptimized(const Dispatch< 2 > &, const ContinuousIndexType &index) const
void MakeZeroInitializer(const TInputImage *const inputImagePtr, VariableLengthVector< RealTypeScalarRealType > &tempZeros) const
A method to generically set all components to zero.
OutputType EvaluateOptimized(const Dispatch< 3 > &, const ContinuousIndexType &index) const
void MakeZeroInitializer(const TInputImage *const, RealTypeScalarRealType &tempZeros) const
virtual OutputType EvaluateUnoptimized(const ContinuousIndexType &index) const
OutputType EvaluateOptimized(const Dispatch< 0 > &, const ContinuousIndexType &) const
typename ContinuousIndexType::ValueType InternalComputationType
OutputType EvaluateOptimized(const DispatchBase &, const ContinuousIndexType &index) const
OutputType EvaluateOptimized(const Dispatch< 1 > &, const ContinuousIndexType &index) const
OutputType EvaluateAtContinuousIndex(const ContinuousIndexType &index) const override
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....