Input description: A sequence of n real or complex values , , sampled at uniform intervals from a function h.
Problem description: The discrete Fourier transform H of h, , for .
Discussion: Although computer scientists tend to be relatively unfamiliar with Fourier transforms, electrical engineers and signal processors eat them for breakfast. Functionally, Fourier transforms provide a way to convert samples of a standard time-series into the ``frequency domain''. This provides a ``dual'' representation of the function, in which certain operations become easier than in the time domain. Applications of Fourier transforms include:
and can be easily computed using Fourier transforms. Note that if the two functions are similar in shape but one is shifted relative to the other (such as and ), the value of will be large at this shift offset . As an application, suppose that we want to detect whether there are any funny periodicities in our random number generator. We can generate a large series of random numbers, turn it into a time series (the ith number at time i), and take the Fourier transform of this series. Any funny spikes will correspond to potential periodicities.
The discrete Fourier transform takes as input n complex numbers , , corresponding to equally spaced points in a time series, and outputs n complex numbers , , each describing a sine function of given frequency. The discrete Fourier transform is defined by
and the inverse Fourier transform is defined by
which enables us move easily between h and H.
Since the output of the discrete Fourier transform consists of n numbers, each of which can be computed using a formula on n numbers, they can be computed in time. The fast Fourier transform (FFT) is an algorithm that computes the discrete Fourier transform in . This is arguably the most important algorithm known, as measured by practical impact, for it opened the door to modern image processing. There are several different algorithms that call themselves FFTs, all of which are based on a divide-and-conquer approach. Essentially, the problem of computing the discrete Fourier transform on n points is reduced to computing two transforms on n/2 points each and is then applied recursively
The FFT usually assumes that n is a power of two. If this is not the case for your data, you are usually better off padding your data with zeros to create elements rather than hunting for a more general code.
Since many image processing systems have strong real-time constraints, FFTs are often implemented in hardware, or at least in assembly language tuned to the particular machine. Be aware of this possibility if the codes below prove too slow.
Implementations: FFTPACK is a package of Fortran subprograms for the fast Fourier transform of periodic and other symmetric sequences, written by P. Swartzrauber. It includes complex, real, sine, cosine, and quarter-wave transforms. A C language translation of the main routines is also provided. FFTPACK resides on Netlib (see Section ) at http://www.netlib.org/fftpack.
Algorithm 545 [Fra79] of the Collected Algorithms of the ACM is an implementation of the fast Fourier transform optimizing virtual memory performance and written in Fortran. See Section for further information.
XTango (see Section ) is an algorithm animation system for UNIX and X-windows, which includes an interesting animation of the fast Fourier transform.
A Pascal implementation of the fast Fourier transform for points appears in [MS91]. For more details, see Section . Sedgewick [Sed92] provides a bare bones implementation of the fast Fourier transform in C++. See Section for details.
Notes: Brigham [Bri74] is an excellent introduction to Fourier transforms and the FFT and is strongly recommended, as is the exposition in [PFTV86]. Expositions in algorithms texts on the fast Fourier transform include [AHU74, Baa88, CLR90, Man89].
Credit for inventing the fast Fourier transform is usually given to Cooley and Tukey [CT65], although it is not completely deserved. See [Bri74] for a complete history.
An interesting divide-and-conquer algorithm for polynomial multiplication [KO63] does the job in time and is discussed in [AHU74, Man89]. An FFT-based algorithm that multiplies two n-bit numbers in time is due to Schönhage and Strassen [SS71] and is presented in [AHU74].
In recent years, wavelets have been proposed to replace Fourier transforms in filtering. See [Dau92] for an introduction to wavelets.
Related Problems: Data compression (see page ), high-precision arithmetic (see page ).