Chapter 5: Discrete Fourier Transform &

Fast Fourier Transform

 Discrete Fourier Transform (DFT)

 Introduction
 DFT and its Properties
 Circular Convolution
 Linear Convolution Using DFT
 Fast Fourier Transform (FFT)
 Decimation-in-Time
 Decimation-in-Frequency

1
Fourier Transform
• The main use of Fourier transform is to find a signal
representation in the frequency domain.
• Signal can be converted from time to frequency domain
and vice versa.
• In the frequency domain, we can see the frequency content
of the signal (as presented by the Fourier coefficients).

2
Signal (time domain) Spectrum (frequency domain)
2.5 60

2
50
1.5

1
40

0.5 Fourier
30
0

-0.5
Transform 20
-1

-1.5 10

-2
0
-2.5 0 5 10 15 20 25 30 35 40 45 50
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

5 15 30 Hz

s(t) = sin(2π(5t)) + Plot of the magnitude of

sin(2π(15t)) + Fourier coeffients vs
sin(2π(30t)) frequency (Hz)

3
 Four Types of Fourier Analysis

1 1 N −1
ck = T xa (t )e − j 2πkF0t dt 1
Tp p
F0 =
Tp ck =
N
 x(n)e− j (2π / N )kn
FS n =0 DFS N −1
∞
xa (t ) =  ck e j 2πkF0t x ( n) =  ck e j (2π / N )kn
k = −∞ k =0

X(Ω)

Ω
∞ ∞
X a (F ) =   x ( n ) e − jΩ n
− j 2πFt
x (t )e
−∞ a
dt X ( e jΩ ) =
n = −∞
CTFT DTFT
∞ 1
2π X ( e
j 2πFt
xa (t ) =  X a ( F )e dF x(n ) = jΩ
) e jΩ n d Ω
−∞ 2π
 Fourier Series

1
ak = T x ( t ) e − j 2 π kf 0 t dt 1 xN (t ) =  kN=− N ak e j 2πkf0t
T0 p f0 = Hz
T0

Major components of Fourier analysis and synthesis

Fourier transform can be derived from Fourier series by making the
period very large (infinite). The discrete frequency will become very
close and approach a continuous spectrum. 5
Continuous Time Fourier Transform
• Continuous Time Fourier transform is used to analyse a non-
periodic, continuous time signal.
• The frequency spectrum is continuous and non-periodic.
• In the frequency domain, frequency F is a continuous value (e.g.,
we can have frequency coefficients Xa(1.99), F=1.99)

∞
X a (F ) =  xa (t )e − j 2πFt dt
−∞ Synthesis equation
Analysis equation ∞
xa (t ) =  X a ( F )e j 2πFt dF
−∞
Discrete Time Fourier Transform

• Discrete Time Fourier Transform (DTFT) is used to analyse a

non-periodic, discrete time signal.
• The frequency spectrum is continuous and periodic.

X(Ω)

∞
X (e jΩ
)=  x ( n ) e − jΩ n Synthesis equation
n = −∞ 1
x(n ) =
2π 2π X (e jΩ )e jΩn dΩ
Analysis equation
2π f
Ω= 7
fs
Limitations of DTFT
• It turns out that an infinite number of samples x[n] are
required to compute DTFT. This would require
infinite memory.
• DTFT cannot be performed by digital processor
because frequency in DTFT is a continuous value.
• Computer only understands and processes discrete
values.
• Discrete Fourier Transform (DFT) is designed to
overcome this problem (by sampling the frequency
signal).

(
X [k ] = X e jΩ k ) Ωk =
2πk k = 0 ,1, K , N − 1
N
8
Relation Between DTFT and DFT

DTFT DFT
∞ N −1 2π
( ) =  x[n]e
X e jΩ − jΩn
X [k ] =  x[n] e
−j
N
kn

n = −∞ n =0

X(Ω) 2π
δΩ =
N

X(kδΩ)

Ω
kδΩ

9
Discrete Fourier Transform
• The DFT of an arbitrary periodic digital signal x[n] (of N samples
in a period) is
N −1 2π
−j kn
X [k ] =  x[n] e N k = 0, 1, ..., N − 1
n =0
Usually the notation is used : WN = e − j 2π N

Thus we can also write X[k] as follows :

N −1
X [k ] =  x [n ] WN
kn
k = 0,1,K, N − 1
n =0
• The inverse DFT is defined as
N −1
1

x[n] is periodic
with period N. x[n ] = X [k ] W − kn
N , n = 0,1,K, N − 1
N k =0
10
10
 N −1
1
x[n ] =
N
 X [k ] W − kn
N , n = 0,1,K, N − 1
k =0

• The frequency spectrum with discrete frequency samples X[k]

(also known as DFT coefficients and obtained via N-point DFT)
corresponds to periodic signal x[n], with period N.
• x[n] = x[n+pN]. The time domain signal (non-periodic) that
corresponds to the DTFT can be recovered by extracting 1 period
of N samples from x[n]. However the length of the non-periodic
signal must contain number of samples ≤ N.
x[n]

L L
−∞ ∞
n
0 1 ……. N-1 11
Computing Discrete Fourier Transform

Example: Consider a sequence x[n] = [1 0 0 1], obtain its

4-point DFT X[k] and plot its magnitude.
N −1
Solution
X [k ] =  x[n] e − j 2 πkn N k = 0,1, K , N − 1 N=4
n =0

for k = 0
3 3
X [0] =  x[ n]e − j 0 =  x[ n] = 1 + 0 + 0 + 1 = 2
n =0 n =0

for k = 1
3
X [1] =  x[n]e − j 2πn / N = 1 + 0 + 0 + 1e − j 2π 3 / 4
n =0

− j 3π / 2 3π 3π
X [1] = 1 + e = 1 + cos( ) − j sin( ) = 1 + j = 2∠45 o
2 2
12
 N −1
X [k ] =  x[n] e − j 2 πkn N k = 0,1,K , N − 1
n =0

for k = 2
3 3
X [ 2] =  x[ n]e − j 2π 2 n / N
=  x[n]e − j 4πn / N
n =0 n =0

X [ 2] = 1 + 0 + 0 + 1e − j 4π 3 / 4 = 1 + 0 + 0 + 1e − j 3π = 1 − 1 = 0

for k = 3
3
X [3] =  x[ n]e − j 2π 3n / N = 1 + 0 + 0 + 1e − j 9π / 2 = 1 − j = 2∠ − 45 o
n =0

X[k] = [2, 1+j, 0, 1-j]

13
 X[k] = [2, 1+j, 0, 1-j]
jw
|X(e )| and |X[k]|
2

1.5

0.5

0
0 0.2 0.4 0.6 0.8 1 (w)
0 1 2 3 (k)

14
MATLAB/Octave Function to Compute DFT
DFT Magnitude
2

1.8

x=[1 0 0 1]; 1.6

1.4

% Do 4-point DFT 1.2

Xk= fft(x,4) 0.8

0.6

0.4

Xk_mag= abs(Xk) ; 0.2

0
0 0 .5 1 1 .5 2 2.5 3

%Adjust frequency to start at 0 F re q u e n c y in d e x k

xvec=0:3; Xk =
stem(xvec,Xk_mag); Columns 1 through 2
xlabel('Frequency index k'); 2.0000 1.0000 +
1.0000i

% Get inverse DFT Columns 3 through 4

x_recon = ifft(Xk); 0 1.0000 - 1.0000i
15
A faster way to compute Discrete Fourier Transform is via the
matrix form.
N −1 2π
−j kn
X [k ] =  x[n] e N k = 0, 1, ..., N − 1
n =0

 DFT in matrix format

16
The complex exponential multiplier WNkn is periodic.

WN = e − j 2π N

WNkn = e − j 2πkn N

For example, for N = 8-point DFT, there are only 8 unique

values of WNkn .

17
Example: Compute the DFT of x(n) = [0 1 2 3].

Solution
The DFT matrix is given by
W40 W40 W40 W40  W40 W40 W40 W40  1 1 1 1 
 0 1 2 3  0 1 2 3  
W4 W4 W4 W4  W4 W4 W4 W4  1 − j − 1 j 
= =
W 0 2 4
W4 W4 W4 6 W 0 2 0
W4 W4 W4 2 1 − 1 1 − 1
 40   4
1  
W4 W43 W46 W49  W40 3 2
W4 W4 W4   1 j − 1 − j 
 X (0)  1 1 1 1   0  6 
 X (1)  1 − j − 1 j  1  − 2 + 2 j 
 =   =  
 X (2) 1 − 1 1 − 1 2  − 2 
      
 X (3)  1 j − 1 − j   3  − 2 − 2 j 
18
Independent Variables of the Frequency Domain
• The horizontal axis of the frequency domain can be
referred to in 3 different ways.
- Frequency index, k (unitless)
- Normalized frequency, Ω (in rad)
2πk
Ωk =
N
- Analog frequency, f (in Hz)
f sk
fk =
N
fs
Frequency resolution (in Hz) is given by .
N
19
Consider an input signal x(n) where
x(n) = 1 for n = 0,1, …7 x(n) = 0 for n > 7
The 32-point DFT is shown below.
9

|X(k)| 5

0
0 5 10 15 20 25 30

Frequency index, k 20
 8

4
|X(k)|
3

0
0 1 2 3 4 5 6 7
Normalized Frequency w 0-2pi
Ω0 Ω1 Normalized frequency, Ω (rad)
21
 Information on sampling frequency is required for this type of plot.
Suppose sampling frequency = 1000Hz
8

5
|X(k)|
4

0
0 100 200 300 400 500 600 700 800 900 1000

Analog frequency, f (Hz)

f0 f1 22
General Properties of DFT

Recall that x[n] obtained via IDFT is periodic with N samples

per period. Circular time shifting will need to be used if the
shifting involves periodic signal.
23
Circular Shifting or Periodic Shifting

• Circular shifting is applicable for x[n] since it is

periodic with period N .
x[n]

L L
−∞ ∞
n
0 1 ……. N-1

• Shifting to the right by one step (delay by one

sample) will cause the sample at time 0 to move to
time 1.
• Shifting within one period will be circular.
24
In the computation of N-point DFT, be aware that x[n] is periodic. N
-point DFT (N = 7) will correspond to the following x[n].
x[n] = [-1 -1 1 2 1 -2 3 -1 -1 1 2 1]
Period, N = 7
n = -3 n=0 n=6 n=7

Let x[n] = [2 1 -2 3 -1 -1 1]. Examples of circular shift of x[n]

are given below:
x[0] = 2
x[n] = [2 1 -2 3 -1 -1 1]
x[6] = 1 x[1] = 1
x[-n7] = [2 1 -1 -1 3 -2 1]
x[5] = -1 x[1- n7] = [1 2 1 -1 -1 3 -2]
x[2] = -2
x[n - 17] = [1 2 1 -2 3 -1 -1]
x[4] = -1 x[3] = 3
25
 If xc[n] is the shifted sequence of x[n] that has been shifted
with m delays via circular shift, it is denoted by
xc[n] = x[<n-m>N] (period = N)
 <n>N = n mod N -> remainder of n/N .
 Example: xc[n] = x[<n-1>5]
xc[0] = x[<0-1>5] = x[4]
x[n] x[<n-1>5] x[<n-3>5]=x[<n+2>5]

0 1 2 3 4 0 1 2 3 4 0 1 2 3 4

Illustration of circular shifting of a finite-length sequence

26
Circular Time Shifting DFT Property

• N-point DFT applied to sequence that has been

circularly shifted is denoted by DFTN{ x[<n-n0>N] }
2π
− j( k ) no
N
DFTN{ x[<n-n0>N] } = X[k] e

• The property states that circular time shift by n0

samples will cause the DFT samples X[k] to

2π
undergo a phase shift of − kno
N

27
Symmetric Properties
• If x[n] is real-valued, then X(k) is conjugate
symmetric.
• Practical signal is real.
• Stated in different ways
X(k) = X*(-k) = X*(N-k)
X(N-k) = X*(k)

• If x[n] is imaginary-valued, X(k) is conjugate

antisymmetric.
• Stated in different ways
X(k) = -X*(-k) = -X*(N-k)
X(N-k) = -X*(k)
28
Example: For x[n] = [2 1 -2 3 -1 -1 1], the 7-point DFT is
partially given as below. Complete the table for X[k], where k = 4, 5
and 6.

X[k]
k Real Imaginary
0 3.0000 0
1 2.1126 -0.7606
2 5.5048 2.6934
3 -2.1174 -4.6815
4
5
6

29
 Solution
X[k]
k Real Imaginary
0 3.0000 0
1 2.1126 -0.7606
2 5.5048 2.6934
3 -2.1174 -4.6815
4 -2.1174 4.6815
5 5.5048 -2.6934
6 2.1126 0.7606

If N is odd, only (N-1)/2 + 1 values are sufficient to fully represent

X[k]. The other DFT coefficients can be obtained from the
symmetric property.
If N is even, (N-2)/2 + 2 values are sufficient to fully represent
X[k]. 30
Circular Convolution
 Linear convolution between two N-length sequences, g[n] and h[n]
results in a N+N-1 length sequence, yL[n]:
N −1
y L [n] =  g [m]h[n − m], n = 0,1, K ,2 N − 2
m =0

 Circular convolution needs to be used if periodic signal is

involved. The following equation states the circular convolution
between sequences g[n] and h[n], both periodic with period = N

N −1
yC [n] =  g [m]h[ n − m N
] = g[n] N h[n]
m =0

 Like linear convolution, circular convolution is commutative :

g[n] N h[n] = h[n] N g[n]
31
 Comparison of Circular and Linear Convolution
Example : Perform linear and circular convolution between g[n] =
[1 2 0 1] and h[n] = [2 2 1 1].
Linear convolution yL[n]:
[1 2 0 1] [1 2 0 1] [1 2 0 1] [1 2 0 1]
[1 1 2 2] [1 1 2 2] [1 1 2 2] [1 1 2 2]
[2 0 0 0] 2 [2 4 0 0] 6 [1 4 0 0] 5 [1 2 0 2] 5

[1 2 0 1] [1 2 0 1] [1 2 0 1]
[1 1 2 2] [1 1 2 2] [1 1 2 2]
[0 2 0 2] 4 [0 0 0 1] 1 [0 0 0 1] 1

Circular convolution yC[n]:

[1 2 0 1] [1 2 0 1] [1 2 0 1] [1 2 0 1]
h[-m4] [ 2 1 1 2 ] h[1-m4] [ 2 2 1 1 ] h[2-m4] [ 1 2 2 1 ] h[3-m4] [ 1 1 2 2 ]
[2 2 0 2] 6 [2 4 0 1] 7 [1 4 0 1] 6 [1 2 0 2] 5
32
The results of linear and circular convolution between
g[n] = [1 2 0 1] and h[n] = [2 2 1 1] are not the same.

33
Making Linear Convolution Equal to Circular
Convolution
• If the length of g[n] is N1 and the length of h[n] is N2 then
linear convolution between g[n] and h[n] will produce a
sequence of length N1 + N2 - 1.

y L [n] = g[n] * h[n] n = 0,1, K , N1 + N 2 − 2

• Circular convolution can be made the same as linear

convolution in one period. The method is to extend the
length of g[n] and h[n] to M = N1 + N2 - 1.
• The length of g[n] and h[n] is extended by zero padding
until length = M = N1 + N2 - 1.
 g [n], 0 ≤ n ≤ N1 − 1 h[n], 0 ≤ n ≤ N 2 − 1
ge[n] =  he [n] = 
0, N1 ≤ n ≤ M − 1 0, N2 ≤ n ≤ M −1

34
 Example : Extend the sequences g[n] and h[n] to length 7 as follows

 g [n], 0 ≤ n ≤ 3 h[n], 0 ≤ n ≤ 3
g e [n] =  he [n] = 
0, 4≤n≤6 0, 4≤n≤6

[1 2 0 1 0 0 0] [1 2 0 1 0 0 0] [1 2 0 1 0 0 0]
he[-m7] [ 2 0 0 0 1 1 2 ] he[1-m7] [ 2 2 0 0 0 1 1 ] he[2-m7] [ 1 2 2 0 0 0 1 ]

[2 0 0 0 0 0 0] 2 [2 4 0 0 0 0 0] 6 [1 4 0 0 0 0 0] 5

[1 2 0 1 0 0 0] [1 2 0 1 0 0 0]
he[3-m7] [ 1 1 2 2 0 0 0 ] he[4-m7] [ 0 1 1 2 2 0 0 ]

[1 2 0 2 0 0 0] 5 [0 2 0 2 0 0 0] 4

[1 2 0 1 0 0 0] [1 2 0 1 0 0 0]
he [5-m7] [ 0 0 1 1 2 2 0 ] he[6-m7] [ 0 0 0 1 1 2 2 ]

[0 0 0 1 0 0 0] 1 [0 0 0 1 0 0 0] 1

The result is the same as linear convolution! 35

• Efficient algorithms are available for computing the DFT of a
finite-duration sequence. They are known collectively as fast
Fourier transform (FFT) algorithms.
• Because these algorithms are available, it is computationally
efficient to implement a convolution of two sequences by the
following procedure:
a) Let g[n] and h[n] be finite-length sequences of length N1
and N2. Denote M = N1 + N2 - 1. Extend the sequences by
zero padding:
 g [n], 0 ≤ n ≤ N1 − 1 h[n], 0 ≤ n ≤ N 2 − 1
ge[n] =  he [n] = 
0, N1 ≤ n ≤ M − 1 0, N2 ≤ n ≤ M −1

b) Compute the M-point DFTs Ge[k] and He[k] of the two

sequences ge[n] and he[n].

36
c) Compute the product Y[k] = Ge[k] He[k] for 0 ≤ k ≤ M - 1.
d) Obtain y[n] by computing the inverse DFT of Y[k].

g[n] Zero- ge[n] (N1+ N2-1)- Ge [k]

padding point DFT
M-point yL[n] =
h[n] he[n] IDFT g[n]∗h[n]
Zero- (N1+ N2-1)-
padding point DFT He[k]

DFT based implementation of the linear convolution

of two finite-length sequences

37
Example: Linear Convolution Using DFT
Given
g[n] = [1 2 0 1] (N1 = 4) and h[n] = [2 2 1 1] (N2 = 4)
Find the linear convolution between g[n] and h[n] through the
DFT.
Solution

a) Apply zero padding

ge[n] = [1 2 0 1 0 0 0]
he[n] = [2 2 1 1 0 0 0]

38
b) Compute Ge[k] and
Ge[k] = [4.0 1.3 – 2i 1.2 - 1.2i -1.0 - 1.8i
-1.0 + 1.8i 1.2 + 1.2i 1.3 + 2i ]
He[k] = [6.0 2.1- 3i 1.3 - 0.7i 0.6 - 1.0i
0.6 + 1i 1.3 + 0.7i 2.1 + 3i ]

c)
Y[k] = Ge[k] He[k]
[24.0000 -3.0794 - 8.2428i 0.6479 - 2.3573i -2.5685 - 0.0170i
-2.5685 + 0.0170i 0.6479 + 2.3573i -3.0794 + 8.2428i]

d) Inverse DFT

yL[n] = [2 6 5 5 4 1 1]
39
Relation Between DTFT and DFT
DTFT DFT
∞ N −1 2π
Xe( )=
jΩ
 x[n]e − jΩn
X [k ] =  x[n] e
−j
N
kn

n = −∞ n =0

• From the above relations, it is obvious that the DFT is obtained by

sampling the DTFT X(ejΩ) on the frequency Ω-axis between
0 ≤ Ω ≤ 2π at Ωk = 2πk/N, k = 0, 1, …, N – 1.

(
X [k ] = X e jΩ k ) Ωk =
2 πk k = 0 ,1, K , N − 1
N

• Given a length-N, x[n] sequence, it is possible to determine a

‘denser’ version of the length-N DFT of the sequence by zero-
padding which approaches X(ejΩ).

40
• Let X(ejΩ) be the DTFT of a length-N sequence x[n]. We wish to
evaluate X(ejΩ) at Ωk = 2πk/M, k = 0, 1, …, M-1 where M >> N.
N −1 N −1
X (Ω k ) =  x[n]e − jΩ k n
=  x[n ]e − j 2πkn M

n =0 n =0

Define a new sequence :

 x[n], 0 ≤ n ≤ N − 1
xe [n] = 
0, N ≤ n ≤ M −1
We thus obtain : M −1
X (Ω k ) =  e
x [n ]e − j 2πkn M

n =0

which is an M-point DFT of xe[n] i.e. Xe[k].

• The DFT Xe[k] can be computed very efficiently using the FFT
algorithm if M is an integer power of 2.
41
The graph below shows the 16-point DFT and 512-point DFT
of 16-length x[n], where x[n] = cos(2πn3/16).

9
X[k]
8

7
X[Ω]
Magnitude

6
M agnitude

4
DFT

0
0 0.2 0.4 0.6 0.8 1
Norm alis ed angular frequency

42
• As we have seen previously, the DTFT is usually approximated by
zero-padding a sequence x[n] before applying DFT.
• Note that zero-padding gives us a high-density spectrum and
provides a better displayed version for plotting. But it DOES NOT
give us a high-resolution spectrum because no new information is
added to the signal; only additional zeros are added in the data.
• To get a high-resolution spectrum, one has to obtain more data
from the experiment or observations.

Example : To illustrate the difference between the high-density

spectrum and the high-resolution spectrum, consider the sequence:
x[n] = cos(0.48πn) + cos(0.52πn)
The following slides illustrate the approximation of the DTFT using
DFT for different observations of x[n].
43
 signal x(n), 0 <= n <= 9

2
N = R = 10
1

x(n) 0

-1

-2

0 1 2 3 4 5 6 7 8 9 10
n
Samples of DTFT Magnitude
10
DFT Magnitude

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
frequency in pi units

R = Length of DFT; N = Length of original data 44

 signal x(n), 0 <= n <= 9 + 40 zeros

1
N = 10, R = 50
x(n) 0

-1

-2

0 5 10 15 20 25 30 35 40 45 50
n
DTFT Magnitude
10
DFT Magnitude

8
Cannot detect the 2
6 frequencies in x[n]
4

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
frequency in pi units

R = Length of DFT; N = Length of original data 45

 Zero padding does not improve frequency resolution

signal x(n), 0 <= n <= 9 + 90 zeros

x(n) 1 N = 10, R = 100

0

-1

-2

0 10 20 30 40 50 60 70 80 90 100
n
DTFT Magnitude
10
Still cannot detect the
DFT Magnitude

8
2 frequencies in x[n]
6

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
frequency in pi units

R = Length of DFT; N = Length of original data 46

 signal x(n), 0 <= n <= 99

1
x(n) 0

-1

-2

0 10 20 30 40 50 60 70 80 90 100
n
DTFT Magnitude
60
DFT Magnitude

N = R = 100
40
Resolution is
sufficient to detect the
20 2 frequencies in x[n]

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
frequency in pi units

R = Length of DFT; N = Length of original data 47

 signal x(n), 0 <= n <= 99 + 100 zeros

1
N = 100, R = 200
x(n) 0

-1

-2

0 20 40 60 80 100 120 140 160 180 200

n
DTFT Magnitude
60
DFT Magnitude

40
Higher density plot

20

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
frequency in pi units

R = Length of DFT; N = Length of original data 48

 signal x(n), 0 <= n <= 99 + 400 zeros

2
N = 100, R = 500
1
x(n) 0

-1

-2

0 50 100 150 200 250 300 350 400 450 500

n
DTFT Magnitude
60
DFT Magnitude

40

20

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
frequency in pi units

R = Length of DFT; N = Length of original data 49

 Spectrum Estimation from Truncated DT Samples

Ideally, the spectrum is computed by using as many samples of signal x[n] as possible
but due to memory limitation, the samples need to be truncated. This is done by
multiplying x[n] with the window function w[n].
This produces a spectrum that is distorted due to convolution of the sinc function
with the impulse function located at the frequencies present in the signal x[n].

Distorted
spectrum, DTFT
of windowed
Ideal signal. DFT is
spectrum obtained by
sampling this
DTFT.
50
 Fast Fourier Transform (FFT)

 DFT plays an important role in analysis, design and

implementation of discrete-time signal processing algorithms and
systems.
 One major problem with the DFT is its extremely high
computational requirement. Thus we need an efficient algorithm
to compute DFT.
 One efficient and fast algorithm used to compute Discrete Fourier
Transform is known as fast Fourier transform (FFT) algorithm.
 We will see later that to achieve the highest efficiency the FFT
algorithms must compute all N values of the DFT.

51
 There are many ways to measure the complexity of an algorithm.
 Computational complexity can be measured by the number of
arithmetic multiplications and additions.
 In terms of multiplications and additions, FFT algorithms can be
of orders of magnitude more efficient than competing algorithms.
 The efficiency of FFT is so high that in many cases the most
efficient way of implementing convolution is by transforming the
sequences involved, multiplying the transform results and inverse-
transforming the product back to the time-domain.
.

52
 N −1
 DFT of length N X [k ] =  x[n]WNkn , k = 0,1,K , N − 1
n =0

 N complex multiplications and (N-1) complex additions for each

value of X[k].
 Thus a total of N2 complex multiplications and N(N-1) ≈ N2
complex additions are required.
 Similar complexity for inverse DFT:
N −1
1
x[n] =
N
 X [k ]W − kn
N , n = 0,1, K , N − 1
k =0

Number of complex additions = N2

Number of complex multiplications = N2
 The number of arithmetic operations required becomes very large for
large values of N, e.g., if N = 1000, then N2 = 106.

53
 Most approaches to improving the efficiency of DFT computation
exploit the symmetry and periodicity properties of WNkn :

WNk + N / 2 = WNkWNN / 2 = −WNk (symmetry property)

WNkn = WNk (n + N ) = WN(k + N )n (periodici ty in n and k )

WN = e − j 2π N

 By using these two properties, it is possible to reduce the number of

operations using algorithms known as the Fast Fourier Transform (FFT).

Number of complex additions = Nlog2N

Number of complex multiplications = Nlog2N
 FFT algorithms are based on the fundamental principle of decomposing
the computation of DFT of a length-N sequence into successively
smaller DFTs.
54
 The radix of the FFT is determined by the size of the smallest
transform. In a radix-2 FFT, the N must be a power of 2 and the
smallest transform is the 2-point DFT.
 Initially, an N-point DFT can be decomposed into two (N/2)-point DFTs.
 Each (N/2)-point DFT can be further decomposed into two (N/4)-point
DFTs, and so on, until a set of 2-point DFTs is obtained.

1
2 - point DFT : X [k ] =  x [n ] e − j 2πkn 2
k = 0,1
n =0

x(0)
x(1)

 Two basic classes of FFT algorithms :

 Decimation-in-time
 Decimation-in-frequency
55
FFT Algorithm : Decimation-in-Time
 In DFT computation we can exploit the symmetry and periodicity
properties of WNkn = e-j(2π/N)kn by decomposing the sequence x[n]
into successively smaller subsequences and computing their DFTs:
decimation-in-time algorithms.
 Consider the case of N = 2v. We can compute X[k] by separating
x[n] into two (N/2)-point sequences consisting of the even-
numbered points in x[n] and odd-numbered points in x[n] :

X [k ] =  x[n ]
W nk
N +  x[n ]
W nk
N
n even n odd

substituting n = 2r for n even and n = 2r+1 for n odd :

( N 2 )−1 ( N 2 )−1
X [k ] =  x[2r ]( WN )
2 rk
+ WNk  x[2r + 1]( )
2 rk
WN
r =0 r =0

56
 WN2 = e −2 j (2π N)
= e − j 2π ( N 2 ) = WN 2

( N 2 )−1 ( N 2 )−1
X [k ] =  x[2r ]WNrk 2 + WNk  x[2r + 1]WNrk 2
=0 42443
1r 4 =0 442444
1r 4 3
G [k ] H [k ]

= G[k ] + WNk H [k ], k = 0,1, K, N − 1

Each of the sums above is recognized as an (N/2)-point DFT,

G[k] being the (N/2)-point DFT of the even-numbered points of
the original sequence and H[k] being the (N/2)-point DFT of the
odd-numbered points of the original sequence.
 Although k = 0, 1, …, N-1, each of the sums must be computed
only for k between 0 and (N/2)-1. Note that G[k] and H[k] are
each periodic in k with period N/2.
57
Example: Decomposing 8-point DFT
X [k ] = G[k ] + WNk H [k ], k = 0,1,K , N − 1
= G[k ] + W8k H [k ], k = 0,1,...7

k = 0, X [0] = G[0] + W80 H [0]

k = 1, X [1] = G[1] + W81H [1]
G[k+4]=G[k]
k = 2, X [2] = G[2] + W82 H [2]
k = 3, X [3] = G[3] + W83 H [3]
H[k+4]=H[k]

k = 4, X [4] = G[4] + W84 H [4]

= G[0] + W84 H [0]
k = 5, X [5] = G[5] + W85 H [5]
= G[1] + W85 H [1]
Try to find X[6] and X[7].
58
 G[0]
x[0] X[0]
G[1] WN0
Even x[2] (N/2)-point X[1]
G[2] WN1
Samples x[4] DFT X[2]
G[3] WN2
x[6] X[3]
WN3

H[0]
x[1] WN4
X[4]
H[1]
x[3] (N/2)-point WN5
X[5]
Odd H[2]
samples
x[5] DFT X[6]
WN6
H[3]
x[7] WN7
X[7]

Flow graph of the decimation-in-time decomposition of an N-point

DFT computation into two (N/2)-point DFT computations (N=8)
The branches entering a node are summed to produce the node variable. When
no coefficient is indicated, the branch transmittance is assumed to be unity. For
other branches, the transmittance of a branch is an integer power of WN.
59
 X [k ] = G[k ] + WNk H [k ], k = 0,1, K , N − 1

 The two (N/2)-point DFTs require approximately 2(N/2)2 complex

multiplications and additions if direct DFT computation is used.
The DFTs will then have to be combined requiring : N complex
multiplications (multiplying the second sum by WNk) and N
complex additions (adding the product obtained to the first sum).
 Therefore, we require at most N+2(N/2)2 = N+(N2/2) complex
multiplications and additions. It can be shown that N+(N2/2) is
less than N2 (the number of multiplications and additions for
direct method).
 In the above example, the original N-point DFT computation is
broken into two (N/2)-point direct DFT computations. If N/2 is
even, we can consider computing each of the (N/2)-point DFTs by
breaking each of them into two (N/4)-point DFTs.
60
 The previous equation for G[k] (replace x[r] with g[r]) thus becomes
( N 2 ) −1 ( N 4 ) −1 ( N 4 ) −1
G [k ] = 
r =0
g [ r ]W rk
N 2 = 
l =0
g [ 2 l ]W 2 lk
N 2 + 
l =0
g [ 2 l + 1]W N( 22l +1)k
( N 4 ) −1 ( N 4 ) −1
= 
l =0
g [ 2 l ]W Nlk 4 + W Nk 2 
l =0
g [ 2 l + 1]W Nlk 4

( N 4 )−1 ( N 4 )−1
Similarly, H [k ] =  h[
l =0
2l ]W lk
N 4 + W k
N 2  h[
l =0
2l + 1]W lk
N 4

 Consequently, we can obtain the (N/2)-point DFTs G[k] by

combining the (N/4)-point DFTs of g[2l] and g[2l+1]. The same
approach can be applied to H[k].
x[0] G[0]
(N/4)-point WN/20
x[4] DFT G[1]
WN/2 1
x[2] G[2]
(N/4)-point WN/2 2
x[6] DFT G[3]
WN/23
61
Flow graph of the decimation-in-time decomposition of an N-point DFT
computation into four (N/4)-point DFT computations (N=8)

x[0] (N/4)-point X[0]

WN0 WN0
x[4] DFT X[1]
WN2 WN1
x[2] WN4 X[2]
(N/4)-point WN2
x[6] DFT X[3]
WN6 WN3

x[1] (N/4)-point X[4]

WN0 WN4
x[5] DFT X[5]
WN2 WN5
x[3] X[6]
(N/4)-point WN4 WN6
x[7] DFT WN7 X[7]
WN6

x[3] Replacing the module into the

WN 0 =1

x[7] (N/4)-point DFT module, we get

W21 =WNN/2 =-1 as shown in the next slide.
62
Complete flow graph of the decimation-in-time decomposition of an N-point
DFT computation into four (N/4)-point DFT computations (N=8)

x[0] X[0]
WN0 WN0 WN0
x[4] X[1]
WN4 WN2 WN1
x[2] WN4 X[2]
WN0 WN2
x[6] X[3]
WN4 WN6 WN3

x[1] WN4 X[4]

WN0 WN0
x[5] X[5]
WN4 WN2 WN5
x[3] WN6 X[6]
WN0 WN4
x[7] X[7]
WN4 WN6 WN7

(m-1)st WN r mth (m-1)st mth

stage stage ≡ stage stage
WN(r+N/2) WN r -1

63
Complete flow graph of the decimation-in-time decomposition of an N-point
DFT computation into four (N/4)-point DFT computations (N=8) after
consideration of the following property :

WNN 2 = e − j ( 2π N )N 2
= e − jπ = −1
WNr + N 2 = WNN 2WNr = −WNr

x[0] X[0]
WN0
x[4] X[1]
-1
WN0
x[2] -1
X[2]
WN0 WN2
x[6] X[3]
-1 -1

WN0
x[1] X[4]
-1
WN0 WN1
x[5] X[5]
-1 WN0 WN2 -1
x[3] -1
X[6]
WN0 WN2 -1 WN3
x[7] X[7]
-1 -1 -1
x(0) 2-point
x(4) DFT Combine X(0)
2-point X(1)
x(2) 2-point DFT's X(2)
x(6) DFT Combine
X(3)
x(1) 2-point 4-point X(4)
x(5) DFT Combine DFT's X(5)
2-point
x(3) 2-point X(6)
x(7) DFT DFT's X(7)

Three stages in the computation of an N=8-point DFT

65
 For the general case of N being power of 2, we can apply the FFT
technique by decomposing the N-point DFT into (N/2)-point
DFTs and then into (N/4)-point DFTs and continue until we are
left with 2-point DFTs.
 The computational saving resulting from the decomposition is as
follows:
N-point (N/2)-point (N/4)-point 2-point
DFT DFT DFT DFT
N2 N + 2(N/2)2 2N + 4(N/4)2 N log2N
No. of complex
Since N = 2v, there can only be v = log2N no. of times
multiplications
decomposition can be carried out. Thus we get total of
& additions
N.v = N log2N complex multiplications & additions.
For large N this represents a significant reduction in
computation compared to the N2 operations required for
the DFT. For example if N=1024, the computational
reduction is 102 times.
66
FFT Input/Output Data Index Bit Reversal
 Let’s look into some of the special properties of the FFT that are
important to FFT software developers and FFT hardware
designers.
 The decimation-in-time phrase refers to how we divided the DFT
input samples into odd and even parts in the derivation of FFT.
 This time decimation leads to the scrambled order of the input
data’s index n.
 The pattern of this shuffled order can be understood with the help
of the table shown in next slide.
 The shuffling of the input data is known as bit reversal because
the scrambled order of the input data index can be obtained by
reversing the bits of the binary representation of the normal input
data index order.
67
 Notice the normal index order in the left column of the above
table and the scrambled order in the right column that
corresponds to the final decimated input index order in
decimation-in-time FFT.
68
Example: Compute the DFT of x(n) = [1 2 3 4]
using FFT.

Solution

Normal order Binary bits Reversed bits Reversed order

0 00 00 0
1 01 10 2
2 10 01 1
3 11 11 3

69
 1+3=4 4 X[0]=4+6=10
x[0]=1

1-3=-2 -2 X[1]=-2+2j
x[2]=3

2+4=6 6 X[2]=4-6=-2
x[1]=2

2-4=-2 2j X[3]=-2-2j
x[3]=4

W40 = 1 W41 = e − j 2π 4 = − j

70
FFT Algorithm : Decimation-in-Frequency
 Besides decimation-in-time, we can also divide X[k] into smaller
and smaller subsequences in the same manner as before:
decimation-in-frequency algorithms.
 We can express X[k] as follows:

( N / 2) −1
X [2r ] =  ( x[n] + x[ n + ( N / 2)])WNrn/ 2 , r = 0, 1, ..., ( N / 2) - 1
n =0
( N / 2) −1
X [2r + 1] =  ( x[n] − x[ n + ( N / 2)])WNnWNrn/ 2 , r = 0, 1, ..., ( N / 2) - 1
n =0

 The complexity is the same as the decimation-in-time FFT.

71
72