Fast Fourier Transform (FFT) - Decomposing Signals with Divide and Conquer Speed

February 6, 2023 — #DSA #Divide and Conquer

Hey everyone, and welcome back to the blog! Ever wondered how your phone processes audio for calls, how JPEG images compress so well, or how engineers analyze complex signals in fields ranging from radio astronomy to medical imaging? A fundamental mathematical tool often at the heart of these processes is the Discrete Fourier Transform (DFT), and the algorithm that makes its computation feasible for large datasets is the Fast Fourier Transform (FFT).

The FFT is a cornerstone algorithm, widely used in signal processing, data compression, solving partial differential equations, and many other scientific and engineering disciplines. It's a brilliant example of how the "Divide and Conquer" strategy can lead to massive performance gains. Let's tune in!

The Problem: Understanding Frequency Components with DFT

The Discrete Fourier Transform (DFT) is a mathematical transformation that converts a finite sequence of equally-spaced samples of a function (often a time-domain signal) into a same-length sequence of equally-spaced samples of its frequency-domain representation. In simpler terms, it decomposes a signal into the different frequencies that make it up.

Given a sequence of $N$ complex numbers $x_0, x_1, \ldots, x_{N-1}$ , its DFT, denoted by $X_k$ , is defined as:

X_k = \sum_{n=0}^{N-1} x_n \cdot e^{-i 2 \pi k n / N}

for $k = 0, 1, \ldots, N-1$ . Here, $e$ is the base of the natural logarithm, $i$ is the imaginary unit ( $i^2 = -1$ ), and $\pi$ is Pi. Each $X_k$ represents the amplitude and phase of a specific frequency component in the original signal.

Conversely, the Inverse Discrete Fourier Transform (IDFT) allows us to reconstruct the original sequence from its frequency components:

x_n = \frac{1}{N} \sum_{k=0}^{N-1} X_k \cdot e^{i 2 \pi k n / N}

A naive, direct computation of the DFT involves $N$ summations for each of the $N$ output values $X_k$ , leading to a time complexity of O(N²). For large $N$ (e.g., processing a few seconds of audio sampled at 44.1 kHz), this becomes prohibitively slow. This is where the FFT shines.

The Solution: Fast Fourier Transform (FFT) - Divide and Conquer!

The Fast Fourier Transform (FFT) is not a different transform but rather a highly efficient algorithm (or family of algorithms) to compute the DFT (and its inverse). The most common FFT algorithms, like the Cooley-Tukey algorithm, dramatically reduce the time complexity from O(N²) to O(N log N).

The magic of the FFT lies in its divide and conquer approach:

Divide: If the length of the input sequence $N$ is not 1 (and typically, for the simplest versions, $N$ is a power of 2), the algorithm splits the input sequence $x_n$ into two smaller subsequences:
- One subsequence containing the even-indexed elements: $x_0, x_2, \ldots, x_{N-2}$ .
- Another subsequence containing the odd-indexed elements: $x_1, x_3, \ldots, x_{N-1}$ .
Each of these subsequences will have $N/2$ elements.
Conquer: Recursively compute the DFT of these two $N/2$ -point subsequences.
- Let $X_{even}$ be the DFT of the even-indexed subsequence.
- Let $X_{odd}$ be the DFT of the odd-indexed subsequence.
Combine: Combine the results from the DFTs of the even and odd subsequences to compute the DFT of the original $N$ -point sequence. This combination step is where the cleverness lies. For $k = 0, 1, \ldots, N/2 - 1$ , the DFT components $X_k$ are calculated as:
- $X_k = X_{even}[k] + W_N^k \cdot X_{odd}[k]$
- $X_{k + N/2} = X_{even}[k] - W_N^k \cdot X_{odd}[k]$
where $W_N^k = e^{-i 2 \pi k / N}$ is called a "twiddle factor" (a complex root of unity). Notice how results for the first half ( $X_k$ ) and the second half ( $X_{k+N/2}$ ) are related, reusing computations.

The base case for the recursion is when the length of the input sequence is 1. The DFT of a single-element sequence $x_0$ is just $X_0 = x_0$ .

Pseudocode for the Recursive FFT (Cooley-Tukey Radix-2 Decimation-in-Time)

Here's a simplified pseudocode for a common recursive FFT algorithm, assuming $N$ is a power of 2:

function fft(x_sequence):
    N = length(x_sequence)

    // Base case: If the sequence has only one element
    if N == 1:
        return x_sequence
    else:
        // Split into even and odd indexed subsequences
        x_even_subsequence = array_of_even_indexed_elements(x_sequence) // e.g., [x[0], x[2], ..., x[N-2]]
        x_odd_subsequence  = array_of_odd_indexed_elements(x_sequence)  // e.g., [x[1], x[3], ..., x[N-1]]

        // Recursively compute FFT for subsequences
        X_even = fft(x_even_subsequence)
        X_odd  = fft(x_odd_subsequence)

        // Combine results
        X_combined = new complex_array of length N

        for k = 0 to (N/2) - 1:
            // Twiddle factor: exp(-i * 2 * pi * k / N)
            twiddle_factor_k = complex_exponential(-2 * PI * k / N)
           
            term_odd_multiplied = twiddle_factor_k * X_odd[k]
           
            X_combined[k]         = X_even[k] + term_odd_multiplied
            X_combined[k + N/2]   = X_even[k] - term_odd_multiplied
           
        return X_combined

Explanation of the Combination Step

The key to the FFT's efficiency is this combination step. The DFT definition can be rewritten by splitting the sum over even and odd indices: $X_k = \sum_{m=0}^{N/2-1} x_{2m} e^{-i 2 \pi k (2m) / N} + \sum_{m=0}^{N/2-1} x_{2m+1} e^{-i 2 \pi k (2m+1) / N}$ $X_k = \sum_{m=0}^{N/2-1} x_{2m} e^{-i 2 \pi k m / (N/2)} + e^{-i 2 \pi k / N} \sum_{m=0}^{N/2-1} x_{2m+1} e^{-i 2 \pi k m / (N/2)}$

Recognizing that the sums are DFTs of the $N/2$ -point even and odd subsequences ( $X_{\text{even}}[k]$ and $X_{\text{odd}}[k]$ respectively, evaluated at frequency $k$ ), we get: $X_k = X_{\text{even}}[k] + W_N^k \cdot X_{\text{odd}}[k]$

For $k \ge N/2$ , we can substitute $k' = k - N/2$ . Using the properties of complex exponentials ( $W_N^{k+N/2} = -W_N^k$ and $X_{\text{even/odd}}[k+N/2] = X_{\text{even/odd}}[k]$ due to periodicity of $N/2$ -point DFTs), we arrive at the formula for $X_{k+N/2}$ . This reuse of computations is what dramatically reduces the operation count.

Correctness of the FFT Algorithm

The correctness of the Fast Fourier Transform algorithm can be rigorously proven using mathematical induction and the fundamental properties of the Discrete Fourier Transform.

Base Case: The DFT of a single-element sequence is the element itself, which the recursive algorithm correctly returns.
Inductive Step: Assuming the recursive calls correctly compute the DFTs for the $N/2$ -point even and odd subsequences, the mathematical derivation of the combination step (as sketched above) shows that these $N/2$ -point DFTs are correctly combined to produce the $N$ -point DFT of the original sequence.

The properties of the complex exponential function (the twiddle factors) are central to this derivation.

Complexity Analysis

Time Complexity: The FFT algorithm, through its divide and conquer approach, achieves a time complexity of O(N log N), where $N$ is the length of the input sequence. The recurrence relation for the time complexity is $T(N) = 2T(N/2) + O(N)$ (where O(N) is for the splitting and combining steps). This recurrence solves to O(N log N). This is a monumental improvement over the O(N²) complexity of the naive DFT computation, especially for large $N$ .
Space Complexity:
- A recursive implementation of FFT typically requires O(N log N) space due to the recursion call stack depth (log N) and the storage of intermediate results at each level (N).
- Iterative (in-place) versions of FFT algorithms can achieve O(N) space complexity, or even O(1) auxiliary space if overwriting the input array is allowed and bit-reversal permutation is handled carefully.

Applications of FFT: Beyond the Math

The Fast Fourier Transform is not just an elegant algorithm; it's a workhorse with a vast range of practical applications across numerous fields:

Signal Processing:
- Audio Signal Processing: Analyzing frequency content of sound, filtering, equalization, audio effects (reverb, chorus), noise reduction, audio compression (like MP3).
- Image Processing: Image filtering, compression (like JPEG), edge detection, pattern recognition.
- Video Compression: Many video codecs use DFT/FFT principles.
Data Compression: As mentioned, fundamental to lossy compression schemes for audio and images.
Solving Partial Differential Equations (PDEs): Spectral methods using FFTs are efficient for solving certain types of PDEs.
Cryptography: Used in some algorithms for large number multiplication (e.g., Schönhage–Strassen algorithm) which can be relevant in cryptographic computations.
Radar and Sonar Signal Processing: Detecting objects, measuring their speed and distance by analyzing the frequency spectrum of reflected signals.
Medical Imaging: Reconstructing images in MRI (Magnetic Resonance Imaging) and CT (Computed Tomography) scans.

Seismic Data Analysis: Analyzing earthquake waves and subsurface structures in geophysics.

Telecommunications: OFDM (Orthogonal Frequency-Division Multiplexing) used in Wi-Fi, LTE, 5G, and digital broadcasting relies heavily on FFT/IFFT for modulating and demodulating signals.

Key Takeaways

The Discrete Fourier Transform (DFT) converts a time-domain signal into its frequency components.
The Fast Fourier Transform (FFT) is a highly efficient algorithm for computing the DFT, reducing time complexity from O(N²) to O(N log N) using a divide and conquer strategy.
The core idea involves recursively breaking down the problem into smaller DFTs of even and odd indexed samples and then cleverly combining the results using "twiddle factors."
FFT is a fundamental algorithm with wide-ranging applications in signal processing, image/audio compression, scientific computing, telecommunications, and many other fields.

The FFT is a beautiful example of algorithmic ingenuity that has had a profound impact on technology, enabling us to process and understand signals and data in ways that would be computationally infeasible otherwise.