Next: , Previous: Complex One-Dimensional DFTs, Up: Tutorial


2.2 Complex Multi-Dimensional DFTs

Multi-dimensional transforms work much the same way as one-dimensional transforms: you allocate arrays of fftw_complex (preferably using fftw_malloc), create an fftw_plan, execute it as many times as you want with fftw_execute(plan), and clean up with fftw_destroy_plan(plan) (and fftw_free).

FFTW provides two routines for creating plans for 2d and 3d transforms, and one routine for creating plans of arbitrary dimensionality. The 2d and 3d routines have the following signature:

     fftw_plan fftw_plan_dft_2d(int n0, int n1,
                                fftw_complex *in, fftw_complex *out,
                                int sign, unsigned flags);
     fftw_plan fftw_plan_dft_3d(int n0, int n1, int n2,
                                fftw_complex *in, fftw_complex *out,
                                int sign, unsigned flags);

These routines create plans for n0 by n1 two-dimensional (2d) transforms and n0 by n1 by n2 3d transforms, respectively. All of these transforms operate on contiguous arrays in the C-standard row-major order, so that the last dimension has the fastest-varying index in the array. This layout is described further in Multi-dimensional Array Format.

FFTW can also compute transforms of higher dimensionality. In order to avoid confusion between the various meanings of the the word “dimension”, we use the term rank to denote the number of independent indices in an array.1 For example, we say that a 2d transform has rank 2, a 3d transform has rank 3, and so on. You can plan transforms of arbitrary rank by means of the following function:

     fftw_plan fftw_plan_dft(int rank, const int *n,
                             fftw_complex *in, fftw_complex *out,
                             int sign, unsigned flags);

Here, n is a pointer to an array n[rank] denoting an n[0] by n[1] by ... by n[rank-1] transform. Thus, for example, the call

     fftw_plan_dft_2d(n0, n1, in, out, sign, flags);

is equivalent to the following code fragment:

     int n[2];
     n[0] = n0;
     n[1] = n1;
     fftw_plan_dft(2, n, in, out, sign, flags);

fftw_plan_dft is not restricted to 2d and 3d transforms, however, but it can plan transforms of arbitrary rank.

You may have noticed that all the planner routines described so far have overlapping functionality. For example, you can plan a 1d or 2d transform by using fftw_plan_dft with a rank of 1 or 2, or even by calling fftw_plan_dft_3d with n0 and/or n1 equal to 1 (with no loss in efficiency). This pattern continues, and FFTW's planning routines in general form a “partial order,” sequences of interfaces with strictly increasing generality but correspondingly greater complexity.

fftw_plan_dft is the most general complex-DFT routine that we describe in this tutorial, but there are also the advanced and guru interfaces, which allow one to efficiently combine multiple/strided transforms into a single FFTW plan, transform a subset of a larger multi-dimensional array, and/or to handle more general complex-number formats. For more information, see FFTW Reference.


Footnotes

[1] The term “rank” is commonly used in the APL, FORTRAN, and Common Lisp traditions, although it is not so common in the C world.