![FFT of Cosine Summation Function](https://upload.wikimedia.org/wikipedia/commons/3/30/FFT_of_Cosine_Summation_Function.svg)
A fast Fourier transform (FFT) is an
algorithm
In mathematics and computer science, an algorithm () is a finite sequence of rigorous instructions, typically used to solve a class of specific problems or to perform a computation. Algorithms are used as specifications for performing ...
that computes the
discrete Fourier transform
In mathematics, the discrete Fourier transform (DFT) converts a finite sequence of equally-spaced samples of a function into a same-length sequence of equally-spaced samples of the discrete-time Fourier transform (DTFT), which is a comple ...
(DFT) of a sequence, or its inverse (IDFT).
Fourier analysis converts a signal from its original domain (often time or space) to a representation in the
frequency domain
In physics, electronics, control systems engineering, and statistics, the frequency domain refers to the analysis of mathematical functions or signals with respect to frequency, rather than time. Put simply, a time-domain graph shows how a s ...
and vice versa. The DFT is obtained by decomposing a
sequence
In mathematics, a sequence is an enumerated collection of objects in which repetitions are allowed and order matters. Like a set, it contains members (also called ''elements'', or ''terms''). The number of elements (possibly infinite) is calle ...
of values into components of different frequencies.
This operation is useful in many fields, but computing it directly from the definition is often too slow to be practical. An FFT rapidly computes such transformations by
factorizing the
DFT matrix
In applied mathematics, a DFT matrix is an expression of a discrete Fourier transform (DFT) as a transformation matrix, which can be applied to a signal through matrix multiplication.
Definition
An ''N''-point DFT is expressed as the multiplicati ...
into a product of
sparse (mostly zero) factors.
As a result, it manages to reduce the
complexity of computing the DFT from
, which arises if one simply applies the definition of DFT, to
, where
is the data size. The difference in speed can be enormous, especially for long data sets where ''N'' may be in the thousands or millions. In the presence of
round-off error, many FFT algorithms are much more accurate than evaluating the DFT definition directly or indirectly. There are many different FFT algorithms based on a wide range of published theories, from simple
complex-number arithmetic to
group theory
In abstract algebra, group theory studies the algebraic structures known as group (mathematics), groups.
The concept of a group is central to abstract algebra: other well-known algebraic structures, such as ring (mathematics), rings, field ...
and
number theory
Number theory (or arithmetic or higher arithmetic in older usage) is a branch of pure mathematics devoted primarily to the study of the integers and integer-valued functions. German mathematician Carl Friedrich Gauss (1777–1855) said, "Mat ...
.
![Simple time domain vs frequency domain](https://upload.wikimedia.org/wikipedia/commons/0/07/Simple_time_domain_vs_frequency_domain.svg)
Fast Fourier transforms are widely used for
applications
Application may refer to:
Mathematics and computing
* Application software, computer software designed to help the user to perform specific tasks
** Application layer, an abstraction layer that specifies protocols and interface methods used in a c ...
in engineering, music, science, and mathematics. The basic ideas were popularized in 1965, but some algorithms had been derived as early as 1805.
In 1994,
Gilbert Strang
William Gilbert Strang (born November 27, 1934), usually known as simply Gilbert Strang or Gil Strang, is an American mathematician, with contributions to finite element theory, the calculus of variations, wavelet analysis and linear algebr ...
described the FFT as "the most important
numerical algorithm
Numerical analysis is the study of algorithms that use numerical approximation (as opposed to symbolic manipulations) for the problems of mathematical analysis (as distinguished from discrete mathematics). It is the study of numerical methods ...
of our lifetime",
and it was included in Top 10 Algorithms of 20th Century by the
IEEE
The Institute of Electrical and Electronics Engineers (IEEE) is a 501(c)(3) professional association for electronic engineering and electrical engineering (and associated disciplines) with its corporate office in New York City and its operat ...
magazine ''Computing in Science & Engineering''.
The best-known FFT algorithms depend upon the
factorization
In mathematics, factorization (or factorisation, see English spelling differences) or factoring consists of writing a number or another mathematical object as a product of several ''factors'', usually smaller or simpler objects of the same kind ...
of ''N'', but there are FFTs with
O(''N'' log ''N'') complexity for all ''N'', even for
prime
A prime number (or a prime) is a natural number greater than 1 that is not a product of two smaller natural numbers. A natural number greater than 1 that is not prime is called a composite number. For example, 5 is prime because the only ways ...
''N''. Many FFT algorithms depend only on the fact that
is an ''N''-th
primitive root of unity, and thus can be applied to analogous transforms over any
finite field
In mathematics, a finite field or Galois field (so-named in honor of Évariste Galois) is a field that contains a finite number of elements. As with any field, a finite field is a set on which the operations of multiplication, addition, subtr ...
, such as
number-theoretic transforms. Since the inverse DFT is the same as the DFT, but with the opposite sign in the exponent and a 1/''N'' factor, any FFT algorithm can easily be adapted for it.
History
The development of fast algorithms for DFT can be traced to
Carl Friedrich Gauss
Johann Carl Friedrich Gauss (; german: Gauß ; la, Carolus Fridericus Gauss; 30 April 177723 February 1855) was a German mathematician and physicist who made significant contributions to many fields in mathematics and science. Sometimes refer ...
's unpublished work in 1805 when he needed it to interpolate the orbit of asteroids
Pallas
Pallas may refer to:
Astronomy
* 2 Pallas asteroid
** Pallas family, a group of asteroids that includes 2 Pallas
* Pallas (crater), a crater on Earth's moon
Mythology
* Pallas (Giant), a son of Uranus and Gaia, killed and flayed by Athena
* Pa ...
and
Juno from sample observations.
His method was very similar to the one published in 1965 by
James Cooley and
John Tukey, who are generally credited for the invention of the modern generic FFT algorithm. While Gauss's work predated even
Joseph Fourier's results in 1822, he did not analyze the computation time and eventually used other methods to achieve his goal.
Between 1805 and 1965, some versions of FFT were published by other authors.
Frank Yates in 1932 published his version called ''interaction algorithm'', which provided
efficient computation of Hadamard and Walsh transforms.
Yates' algorithm is still used in the field of statistical design and analysis of experiments. In 1942,
G. C. Danielson Gordon Charles Danielson (October 28, 1912 - September 30, 1983) was a Distinguished Professor in Sciences and Humanities in 1964 at Iowa State University at Ames, Iowa.
His name was added to the Distinguished Professor Award Wall in Beardsher H ...
and
Cornelius Lanczos
__NOTOC__
Cornelius (Cornel) Lanczos ( hu, Lánczos Kornél, ; born as Kornél Lőwy, until 1906: ''Löwy (Lőwy) Kornél''; February 2, 1893 – June 25, 1974) was a Hungarian-American and later Hungarian-Irish mathematician and physicist. Accor ...
published their version to compute DFT for
x-ray crystallography
X-ray crystallography is the experimental science determining the atomic and molecular structure of a crystal, in which the crystalline structure causes a beam of incident X-rays to diffract into many specific directions. By measuring the angles ...
, a field where calculation of Fourier transforms presented a formidable bottleneck.
While many methods in the past had focused on reducing the constant factor for
computation by taking advantage of "symmetries", Danielson and Lanczos realized that one could use the "periodicity" and apply a "doubling trick" to "double
/nowiki> with only slightly more than double the labor", though like Gauss they did not analyze that this led to scaling.
James Cooley and John Tukey independently rediscovered these earlier algorithms and published a more general FFT in 1965 that is applicable when ''N'' is composite and not necessarily a power of 2, as well as analyzing the scaling. Tukey came up with the idea during a meeting of President Kennedy
John Fitzgerald Kennedy (May 29, 1917 – November 22, 1963), often referred to by his initials JFK and the nickname Jack, was an American politician who served as the 35th president of the United States from 1961 until assassination of Joh ...
's Science Advisory Committee where a discussion topic involved detecting nuclear tests by the Soviet Union by setting up sensors to surround the country from outside. To analyze the output of these sensors, an FFT algorithm would be needed. In discussion with Tukey, Richard Garwin
Richard Lawrence Garwin (born April 19, 1928) is an American physicist, best known as the author of the first hydrogen bomb design.
In 1978, Garwin was elected a member of the National Academy of Engineering for contributing to the application ...
recognized the general applicability of the algorithm not just to national security problems, but also to a wide range of problems including one of immediate interest to him, determining the periodicities of the spin orientations in a 3-D crystal of Helium-3. Garwin gave Tukey's idea to Cooley (both worked at IBM's Watson labs) for implementation. Cooley and Tukey published the paper in a relatively short time of six months. As Tukey did not work at IBM, the patentability of the idea was doubted and the algorithm went into the public domain, which, through the computing revolution of the next decade, made FFT one of the indispensable algorithms in digital signal processing.
Definition
Let , …, be complex number
In mathematics, a complex number is an element of a number system that extends the real numbers with a specific element denoted , called the imaginary unit and satisfying the equation i^= -1; every complex number can be expressed in the fo ...
s. The DFT is defined by the formula
:
where is a primitive th root of 1.
Evaluating this definition directly requires operations: there are ''N'' outputs ''X''''k'', and each output requires a sum of ''N'' terms. An FFT is any method to compute the same results in operations. All known FFT algorithms require operations, although there is no known proof that lower complexity is impossible.
To illustrate the savings of an FFT, consider the count of complex multiplications and additions for data points. Evaluating the DFT's sums directly involves complex multiplications and complex additions, of which operations can be saved by eliminating trivial operations such as multiplications by 1, leaving about 30 million operations. In contrast, the radix-2 Cooley–Tukey algorithm, for ''N'' a power of 2, can compute the same result with only complex multiplications (again, ignoring simplifications of multiplications by 1 and similar) and complex additions, in total about 30,000 operations — a thousand times less than with direct evaluation. In practice, actual performance on modern computers is usually dominated by factors other than the speed of arithmetic operations and the analysis is a complicated subject (for example, see Frigo & Johnson
Johnson is a surname of Anglo-Norman origin meaning "Son of John". It is the second most common in the United States and 154th most common in the world. As a common family name in Scotland, Johnson is occasionally a variation of ''Johnston'', a ...
, 2005), but the overall improvement from to remains.
Algorithms
Cooley–Tukey algorithm
By far the most commonly used FFT is the Cooley–Tukey algorithm. This is a divide-and-conquer algorithm
In computer science, divide and conquer is an algorithm design paradigm. A divide-and-conquer algorithm recursively breaks down a problem into two or more sub-problems of the same or related type, until these become simple enough to be solved dire ...
that recursively
Recursion (adjective: ''recursive'') occurs when a thing is defined in terms of itself or of its type. Recursion is used in a variety of disciplines ranging from linguistics to logic. The most common application of recursion is in mathematics ...
breaks down a DFT of any composite
Composite or compositing may refer to:
Materials
* Composite material, a material that is made from several different substances
** Metal matrix composite, composed of metal and other parts
** Cermet, a composite of ceramic and metallic materials
...
size into many smaller DFTs of sizes and , along with multiplications by complex roots of unity
In mathematics, a root of unity, occasionally called a de Moivre number, is any complex number that yields 1 when raised to some positive integer power . Roots of unity are used in many branches of mathematics, and are especially important in ...
traditionally called twiddle factor
A twiddle factor, in fast Fourier transform (FFT) algorithms, is any of the trigonometric constant coefficients that are multiplied by the data in the course of the algorithm. This term was apparently coined by Gentleman & Sande in 1966, and has ...
s (after Gentleman and Sande, 1966).
This method (and the general idea of an FFT) was popularized by a publication of Cooley and Tukey in 1965, but it was later discovered that those two authors had independently re-invented an algorithm known to Carl Friedrich Gauss
Johann Carl Friedrich Gauss (; german: Gauß ; la, Carolus Fridericus Gauss; 30 April 177723 February 1855) was a German mathematician and physicist who made significant contributions to many fields in mathematics and science. Sometimes refer ...
around 1805 (and subsequently rediscovered several times in limited forms).
The best known use of the Cooley–Tukey algorithm is to divide the transform into two pieces of size ''N''/2 at each step, and is therefore limited to power-of-two sizes, but any factorization can be used in general (as was known to both Gauss and Cooley/Tukey). These are called the ''radix-2'' and ''mixed-radix'' cases, respectively (and other variants such as the split-radix FFT The split-radix FFT is a fast Fourier transform (FFT) algorithm for computing the discrete Fourier transform (DFT), and was first described in an initially little-appreciated paper by R. Yavne (1968) and subsequently rediscovered simultaneously by ...
have their own names as well). Although the basic idea is recursive, most traditional implementations rearrange the algorithm to avoid explicit recursion. Also, because the Cooley–Tukey algorithm breaks the DFT into smaller DFTs, it can be combined arbitrarily with any other algorithm for the DFT, such as those described below.
Other FFT algorithms
There are FFT algorithms other than Cooley–Tukey.
For ''N'' = ''N''1''N''2 with coprime
In mathematics, two integers and are coprime, relatively prime or mutually prime if the only positive integer that is a divisor of both of them is 1. Consequently, any prime number that divides does not divide , and vice versa. This is equivale ...
''N''1 and ''N''2, one can use the prime-factor (Good–Thomas) algorithm (PFA), based on the Chinese remainder theorem, to factorize the DFT similarly to Cooley–Tukey but without the twiddle factors. The Rader–Brenner algorithm (1976) is a Cooley–Tukey-like factorization but with purely imaginary twiddle factors, reducing multiplications at the cost of increased additions and reduced numerical stability; it was later superseded by the split-radix variant of Cooley–Tukey (which achieves the same multiplication count but with fewer additions and without sacrificing accuracy). Algorithms that recursively factorize the DFT into smaller operations other than DFTs include the Bruun and QFT algorithms. (The Rader–Brenner and QFT algorithms were proposed for power-of-two sizes, but it is possible that they could be adapted to general composite ''N''. Bruun's algorithm applies to arbitrary even composite sizes.) Bruun's algorithm, in particular, is based on interpreting the FFT as a recursive factorization of the polynomial
In mathematics, a polynomial is an expression consisting of indeterminates (also called variables) and coefficients, that involves only the operations of addition, subtraction, multiplication, and positive-integer powers of variables. An example ...
''z''''N'' − 1, here into real-coefficient polynomials of the form ''z''''M'' − 1 and ''z''2''M'' + ''az''''M'' + 1.
Another polynomial viewpoint is exploited by the Winograd FFT algorithm, which factorizes ''z''''N'' − 1 into cyclotomic polynomials—these often have coefficients of 1, 0, or −1, and therefore require few (if any) multiplications, so Winograd can be used to obtain minimal-multiplication FFTs and is often used to find efficient algorithms for small factors. Indeed, Winograd showed that the DFT can be computed with only O(''N'') irrational multiplications, leading to a proven achievable lower bound on the number of multiplications for power-of-two sizes; unfortunately, this comes at the cost of many more additions, a tradeoff no longer favorable on modern processors
A central processing unit (CPU), also called a central processor, main processor or just processor, is the electronic circuitry that executes instructions comprising a computer program. The CPU performs basic arithmetic, logic, controlling, ...
with hardware multipliers. In particular, Winograd also makes use of the PFA as well as an algorithm by Rader for FFTs of ''prime'' sizes.
Rader's algorithm, exploiting the existence of a generator for the multiplicative group
A group is a number of persons or things that are located, gathered, or classed together.
Groups of people
* Cultural group, a group whose members share the same cultural identity
* Ethnic group, a group whose members share the same ethnic ide ...
modulo prime ''N'', expresses a DFT of prime size ''N'' as a cyclic convolution
In mathematics (in particular, functional analysis), convolution is a mathematical operation on two functions ( and ) that produces a third function (f*g) that expresses how the shape of one is modified by the other. The term ''convolution'' ...
of (composite) size ''N'' − 1, which can then be computed by a pair of ordinary FFTs via the convolution theorem
In mathematics, the convolution theorem states that under suitable conditions the Fourier transform of a convolution of two functions (or signals) is the pointwise product of their Fourier transforms. More generally, convolution in one domain (e.g ...
(although Winograd uses other convolution methods). Another prime-size FFT is due to L. I. Bluestein, and is sometimes called the chirp-z algorithm
The chirp Z-transform (CZT) is a generalization of the discrete Fourier transform (DFT). While the DFT samples the Z plane at uniformly-spaced points along the unit circle, the chirp Z-transform samples along spiral arcs in the Z-plane, correspon ...
; it also re-expresses a DFT as a convolution, but this time of the ''same'' size (which can be zero-padded to a power of two
A power of two is a number of the form where is an integer, that is, the result of exponentiation with number two as the base and integer as the exponent.
In a context where only integers are considered, is restricted to non-negativ ...
and evaluated by radix-2 Cooley–Tukey FFTs, for example), via the identity
:
Hexagonal fast Fourier transform The fast Fourier transform (FFT) is an important tool in the fields of image and signal processing. The hexagonal fast Fourier transform (HFFT) uses existing FFT routines to compute the discrete Fourier transform (DFT) of images that have been cap ...
(HFFT) aims at computing an efficient FFT for the hexagonally-sampled data by using a new addressing scheme for hexagonal grids, called Array Set Addressing (ASA).
FFT algorithms specialized for real or symmetric data
In many applications, the input data for the DFT are purely real, in which case the outputs satisfy the symmetry
:
and efficient FFT algorithms have been designed for this situation (see e.g. Sorensen, 1987). One approach consists of taking an ordinary algorithm (e.g. Cooley–Tukey) and removing the redundant parts of the computation, saving roughly a factor of two in time and memory. Alternatively, it is possible to express an ''even''-length real-input DFT as a complex DFT of half the length (whose real and imaginary parts are the even/odd elements of the original real data), followed by O(''N'') post-processing operations.
It was once believed that real-input DFTs could be more efficiently computed by means of the discrete Hartley transform
A discrete Hartley transform (DHT) is a Fourier-related transform of discrete, periodic data similar to the discrete Fourier transform (DFT), with analogous applications in signal processing and related fields. Its main distinction from the DFT is ...
(DHT), but it was subsequently argued that a specialized real-input DFT algorithm (FFT) can typically be found that requires fewer operations than the corresponding DHT algorithm (FHT) for the same number of inputs. Bruun's algorithm (above) is another method that was initially proposed to take advantage of real inputs, but it has not proved popular.
There are further FFT specializations for the cases of real data that have even/odd symmetry, in which case one can gain another factor of roughly two in time and memory and the DFT becomes the discrete cosine/ sine transform(s) ( DCT/DST
Daylight saving time (DST), also referred to as daylight savings time or simply daylight time (United States, Canada, and Australia), and summer time (United Kingdom, European Union, and others), is the practice of advancing clocks (typicall ...
). Instead of directly modifying an FFT algorithm for these cases, DCTs/DSTs can also be computed via FFTs of real data combined with O(''N'') pre- and post-processing.
Computational issues
Bounds on complexity and operation counts
A fundamental question of longstanding theoretical interest is to prove lower bounds on the complexity and exact operation counts of fast Fourier transforms, and many open problems remain. It is not rigorously proved whether DFTs truly require (i.e., order '''' or greater) operations, even for the simple case of power of two
A power of two is a number of the form where is an integer, that is, the result of exponentiation with number two as the base and integer as the exponent.
In a context where only integers are considered, is restricted to non-negativ ...
sizes, although no algorithms with lower complexity are known. In particular, the count of arithmetic operations is usually the focus of such questions, although actual performance on modern-day computers is determined by many other factors such as cache
Cache, caching, or caché may refer to:
Places United States
* Cache, Idaho, an unincorporated community
* Cache, Illinois, an unincorporated community
* Cache, Oklahoma, a city in Comanche County
* Cache, Utah, Cache County, Utah
* Cache County ...
or CPU pipeline
In computing, a pipeline, also known as a data pipeline, is a set of data processing elements connected in series, where the output of one element is the input of the next one. The elements of a pipeline are often executed in parallel or in time- ...
optimization.
Following work by Shmuel Winograd
__NOTOC__
Shmuel Winograd ( he, שמואל וינוגרד; January 4, 1936 – March 25, 2019) was an Israeli-American computer scientist, noted for his contributions to computational complexity. He has proved several major results regarding the ...
(1978), a tight Θ(''N'') lower bound is known for the number of real multiplications required by an FFT. It can be shown that only irrational real multiplications are required to compute a DFT of power-of-two length . Moreover, explicit algorithms that achieve this count are known (Heideman & Burrus, 1986; Duhamel, 1990). However, these algorithms require too many additions to be practical, at least on modern computers with hardware multipliers (Duhamel, 1990; Frigo & Johnson
Johnson is a surname of Anglo-Norman origin meaning "Son of John". It is the second most common in the United States and 154th most common in the world. As a common family name in Scotland, Johnson is occasionally a variation of ''Johnston'', a ...
, 2005).
A tight lower bound is not known on the number of required additions, although lower bounds have been proved under some restrictive assumptions on the algorithms. In 1973, Morgenstern proved an lower bound on the addition count for algorithms where the multiplicative constants have bounded magnitudes (which is true for most but not all FFT algorithms). Pan (1986) proved an lower bound assuming a bound on a measure of the FFT algorithm's "asynchronicity", but the generality of this assumption is unclear. For the case of power-of-two ''N'', Papadimitriou (1979) argued that the number of complex-number additions achieved by Cooley–Tukey algorithms is ''optimal'' under certain assumptions on the graph
Graph may refer to:
Mathematics
*Graph (discrete mathematics), a structure made of vertices and edges
**Graph theory, the study of such graphs and their properties
*Graph (topology), a topological space resembling a graph in the sense of discre ...
of the algorithm (his assumptions imply, among other things, that no additive identities in the roots of unity are exploited). (This argument would imply that at least real additions are required, although this is not a tight bound because extra additions are required as part of complex-number multiplications.) Thus far, no published FFT algorithm has achieved fewer than complex-number additions (or their equivalent) for power-of-two ''N''.
A third problem is to minimize the ''total'' number of real multiplications and additions, sometimes called the "arithmetic complexity" (although in this context it is the exact count and not the asymptotic complexity that is being considered). Again, no tight lower bound has been proven. Since 1968, however, the lowest published count for power-of-two ''N'' was long achieved by the split-radix FFT algorithm The split-radix FFT is a fast Fourier transform (FFT) algorithm for computing the discrete Fourier transform (DFT), and was first described in an initially little-appreciated paper by R. Yavne (1968) and subsequently rediscovered simultaneously by ...
, which requires real multiplications and additions for ''N'' > 1. This was recently reduced to (Johnson and Frigo, 2007; Lundy and Van Buskirk, 2007). A slightly larger count (but still better than split radix for ''N'' ≥ 256) was shown to be provably optimal for ''N'' ≤ 512 under additional restrictions on the possible algorithms (split-radix-like flowgraphs with unit-modulus multiplicative factors), by reduction to a satisfiability modulo theories
In computer science and mathematical logic, satisfiability modulo theories (SMT) is the problem of determining whether a mathematical formula is satisfiable. It generalizes the Boolean satisfiability problem (SAT) to more complex formulas involvi ...
problem solvable by brute force (Haynal & Haynal, 2011).
Most of the attempts to lower or prove the complexity of FFT algorithms have focused on the ordinary complex-data case, because it is the simplest. However, complex-data FFTs are so closely related to algorithms for related problems such as real-data FFTs, discrete cosine transforms, discrete Hartley transform
A discrete Hartley transform (DHT) is a Fourier-related transform of discrete, periodic data similar to the discrete Fourier transform (DFT), with analogous applications in signal processing and related fields. Its main distinction from the DFT is ...
s, and so on, that any improvement in one of these would immediately lead to improvements in the others (Duhamel & Vetterli, 1990).
Approximations
All of the FFT algorithms discussed above compute the DFT exactly (i.e. neglecting floating-point errors). A few "FFT" algorithms have been proposed, however, that compute the DFT ''approximately'', with an error that can be made arbitrarily small at the expense of increased computations. Such algorithms trade the approximation error for increased speed or other properties. For example, an approximate FFT algorithm by Edelman et al. (1999) achieves lower communication requirements for parallel computing with the help of a fast multipole method __NOTOC__
The fast multipole method (FMM) is a numerical technique that was developed to speed up the calculation of long-ranged forces in the ''n''-body problem. It does this by expanding the system Green's function using a multipole expansion, w ...
. A wavelet
A wavelet is a wave-like oscillation with an amplitude that begins at zero, increases or decreases, and then returns to zero one or more times. Wavelets are termed a "brief oscillation". A taxonomy of wavelets has been established, based on the num ...
-based approximate FFT by Guo and Burrus (1996) takes sparse inputs/outputs (time/frequency localization) into account more efficiently than is possible with an exact FFT. Another algorithm for approximate computation of a subset of the DFT outputs is due to Shentov et al. (1995). The Edelman algorithm works equally well for sparse and non-sparse data, since it is based on the compressibility (rank deficiency) of the Fourier matrix itself rather than the compressibility (sparsity) of the data. Conversely, if the data are sparse—that is, if only ''K'' out of ''N'' Fourier coefficients are nonzero—then the complexity can be reduced to O(''K''log(''N'')log(''N''/''K'')), and this has been demonstrated to lead to practical speedups compared to an ordinary FFT for ''N''/''K'' > 32 in a large-''N'' example (''N'' = 222) using a probabilistic approximate algorithm (which estimates the largest ''K'' coefficients to several decimal places).
Accuracy
FFT algorithms have errors when finite-precision floating-point arithmetic is used, but these errors are typically quite small; most FFT algorithms, e.g. Cooley–Tukey, have excellent numerical properties as a consequence of the pairwise summation In numerical analysis, pairwise summation, also called cascade summation, is a technique to sum a sequence of finite-precision floating-point numbers that substantially reduces the accumulated round-off error compared to naively accumulating the s ...
structure of the algorithms. The upper bound on the relative error
The approximation error in a data value is the discrepancy between an exact value and some '' approximation'' to it. This error can be expressed as an absolute error (the numerical amount of the discrepancy) or as a relative error (the absolute e ...
for the Cooley–Tukey algorithm is , compared to for the naïve DFT formula, where ε is the machine floating-point relative precision. In fact, the root mean square (rms) errors are much better than these upper bounds, being only for Cooley–Tukey and for the naïve DFT (Schatzman, 1996). These results, however, are very sensitive to the accuracy of the twiddle factors used in the FFT (i.e. the trigonometric function
In mathematics, the trigonometric functions (also called circular functions, angle functions or goniometric functions) are real functions which relate an angle of a right-angled triangle to ratios of two side lengths. They are widely used in a ...
values), and it is not unusual for incautious FFT implementations to have much worse accuracy, e.g. if they use inaccurate trigonometric recurrence formulas. Some FFTs other than Cooley–Tukey, such as the Rader–Brenner algorithm, are intrinsically less stable.
In fixed-point arithmetic, the finite-precision errors accumulated by FFT algorithms are worse, with rms errors growing as for the Cooley–Tukey algorithm (Welch, 1969). Achieving this accuracy requires careful attention to scaling to minimize loss of precision, and fixed-point FFT algorithms involve rescaling at each intermediate stage of decompositions like Cooley–Tukey.
To verify the correctness of an FFT implementation, rigorous guarantees can be obtained in time by a simple procedure checking the linearity, impulse-response, and time-shift properties of the transform on random inputs (Ergün, 1995).
Multidimensional FFTs
As defined in the multidimensional DFT article, the multidimensional DFT
:
transforms an array ''x''n with a ''d''-dimensional vector
Vector most often refers to:
*Euclidean vector, a quantity with a magnitude and a direction
*Vector (epidemiology), an agent that carries and transmits an infectious pathogen into another living organism
Vector may also refer to:
Mathematic ...
of indices by a set of ''d'' nested summations (over for each ''j''), where the division n/N, defined as , is performed element-wise. Equivalently, it is the composition of a sequence of ''d'' sets of one-dimensional DFTs, performed along one dimension at a time (in any order).
This compositional viewpoint immediately provides the simplest and most common multidimensional DFT algorithm, known as the row-column algorithm (after the two-dimensional case, below). That is, one simply performs a sequence of ''d'' one-dimensional FFTs (by any of the above algorithms): first you transform along the ''n''1 dimension, then along the ''n''2 dimension, and so on (or actually, any ordering works). This method is easily shown to have the usual complexity, where is the total number of data points transformed. In particular, there are ''N''/''N''1 transforms of size ''N''1, etcetera, so the complexity of the sequence of FFTs is:
:
In two dimensions, the ''x''k can be viewed as an matrix, and this algorithm corresponds to first performing the FFT of all the rows (resp. columns), grouping the resulting transformed rows (resp. columns) together as another matrix, and then performing the FFT on each of the columns (resp. rows) of this second matrix, and similarly grouping the results into the final result matrix.
In more than two dimensions, it is often advantageous for cache
Cache, caching, or caché may refer to:
Places United States
* Cache, Idaho, an unincorporated community
* Cache, Illinois, an unincorporated community
* Cache, Oklahoma, a city in Comanche County
* Cache, Utah, Cache County, Utah
* Cache County ...
locality to group the dimensions recursively. For example, a three-dimensional FFT might first perform two-dimensional FFTs of each planar "slice" for each fixed ''n''1, and then perform the one-dimensional FFTs along the ''n''1 direction. More generally, an asymptotically optimal
In computer science, an algorithm is said to be asymptotically optimal if, roughly speaking, for large inputs it performs at worst a constant factor (independent of the input size) worse than the best possible algorithm. It is a term commonly en ...
cache-oblivious algorithm
In computing, a cache-oblivious algorithm (or cache-transcendent algorithm) is an algorithm designed to take advantage of a processor cache without having the size of the cache (or the length of the cache lines, etc.) as an explicit parameter. An ...
consists of recursively dividing the dimensions into two groups and that are transformed recursively (rounding if ''d'' is not even) (see Frigo and Johnson, 2005). Still, this remains a straightforward variation of the row-column algorithm that ultimately requires only a one-dimensional FFT algorithm as the base case, and still has O(''N'' log ''N'') complexity. Yet another variation is to perform matrix transpositions in between transforming subsequent dimensions, so that the transforms operate on contiguous data; this is especially important for out-of-core In computing, external memory algorithms or out-of-core algorithms are algorithms that are designed to process data that are too large to fit into a computer's main memory at once. Such algorithms must be optimized to efficiently fetch and access d ...
and distributed memory
In computer science, distributed memory refers to a multiprocessor computer system in which each processor has its own private memory. Computational tasks can only operate on local data, and if remote data are required, the computational task mu ...
situations where accessing non-contiguous data is extremely time-consuming.
There are other multidimensional FFT algorithms that are distinct from the row-column algorithm, although all of them have complexity. Perhaps the simplest non-row-column FFT is the vector-radix FFT algorithm The vector-radix FFT algorithm, is a multidimensional fast Fourier transform (FFT) algorithm, which is a generalization of the ordinary Cooley–Tukey FFT algorithm that divides the transform dimensions by arbitrary radices. It breaks a multidimens ...
, which is a generalization of the ordinary Cooley–Tukey algorithm where one divides the transform dimensions by a vector of radices at each step. (This may also have cache benefits.) The simplest case of vector-radix is where all of the radices are equal (e.g. vector-radix-2 divides ''all'' of the dimensions by two), but this is not necessary. Vector radix with only a single non-unit radix at a time, i.e. , is essentially a row-column algorithm. Other, more complicated, methods include polynomial transform algorithms due to Nussbaumer (1977), which view the transform in terms of convolutions and polynomial products. See Duhamel and Vetterli (1990) for more information and references.
Other generalizations
An generalization to spherical harmonics
In mathematics and physical science, spherical harmonics are special functions defined on the surface of a sphere. They are often employed in solving partial differential equations in many scientific fields.
Since the spherical harmonics form ...
on the sphere ''S''2 with ''N''2 nodes was described by Mohlenkamp, along with an algorithm conjectured (but not proven) to have complexity; Mohlenkamp also provides an implementation in the libftsh library. A spherical-harmonic algorithm with complexity is described by Rokhlin and Tygert.
The fast folding algorithm In signal processing, the fast folding algorithm (Staelin, 1969) is an efficient algorithm for the detection of approximately- periodic events within time series data. It computes superpositions of the signal modulo various window sizes simultaneou ...
is analogous to the FFT, except that it operates on a series of binned waveforms rather than a series of real or complex scalar values. Rotation (which in the FFT is multiplication by a complex phasor) is a circular shift of the component waveform.
Various groups have also published "FFT" algorithms for non-equispaced data, as reviewed in Potts ''et al.'' (2001). Such algorithms do not strictly compute the DFT (which is only defined for equispaced data), but rather some approximation thereof (a non-uniform discrete Fourier transform, or NDFT, which itself is often computed only approximately). More generally there are various other methods of spectral estimation.
Applications
The FFT is used in digital recording, sampling, additive synthesis and pitch correction
Pitch correction is an electronic effects unit or audio software that changes the intonation (highness or lowness in pitch) of an audio signal so that all pitches will be notes from the equally tempered system (i.e., like the pitches on a piano) ...
software.
The FFT's importance derives from the fact that it has made working in the frequency domain equally computationally feasible as working in the temporal or spatial domain. Some of the important applications of the FFT include:
* fast large-integer multiplication algorithm
A multiplication algorithm is an algorithm (or method) to multiply two numbers. Depending on the size of the numbers, different algorithms are more efficient than others. Efficient multiplication algorithms have existed since the advent of the d ...
s and polynomial multiplication,
* efficient matrix–vector multiplication for Toeplitz, circulant and other structured matrices,
* filtering algorithms (see overlap–add and overlap–save methods),
* fast algorithms for discrete cosine or sine transforms (e.g. fast DCT used for JPEG and MPEG
The Moving Picture Experts Group (MPEG) is an alliance of working groups established jointly by ISO and IEC that sets standards for media coding, including compression coding of audio, video, graphics, and genomic data; and transmission and f ...
/MP3
MP3 (formally MPEG-1 Audio Layer III or MPEG-2 Audio Layer III) is a coding format for digital audio developed largely by the Fraunhofer Society in Germany, with support from other digital scientists in the United States and elsewhere. Origin ...
encoding and decoding),
* fast Chebyshev approximation,
* solving difference equations,
* computation of isotopic distributions.
* modulation and demodulation of complex data symbols using orthogonal frequency division multiplexing (OFDM) for 5G, LTE, Wi-Fi, DSL, and other modern communication systems.
Research areas
; Big FFTs: With the explosion of big data in fields such as astronomy, the need for 512K FFTs has arisen for certain interferometry calculations. The data collected by projects such as WMAP
The Wilkinson Microwave Anisotropy Probe (WMAP), originally known as the Microwave Anisotropy Probe (MAP and Explorer 80), was a NASA spacecraft operating from 2001 to 2010 which measured temperature differences across the sky in the cosmic mic ...
and LIGO
The Laser Interferometer Gravitational-Wave Observatory (LIGO) is a large-scale physics experiment and observatory designed to detect cosmic gravitational waves and to develop gravitational-wave observations as an astronomical tool. Two large ...
require FFTs of tens of billions of points. As this size does not fit into main memory, so called out-of-core FFTs are an active area of research.
; Approximate FFTs: For applications such as MRI, it is necessary to compute DFTs for nonuniformly spaced grid points and/or frequencies. Multipole based approaches can compute approximate quantities with factor of runtime increase.
; Group FFTs: The FFT may also be explained and interpreted using group representation theory
In the mathematical field of representation theory, group representations describe abstract groups in terms of bijective linear transformations of a vector space to itself (i.e. vector space automorphisms); in particular, they can be used to re ...
allowing for further generalization. A function on any compact group, including non-cyclic, has an expansion in terms of a basis of irreducible matrix elements. It remains active area of research to find efficient algorithm for performing this change of basis. Applications including efficient spherical harmonic
In mathematics and physical science, spherical harmonics are special functions defined on the surface of a sphere. They are often employed in solving partial differential equations in many scientific fields.
Since the spherical harmonics form ...
expansion, analyzing certain Markov process
A Markov chain or Markov process is a stochastic model describing a sequence of possible events in which the probability of each event depends only on the state attained in the previous event. Informally, this may be thought of as, "What happe ...
es, robotics etc.
; Quantum FFTs: Shor's fast algorithm for integer factorization on a quantum computer has a subroutine to compute DFT of a binary vector. This is implemented as sequence of 1- or 2-bit quantum gates now known as quantum FFT, which is effectively the Cooley–Tukey FFT realized as a particular factorization of the Fourier matrix. Extension to these ideas is currently being explored.
Language reference
See also
FFT-related algorithms:
* Goertzel algorithm
The Goertzel algorithm is a technique in digital signal processing (DSP) for efficient evaluation of the individual terms of the discrete Fourier transform (DFT). It is useful in certain practical applications, such as recognition of dual-tone mult ...
– computes individual terms of discrete Fourier transform
FFT implementations:
* ALGLIB
ALGLIB is a cross-platform open source numerical analysis and data processing library. It can be used from several programming languages (C++, C#, VB.NET, Python, Delphi).
ALGLIB started in 1999 and has a long history of steady development wi ...
– a dual/GPL-licensed C++ and C# library (also supporting other languages), with real/complex FFT implementation
* FFTPACK
FFTPACK is a package of Fortran subroutines for the fast Fourier transform. It includes complex, real, sine, cosine, and quarter-wave transforms. It was developed by Paul Swarztrauber of the National Center for Atmospheric Research, and is i ...
– another Fortran FFT library (public domain)
* Architecture-specific:
** Arm Performance Libraries
** Intel Integrated Performance Primitives
Intel Integrated Performance Primitives (Intel IPP) is a multi-threaded software library of functions for multimedia and data processing applications, produced by Intel.
The library supports Intel and compatible processors and is available for Li ...
** Intel Math Kernel Library
Intel oneAPI Math Kernel Library (Intel oneMKL; formerly Intel Math Kernel Library or Intel MKL) is a library of optimized math routines for science, engineering, and financial applications. Core math functions include BLAS, LAPACK, ScaLAPACK, ...
* Many more implementations are available, for CPUs and GPUs, such as PocketFFT for C++
Other links:
* Odlyzko–Schönhage algorithm In mathematics, the Odlyzko–Schönhage algorithm is a fast algorithm for evaluating the Riemann zeta function at many points, introduced by . The main point is the use of the fast Fourier transform to speed up the evaluation of a finite Dirichle ...
applies the FFT to finite Dirichlet series
In mathematics, a Dirichlet series is any series of the form
\sum_^\infty \frac,
where ''s'' is complex, and a_n is a complex sequence. It is a special case of general Dirichlet series.
Dirichlet series play a variety of important roles in analy ...
* Schönhage–Strassen algorithm
The Schönhage–Strassen algorithm is an asymptotically fast multiplication algorithm for large integers. It was developed by Arnold Schönhage and Volker Strassen in 1971.A. Schönhage and V. Strassen,Schnelle Multiplikation großer Zahlen, '' ...
– asymptotically fast multiplication algorithm for large integers
* Butterfly diagram – a diagram used to describe FFTs
* Spectral music
Spectral music uses the acoustic properties of sound – or sound spectra – as a basis for composition.
Definition
Defined in technical language, spectral music is an acoustic musical practice where compositional decisions are often inform ...
(involves application of DFT analysis to musical composition)
* Spectrum analyzer
A spectrum analyzer measures the magnitude of an input signal versus frequency within the full frequency range of the instrument. The primary use is to measure the power of the spectrum of known and unknown signals. The input signal that most co ...
– any of several devices that perform spectrum analysis, often via a DFT
* Time series
In mathematics, a time series is a series of data points indexed (or listed or graphed) in time order. Most commonly, a time series is a sequence taken at successive equally spaced points in time. Thus it is a sequence of discrete-time data. Ex ...
* Fast Walsh–Hadamard transform
* Generalized distributive law The generalized distributive law (GDL) is a generalization of the distributive property which gives rise to a general message passing algorithm. It is a synthesis of the work of many authors in the information theory, digital communications, signal ...
* Least-squares spectral analysis
Least-squares spectral analysis (LSSA) is a method of estimating a frequency spectrum, based on a least squares fit of sinusoids to data samples, similar to Fourier analysis. Fourier analysis, the most used spectral method in science, generally ...
* Multidimensional transform
In mathematical analysis and applications, multidimensional transforms are used to analyze the frequency content of signals in a domain of two or more dimensions.
Multidimensional Fourier transform
One of the more popular multidimensional tran ...
* Multidimensional discrete convolution
* Fast Fourier Transform Telescope
References
Further reading
*
*
*
*
*
*
* (NB. Contains extensive bibliography.)
* Elena Prestini: "The Evolution of Applied Harmonic Analysis", Springer, ISBN 978-0-8176-4125-2 (2004), Sec.3.10 'Gauss and the asteroids: history of the FFT'.
External links
Fast Fourier Transform for Polynomial Multiplication
fast Fourier algorithm
*
Fast Fourier Transforms
', Connexions online book edited by Charles Sidney Burrus, with chapters by Charles Sidney Burrus, Ivan Selesnick, Markus Pueschel, Matteo Frigo, and Steven G. Johnson (2008)
Fast Fourier transform — FFT
FFT programming in C++ the Cooley–Tukey algorithm
Online documentation, links, book, and code
* Sri Welaratna,
Thirty years of FFT analyzers
", ''Sound and Vibration'' (January 1997, 30th anniversary issue) a historical review of hardware FFT devices
ALGLIB FFT Code
a dual/GPL-licensed multilanguage (VBA, C++, Pascal, etc.) numerical analysis and data processing library
SFFT: Sparse Fast Fourier Transform
MIT's sparse (sub-linear time) FFT algorithm, sFFT, and implementation
VB6 FFT
a VB6 optimized library implementation with source code
a visual interactive intro to Fourier transforms and FFT methods
{{Authority control
FFT algorithms
Digital signal processing
Discrete transforms