siware.dev

# Discrete Fourier Transform

The idea of Fourier is that a function (signal) can be decomposed as infinite number of sine/cosine waves (harmonics). We can consider $sin$ and $cos$ as the basis function and they form an orthogonal basis.

We can lookup on the Internet regarding the topic of Fourier and we will find there are extensive knowledge. Studying the terminology is important to avoid confusion. For starter, the concept of Fourier manifested into many variants:

• Continuous Time Discrete Fourier Series (CT-DFS)
• Discrete Time Discrete Fourier Series (DT-DFS)
• Continuous Time Fourier Transform (CTFT)
• Discrete Time Fourier Transform (DTFT)
• Discrete Fourier Transform (DFT)
• Fast Fourier Transform (FFT)

The differences can be summarized in this table:

Series Input Periodicity Input Resolution Output Periodicity Output Resolution
CT-DFS Infinite Periodic Continuous Aperiodic Discrete
DT-DFS Infinite Periodic Discrete Periodic Discrete
CTFT Infinite Aperiodic Continuous Aperiodic Continuous
DTFT Infinite Aperiodic Discrete Periodic Continuous
DFT Finite Aperiodic Discrete Periodic Discrete

Note that FFT is an implmentation of DFT and practically what we use for computation. For more details on the differences, consult http://complextoreal.com/wp-content/uploads/2013/04/fft5.pdf.

## Fourier Transform Pair

Forward Fourier Transform $\mathcal{F} \{f(x)\}$ is mathematical operation that transforms a function of time (or space) into a function of frequency. Backward (Inverse) Fourier Transform $\mathcal{F}^{-1}\{\hat{f}(\omega)\}$ reconstructs the function of time from function of frequency. Together they form Fourier Transform pair.

Fourier Transform is commonly expressed in terms of complex exponentials (instead of sine/cosine). Recall Euler’s Formula is $e^{ix} = cos(x) + i \; sin(x)$.

$\quad \langle e^{imx}, e^{-inx} dx \rangle = \displaystyle\int_{-\pi}^{\pi} e^{imx} e^{-inx} dx = \begin{cases} 0 & \quad m \neq n \\ 2\pi & \quad m = n \end{cases}$

$\quad \langle e^{i 2\pi mx}, e^{-i 2\pi nx} dx \rangle = \displaystyle\int_{-\pi}^{\pi} e^{i 2\pi mx} e^{-i 2\pi nx} dx = \begin{cases} 0 & \quad m \neq n \\ 1 & \quad m = n \end{cases}$

This makes chosing $e^{i 2\pi \xi x}$ as basis function convenient. Note that $\xi$ is used to denote an ordinary frequency (Hz) whereas $\omega = 2 \pi \xi$ is for angular frequency (rad/s).

Continuous Time Fourier Transform (CTFT) pair is given as:

• Forward Transform: $\hat{f}(\xi) = \displaystyle\int_{-\infty}^\infty f(x) \; e^{-i 2\pi \xi x} dx$
• Backward Transform: $f(x) = \displaystyle\int_{-\infty}^\infty \hat{f}(\xi) \; e^{i 2\pi \xi x} d\xi$

In Discrete Time Fourier Transform (DTFT), the function $f(x)$ is sampled with period $T_s$ (freq $\frac{1}{\xi_s}$) and this can be expressed as:

$\quad f_s(x) = \displaystyle\sum_{n=-\infty}^{\infty} f(nT_s) \; \delta(x - nT_s)$

plugging $f_s$ to CTFT forward transform yields:

$\quad \hat{f}(\xi) = \displaystyle\sum_{n=-\infty}^\infty f(nT_s) \; e^{-i 2\pi \xi nT_s}$

and let the n-th sample $f_n = f(nT_s)$ and $\hat{\xi} = \xi \; T_s = \xi / \xi_s$, this simplifies to:

$\quad \hat{f}(\hat{\xi}) = \displaystyle\sum_{n=-\infty}^{\infty} f_n \; e^{-i 2 \pi \hat{\xi} n}$

Note that the output of CTFT and DTFT are continuous function of frequency. In addition DTFT requires infinite sum which is not suitable for computation. Discrete Fourier Transform (DFT) approximates the result by sampling the output function.

Given $L$ samples of a function $f$ and $\hat{f}$ is discretized by $N$, DFT can be expressed as:

$\quad \hat{f_k} = \displaystyle\sum_{n=0}^{L-1} f_n \; e^{-i \frac{2\pi}{N} k n}$

Note: in most implementation $L = N$.

To summarize, Fourier Transform pairs are:

Forward Transform Backward (Inverse) Transform
CTFT $\hat{f}(\xi) = \displaystyle\int_{-\infty}^\infty f(x) \; e^{-i 2\pi \xi x} dx$ $f(x) = \displaystyle\int_{-\infty}^\infty \hat{f}(\xi) \; e^{i 2\pi \xi x} d\xi$
DFT $\hat{f_k} = \displaystyle\sum_{n=0}^{N-1} f_n \; e^{-i \frac{2\pi}{N} k n}$ $f_n = \frac{1}{N} \displaystyle\sum_{k=0}^{N-1} \hat{f_k} \; e^{i \frac{2\pi}{N} k n}$

## Properties of Fourier Transform

#### Spectral Derivative

The benefit of doing computation in Fourier space is that it is easy to compute derivative (a.k.a Spectral derivative).

$\mathcal{F}\{\frac{d}{dx} f(x)\} = \displaystyle\int_{-\infty}^\infty \frac{df}{dx} e^{-i \omega x} dx$

$= -\displaystyle\int_{-\infty}^\infty f(x) (-i \omega e^{-i \omega x})dx + \Bigl[ f(x) e^{-i \omega x}\Bigr]_{-\infty}^\infty$

$= i \omega \displaystyle\int_{-\infty}^\infty f(x) e^{-i \omega x} dx$

$= i \omega \mathcal{F}\{f(x)\}$

It turns out that computing FFT-based spectral derivative isn’t as straightforward as above. For complete discussion on this topic, please consult https://math.mit.edu/~stevenj/fft-deriv.pdf.

In summary, algorithm to compute first derivative:

1. Given $f_n$ for $0 \le n \lt N$, use an FFT to compute $\hat{f_k}$ for $0 \le k \lt N$
2. Obtain $\hat{f_k}^\prime$ by multiplying $\begin{cases} i \frac{2 \pi}{N} k & \quad k \lt N/2 \\ i \frac{2 \pi}{N} (k - N) & \quad k \gt N/2 \\ 0 & \quad k = N/2 \; (N \; is \; even) \end{cases}$
3. Compute $f_n^\prime$ from $\hat{f_k}^\prime$ via inverse FFT

To compute second derivative:

1. Given $f_n$ for $0 \le n \lt N$, use an FFT to compute $\hat{f_k}$ for $0 \le k \lt N$
2. Obtain $\hat{f_k}^{\prime\prime}$ by multiplying $\begin{cases} -(\frac{2 \pi}{N} k)^2 & \quad k \le N/2 \\ -(\frac{2 \pi}{N} (k-N))^2 & \quad k \gt N/2 \end{cases}$
3. Compute $f_n^{\prime\prime}$ from $\hat{f_k}^{\prime\prime}$ via inverse FFT

#### Convolution

Convolution is defined as $(f * g) = \displaystyle\int_{-\infty}^{\infty} f(x - \xi) g(\xi) d\xi$

Fourier Transform property for convolution is that:

$\mathcal{F}\{(f * g)\} = \mathcal{F}\{f\} \mathcal{F}\{g\} = \hat{f} \hat{g}$

This basically means convolution equals to multiplication in Fourier space. This is another benefit of Fourier Transform.

#### Parseval’s Theorem

The energy in $f$ is proportional with energy in its fourier transform $\hat{f}$

$\displaystyle\int_{-\infty}^{\infty} \lvert \hat{f(\omega)} \rvert^2 d\omega = 2 \pi \displaystyle\int_{-\infty}^{\infty} \lvert f(x) \rvert^2 dx$

This theorem is useful because it verifies that if we truncate the small Fourier coefficients, we are still capturing most of function $f$.

#### Linear Operator

It’s useful to know that Fourier Transform is a linear operator:

$\mathcal{F}\{\alpha f(x) + \beta g(x)\} = \alpha \mathcal{F}\{f\} + \beta \mathcal{F}\{g\}$