For a length N finite sequence {x[n]}0n−1, the Discrete Fourier Transform of the signal is a length N finite sequence {X[k]}0n−1 where
X[k]=∑n=0N−1x[n]e−jN2πkn
One way to interpret the DFT is in terms of the Fourier series for a disrete periodic signal x~[n]=x[((n))N] where the ((n))N=nmodN. Recall that the coefficient of the kth term of the Fourier Series is
ak=N1∑n=0N−1x[n]e−jN2πkn
Notice that the ak of the Fourier Series are the DFT values except scaled by a factor of N. In other words, if we extend a finite signal periodically, then the DFT and the DTFS are the same up to a constant scale factor. This gives an intuitive inverse DFT.
For a length N finite sequence {X[k]}0N−1 representing the DFT of a finite perioidc signal {x[n]}0N−1, the inverse DFT is given by
x[n]=N1∑k=0N−1X[k]ejN2πkn
Notice that the DFT and the IDFT are very similar in form. It turns out that the IDFT can be expressed as a DFT of X∗[k]. Namely
IDFT{X[k]}=N1DFT{X⋆[k]}⋆
Further intuition for the DFT comes by relating it to the DTFT. Suppose we have a finite signal x[n] which is 0 for n<0 and n>N−1. The DTFT of this signal is
X(ω)=∑n=−∞∞x[n]e−jωn=∑n=0N−1x[n]e−jωn
Suppose we sample the DTFT at intervals of N2πk, then the kth sample is given by
X[k]=X(N2πk)=∑n=0N−1x[n]e−jN2πkn
Thus we can think of the DFT as a N evenly spaced samples of the DTFT. One important point to notice is that while the DTFT is often centered around 0 (meaning it is plotted over a range from −π to π), because we are summing from 0 to N-1 in the DFT, the DFT coefficients are centered around π, and thus they are plotted on a range fo [0,2π−N2π]
x[n]⊛y[n]↔X[k]Y[k]
x[n]⊛y[n]=∑m=0N−1x[m]y[((n−m))N]
Because multiplying DFT coefficients performs a specific case of convolution, we can compute a linear convolution using the circular convolution. Suppose we have two finite signals {x[n]}0L−1 and {h[n]}0P−1 The linear convolution of these two signals will be length L+P−1, so in order to take an IDFT and get L+P−1 samples, we need to take at least N≤L+P−1 points.
Pad each vector to length L+P−1
Compute X[k]H[k]
If N is smaller than L−P+1, the result is akin to aliasing in the time domain. To see why, consider that the DFT coefficients are essentially the DTFS coefficients of the periodic extension of x[n] (denote x~[n]).
x~[n]=∑r=−∞∞x[n−rN]
Y(ejω)=X(ejω)H(ejω)
y~[n]=∑r=−∞∞y[n−rN].
Notice that if N is not large enough, then these copies will be overlapping (a.k.a aliasing). Since the DFT is just sampling the DTFT, the circular convolution will represent the true convolution so long as the copies don’t overlap.
In a discrete time system, the input signal might have a very long length, making it impractical to be stored in a computer’s memory or to compute the DFT of it all at once (especially if we have a real-time system). Thus to compute the output of a digital filter (with impulse response of length P), we need to compute the DFT in blocks shorter than the signal.
Decompose x[n] into nonoverlapping segments of length L
x[n]=∑rxr[n]xr[n]={x[n]0rL≤n≤(r+1)Lelse.
y[n]=x[n]∗h[n]=∑rxr[n]∗h[n].
Zero pad xr[n] and h[n] to length N≥L+P−1 to prevent time-domain aliasing
The neighboring outputs overlap in the last P−1 points, so add the overlapping sections together to get the final output
Divide x[n] into sections of length L such that each section overlaps the previous by P−1 points
xr[n]=x[n+r(L−P+1)−P+1]0≤n≤L−1
Zero pad xr[n] and h[n] to length N≥L+P−1 to prevent time domain aliasing.
The first P−1 samples of the output will be incorrect, so we can discard them.
y[n]=∑r=0∞yr[n−r(L−P+1)+P−1]yr[n]={xr[n]∗h[n]0P−1≤n≤L−1 else
The DFT gives us an easy way to do convolutions. Unfortunately naiively, computing it is an O(N2) operation because we must sum together N elements to compute N different coefficients. Thankfully, there is a fast algorithm which can compute the DFT in O(NlogN) time so we can compute convolutions quickly.
The fast fourier transform is an algorithm which computes the DFT efficiently in O(NlogN)time.
The Nth roots of unity are the complex roots of WNN=1.
WNk=e−jN2πk
WNN−kn=WN−kn=(WNkn)⋆
Wkn=WNk(n+N)=WN(k+N)n
When squared, the Nth roots of unity become the 2Nth roots of unity.
WN2=W2N
Using Theorem 1, Theorem 2, Theorem 3, we can take two approaches to the FFT: decimation in time, which splits x[n] into smaller subsequences, and decimation in frequency which splits X[n] into smaller subsequences.
The idea here is too break x[n] into smaller subsequences. We assume that N is a power of 2 for simplicity.
X[k]=∑n=0N−1x[n]WNkn=∑n evenx[n]WNkn+∑n oddx[n]WNkn
We let n=2r and n=2r+1.
X[k]=r=0∑2N−1x[2r]WN2rk+r=0∑2N−1x[2r+1]WNk(2r+1)=r=0∑2N−1x[2r]W2Nrk+WNkr=0∑2N−1x[2r+1]W2Nkr
∴X[k]=G[k]+WNkH[k]
Both G and H are 2N periodic, and notice that
WNk+2N=e−jN2π(k+2N)=−WNk.
This means once we compute G[k] and H[k] we can compute X[k] easily because
X[k]=G[k]+WNkH[k]X[k+2N]=G[k]−WNkH[k]for k∈[0,2N)
We can continue this relationship recursively downwards. Once we get too N=2, we can represet this as a simple butterfly operation:
X[0]=x[0]+x[1]X[1]=x[0]−x[1].
X[2r]X[2r+1]=n=0∑2N−1x[n]WN2rn+n=0∑2N−1x[n+2N]WN2r(n+2N)=W2Nrnn=0∑2N−1(x[n]+x[n+2N])=W2Nrnn=0∑2N−1(x[n]−x[n+2N])