Chris@10: Next: Complex Multi-Dimensional DFTs, Chris@10: Previous: Tutorial, Chris@10: Up: Tutorial Chris@10:
Chris@10: Plan: To bother about the best method of accomplishing an accidental result. Chris@10: [Ambrose Bierce, The Enlarged Devil's Dictionary.] Chris@10:Chris@10: Chris@10:
The basic usage of FFTW to compute a one-dimensional DFT of size
Chris@10: N is simple, and it typically looks something like this code:
Chris@10:
Chris@10:
#include <fftw3.h>
Chris@10: ...
Chris@10: {
Chris@10: fftw_complex *in, *out;
Chris@10: fftw_plan p;
Chris@10: ...
Chris@10: in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N);
Chris@10: out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N);
Chris@10: p = fftw_plan_dft_1d(N, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
Chris@10: ...
Chris@10: fftw_execute(p); /* repeat as needed */
Chris@10: ...
Chris@10: fftw_destroy_plan(p);
Chris@10: fftw_free(in); fftw_free(out);
Chris@10: }
Chris@10:
Chris@10: You must link this code with the fftw3 library. On Unix systems,
Chris@10: link with -lfftw3 -lm.
Chris@10:
Chris@10:
The example code first allocates the input and output arrays. You can
Chris@10: allocate them in any way that you like, but we recommend using
Chris@10: fftw_malloc, which behaves like
Chris@10: malloc except that it properly aligns the array when SIMD
Chris@10: instructions (such as SSE and Altivec) are available (see SIMD alignment and fftw_malloc). [Alternatively, we provide a convenient wrapper function fftw_alloc_complex(N) which has the same effect.]
Chris@10:
Chris@10:
Chris@10:
The data is an array of type fftw_complex, which is by default a
Chris@10: double[2] composed of the real (in[i][0]) and imaginary
Chris@10: (in[i][1]) parts of a complex number.
Chris@10:
Chris@10: The next step is to create a plan, which is an object
Chris@10: that contains all the data that FFTW needs to compute the FFT.
Chris@10: This function creates the plan:
Chris@10:
Chris@10:
fftw_plan fftw_plan_dft_1d(int n, fftw_complex *in, fftw_complex *out, Chris@10: int sign, unsigned flags); Chris@10:Chris@10:
Chris@10: The first argument, n, is the size of the transform you are
Chris@10: trying to compute. The size n can be any positive integer, but
Chris@10: sizes that are products of small factors are transformed most
Chris@10: efficiently (although prime sizes still use an O(n log n) algorithm).
Chris@10:
Chris@10:
The next two arguments are pointers to the input and output arrays of Chris@10: the transform. These pointers can be equal, indicating an Chris@10: in-place transform. Chris@10: Chris@10: Chris@10:
The fourth argument, sign, can be either FFTW_FORWARD
Chris@10: (-1) or FFTW_BACKWARD (+1),
Chris@10: and indicates the direction of the transform you are interested in;
Chris@10: technically, it is the sign of the exponent in the transform.
Chris@10:
Chris@10:
The flags argument is usually either FFTW_MEASURE or
Chris@10: FFTW_ESTIMATE. FFTW_MEASURE instructs FFTW to run
Chris@10: and measure the execution time of several FFTs in order to find the
Chris@10: best way to compute the transform of size n. This process takes
Chris@10: some time (usually a few seconds), depending on your machine and on
Chris@10: the size of the transform. FFTW_ESTIMATE, on the contrary,
Chris@10: does not run any computation and just builds a
Chris@10: reasonable plan that is probably sub-optimal. In short, if your
Chris@10: program performs many transforms of the same size and initialization
Chris@10: time is not important, use FFTW_MEASURE; otherwise use the
Chris@10: estimate.
Chris@10:
Chris@10:
You must create the plan before initializing the input, because
Chris@10: FFTW_MEASURE overwrites the in/out arrays.
Chris@10: (Technically, FFTW_ESTIMATE does not touch your arrays, but you
Chris@10: should always create plans first just to be sure.)
Chris@10:
Chris@10:
Once the plan has been created, you can use it as many times as you
Chris@10: like for transforms on the specified in/out arrays,
Chris@10: computing the actual transforms via fftw_execute(plan):
Chris@10:
void fftw_execute(const fftw_plan plan); Chris@10:Chris@10:
Chris@10: The DFT results are stored in-order in the array out, with the
Chris@10: zero-frequency (DC) component in out[0].
Chris@10: If in != out, the transform is out-of-place and the input
Chris@10: array in is not modified. Otherwise, the input array is
Chris@10: overwritten with the transform.
Chris@10:
Chris@10:
If you want to transform a different array of the same size, you
Chris@10: can create a new plan with fftw_plan_dft_1d and FFTW
Chris@10: automatically reuses the information from the previous plan, if
Chris@10: possible. Alternatively, with the “guru” interface you can apply a
Chris@10: given plan to a different array, if you are careful.
Chris@10: See FFTW Reference.
Chris@10:
Chris@10:
When you are done with the plan, you deallocate it by calling
Chris@10: fftw_destroy_plan(plan):
Chris@10:
void fftw_destroy_plan(fftw_plan plan); Chris@10:Chris@10:
If you allocate an array with fftw_malloc() you must deallocate
Chris@10: it with fftw_free(). Do not use free() or, heaven
Chris@10: forbid, delete.
Chris@10:
Chris@10: FFTW computes an unnormalized DFT. Thus, computing a forward
Chris@10: followed by a backward transform (or vice versa) results in the original
Chris@10: array scaled by n. For the definition of the DFT, see What FFTW Really Computes.
Chris@10:
Chris@10:
Chris@10:
If you have a C compiler, such as gcc, that supports the
Chris@10: C99 standard, and you #include <complex.h> before
Chris@10: <fftw3.h>, then fftw_complex is the native
Chris@10: double-precision complex type and you can manipulate it with ordinary
Chris@10: arithmetic. Otherwise, FFTW defines its own complex type, which is
Chris@10: bit-compatible with the C99 complex type. See Complex numbers.
Chris@10: (The C++ <complex> template class may also be usable via a
Chris@10: typecast.)
Chris@10:
Chris@10: To use single or long-double precision versions of FFTW, replace the
Chris@10: fftw_ prefix by fftwf_ or fftwl_ and link with
Chris@10: -lfftw3f or -lfftw3l, but use the same
Chris@10: <fftw3.h> header file.
Chris@10:
Chris@10:
Chris@10:
Many more flags exist besides FFTW_MEASURE and
Chris@10: FFTW_ESTIMATE. For example, use FFTW_PATIENT if you're
Chris@10: willing to wait even longer for a possibly even faster plan (see FFTW Reference).
Chris@10: You can also save plans for future use, as described by Words of Wisdom-Saving Plans.
Chris@10:
Chris@10:
Chris@10:
Chris@10: