Lesson 21

Fast Fourier Transform

The work of Cooley and Tukey showed how to calculate the DFT with complexity $N\log N$ (called the Fast Fourier Transform) instead of complexity $N^2$ using the direct algorithm. The fft command that you use in MATLAB implements a Fast Fourier Transform.

Examine:


\begin{displaymath}X[k] = \sum_{n=0}^{N-1} x[n] W^{kn} \end{displaymath}

There are approximately $N^2$ complex multiplications and additions required to implement this ($N$ for each value of $X[k]$).

If $N=2^{10}=1024$, then $N^2 = 2^{20} = 10^6$, a very large number!

However, the FFT would only require about 5000, a substantial savings in complexity (the actual calculation is $\frac{N}{2} \log_2N$).








There are a number of different FFT algorithms that exist including decimation-in-time and decimation-in-frequency.

The primary idea is to split up the size-$N$ DFT into $\frac{N}{2}$ DFTs of length 2 each.

You split the sum into 2 subsequences of length $\frac{N}{2}$ and continue all the way down until you have $\frac{N}{2}$ subsequences of length 2.

First break $x[n]$ into even and odd subsequences:


\begin{displaymath}X[k] = \sum_{n {\rm even}} x[n] W^{kn} + \sum_{n {\rm odd}} x[n] W^{kn} \end{displaymath}

Now let $n = 2m$ for even numbers and $n=2m+1$ for odd numbers:


\begin{displaymath}X[k] = \sum_{m=0}^{\frac{N}{2} -1} x[2m] W^{2mk} + \sum_{m=0}^{\frac{N}{2} -1}x[2m+1] W^{k(2m+1)} = \end{displaymath}


\begin{displaymath}\sum_{m=0}^{\frac{N}{2} -1} x[2m] (W^2)^{mk} + W^k \sum_{m=0}^{\frac{N}{2} -1}x[2m+1] (W^2)^{mk} =\end{displaymath}


\begin{displaymath}X_e[k] + W^k X_o[k] \end{displaymath}

$X_e[k]$ and $X_o[k]$ are both the DFT of a $\frac{N}{2}$ point sequence.

$W^k$ is often referred to as the ``twiddle factor.''

Now break up the size $\frac{N}{2}$ subsequences in half by letting $m=2p$:


\begin{displaymath}X_e[k] = \sum_{p=0}^{\frac{N}{4} -1} x[4p] (W^{4})^{kp} + W^{2k} \sum_{p=0}^{\frac{N}{4} -1}x[4p+2] (W^{4})^{kp} = \end{displaymath}

The first subsequence here is the terms $x[0], x[4], \dots$ and the second is $x[2], x[6], \dots.$

Also, we have that:


\begin{displaymath}W_N^{\frac{N}{2}} = -1\end{displaymath}



\begin{displaymath}Y[k] = \sum_{n=0}^{1} y[n] W_2^{kn} = y[0] + W_2^k y[1] \end{displaymath}


\begin{displaymath}W_2 = e^{-\frac{j2\pi}{2}} = e^{-j\pi} = -1\end{displaymath}

So we get,

\begin{displaymath}Y[k] = y[0] + (-1)^k y[1]\end{displaymath}

and:


\begin{displaymath}Y[0] = y[0] + y[1] \end{displaymath}


\begin{displaymath}Y[1] = y[0] - y[1] \end{displaymath}

\fbox{Ex.}

This was a problem I had on a DSP final exam in 1984:

Express the DFT of the 9-point sequence $\{ x[0], x[1], \dots, x[8] \}$ in terms of the DFTs of 3-point sequences $\{x[0],x[3],x[6]\}, \{x[1],x[4],x[7]\},$ and $\{x[2],x[5],x[8]\}$

We start with:


\begin{displaymath}X[k] = \sum_{m=0}^2 x[3m] W_9^{3mk} +\sum_{m=0}^2 x[3m+1] W_9^{(3m+1)k} + \sum_{m=0}^2 x[3m+2] W_9^{(3m+2)k} \end{displaymath}