CUDA Cubic B-Spline Interpolation (CI)

Version 1.2

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). Read this, if you are using CUDA 5. This software has been released under a revised BSD style license.

Getting started

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

1D textures:
cubicTex1D(texture tex, float x)

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 filtered 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 1D and 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 filtered 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.


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 filtering 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 prefilter, and using it is rather simple. The interface for the 2D case is:
CubicBSplinePrefilter2D(float* image, uint pitch, uint width, uint height);

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

The image and volume should point to the GPU memory containing the data samples. Note that the sample values will be replaced by the cubic B-spline coefficients. The pitch variable describes the width of a row in the image in bytes. width, height, and depth describe the extent of the data in pixels/voxels. Instead of float*, it is also possible to pass float2*, float3*, or float4* data (e.g., for RGB color data).

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:
cudaPitchedPtr 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, cudaPitchedPtr 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 cudaPitchedPtr 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 output of the pre-filter functions into a texture:
void CreateTextureFromVolume(texture* tex, cudaArray** texArray, const cudaPitchedPtr volume, cudaExtent extent, bool onDevice);

Using layered textures

It is possible to use the cubic interpolation texture lookups with layered textures. In order to do this, you need to add:
#define cudaTextureType2DLayered
to your code, and replace the calls to:
tex2DLayered(texture tex, float x, float y, int layer)
cubicTex2DLayered(texture tex, float x, float y, int layer)
For 1D it works similarly.

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. 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.


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 GPU implementation is described in [2]. The fast cubic B-spline interpolation is an adapted version of the method introduced by Sigg and Hadwiger [3]. A description of the adapted algorithm, its merits and its drawbacks is given in [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. Daniel Ruijters and Philippe Thévenaz, "GPU Prefilter for Accurate Cubic B-Spline Interpolation," The Computer Journal, vol. 55, no. 1, pp. 15-20, January 2012.

  3. 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.

  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-2013 Danny Ruijters