Berkeley Notes
  • Introduction
  • EE120
    • Introduction to Signals and Systems
    • The Fourier Series
    • The Fourier Transform
    • Generalized transforms
    • Linear Time-Invariant Systems
    • Feedback Control
    • Sampling
    • Appendix
  • EE123
    • The DFT
    • Spectral Analysis
    • Sampling
    • Filtering
  • EECS126
    • Introduction to Probability
    • Random Variables and their Distributions
    • Concentration
    • Information Theory
    • Random Processes
    • Random Graphs
    • Statistical Inference
    • Estimation
  • EECS127
    • Linear Algebra
    • Fundamentals of Optimization
    • Linear Algebraic Optimization
    • Convex Optimization
    • Duality
  • EE128
    • Introduction to Control
    • Modeling Systems
    • System Performance
    • Design Tools
    • Cascade Compensation
    • State-Space Control
    • Digital Control Systems
    • Cayley-Hamilton
  • EECS225A
    • Hilbert Space Theory
    • Linear Estimation
    • Discrete Time Random Processes
    • Filtering
  • EE222
    • Real Analysis
    • Differential Geometry
    • Nonlinear System Dynamics
    • Stability of Nonlinear Systems
    • Nonlinear Feedback Control
Powered by GitBook
On this page
  • Convolution and the DFT
  • Circular Convolution
  • Linear Convolution with the DFT
  • Block Convolutions
  • FFT
  • Decimation in Time
  • Decimation in Frequency

Was this helpful?

  1. EE123

The DFT

Whereas the CTFT takes a continuous signal and outputs a continuous frequency spectrum and the DTFT takes a discrete signal and outputs a continuous, periodic frequecy spectrum, the Discrete Fourier Transform takes a discrete finite signal and outputs a discrete frequency spectrum. This is useful for signal processing because we cannot store infinite signals in a computer’s memory.

Definition 1

For a length NNN finite sequence {x[n]}0n−1\{x[n]\}^{n-1}_{0}{x[n]}0n−1​, the Discrete Fourier Transform of the signal is a length N finite sequence {X[k]}0n−1\{X[k]\}^{n-1}_{0}{X[k]}0n−1​ where

X[k]=∑n=0N−1x[n]e−j2πNknX[k] = \sum_{n=0}^{N-1}{x[n]e^{-j\frac{2\pi}{N}kn}}X[k]=∑n=0N−1​x[n]e−jN2π​kn

One way to interpret the DFT is in terms of the Fourier series for a disrete periodic signal x~[n]=x[((n))N]\tilde{x}[n]=x[((n))_N]x~[n]=x[((n))N​] where the ((n))N=nmod  N((n))_N = n \mod N((n))N​=nmodN. Recall that the coefficient of the kth term of the Fourier Series is

ak=1N∑n=0N−1x[n]e−j2πNkna_k = \frac{1}{N}\sum_{n=0}^{N-1}{x[n]e^{-j\frac{2\pi}{N}kn}}ak​=N1​∑n=0N−1​x[n]e−jN2π​kn

Notice that the aka_kak​ of the Fourier Series are the DFT values except scaled by a factor of NNN. In other words, if we extend a finite signal periodically, then the DFT and the DTFS are the same up to a constant scale factor. This gives an intuitive inverse DFT.

Definition 2

For a length N finite sequence {X[k]}0N−1\{X[k]\}^{N-1}_{0}{X[k]}0N−1​ representing the DFT of a finite perioidc signal {x[n]}0N−1\{x[n]\}^{N-1}_{0}{x[n]}0N−1​, the inverse DFT is given by

x[n]=1N∑k=0N−1X[k]ej2πNknx[n] = \frac{1}{N}\sum_{k=0}^{N-1}{X[k]e^{j\frac{2\pi}{N}kn}}x[n]=N1​∑k=0N−1​X[k]ejN2π​kn

Notice that the DFT and the IDFT are very similar in form. It turns out that the IDFT can be expressed as a DFT of X∗[k]X^*[k]X∗[k]. Namely

IDFT{X[k]}=1NDFT{X⋆[k]}⋆IDFT\{X[k]\} = \frac{1}{N}DFT\{X^\star[k]\}^\starIDFT{X[k]}=N1​DFT{X⋆[k]}⋆

Further intuition for the DFT comes by relating it to the DTFT. Suppose we have a finite signal x[n]x[n]x[n] which is 000 for n<0n < 0n<0 and n>N−1n > N-1n>N−1. The DTFT of this signal is

X(ω)=∑n=−∞∞x[n]e−jωn=∑n=0N−1x[n]e−jωnX(\omega) = \sum_{n=-\infty}^{\infty}{x[n]e^{-j\omega n}} = \sum_{n=0}^{N-1}{x[n]e^{-j\omega n}}X(ω)=∑n=−∞∞​x[n]e−jωn=∑n=0N−1​x[n]e−jωn

Suppose we sample the DTFT at intervals of 2πNk\frac{2\pi}{N}kN2π​k, then the kth sample is given by

X[k]=X(2πNk)=∑n=0N−1x[n]e−j2πNknX[k] = X\left(\frac{2\pi}{N}k\right) = \sum_{n=0}^{N-1}{x[n]e^{-j\frac{2\pi}{N}k n}}X[k]=X(N2π​k)=∑n=0N−1​x[n]e−jN2π​kn

Thus we can think of the DFT as a NNN evenly spaced samples of the DTFT. One important point to notice is that while the DTFT is often centered around 0 (meaning it is plotted over a range from −π-\pi−π to π\piπ), because we are summing from 0 to N-1 in the DFT, the DFT coefficients are centered around π\piπ, and thus they are plotted on a range fo [0,2π−2πN][0, 2\pi - \frac{2\pi}{N}][0,2π−N2π​]

Convolution and the DFT

Circular Convolution

When the DFT coefficients of two signals are multiplied, the resulting coefficients describe a circular convolution of the original two signals.

x[n]⊛y[n]↔X[k]Y[k]x[n]\circledast y[n] \leftrightarrow X[k]Y[k]x[n]⊛y[n]↔X[k]Y[k]

Definition 3

A circular convolution between two finite sequences is given by

x[n]⊛y[n]=∑m=0N−1x[m]y[((n−m))N]x[n]\circledast y[n] = \sum_{m=0}^{N-1}x[m]y[((n-m))_N]x[n]⊛y[n]=∑m=0N−1​x[m]y[((n−m))N​]

The mechanics of the circular convolution are the same as that of the regular convolution, except the signal is circularly shifted as shown in Figure 1.

A circular convolution is equivalent to a periodic convolution over a single period.

Linear Convolution with the DFT

Because multiplying DFT coefficients performs a specific case of convolution, we can compute a linear convolution using the circular convolution. Suppose we have two finite signals {x[n]}0L−1\{x[n]\}_0^{L-1}{x[n]}0L−1​ and {h[n]}0P−1\{h[n]\}_0^{P-1}{h[n]}0P−1​ The linear convolution of these two signals will be length L+P−1L+P-1L+P−1, so in order to take an IDFT and get L+P−1L+P-1L+P−1 samples, we need to take at least N≤L+P−1N\le L+P-1N≤L+P−1 points.

  1. Pad each vector to length L+P−1L+P-1L+P−1

  2. Compute X[k]H[k]X[k]H[k]X[k]H[k]

  3. Take the Inverse DFT

If NNN is smaller than L−P+1L-P+1L−P+1, the result is akin to aliasing in the time domain. To see why, consider that the DFT coefficients are essentially the DTFS coefficients of the periodic extension of x[n]x[n]x[n] (denote x~[n]\tilde{x}[n]x~[n]).

x~[n]=∑r=−∞∞x[n−rN]\tilde{x}[n]=\sum_{r=-\infty}^{\infty}x[n-rN]x~[n]=∑r=−∞∞​x[n−rN]

If we compute the DTFT of each periodic extension, then

Y(ejω)=X(ejω)H(ejω)Y(e^{j\omega})=X(e^{j\omega})H(e^{j\omega})Y(ejω)=X(ejω)H(ejω)

and the IDTFT of this will be

y~[n]=∑r=−∞∞y[n−rN].\tilde{y}[n] = \sum_{r=-\infty}^{\infty}y[n-rN].y~​[n]=∑r=−∞∞​y[n−rN].

Notice that if NNN is not large enough, then these copies will be overlapping (a.k.a aliasing). Since the DFT is just sampling the DTFT, the circular convolution will represent the true convolution so long as the copies don’t overlap.

Block Convolutions

In a discrete time system, the input signal might have a very long length, making it impractical to be stored in a computer’s memory or to compute the DFT of it all at once (especially if we have a real-time system). Thus to compute the output of a digital filter (with impulse response of length PPP), we need to compute the DFT in blocks shorter than the signal.

The first method of block convolution is the overlap-add method.

  1. Decompose x[n]x[n]x[n] into nonoverlapping segments of length LLL

    x[n]=∑rxr[n]xr[n]={x[n]rL≤n≤(r+1)L0else.x[n] = \sum_{r}x_r[n] \qquad x_r[n] = \begin{cases} x[n] & rL \le n \le (r+1)L\\ 0 & \text{else}. \end{cases}x[n]=∑r​xr​[n]xr​[n]={x[n]0​rL≤n≤(r+1)Lelse.​

  2. Since convolution is linear,

    y[n]=x[n]∗h[n]=∑rxr[n]∗h[n].y[n] = x[n]*h[n]=\sum_r{x_r[n]*h[n]}.y[n]=x[n]∗h[n]=∑r​xr​[n]∗h[n].

  3. Zero pad xr[n]x_r[n]xr​[n] and h[n]h[n]h[n] to length N≥L+P−1N\ge L+P-1N≥L+P−1 to prevent time-domain aliasing

  4. Compute the DFTs, multiply them, and take the inverse DFT.

  5. The neighboring outputs overlap in the last P−1P-1P−1 points, so add the overlapping sections together to get the final output

The other method of block convolution is the overlap-save method.

  1. Divide x[n]x[n]x[n] into sections of length LLL such that each section overlaps the previous by P−1P-1P−1 points

    xr[n]=x[n+r(L−P+1)−P+1]0≤n≤L−1x_r[n]=x[n+r(L-P+1)-P+1] \qquad 0 \le n \le L-1xr​[n]=x[n+r(L−P+1)−P+1]0≤n≤L−1

  2. Zero pad xr[n]x_r[n]xr​[n] and h[n]h[n]h[n] to length N≥L+P−1N\ge L+P-1N≥L+P−1 to prevent time domain aliasing.

  3. Compute the DFTs, multiply the coefficients, and compute the inverse DFT.

  4. The first P−1P-1P−1 samples of the output will be incorrect, so we can discard them.

    y[n]=∑r=0∞yr[n−r(L−P+1)+P−1]yr[n]={xr[n]∗h[n]P−1≤n≤L−10 elsey[n]=\sum_{r=0}^{\infty}y_r[n-r(L-P+1)+P-1] \qquad y_r[n]= \begin{cases} x_r[n]*h[n] & P-1\le n \le L-1\\ 0 & \text{ else} \end{cases}y[n]=∑r=0∞​yr​[n−r(L−P+1)+P−1]yr​[n]={xr​[n]∗h[n]0​P−1≤n≤L−1 else​

FFT

The DFT gives us an easy way to do convolutions. Unfortunately naiively, computing it is an O(N2)O(N^2)O(N2) operation because we must sum together NNN elements to compute NNN different coefficients. Thankfully, there is a fast algorithm which can compute the DFT in O(Nlog⁡N)O(N\log N)O(NlogN) time so we can compute convolutions quickly.

Definition 4

The fast fourier transform is an algorithm which computes the DFT efficiently in O(Nlog⁡N)O(N\log N)O(NlogN)time.

It works by exploiting properties of the Nth roots of unity.

Definition 5

The Nth roots of unity are the complex roots of WNN=1W_N^N=1WNN​=1.

WNk=e−j2πkNW_N^k=e^{-j\frac{2\pi k}{N}}WNk​=e−jN2πk​

The roots of unity have the following properties.

Theorem 1

The Nth roots of unity are conjugate symmetric.

WNN−kn=WN−kn=(WNkn)⋆W_N^{N-kn} = W_N^{-kn} = (W_N^{kn})^\starWNN−kn​=WN−kn​=(WNkn​)⋆

Theorem 2

The Nth roots of unity are periodic in N.

Wkn=WNk(n+N)=WN(k+N)nW_{kn} = W_N^{k(n+N)} = W_N^{(k+N)n}Wkn​=WNk(n+N)​=WN(k+N)n​

Theorem 3

When squared, the Nth roots of unity become the N2\frac{N}{2}2N​th roots of unity.

WN2=WN2W_N^2 = W_\frac{N}{2}WN2​=W2N​​

Using Theorem 1, Theorem 2, Theorem 3, we can take two approaches to the FFT: decimation in time, which splits x[n]x[n]x[n] into smaller subsequences, and decimation in frequency which splits X[n]X[n]X[n] into smaller subsequences.

Decimation in Time

The idea here is too break x[n]x[n]x[n] into smaller subsequences. We assume that NNN is a power of 2 for simplicity.

X[k]=∑n=0N−1x[n]WNkn=∑n evenx[n]WNkn+∑n oddx[n]WNknX[k]=\sum_{n=0}^{N-1}x[n]W_N^{kn} = \sum_{\text{n even}}x[n]W_N^{kn}+\sum_{\text{n odd}}x[n]W_N^{kn}X[k]=∑n=0N−1​x[n]WNkn​=∑n even​x[n]WNkn​+∑n odd​x[n]WNkn​

We let n=2rn=2rn=2r and n=2r+1n=2r+1n=2r+1.

X[k]=∑r=0N2−1x[2r]WN2rk+∑r=0N2−1x[2r+1]WNk(2r+1)=∑r=0N2−1x[2r]WN2rk+WNk∑r=0N2−1x[2r+1]WN2kr\begin{aligned} X[k] &= \sum_{r=0}^{\frac{N}{2}-1}x[2r]W_N^{2rk}+\sum_{r=0}^{\frac{N}{2}-1}x[2r+1]W_N^{k(2r+1)}\\ &= \sum_{r=0}^{\frac{N}{2}-1}x[2r]W_{\frac{N}{2}}^{rk}+W_N^k\sum_{r=0}^{\frac{N}{2}-1}x[2r+1]W_{\frac{N}{2}}^{kr}\\\end{aligned}X[k]​=r=0∑2N​−1​x[2r]WN2rk​+r=0∑2N​−1​x[2r+1]WNk(2r+1)​=r=0∑2N​−1​x[2r]W2N​rk​+WNk​r=0∑2N​−1​x[2r+1]W2N​kr​​

These are just the DFTs of the even and odd elements of the signal!

∴X[k]=G[k]+WNkH[k]\begin{aligned} \therefore X[k] = G[k] + W_N^kH[k]\end{aligned}∴X[k]=G[k]+WNk​H[k]​

Both GGG and HHH are N2\frac{N}{2}2N​ periodic, and notice that

WNk+N2=e−j2πN(k+N2)=−WNk.W_N^{k+\frac{N}{2}}=e^{-j\frac{2\pi}{N}(k+\frac{N}{2})}= -W_N^k.WNk+2N​​=e−jN2π​(k+2N​)=−WNk​.

This means once we compute G[k]G[k]G[k] and H[k]H[k]H[k] we can compute X[k]X[k]X[k] easily because

X[k]=G[k]+WNkH[k]X[k+N2]=G[k]−WNkH[k]for k∈[0,N2)X[k] = G[k]+W_N^kH[k]\qquad X\left[k+\frac{N}{2}\right]=G[k]-W_N^kH[k] \qquad \text{for }k\in\left[0, \frac{N}{2}\right)X[k]=G[k]+WNk​H[k]X[k+2N​]=G[k]−WNk​H[k]for k∈[0,2N​)

We can continue this relationship recursively downwards. Once we get too N=2N=2N=2, we can represet this as a simple butterfly operation:

X[0]=x[0]+x[1]X[1]=x[0]−x[1].X[0] = x[0]+x[1] \qquad X[1] = x[0]-x[1].X[0]=x[0]+x[1]X[1]=x[0]−x[1].

Decimation in Frequency

The decimation in frequency approach is very similar to the decimation in time approach except instead we split the frequency components

X[2r]=∑n=0N2−1x[n]WN2rn+∑n=0N2−1x[n+N2]WN2r(n+N2)=WN2rn∑n=0N2−1(x[n]+x[n+N2])X[2r+1]=WN2rn∑n=0N2−1(x[n]−x[n+N2])\begin{aligned} X[2r] &= \sum_{n=0}^{\frac{N}{2}-1}x[n]W_N^{2rn}+\sum_{n=0}^{\frac{N}{2}-1}x\left[n+\frac{N}{2}\right]W_N^{2r\left(n+\frac{N}{2}\right)}=W_{\frac{N}{2}}^{rn}\sum_{n=0}^{\frac{N}{2}-1}\left(x[n]+x\left[n+\frac{N}{2}\right]\right)\\ X[2r+1] &= W_{\frac{N}{2}}^{rn}\sum_{n=0}^{\frac{N}{2}-1}\left(x[n]-x\left[n+\frac{N}{2}\right]\right)\end{aligned}X[2r]X[2r+1]​=n=0∑2N​−1​x[n]WN2rn​+n=0∑2N​−1​x[n+2N​]WN2r(n+2N​)​=W2N​rn​n=0∑2N​−1​(x[n]+x[n+2N​])=W2N​rn​n=0∑2N​−1​(x[n]−x[n+2N​])​

PreviousEE123NextSpectral Analysis

Last updated 3 years ago

Was this helpful?

Figure 1: A circular shift