Go to the first, previous, next, last section, table of contents.

FFTW Reference

This chapter provides a complete reference for all FFTW functions. Programs using FFTW should be linked with -lfftw -lm.

Complex Numbers

By including <fftw.h>, you will have access to the following definitions:

typedef double FFTW_REAL;

typedef struct {
     FFTW_REAL re, im;
} FFTW_COMPLEX;

#define c_re(c)  ((c).re)
#define c_im(c)  ((c).im)

All FFTW operations are performed on the FFTW_COMPLEX data type. The two macros c_re and c_im retrieve, respectively, the real and imaginary parts of a complex number.

Users who wish to work in single precision rather than double precision merely need to change the definition of FFTW_REAL from double to float in fftw.h and then recompile the library.

Plan Creation for One-dimensional Transforms

#include <fftw.h>

fftw_plan fftw_create_plan(int n, fftw_direction dir,
                           int flags);

The function fftw_create_plan creates a plan, which is a data structure containing all the information that fftw needs in order to compute the 1D Fourier Transform. You can create as many plans as you need, but only one plan for a given array size is required (a plan can be reused many times).

fftw_create_plan returns a valid plan, or NULL if, for some reason, the plan can't be created. In the default installation, this can't happen, but it is possible to configure FFTW in such a way that some input sizes are forbidden, and FFTW cannot create a plan.

Arguments

Plan Creation for Multi-dimensional Transforms

#include <fftw.h>

fftwnd_plan fftwnd_create_plan(int rank, const int *n,
                               fftw_direction dir, int flags);

fftwnd_plan fftw2d_create_plan(int nx, int ny,
                               fftw_direction dir, int flags);

fftwnd_plan fftw3d_create_plan(int nx, int ny, int nz,
                               fftw_direction dir, int flags);

The function fftwnd_create_plan creates a plan, which is a data structure containing all the information that fftwnd needs in order to compute a multi-dimensional Fourier Transform. You can create as many plans as you need, but only one plan for a given array size is required (a plan can be reused many times). The functions fftw2d_create_plan and fftw3d_create_plan are optional, alternative interfaces to fftwnd_create_plan for two and three dimensions, respectively.

fftwnd_create_plan returns a valid plan, or NULL if, for some reason, the plan can't be created. This can happen if memory runs out or if the arguments are invalid in some way (e.g. if rank < 0).

Arguments

Computing the one-dimensional transform

#include <fftw.h>

void fftw(fftw_plan plan, int howmany,
          FFTW_COMPLEX *in, int istride, int idist,
          FFTW_COMPLEX *out, int ostride, int odist);

The function fftw computes the one-dimensional Fourier Transform, using a plan created by fftw_create_plan (See section Plan Creation for One-dimensional Transforms.)

Arguments

Computing the multi-dimensional transform

#include <fftw.h>

void fftwnd(fftwnd_plan plan, int howmany,
            FFTW_COMPLEX *in, int istride, int idist,
            FFTW_COMPLEX *out, int ostride, int odist);

The function fftwnd computes the multi-dimensional Fourier Transform, using a plan created by fftwnd_create_plan (See section Plan Creation for Multi-dimensional Transforms.) (Note that the plan determines the rank and dimensions of the array to be transformed.)

Arguments

Destroying a one-dimensional plan

#include <fftw.h>

void fftw_destroy_plan(fftw_plan plan);

The function fftw_destroy_plan frees the plan plan, and releases all the memory associated with it. After destruction, a plan is no longer valid.

Destroying a multi-dimensional plan

#include <fftw.h>

void fftwnd_destroy_plan(fftwnd_plan plan);

The function fftwnd_destroy_plan frees the plan plan, and releases all the memory associated with it. After destruction, a plan is no longer valid.

How to install your own memory allocator

#include <fftw.h>

extern void *(*fftw_malloc_hook) (size_t n);
extern void (*fftw_free_hook) (void *p);

FFTW needs some memory for internal use. It ordinarily calls malloc and free to allocate and release memory; if malloc fails, FFTW prints an error message and exits. This behavior may be undesirable in some applications; consequently, FFTW provide means by which you can install your own memory allocator and take whatever error-correcting action you find appropriate. The variables fftw_malloc_hook and fftw_free_hook are pointers to functions, and they are normally NULL. If you set those variables to point to another function, FFTW will use your function instead of its own internal malloc and free. fftw_malloc_hook must point to a malloc-like function, and fftw_free_hook must point to a free-like function.


Go to the first, previous, next, last section, table of contents.