Allowing generic output type in VectorLinearInterpolateImageFunction::EvaluateAtContinuousIndex

[mathdugre]

Description

The current design restricts the output type of the interpolation to double. As shown here NumericTraits<typename TInputImage::PixelType>::RealType evaluates to double regardless of the input type.

We found this was an issue while profiling ANTs registrationSyN. Using ITK templating, the pipeline offers the option for computation to be done in either float or double. Unexpectedly, the execution time of the pipeline is significantly longer when using float compared to double (See figure below). Our observation shows that the time difference is due to the cast from float to double, and possible back to float in ANTs.

Caption: Mean time per iteration for the interpolation stages and levels in ANTs registrationSyN.

Looking at gdb output we see that the input data is an 3D image of type float

Thread 1 "antsRegistratio" hit Breakpoint 21, itk::VectorLinearInterpolateImageFunction<itk::Image<itk::Vector<float, 3u>, 3u>, float>::EvaluateAtContinuousIndex (this=0x55be27bacc00, index=...)                                                             
    at /tmp/ants/build/staging/include/ITK-5.3/itkVectorLinearInterpolateImageFunction.hxx:48                                                                                                                                                                  
48      VectorLinearInterpolateImageFunction<TInputImage, TCoordRep>::EvaluateAtContinuousIndex( 

but the output is of type double.

itk::FixedArray<double, 3u>::Fill (value=<optimized out>, this=<optimized out>) at /tmp/ants/build/staging/include/ITK-5.3/itkVectorLinearInterpolateImageFunction.hxx:70
70        output.Fill(0.0);

Impact analysis

Currently, the downside is that only double is accepted as an output type for the interpolation.

One the one hand, allowing generic type for the output would prevent confusion and potential slowdowns, as with the ANTs case above. Alternatively, the same type as the input data could be use, with similar benefits.
On the other hand, using different output type could affect the precision of the calculation if lower-precision data type are used, or increase the runtime if higher-precision data type are used. i.e. trade-off between performance and precision.

Due to the potential change in precision incurred by this change, this may affected substantially some tools depending on this function.

Expected behavior

Actual behavior

Versions

Additional Information

[github-actions]

Thank you for contributing an issue! 🙏

Welcome to the ITK community! 🤗👋☀️

We are glad you are here and appreciate your contribution. Please keep in mind our community participation guidelines. 📜
Also, please check existing open issues and consider discussion on the ITK Discourse. 📖

This is an automatic message. Allow for time for the ITK community to be able to read the issue and comment on it.

[thewtex]

Thank you for the nice investigation and report.

Adding the capability to output double does seem reasonable.

[mathdugre]

Thank you for the quick reply.

Adding the capability to output double does seem reasonable.

I might have miss expressed myself. Currently, the function will always output as a double. My suggestion was to allow other types to be returned.
For example, in the ANTs case above, returning data in float would prevent the casting overhead and increase the overall performance of the pipeline.

I would be happy to work on a PR, but I might need some directions as I'm not super familiar with the ITK codebase nor C.

[thewtex]

Ah, yes, returning float make sense, too.

A patch is welcome and appreciated. We are happy to provide guidance on implementing the desired functionality in ITK and ANTs. Our contributing guide provides general direction, but feel free to reach out if you get stuck on specifics.