# Resampling to Isotropic: Signal Processing Theory
**URL:** https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403
**Category:** Algorithms
**Created:** [November 7, 2018, 8:04pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403 "2018-11-07T20:04:33Z")
**Posts on this page:** 19
**Page:** 1
### Author:  [@Besler](https://discourse.itk.org/u/Besler)
#### Post date: [November 7, 2018, 8:04pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/1 "2018-11-07T20:04:33Z")
Hi ITK Community!
Sanity check here about Nyquist-Shannon sampling theorem and resampling images to be isotropic. I have read through the online example [1]. I have also read the sinc approximation paper [2].
**Setup** I have CT data that is captured at 0.5 x 0.5 x 2.0 mm resolution. I want my data in 1.0 x 1.0 x 1.0 mm resolution for quantification purposes (FEA, actually).
**Theory** Following classsic DSP theory, I would do the following.
I apply a low pass filter (say Gaussian) in the X and Y directions at a spatial frequency of 0.5 mm^-1 (Sampling frequency is 1 mm^-1 so Nyquist frequency is 0.5 mm^-1).
Then, I would drop samples in X and Y and zero order hold in Z.
Then, I would apply a low pass filter in Z at 0.5 mm^-1.
This all makes sense from a classic DSP approach. I avoid aliasing in X and Y with prefiltering and I get a reconstructed filter with postfiltering in Z.
**Interpolation Filter** The problem arrises with the concept of interpolation in the ITK Resample filter.
In the online example [1], X and Y are blurred with a Gaussian filter and a linear interpolator is used for resampling. Let’s say I use a Hamming Windowed Sinc interpolator. I map a point from the isotropic image to the fixed image, sample by convolving the Hamming Windowed Sinc, and store that value. This effectively causes blurring in my input image, similar to the effect of prefiltering.
**Question** If I use a proper interpolation filter (3rd order BSpline, windowed sinc) I do not have to prefilter my image. I will have minimal aliasing because I am effectively prefiltering my image in the interpolator. Is this correct?
[1] [https://itk.org/Doxygen/html/Examples\_2Filtering\_2ResampleVolumesToBeIsotropic\_8cxx-example.html](https://itk.org/Doxygen/html/Examples_2Filtering_2ResampleVolumesToBeIsotropic_8cxx-example.html)
[2] [https://www.researchgate.net/profile/WJ\_Niessen/publication/221401597\_Quantitative\_Comparison\_of\_Sinc-Approximating\_Kernels\_for\_Medical\_Image\_Interpolation/links/0deec525e649cb45a6000000/Quantitative-Comparison-of-Sinc-Approximating-Kernels-for-Medical-Image-Interpolation.pdf](https://www.researchgate.net/profile/WJ_Niessen/publication/221401597_Quantitative_Comparison_of_Sinc-Approximating_Kernels_for_Medical_Image_Interpolation/links/0deec525e649cb45a6000000/Quantitative-Comparison-of-Sinc-Approximating-Kernels-for-Medical-Image-Interpolation.pdf)
---
### Author:  [@phcerdan](https://discourse.itk.org/u/phcerdan)
#### Post date: [November 7, 2018, 8:26pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/2 "2018-11-07T20:26:47Z")
> [@Besler](#):
>
> **Question** If I use a proper interpolation filter (3rd order BSpline, windowed sinc) I do not have to prefilter my image
I think so.
Have a look at: [LesionSizingToolkit/include/itkIsotropicResamplerImageFilter.hxx at master · InsightSoftwareConsortium/LesionSizingToolkit · GitHub](https://github.com/InsightSoftwareConsortium/LesionSizingToolkit/blob/master/include/itkIsotropicResamplerImageFilter.hxx)
---
### Author:  [@Besler](https://discourse.itk.org/u/Besler)
#### Post date: [November 7, 2018, 8:31pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/3 "2018-11-07T20:31:37Z")
Thank you for the example! Here they’re doing 3rd order BSpline interpolation, no prefiltering. This agrees with my understanding.
Thank You!
---
### Author:  [@warfield](https://discourse.itk.org/u/warfield)
#### Post date: [November 7, 2018, 11:10pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/4 "2018-11-07T23:10:14Z")
A good explanation of the difference between prefiltering and interpolating, and what has been called ‘generalized interpolation’ is available here:
[http://bigwww.epfl.ch/publications/thevenaz9901.pdf](http://bigwww.epfl.ch/publications/thevenaz9901.pdf)
ITK implements a version of the code that is here:
[http://bigwww.epfl.ch/thevenaz/interpolation/](http://bigwww.epfl.ch/thevenaz/interpolation/)
and is also explained in this paper:
Unser et al. [https://doi-org/10.1109/78.193220](https://doi-org/10.1109/78.193220)
---
### Author:  [@phcerdan](https://discourse.itk.org/u/phcerdan)
#### Post date: [November 8, 2018, 12:02am UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/5 "2018-11-08T00:02:19Z")
What a good and complete description [http://bigwww.epfl.ch/publications/thevenaz9901.pdf](http://bigwww.epfl.ch/publications/thevenaz9901.pdf)
Thanks @warfield.
We should create a sticky post holding high-quality papers/reviews like that on important topics.
---
### Author:  [@Besler](https://discourse.itk.org/u/Besler)
#### Post date: [November 8, 2018, 5:10am UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/6 "2018-11-08T05:10:41Z")
Wow, thank you @warfield. That was a journey.
I’m going to check my understanding here.
**Interpolation** Represent the continuous signal as the convolution of the discrete signal with an interpolation function.
**Generalized Interpolation** Represent the continuous signal as the convolution of a coefficient image and basis functions. The distinction is to calculate the coefficients from the image _before_ interpolating (don’t assume the coefficients are the discrete signal).
**Aliasing**
From what I’m reading, antialias depends on your implementation. If decimation is implemented (Unser et al. [https://doi-org/10.1109/78.193220](https://doi-org/10.1109/78.193220), section D. Least Squares Fitting), then this acts like the antialiasing filter.
So, in general, antialias filtering is needed. Interpolation provides a continuous representation for the discrete signal. Since you want to store that signal with fewer samples, you need to remove frequency content above the Nyquest sample frequency before building a representation (i.e. before resampling).
For the implementation in ITK, antialiasing filter is needed as per the example.
This raises concerns about how I’ve been doing my registrations! I suppose you should always register to a high resolution image so you don’t need antialias filtering, among other reasons (fidelity, etc.)
P.S. I was mixing up the terms prefilteringand antialias filtering. Antialias filtering means to remove frequency content for discrete representation. Prefiltering, in the context of interpolation, means to determine the coefficients.
---
### Author:  [@grothausmann.roman](https://discourse.itk.org/u/grothausmann.roman)
#### Post date: [December 13, 2019, 12:46pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/7 "2019-12-13T12:46:25Z")
I was about to post a very similar topic when I found this discussion.
Many thanks for all the informative feedback and the pointers for further reading.
This leaves me with 3 remaining questions:
Regarding
> For the implementation in ITK, antialiasing filter is needed as per the example.
does that mean that even when using [IsotropicResamplerImageFilter](https://itk.org/Doxygen50/html/classitk_1_1IsotropicResamplerImageFilter.html) instead of the [ResampleImageFilter](https://itk.org/Doxygen50/html/classitk_1_1ResampleImageFilter.html) (as used in [ResampleVolumesToBeIsotropic](https://itk.org/Doxygen/html/Examples_2Filtering_2ResampleVolumesToBeIsotropic_8cxx-example.html)) it is still essential to apply low pass filtering (“remove frequency content above the Nyquest sample frequency”) for the spatial directions to be down-sampled? As far as I understand the code of [IsotropicResamplerImageFilter](https://github.com/InsightSoftwareConsortium/LesionSizingToolkit/blob/master/include/itkIsotropicResamplerImageFilter.hxx) it mainly only differs by using BSpline interpolation instead of linear interpolation as in [ResampleVolumesToBeIsotropic](https://itk.org/Doxygen/html/Examples_2Filtering_2ResampleVolumesToBeIsotropic_8cxx-example.html). Is that correct or am I missing something here?
So in other words it would NOT be appropriate to remove the low pass filtering from [ResampleVolumesToBeIsotropic](https://itk.org/Doxygen/html/Examples_2Filtering_2ResampleVolumesToBeIsotropic_8cxx-example.html) when replacing the linear interpolation with BSpline interpolation?
> ITK implements a version of the code that is here: [http://bigwww.epfl.ch/thevenaz/interpolation/](http://bigwww.epfl.ch/thevenaz/interpolation/)
The ITK implementation referred to here is [IsotropicResamplerImageFilter](https://github.com/InsightSoftwareConsortium/LesionSizingToolkit/blob/master/include/itkIsotropicResamplerImageFilter.hxx) or [ResampleVolumesToBeIsotropic](https://itk.org/Doxygen/html/Examples_2Filtering_2ResampleVolumesToBeIsotropic_8cxx-example.html)?
And my original question:
What would be the correct way to generalize e.g. [ResampleVolumesToBeIsotropic](https://itk.org/Doxygen/html/Examples_2Filtering_2ResampleVolumesToBeIsotropic_8cxx-example.html) to resample arbitrary voxel spacing relations? For example if the final intended voxel spacing in x is close to the input voxel size in x I would think that the blurring should depend on the difference of the two. [IsotropicResamplerImageFilter](https://github.com/InsightSoftwareConsortium/LesionSizingToolkit/blob/master/include/itkIsotropicResamplerImageFilter.hxx) only seems to handle the case in which the desired output spacing equals the input spacing. E.g. in the case of a 1:10 voxel size ratio and a desired iso output spacing equaling to the z-spacing: What low pass filter should be used (Gauss, Hamming, etc) for xy and how to calculate appropriate parameters for the filter (e.g. sigma)?
---
### Author:  [@dzenanz](https://discourse.itk.org/u/dzenanz)
#### Post date: [December 13, 2019, 3:17pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/8 "2019-12-13T15:17:42Z")
> [@grothausmann.roman](#):
>
> ITK implements a version of the code that is here: [Spline Interpolation](http://bigwww.epfl.ch/thevenaz/interpolation/)
As [IsotropicResamplerImageFilter](https://github.com/InsightSoftwareConsortium/LesionSizingToolkit/blob/master/include/itkIsotropicResamplerImageFilter.hxx) is in a remote module, “ITK implements” probably does not refer to it.
Looking at the code of [IsotropicResamplerImageFilter](https://github.com/InsightSoftwareConsortium/LesionSizingToolkit/blob/master/include/itkIsotropicResamplerImageFilter.hxx), it does not include smoothing so it needs to be applied separately.
---
### Author:  [@grothausmann.roman](https://discourse.itk.org/u/grothausmann.roman)
#### Post date: [December 16, 2019, 10:23am UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/9 "2019-12-16T10:23:21Z")
Thanks @dzenanz for your reply.
So I would conclude that “ITK implements” refers to the BSpline interpolation implementation of ITK and not to [ResampleVolumesToBeIsotropic](https://itk.org/Doxygen/html/Examples_2Filtering_2ResampleVolumesToBeIsotropic_8cxx-example.html) (which seems to originate from [@luisibanez from 2005](https://github.com/InsightSoftwareConsortium/ITK/blob/dd76165be51ee20a5eeba90126c3a51924a43e98/Examples/Filtering/ResampleVolumesToBeIsotropic.cxx) and was ever since using linear interpolation).
@warfield is that what you meant?
---
### Author:  [@matt.mccormick](https://discourse.itk.org/u/matt.mccormick)
#### Post date: [December 16, 2019, 4:53pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/10 "2019-12-16T16:53:50Z")
To summarize on a high level:
- An **anti-aliasing filter** , e.g. Gaussian smoothing, is required to remove high frequencies that will cause artifacts when **downsampling**.
- An **accurate interpolator** , e.g. windowed-sinc interpolator, is required to avoid artifacts when **upsampling**.
---
### Author:  [@warfield](https://discourse.itk.org/u/warfield)
#### Post date: [December 16, 2019, 5:07pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/11 "2019-12-16T17:07:28Z")
I meant that the BSpline interpolation uses an algorithmic strategy and code that is based on what Thevanez originally published.
---
### Author:  [@grothausmann.roman](https://discourse.itk.org/u/grothausmann.roman)
#### Post date: [December 17, 2019, 12:36pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/12 "2019-12-17T12:36:12Z")
Many thanks for these clearifications @matt.mccormick and @warfield
> [@matt.mccormick](#):
>
> An **anti-aliasing filter** , e.g. Gaussian smoothing, is required to remove high frequencies that will cause artifacts when **downsampling**.
So what I wonder then is how to choose the correct sigma for Gaussian smoothing when downsampling. [ResampleVolumesToBeIsotropic](https://github.com/InsightSoftwareConsortium/ITK/blob/246c02bdad8c31e79f1ef7c9085f4eace30eafd0/Examples/Filtering/ResampleVolumesToBeIsotropic.cxx) uses the new pixel size directly:
>
However, I would think that is not a very appropriate choice in case the original pixel size is very similar, e.g. for downsampling from 0.9mm to 1mm pixel size. Would it be better to choose the difference of the original and new pixel size for sigma or would that be too small? Or should it be a function that depends on the difference?
> [@matt.mccormick](#):
>
> An **accurate interpolator** , e.g. windowed-sinc interpolator, is required to avoid artifacts when **upsampling**
Despite the argument that sophistication is not appropriate in the case of the example
>
and therefore linear interpolation is used
>
it might be better to change the example to use BSpline or windowed-sinc interpolator to be a better role model in general, as it seems this example is uses as a base for more general cases.
---
### Author:  [@dzenanz](https://discourse.itk.org/u/dzenanz)
#### Post date: [December 17, 2019, 2:11pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/13 "2019-12-17T14:11:42Z")
> [@grothausmann.roman](#):
>
> choose the difference of the original and new pixel size for sigma
This sounds much better to me than the target spacing. But ff some of the spacings are already close to your target, maybe it is better not to resample along that dimension. E.g. instead of resampling 0.9x0.9x1.2 into 1x1x1, maybe resample it into 0.9x0.9x0.9 or 1.2x1.2x1.2, or something of that sort.
Windowed-sinc interpolator is quite compute intensive, and is maybe not appropriate for on-the-fly execution. But it is probably the best choice for one-off resampling of all kinds. And as computation power has grown much more than image resolutions since that example was originally written, perhaps it is time to update it to use BSpline or windowed-sinc interpolator. Do you want to do it @grothausmann.roman?
---
### Author:  [@matt.mccormick](https://discourse.itk.org/u/matt.mccormick)
#### Post date: [December 17, 2019, 2:56pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/14 "2019-12-17T14:56:43Z")
> [@grothausmann.roman](#):
>
> However, I would think that is not a very appropriate choice in case the original pixel size is very similar, e.g. for downsampling from 0.9mm to 1mm pixel size. Would it be better to choose the difference of the original and new pixel size for sigma or would that be too small?
Ideally, sigma should be specified in each direction based on the physical sampling distance in the output image.
---
### Author:  [@grothausmann.roman](https://discourse.itk.org/u/grothausmann.roman)
#### Post date: [February 25, 2020, 10:55am UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/15 "2020-02-25T10:55:20Z")
Thanks for all your feedback, perhaps we can discuss this easier when meeting at WBIR2020.
---
### Author:  [@grothausmann.roman](https://discourse.itk.org/u/grothausmann.roman)
#### Post date: [June 10, 2020, 4:09pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/16 "2020-06-10T16:09:31Z")
Since WBIR2020 is postponed, continuing the discussion here.
I found this article from Cardoso M.J., Modat M., Vercauteren T., Ourselin S. (2015) Scale Factor Point Spread Function Matching: Beyond Aliasing in Image Resampling – MICCAI 2015 [https://doi.org/10.1007/978-3-319-24571-3\_81](https://doi.org/10.1007/978-3-319-24571-3_81) [PDF](http://discovery.ucl.ac.uk/1469251/1/scale-factor-point-5.pdf) which provides the formula `sigma= sqrt((isoSpacing^2 - inputSpacing[0]^2)/(2*sqrt(2*ln(2)))^2)`, which is not a function of the difference but based on difference and also leads to sigma= 0 in case input and output spacing are equal.
@dzenanz @matt.mccormick What do you think? Is its application appropriate in the discussed case, i.e. is its use in my [resample-iso](https://github.com/romangrothausmann/ITK-CLIs/blob/3165433931a5d75d7014567920fad136fd4555a0/resample-iso.cxx#L65-L146) correct?
---
### Author:  [@matt.mccormick](https://discourse.itk.org/u/matt.mccormick)
#### Post date: [June 10, 2020, 5:34pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/17 "2020-06-10T17:34:43Z")
@grothausmann.roman 👍 looks great.
---
### Author:  [@lassoan](https://discourse.itk.org/u/lassoan)
#### Post date: [October 15, 2020, 6:49pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/18 "2020-10-15T18:49:32Z")
David Gobbi implemented an [image resize filter in VTK that automatically applies anti-aliasing and interpolation as needed](https://vtk.org/doc/nightly/html/classvtkImageReslice.html#details) (it computes optimal anti-aliasing filter parameters). It would be nice if ITK would offer a similar feature. Ideally, it would be just a flag in itk::ResampleImageFilter to enable/disable anti-aliasing (to choose between higher-quality vs. faster resize).
---
### Author:  [@dzenanz](https://discourse.itk.org/u/dzenanz)
#### Post date: [October 15, 2020, 10:27pm UTC](https://discourse.itk.org/t/resampling-to-isotropic-signal-processing-theory/1403/19 "2020-10-15T22:27:02Z")
I still have this on my todo list. But it is now likely I will get to this in the next few months.