CUDA Cubic B-Spline Interpolation (CI)

Version 0.9

This read me serves as a quick guide to using the CUDA Cubic B-Spline Interpolation (abbreviated as CI) code. The most recent version of CI and some background information can be found online. To execute and compile CI you need CUDA and the CUDA SDK (2.0 or higher) This software has been released under a revised BSD style license.

Getting started

If you want to simply replace linear interpolation by cubic interpolation, then all you need to do is to include the appropriate header and replace your tex2D and tex3D calls by one of the following:

2D textures:
cubicTex2D(texture tex, float x, float y)
cubicTex2D(texture tex, float2 coord)

3D textures:
cubicTex3D(texture tex, float x, float y, float z)
cubicTex3D(texture tex, float3 coord)

Whereby the texture coordinates coord are expressed in absolute pixel respectively voxel coordinates (thus not in normalized coordinates).

It is also possible to query the cubicly interpolated 1st order derivative of the texture at a coordinate coord; e.g. in the x-direction:
cubicTex3D_1st_derivative_x(texture tex, float3 coord)

The calls for the y- and z-direction are similar, and also the derivatives of 2D textures can be retrieved in the same way. The cost of querying the derivative in a single direction costs the same amount of time as a normal cubicly interpolated lookup in the texture. The gradient of the texture at a specified coordinate can be composed by querying the derivatives in x-, y- and z-direction at that location.

Pre-filtering

When the approach described above is directly applied, it will result in smoothened images. This is caused by the fact that the cubic B-spline interpolation yields a function that does not pass through its coefficients (i.e. texture values). In order to wind up with a cubic B-spline interpolated image that passes through the original samples, we need to pre-filter the texture, as is beautifully described by Philippe Thévenaz et al.

Luckily, CI also provides a CUDA implementation of this filter, and using it is rather simple. The interface for the 2D case is:
CubicBSplinePrefilter2D(float* image, uint width, uint height);

and for the 3D case:
CubicBSplinePrefilter3D(float* volume, uint width, uint height, uint depth);

Note that image and volume should point to a contineous chunk of GPU memory, and that the sample values will be replaced by the cubic B-spline coefficients.

Getting your data there

In order to make your even easier, CI also provides some routines to transfer your data to the GPU memory, and back. Copying your sample values to the GPU memory can be accomplished by this function:
float* CopyVolumeHostToDevice(const float* host, uint width, uint height, uint depth);

The routine allocates GPU memory and copies the sample values from CPU to GPU memory. The pointer to the CPU memory is passed as host, and a pointer to the GPU memory is returned. The counterpart that copies data from the GPU memory to the CPU memory is called:
void CopyVolumeDeviceToHost(float* host, float* device, uint width, uint height, uint depth);

Note that the host destination CPU memory should be allocated before CopyVolumeDeviceToHost is called. This routine will also free the GPU memory.

Often you will have your original data in a different format than float, and for large data sets it costs some time (and memory) to cast everything to float. Therefore CI also provides a number of functions that copy and cast your data immediately to float on the GPU, which is faster and easier. In C++ you can use the following template function:
template extern float* CastVolumeHostToDevice(const T* host, uint width, uint height, uint depth);

The usage of the parameters is the same as for CopyVolumeHostToDevice, except that host can be of the type uchar, schar, ushort or short. Note that the sample values will be normalized, meaning that the maximum value (e.g. 255 for uchar) will be mapped to 1.0.

The following function can be used to copy the continuous chunk of memory used by the pre-filter functions into a texture:
void CreateTextureFromVolume(texture* tex, cudaArray** texArray, const T* volume, cudaExtent extent, bool onDevice);

Example programs

Some example programs are provided along with the CI code, to illustrate the usage of the various routines. In order to compile the example programs, you first need to make sure that the CUDA SDK examples can be compiled (on the Mac you might need to patch the build script first). Then it is simply a matter of opening the Visual Studio solution on Windows machines, or running make in the example folder on the Mac or Linux.

Background

More background information to the CI code is provided online. A comprehensive discussion of uniform B-spline interpolation and the pre-filter can be found in [1]. The fast cubic B-spline interpolation is an adapted version of the method introduced by Sigg and Hadwiger [2]. A description of the adapted algorithm, its merits and its drawbacks is given in [3,4].

  1. Philippe Thévenaz, Thierry Blu, and Michael Unser, "Interpolation Revisited," IEEE Transactions on Medical Imaging, vol. 19, no. 7, pp. 739-758, July 2000.

  2. Christian Sigg and Markus Hadwiger, "Fast Third-Order Texture Filtering," In GPU Gems 2: Programming Techniques for High-Performance Graphics and General-Purpose Computation, Matt Pharr (ed.), Addison-Wesley; chapter 20, pp. 313-329, 2005.

  3. Daniel Ruijters, Bart M. ter Haar Romeny, and Paul Suetens, "Accuracy of GPU-based B-Spline Evaluation," In Proc. Tenth IASTED International Conference on Computer Graphics and Imaging (CGIM), Innsbruck, Austria, pp. 117-122, February 13-15, 2008.

  4. Daniel Ruijters, Bart M. ter Haar Romeny, and Paul Suetens, "Efficient GPU-Based Texture Interpolation using Uniform B-Splines," Journal of Graphics Tools, vol. 13, no. 4, pp. 61-69, 2008.

Copyright 2008-2010 Danny Ruijters