Frequency Analysis from the Basics (9) "Fast Fourier Transform (FFT)"

In the previous article, we explained the derivation of the Discrete Fourier Transform (DFT). The DFT is a fundamental and extremely useful formula when performing Fourier transforms numerically. However, calculating it according to the formula requires an enormous amount of computation and has not been widely used in practice. Since the development of the Fast Fourier Transform (FFT) algorithm by Cooley and Turkey in 1965, it has become very popular and is now one of the most useful technologies in digital signal processing. Recently, even TV broadcasting has gone digital, and FFT technology is used in a digital modulation technique called OFDM.

This time, I will briefly explain the calculation method for FFT (Fast Fourier Transform).

As explained last time, DFT is used to obtain the time signals of N points. x_n Frequency spectrum X at N points (n=0, ..., N-1) _k (k=0, ..., N-1) is obtained. That is;

.................................（1）

Here,

.................................（2）

Normally, the DFT in equation (1) is calculated using complex numbers, and N2 complex multiplications are required to find Xk for N points. However, if N can be expressed as a product of multiple integers, the Fast Firmware Technique (FFT) is a method that performs the DFT calculation quickly by cleverly decomposing the operation in equation (1) into smaller groups and using the periodicity of W (period N) to reduce the number of multiplications as much as possible.

The number of calculation points N (the number of sample points in the time function) in an FFT is usually chosen to be a power of 2 (N=2m) due to the simplicity of the calculation partition and compatibility with the implementation of digital computer programs. For example, N=1024 (m=10), N=2048 (m=11), etc. This is called a two-basis FFT.

In the following explanation, we will use N=8 (m=3). Substituting N=8 into equation (1), we get:

.................................（3）

.................................（4）

This is how it works. Here, W has the function of rotating the unit circle by π/4 (45°) in the negative angular direction (clockwise), so it is called the rotation factor. This rotation factor W has the same value for a full rotation, and the sign is reversed for a half rotation. The specific points on the unit circle when N=8 are shown in Figure 1.

Next, if we divide the time signal x _n into even and odd sequences, then equation (3) becomes:

.................................（5）

Here, let r = 0, ..., 3.

Equation (5) shows that it is the sum of the even-term 4-point DFT and the odd-term 4-point DFT multiplied by the rotation factor ^Wk. Note that ^Wk = -Wk^-4 (Figure 1), and the flowchart for this operation is shown in Figure 2.

In exactly the same way, by decomposing a 4-point DFT into two 2-point DFTs, the final FFT of an 8-point DFT will have a signal flow diagram like the one in Figure 3. This method is called the time decimation method because it obtains the result by dividing the time signal.

If you look closely at the flowchart in Figure 3, you will see that every calculation sequence starts from the same pair of nodes as the previous sequence. In other words, the FFT operation is based on product-sum pairs of the form given by equation (6).

.................................（6）

When this operation is represented as a flowchart, it looks like Figure 4, and because it resembles the shape of a butterfly's wings, it is called the butterfly operation.

Equation (6) means that if we denote a pair of nodal data as xα and xβ, this data is retrieved from the memory of the time signal (complex number), the calculation on the right-hand side is performed, and the result is stored in the same memory of xα and xβ. In this way, the butterfly operation, which is the basic operation of FFT, can be performed with only a pair of nodal data, and is therefore called an in-place memory operation. It has the advantage of being very efficient because it requires less memory and fewer memory accesses.

There is one more point to note regarding the flowchart in Figure 3. The sequence of time signals xn is no longer in the usual sampled order. Therefore, before performing the FFT operation, it is necessary to rearrange it as shown in Figure 3. This can be easily done by converting the ordinal numbers of the sequence into binary.

Table 1 shows an example for the case where N=8. In this way, by rewriting the sequence number of the original sequence into binary and then reversing the bit sequence, it is possible to derive the arrangement shown in Figure 3. This operation is called bit reverse (BTR). This operation can also be easily performed on a standard computer.

Table 1: Sorting by bit reverse order when N=8

Next, by dividing the DFT calculation result X _k (k=0, ..., N-1) into even and odd numbers, we can derive the frequency decimation method in the same way as the time decimation method. Figure 5 shows the signal flow diagram of the 8-point FFT using the frequency decimation method when N=8.

In this case, the butterfly operation is given by equation (7), and the signal flow diagram is shown in Figure 6.

.................................（7）

As shown in Figure 5, in the frequency decimation method, the resulting _Xk is usually not in frequency order, requiring BTR processing. By performing the reverse processing shown in Table 1, the correct ordered sequence can be obtained.

There are various other calculation methods, such as applying BTR to the result using time decimation or applying BTR to the frequency decimation method first. However, these methods are not commonly used because they require BTR processing when accessing the rotation factor W.

Let's consider how much the computation time can be reduced by this FFT operation compared to a normal DFT. Assuming that the butterfly operation in equation (6) performs one complex number multiplication, the N-point FFT performs NN log ₂ 2 butterfly operations. Compared to the number of multiplications N2 in the DFT, this is an astonishing reduction: approximately 1/200 when N=1024 (= ^2¹⁰) and approximately 1/2300 when N=16384 (= ^2¹⁴).
Furthermore, reducing the number of multiplication steps also leads to a reduction in calculation errors due to rounding errors, which is considered a major advantage of the FFT algorithm.

Finally, here's a summary.

1. FFT is an algorithm that performs DFT calculations at high speed.
2. The number of points N used in the FFT calculation is usually a power of 2, such as 1024 or 2048, and this is called a two-basis FFT.
3. There are two main types of FFT calculations: time decimation and frequency decimation.
4. The basic operation of FFT is called the butterfly operation, and since it is an in-place memory operation, it can save memory.
5. Compared to the standard DFT, the FFT can significantly reduce the number of multiplications, so it can be expected to not only speed up calculations but also reduce calculation errors.

【keyword】
DFT (Discrete Fourier Transform), FFT (Fast Fourier Transform), 2-Base FFT, Rotation Factor, Time Decimation, Butterfly Operation, Bit Reverse (BTR), In-Place Memory Operation, Frequency Decimation

[Reference materials]

1. "The Fast Fourier Transform" by E. Oran Brigham, published by Science and Technology Press.
2. "Digital Fourier Analysis (I) - Fundamentals" by Kenichi Kido, published by Corona Publishing Co.
3. "Signal Processing" by Iwao Morishita and Hidefumi Obata, published by Asakura Shoten.

(Excerpt from the email newsletter issued on May 23, 2013)