[Insight-developers] Proposed requirements for convolution filters in ITK

Thu Jun 2 14:05:04 EDT 2011

Hi Dirk,

Thanks for your thoughtful comments. My responses are below.

> "The output image should be the same size as the first input image."
> Perhaps the main reason for this is for easy streaming.  But in my
> implementations, I generally set the output to be the full convolution
> size (N+M-1,N+M-1), where N is the size of image 1, and M is the size of
> image 2.  As long as we have an accessor for this image and perhaps one
> for only the valid region, I think that is fine.  We should keep in mind
> that if the images are the same size, the valid region is only 1 pixel.

Would having a "mode" that determines the output region of the image
obtained via GetOutput() be acceptable? Or were you thinking of having
different methods such as GetValidRegionOutput(), GetFullOuput(),
etc...

> "Treatment of regions beyond the boundaries of the input images should
> be controllable by the user."
> I believe that no handling of boundaries is necessary because the
> algorithm should simply calculate the convolution in the overlap region
> as the images are convolved.  This can be done explicitly in the spatial
> domain, and it can be accomplished in the Fourier domain by padding both
> images first with zeros so that their size is (N+M-1,N+M-1), where N is
> the size of image 1, and M is the size of image 2.

I could have said what I meant more clearly. Consider the process of
calculating the inner product of the kernel with the image at a pixel
that is less than half the kernel radius away from a boundary in some
dimension. Some kernel pixels will fall outside the image boundaries.
What I meant to say was that the user should be able to specify what
values should be used in place of the image pixel values for those
kernel pixels. Zero is one option, but there are others such as those
I cited in my original email that may be more desirable for certain
applications. As an interesting side note, the current
ConvolutionImageFilter is implicitly padding with the ZeroFluxNeumann
boundary condition, but it takes a little digging to discover that.

The convolution filters should also be a "no padding" option to handle
the case where a user takes care of the padding herself with a method
not offered by the convolution filter. In that case, the spatial
domain filter should probably perform a cyclic convolution to produce
the same results as the Fourier domain filter, and an exception should
be thrown if the inputs images are not the same size.

> "Normalization of the kernel image should be an option. Default is off."
> Do you mean normalization to have a sum of one?

Yes.

> I agree that it should
> be off by default.  For example, if the user wants to convolve two
> images, it would be confusing if the algorithm normalizes the second
> image.

That's a good point.

>The user should have to explicitly specify normalization or do
> the normalization himself/herself beforehand.

Now that I've though about it longer, I am in favor of letting the
user take care of normalization himself. It would make the convolution
filter code somewhat less complex (and require less testing). Adding
normalization wouldn't be hard to add in an application.

> "VNL operates only on images whose sizes are powers of two."
> This is not entirely true; VNL requires that the image size be multiples
> of 2s, 3s, and 5s, and this is far less restricting (and much faster).
> It also turns out that, even though FFTW handles any size, it is much
> faster when the images are multiples of 2s, 3s, and 5s.  In my
> implementation, I wrote a function that calculates the next largest
> image size in each dimension that is a multiple of 2s, 3s, and 5s and
> then pads the image with zeros before doing the Fourier domain
> computation.  Thus, the padding size should not be a parameter.

Thanks for the clarification about VNL's restrictions. I agree that
the padding size should all be done internally.

> At the
> end, it crops the correlation map back to the correct size.  Then the
> result can be compared with the spatial domain version since the spatial
> domain doesn't require any padding.

Yes. I was a little confused when I wrote that, but as long as the
boundary conditions are handled the same and the Fourier domain result
is cropped, comparison shouldn't be an issue.

> I have already implemented a filter that addresses most of these
> concerns in the context of Normalized Cross-Correlation in the Fourier
> domain using VNL.  It is based on my "Masked FFT Registration" CVPR 2010
> paper.  Many parts of that code should be relevant for the convolution
> filter, and I can share it with you if that is helpful.

Please do. Also, if you have a chance, take a look at Gaetan's "FFT
based convolution' article in the Insight Journal
(http://hdl.handle.net/10380/3154). Once we arrive at a consensus as
to what the features the convolution filters should support, I am
planning on modifying Gaetan's contribution accordingly and submitting
that as a gerrit patch.

Thanks,
Cory

>
>
> Date: Wed, 1 Jun 2011 11:17:02 -0400
> From: Cory Quammen <cquammen at cs.unc.edu>
> Subject: [Insight-developers] Proposed requirements for convolution
>        filters in ITK
> To: Insight Developers <insight-developers at itk.org>,    Nicholas
> Tustison
>        <ntustison at gmail.com>,  Luis Ibanez <luis.ibanez at kitware.com>
> Message-ID: <BANLkTikVqQuHX7ibA_zW0sXQz7H+0vMvWQ at mail.gmail.com>
> Content-Type: text/plain; charset=ISO-8859-1
>
> Here is what I would like to see in the convolution filters in ITK.
>
> - Convolution filters should have the same interface to make it easy
> to switch between the spatial domain and Fourier domain methods. This
> suggests a class hierarchy for image convolution filters, perhaps
> organized as follows:
>
> --ConvolutionImageFilter
> ---SpatialDomainConvolutionImageFilter
> ---FourierDomainConvolutionImageFilter
>
> - Each algorithm should take two inputs. The first is treated as the
> image to convolve and the second as the convolution kernel.
>
> - The output image should be the same size as the first input image.
>
> - Both the image to convolve and the kernel can have an arbitrary
> size. There should be no assumption that requires the kernel to be
> smaller than the input image in any dimension or to have an even or
> odd size.
>
> - Results from the spatial domain implementation and Fourier domain
> implementation should differ only by a tolerance when the same
> settings are used for each filter and the pixel type of the images is
> compatible with the underlying FFT routines used in the Fourier domain
> implementation.
>
> - The spacing of the two input images should be required to be the
> same (other specialized convolution algorithms may handle other
> cases). *Or* the image with coarser spacing should be resampled to the
> same spacing as the image with finer spacing.
>
> - The user must have the ability to specify which index is treated as
> the center of the second input. This should default to
> floor((2*index+size-1)/2) where index and size are from the largest
> possible region of the kernel image input.
>
> - The treatment of regions beyond the boundaries of the input images
> should be controllable by the user. Typical choices include constant
> padding, periodic, mirror, and zero-flux Neumann boundary conditions.
> Default should be constant padding with a constant of zero.
>
> - Normalization of the kernel image should be an option. Default is off.
>
> Outstanding question:
>
> - VNL operates only on images whose sizes are powers of two. For best
> performance, FFTW requires that the largest prime factor in an image
> size should be 13 or less. Input images should be padded to meet those
> size constraints in the Fourier domain filter, but no such restriction
> exists for padding in the spatial domain filter. Nevertheless, to
> compare outputs of the spatial and Fourier domain filters, users
> should be able to set the padding size in the spatial domain filter to
> be the same as the padding in the Fourier domain filter. Does exposing
> a class method to control this justified if it is used only for
> testing? Or should testing for agreement between the spatial and
> Fourier domain filters be restricted to cases where no special padding
> for VNL or FFTW is required?
>
> Please let me know what you think.
>
> Thanks,
> Cory
>
> _______________________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html
>
> Kitware offers ITK Training Courses, for more information visit:
> http://kitware.com/products/protraining.html
>
> Please keep messages on-topic and check the ITK FAQ at:
> http://www.itk.org/Wiki/ITK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.itk.org/mailman/listinfo/insight-developers
>

-- 
Cory Quammen
Computer Integrated Systems for Microscopy and Manipulation (CISMM)
Department of Computer Science
University of North Carolina at Chapel Hill
http://www.cs.unc.edu/~cquammen