Radix-2 Decimation in Time FFT Algorithm

Contents

Introduction

In the previous sections we considered DFT
and its properties
and also showed
the Fast Fourier Transform (FFT) principle on the basis of the divide and conquer procedure.

In this section we will consider one of the most widespread FFT algorithms: radix-2 decimation in time FFT [1-4].

Let’s give the DFT equation once again:

$$
S(k) = \sum_{n = 0 }^{N-1} s(n) \cdot \exp \left( -j \cdot \frac{2\pi}{N} \cdot n \cdot k\right), \text{ } k = 0 \ldots N-1.
$$

(1)

The divide and conquer procedure was considered in the previous section.

We will operate on twiddle factors:

$$
W_N^{ k} = \exp \left( -j \cdot \frac{2\pi}{N} \cdot k\right),\text{ } k = 0 \ldots N-1,
$$

(2)

as we already did upon considering DFT properties. Then (1) taking into account (2) takes the more simplified form:

$$
S(k) = \sum_{n = 0 }^{N-1} s(n) \cdot W_N^{n \cdot k}, \text{ } k = 0 \ldots N-1.
$$

(3)

Initial Sequence Divide by Decimation in Time

Decimation in time consists of dividing the input sequence $s(n)$, $n = 0 \ldots N-1$,
into two sequences of half duration $s_0(m)$ and $s_1(m)$, $m = 0 \ldots \frac{N}{2}-1$,
therefore, $s_0(m) = s(2 \cdot m)$, and $s_1(m) = s(2 \cdot m + 1)$.
The sequence $s_0(m)$ contains samples with even indexes, and $s_1(m)$ with odd ones.
Decimation in time for $N = 8$ is visually presented in Figure 1.

Figure 1. Decimation in Time for $N=8$

Divide and Conquer Procedure

Consider DFT of the signal which is decimated in time:

$$
S(k) = \sum_{m = 0}^{\frac{N}{2}-1} s( 2 \cdot m) \cdot W_N^{2 \cdot m \cdot k} +
\sum_{m= 0}^{\frac{N}{2}-1} s( 2 \cdot m+1) \cdot W_N^{(2 \cdot m+1) \cdot k} = \\
= \sum_{m= 0}^{\frac{N}{2}-1} s( 2 \cdot m) \cdot W_N^{2 \cdot m \cdot k} +
W_N^k \cdot \sum_{m= 0}^{\frac{N}{2}-1} s( 2 \cdot m+1) \cdot W_N^{2 \cdot m \cdot k} , \text{ } k = 0 \ldots N-1.
$$

(4)

If we consider only the first DFT half $S(k)$ for indexes $k = 0 \ldots \frac{N}{2}-1$,
and also to take into account that

$$
W_N^{2 \cdot m \cdot k} = \exp \left( -j \cdot \frac{2\pi}{N} \cdot 2 \cdot m \cdot k \right) =
\exp \left( -j \cdot \frac{2\pi}{\frac{N}{2}} \cdot m \cdot k \right) = W_{\frac{N}{2}}^{m \cdot k},
$$

(5)

then (4) will be transformed to the form:

$$
S(k)
= \sum_{m= 0}^{\frac{N}{2}-1} s( 2 \cdot m) \cdot W_{\frac{N}{2}}^{m \cdot k} +
W_N^k \cdot \sum_{m= 0}^{\frac{N}{2}-1} s( 2 \cdot m+1) \cdot W_{\frac{N}{2}}^{m \cdot k} =
$$
$$
= S_0(k) + W_N^k \cdot S_1(k), \text{ } k = 0 \ldots \frac{N}{2}-1,
$$

(6)

where

$$
S_0(k) = \sum_{m= 0}^{\frac{N}{2}-1} s_0(m) \cdot W_{\frac{N}{2}}^{m \cdot k}, \text{ } k = 0 \ldots \frac{N}{2}-1;
$$
$$
S_1(k) = \sum_{m= 0}^{\frac{N}{2}-1} s_1(m) \cdot W_{\frac{N}{2}}^{m \cdot k}, \text{ } k = 0 \ldots \frac{N}{2}-1,
$$

(7)

$\frac{N}{2}-$points DFT of decimated sequences $s_0(m)$ and $s_1(m)$, $m = 0 \ldots \frac{N}{2}-1$ respectively.

Thus, decimation in time can be considered as an algorithm of dividing the sequence into two sequences of half duration. The first DFT half is the DFT sum $S_0(k)$ of the “even” sequence $s_0(m)$ and DFT $S_1(k)$ is of the “odd” sequence $s_1(m)$ multiplied by twiddle factors $W_N^k$.

Now analyze the second DFT half $S\left( k+\frac{N}{2} \right)$ for $k = 0 \ldots \frac{N}{2}-1$:

$$
S\left( k+\frac{N}{2} \right) = \sum_{m= 0}^{\frac{N}{2} - 1} s(2\cdot m) W_{N}^{2\cdot m \cdot \left( k+\frac{N}{2}\right)} + \sum_{m= 0}^{\frac{N}{2} - 1} s(2\cdot m+1) W_{N}^{(2\cdot m + 1) \cdot \left( k+\frac{N}{2}\right)}, \text{ } k = 0 \ldots \frac{N}{2} - 1.
$$

(8)

Consider twiddle factors $W_{N}^{2\cdot m \cdot \left( k+\frac{N}{2}\right)} $in more detail:

$$
W_{N}^{2\cdot m \cdot \left( k+\frac{N}{2}\right)} =
\underbrace{W_{N}^{2\cdot m \cdot k}}_{W_{\frac{N}{2}}^{ m \cdot k}}
\cdot W_{N}^{2 \cdot m \cdot \frac{N}{2}} =
W_{\frac{N}{2}}^{ m \cdot k} \cdot \underbrace{W_{N}^{ m \cdot N}}_{1} =
W_{\frac{N}{2}}^{ m \cdot k}.
$$

(9)

It is possible to simplify $W_{N}^{(2\cdot m + 1) \cdot \left( k+\frac{N}{2}\right)}$ in a similar way:

$$
W_{N}^{(2\cdot m + 1) \cdot \left( k+\frac{N}{2}\right)} =
\underbrace{W_{N}^{2\cdot m \cdot k}}_{W_{\frac{N}{2}}^{m \cdot k}}
\underbrace{\cdot W_{N}^{2 \cdot m \cdot \frac{N}{2}}}_{1}
\cdot W_{N}^{k}
\underbrace{\cdot W_{N}^{\frac{N}{2}}}_{-1} = -W_{N}^{k} \cdot W_{\frac{N}{2}}^{m \cdot k}.
$$

(10)

Then (8) taking into account (9) and (10) takes the form:

$$
S\left( k+\frac{N}{2} \right) = \sum_{m= 0}^{\frac{N}{2} - 1} s(2\cdot m) W_{\frac{N}{2}}^{m \cdot k} - W_{N}^{k} \cdot \sum_{m= 0}^{\frac{N}{2} - 1} s(2\cdot m+1) \cdot W_{\frac{N}{2}}^{m \cdot k} = \\
= S_0(k) - W_N^k \cdot S_1(k), \text{ } k = 0 \ldots \frac{N}{2} - 1.
$$

(11)

Using the equation for the first (6) and the second (11) half of DFT, it is definitely possible to write "the conquer procedure" as:

\begin{eqnarray}
S(k) &=& S_0(k) + W_N^k \cdot S_1(k), \text{ } k = 0 \ldots \frac{N}{2} - 1;
\\
S\left(k+\frac{N}{2} \right) &=& S_0(k) - W_N^k \cdot S_1(k), \textrm{ } k = 0 \ldots \frac{N}{2} - 1.
\end{eqnarray}

(12)

Butterfly Diagram

Equation (12) combines two $\frac{N}{2}-$points DFT $S_0(k)$ and $S_1(k)$, $k = 0 \ldots \frac{N}{2}-1$, of decimated signals of half duration $s_0(m) = s(2\cdot m)$ and $s_1(m) = s(2\cdot m+1)$, $m = 0 \ldots \frac{N}{2}-1$, in the resulting $N-$points DFT $S(k)$, $k= 0 \ldots N-1$, of a input signal.

Graphically the divide and conquer process is presented in Figure 2.

Figure 2. Butterfly Diagram

Because of a specific diagram form it takes the name “butterfly”. This divide and conquer procedure is the main one upon creating radix-2 FFT algorithms.

In Figure 3 the FFT algorithm diagram according to (12) is presented for $N = 8$.

Figure 3. Decimation in Time FFT Algorithm Diagram for $N = 8$

Full Diagram of Decimation in Time FFT Algorithm. Bit Reversal

Equation (12) presents the conquer procedure to calculate $N-$points DFT $S(k)$, $k = 0 \ldots N-1$, by means of two $\frac{N}{2}-$points DFT $S_0(k)$ and $S_1(k)$, $k = 0 \ldots \frac{N}{2}-1,$ of even and odd decimated sequences $s(2 \cdot m)$ and $s(2\cdot m+1)$, $m = 0 \ldots \frac{N}{2}-1$.

The same procedure can be applied to calculate each of $\frac{N}{2}-$points DFT $S_0(k)$ and $S_1(k)$ by means of two $\frac{N}{4}-$points DFT. Then for $N = 2^L$ it is possible to perform $L-1$ sequence division step into “even” and “odd” and after that to conquer the spectrum with $L$ steps. As a result we will get the full FFT algorithm diagram. As an example a full decimation in time FFT diagram for $N = 8$ is given in Figure 4.

Figure 4. Decimation in Time FFT Algorithm Diagram for $N=8$

At the first step input signal samples are interchanged and the initial sequence is divided into “even” and “odd” sequences. Then “even” and “odd” sequences, in their turn, are divided into “even” and “odd” sequences again.

This procedure is called bit reversal, in such a way it is possible to renumber samples having rewritten the sample number in the binary system in the reversed direction.

For example, $s(4)$ has an index in the decimal number system $4_{10} = 100_2,$. If $100_2$ is rewritten from right to left, then we will get $001_2$, that is $s(4)$ after division into “even-odd” before the first operation the “butterfly” will take the place of the sample $s(1)$ which in its turn will take the place $s(4)$.

In a similar way all samples will be interchanged, at the same time some of them will save their places, in particular $s(2)$ because if $2_{10} = 010_2$ is rewritten from right to left, then all the same $010_2$ will save its place, similarly $s(0)$, $s(5)$ and $s(7)$.

It is important to note that this renumbering method has to be applied upon putting the number into the binary system consisting from $L$ bits. In the given example $N = 2^L = 8$, so 3 bits of binary number were used, but if $N$ were equal to 16, then it would be necessary to write down the number upon using 4 bits. In this case $2_{10} = 0010_2,$ and after reversal we will get $0100_2 = 4_{10}$ that is when $N=16$, the sample $s(2)$ will not save its place, and will be interchange with $s(4)$.

After the bit reversal we get four 2-points DFT:

\begin{eqnarray}
S_{00}(0) = s(0) + W_2^0 \cdot s(4);\nonumber \\
S_{00}(1) = s(0) - W_2^0 \cdot s(4);\nonumber \\
\nonumber \\
S_{01}(0) = s(2) + W_2^0 \cdot s(6);\nonumber \\
S_{01}(1) = s(2) - W_2^0 \cdot s(6);\nonumber \\
\\
S_{02}(0) = s(1) + W_2^0 \cdot s(5);\nonumber \\
S_{02}(1) = s(1) - W_2^0 \cdot s(5);\nonumber \\
\nonumber \\
S_{03}(0) = s(3) + W_2^0 \cdot s(7);\nonumber \\
S_{03}(1) = s(3) - W_2^0 \cdot s(7).\nonumber
\end{eqnarray}

(13)

Two 4-points DFT are formed on the basis of four 2-points DFT:

\begin{eqnarray}
S_{10}(0) = S_{00}(0) + W_4^0 \cdot S_{01}(0);\nonumber \\
S_{10}(1) = S_{00}(1) + W_4^1 \cdot S_{01}(1);\nonumber \\
S_{10}(2) = S_{00}(0) - W_4^0 \cdot S_{01}(0);\nonumber \\
S_{10}(3) = S_{00}(1) - W_4^1 \cdot S_{01}(1);\nonumber \\
\\
S_{11}(0) = S_{02}(0) + W_4^0 \cdot S_{03}(0);\nonumber \\
S_{11}(1) = S_{02}(1) + W_4^1 \cdot S_{03}(1);\nonumber \\
S_{11}(2) = S_{02}(0) - W_4^0 \cdot S_{03}(0);\nonumber \\
S_{11}(3) = S_{02}(1) - W_4^1 \cdot S_{03}(1).\nonumber
\end{eqnarray}

(14)

And the full input signal DFT is formed at the last level:

\begin{eqnarray}
S(0) = S_{10}(0) + W_8^0 \cdot S_{11}(0);\nonumber \\
S(1) = S_{10}(1) + W_8^1 \cdot S_{11}(1);\nonumber \\
S(2) = S_{10}(2) + W_8^2 \cdot S_{11}(2);\nonumber \\
S(3) = S_{10}(3) + W_8^3 \cdot S_{11}(3);\nonumber \\
S(4) = S_{10}(0) - W_8^0 \cdot S_{11}(0);\nonumber \\
S(5) = S_{10}(1) - W_8^1 \cdot S_{11}(1);\nonumber \\
S(6) = S_{10}(2) - W_8^2 \cdot S_{11}(2);\nonumber \\
S(7) = S_{10}(3) - W_8^3 \cdot S_{11}(3).
\end{eqnarray}

(15)

Conclusions

In this section we considered one of the most widespread fast Fourier transform algorithms: an algorithm with decimation in time.

The divide procedure by decimation in time is shown and the divide and conquer procedure on the basis of the Butterfly Diagram according to (12) is given.

The complete FFT diagram with decimation in time was given for $N = 8$.

In the following section we will consider one more radix-2 FFT algorithm: decimation in frequency FFT algorithm.

You can provide your any feedbacks, questions and wishes in the message board or в guestbook.

Reference

[1]
James W. Cooley and John W. Tukey
An Algorithm for the Machine Calculation of Complex Fourier Series.
Mathematics of Computation, 1965 p. 297–301.

[2]
Richard E. Blahut
Fast Algorithms for Signal Processing.
Cambridge University Press, 2010.

[3]
Nussbaumer Henri J.
Fast Fourier Transform and Convolution Algorithms.
Second Corrected and Updated Edition.
Springer-Verlag, 1982.

[4]
Oppenheim Alan V. and Schafer Ronald W.
Discrete-Time Signal Processing
Second Edition.
Prentice-Hall, New Jersey, 1999.