Transcript Discrete Fourier Transform (DFT)

Efficient Computation of the DFT:
The FFT Algorithm
Hafiz Malik
Dept. of Electrical & Computer Engineering
The University of Michigan-Dearborn
[email protected]
http://www-perosnal-engin.umd.umich.edu/~hafiz
7/17/2015
DFT - Hafiz Malik
1
Acknowledgement
• These lecture slides are based on
– Chapter 5 of
• Digital Signal Processing: Principles, Algorithms and
Applications 3rd e
• John G. Proakis and Dimitris K Manolakis
– Lecture slides of Prof. Richard M. Stern
• Department of Electrical and Computer Engineering,
Carnegie Mellon University
7/17/2015
DFT - Hafiz Malik
2
The Cooley-Tukey decimation-in-time algorithm
• Consider the DFT algorithm for an integer power of 2, N  2
X[k] 
N1

x[n]WN nk 
n0
N1

x[n]e j2nk / N ; WN  e j2 / N
n0
• Create separate sums for even and odd values of n:
 x[n]WN nk   x[n]WN nk
X[k] 
n even
n odd
• Letting n  2r for n even and n  2r 1 for n odd, we obtain
X[k] 
 N / 21

r 0
7/17/2015
 N / 2 1
x[2r]WN 2rk 

r 0
DFT - Hafiz Malik
x[2r  1]WN  2r1 k
3
The Cooley-Tukey DIT algorithm
• Splitting indices in time, we have obtained
X[k] 
 N / 21

r 0
 N / 2 1
x[2r]WN 2rk 

r 0
x[2r  1]WN  2r1 k
• But W 2  e  j 2 2 / N  e  j2 /( N / 2)  W
and W 2rk W k  W k W rk
N
N/2
N
N
N N/2
So …
X[k] 
(N/ 2)1
(N/ 2)1
n0
n0

x[2r]WNrk/ 2  WNk
N/2-point DFT of x[2r]
7/17/2015

x[2r 1]WNrk/ 2
N/2-point DFT of x[2r+1]
DFT - Hafiz Malik
4
Computational Gain so far …
• We have split the DFT computation into two halves:
X[k] 

N 1

x[n]WN nk
k 0
( N/ 2)1
( N/ 2)1
n0
n0

x[2r]WNrk/ 2  WNk

x[2r 1]WNrk/ 2
• Have we gained anything? Consider the nominal number of
multiplications for N  8
– Original form produces 82  64 multiplications
– New form produces 2(42 )  8  40 multiplications
– So we’re already ahead ….. Let’s keep going!!
7/17/2015
DFT - Hafiz Malik
5
Signal flowgraph notation
• In generalizing this formulation, it is most convenient to adopt a graphic
approach …
• Signal flowgraph notation describes the three basic DSP operations:
– Addition
x[n]
– Multiplication by a constant
y[n]
x[n]+y[n]
– Delay
a
x[n]
x[n]
7/17/2015
DFT - Hafiz Malik
ax[n]
z-1
x[n-1]
6
Signal flowgraph representation of 8-point DFT
k
• Recall that the DFT is now of the form X[k]  G[k]  WN H[k]
• The DFT in (partial) flowgraph notation:
7/17/2015
DFT - Hafiz Malik
7
Continuing with the decomposition …
• So why not break up into additional DFTs? Let’s take the
upper 4-point DFT and break it up into two 2-point DFTs:
7/17/2015
DFT - Hafiz Malik
8
The complete decomposition into 2-point DFTs
7/17/2015
DFT - Hafiz Malik
9
Now let’s take a closer look at the 2-point DFT
• The expression for the 2-point DFT is:
X[k] 
1

n0
x[n]W2nk 
1

x[n]e  j 2nk / 2
n0
• Evaluating for k  0,1 we obtain
X[0]  x[0]  x[1]
X[1]  x[0]  e  j 21 / 2 x[1]  x[0]  x[1]
which in signal flowgraph notation looks like ...
This topology is referred to as the
basic butterfly
7/17/2015
DFT - Hafiz Malik
10
The complete 8-point DIT FFT Structure
7/17/2015
DFT - Hafiz Malik
11
Computational Complexity of N-point FFTs
• Let
 where   log (N)
N

2
2 butterflys/column)(2
• (log2(N) columns)(N/2
mults/butterfly)
or ~N log 2 (N)
multiplys
7/17/2015
DFT - Hafiz Malik
12
Additional timesavers: reducing multiplications
in the basic butterfly
• As we derived it, the basic butterfly is of the form
WNr
WNr N / 2
• Since W N / 2  1 we can reducing computation by 2 by
N
premultiplying byWNr
1
WNr
7/17/2015
1
DFT - Hafiz Malik
13
Alternate FFT structures
• We developed the basic decimation-in-time
(DIT) FFT structure, but other forms are
possible simply by rearranging the branches of
the signal flowgraph
• Consider the rearranged signal flow diagrams
on the following panels …..
7/17/2015
DFT - Hafiz Malik
16
Alternate DIT FFT structures
• DIT structure with input bit-reversed, output natural
7/17/2015
DFT - Hafiz Malik
17
Alternate DIT FFT structures
• DIT structure with input natural, output bit-reversed:
7/17/2015
DFT - Hafiz Malik
18
Alternate DIT FFT structures
• DIT structure with both input and output natural:
7/17/2015
DFT - Hafiz Malik
19
Alternate DIT FFT structures
• DIT structure with same structure for each stage:
7/17/2015
DFT - Hafiz Malik
20
Comments on alternate FFT structures
• A method to avoid bit-reversal in filtering
operations is:
– Compute forward transform using natural input, bitreversed output (as in OSB 9.10)
– Multiply DFT coefficients of input and filter response
(both in bit-reversed order)
– Compute inverse transform of product using bitreversed input and natural output (as in OSB 9/14)
• Latter two topologies (as in OSB 9.15 and 9.16)
are now rarely used
7/17/2015
DFT - Hafiz Malik
21
The decimation-in-frequency (DIF)
FFT algorithm
• Introduction: Decimation in frequency is an
alternate way of developing the FFT algorithm
• It is different from decimation in time in its
development, although it leads to a very
similar structure
7/17/2015
DFT - Hafiz Malik
22
FFT: Decimation in Frequency
• Consider the original DFT equation ….
N 1
X (k )   x(n)WNnk
n 0
• Separate the first half and the second half of time samples:
X (k ) 
( N / 2 ) 1

x ( n)W
n 0

( N / 2 ) 1

x ( n)W
n 0

nk
N
nk
N
N 1


n N / 2
W
 x(n)  (1)
( N / 2 ) 1
n 0
k
x( n)WNnk
( N / 2 ) 1
( N / 2) k
N

n 0
x( n  N / 2)WNnk

x ( n  N / 2) WNnk
• Note that these are not N/2-point DFTs
7/17/2015
DFT - Hafiz Malik
23
FFT: Decimation in Frequency
X (k ) 
 x(n)  (1)
( N / 2 ) 1
k
n 0

x(n  ( N / 2)) WNnk
• For k even, let k  2r
 x(n)  (1)
( N / 2 ) 1
X (k ) 
n 0

2r

x(n  ( N / 2)) WNn 2 r
( N / 2 ) 1
 x(n)  x(n  ( N / 2))W
nr
N /2
n 0
• For k odd, let k  2r 1
X (k ) 
 x(n)  (1)
( N / 2 ) 1
n 0

2r

(1) x(n  ( N / 2)) WNn ( 2 r 1)
( N / 2 ) 1
 x(n)  x(n  ( N / 2))W W
n 0
n
N
nr
N /2
• These expressions are the N/2-point DFTs of
n
x
(
n
)

x
(
n

(
N
/
2
))
and
[
x
(
n
)

x
(
n

(
N
/
2
))]
W
N
7/17/2015
DFT - Hafiz Malik
24
These equations describe the
following structure:
7/17/2015
DFT - Hafiz Malik
25
Continuing by decomposing the odd
and even output points results …
7/17/2015
DFT - Hafiz Malik
26
… and replacing the N/4-point DFTs by
butterflys we obtain
7/17/2015
DFT - Hafiz Malik
27
The DIF FFT is the transpose of the DIT FFT
• To obtain flowgraph transposes:
– Reverse direction of flowgraph arrows
– Interchange input(s) and output(s)
•
DIT butterfly:
DIF butterfly:
• Comment:
– We will revisit transposed forms again in our discussion of
filter implementation
7/17/2015
DFT - Hafiz Malik
28
The DIF FFT is the transpose of the DIT FFT
• Comparing DIT and DIF structures:
DIT FFT structure:
DIF FFT structure:
• Alternate forms for DIF FFTs are similar to those of
DIT FFTs
7/17/2015
DFT - Hafiz Malik
29
Alternate DIF FFT structures
• DIF structure with input natural, output bit-reversed:
7/17/2015
DFT - Hafiz Malik
30
Alternate DIF FFT structures
• DIF structure with input bit-reversed, output natural:
7/17/2015
DFT - Hafiz Malik
31
Alternate DIF FFT structures
• DIF structure with both input and output natural:
7/17/2015
DFT - Hafiz Malik
32
Alternate DIF FFT structures
• DIF structure with same structure for each stage:
7/17/2015
DFT - Hafiz Malik
33
FFT structures for other DFT sizes
• Can we do anything when the DFT size N is not
an integer power of 2 (the non-radix 2 case)?
• Yes! Consider a value of N that is not a power
of 2, but that still is highly factorable …
N let
p1 p2 p3 p4 ... p ; q1  N / p1, q2  N / p1 p2 , etc.
• Let
Then
X[k] 
N 1

n0
q1 1
x[n]WNnk
q1 1
q1 1
p1rk
( p1r1)k
( p r 2)k
  x[ p1r]WN   x[ p1r 1]WN
  x[ p1r  2]WN 1
 ...
r0
r0
r0
7/17/2015
DFT - Hafiz Malik
34
Non-radix 2 FFTs (continued)
• An arbitrary term of the sum on the previous
panel is
q1 1

x[ p1r  l]WN( p1rl)k
r0
q1 1
q1 1
p1rk lk
  x[ p1r  l]WN WN  WNlk  x[ p1r  l]Wqrk
1
r0
r0
• This is, of course, a DFT of size of points spaced
by
q1
p1
7/17/2015
DFT - Hafiz Malik
35
Non-radix 2 FFTs (continued)
• In general, for the first decomposition we use
• Comments:
X[k] 
p1 1

l 0
WNlk
q1 1

r0
x[ p1r  l]Wqrk
1
– This procedure can be repeated for subsequent factors of
N
– The amount of computational savings depends on the
extent to which N is “composite”, able to be factored into
small integers
– Generally the smallest factors possible used, with the
exception of some use of radix-4 and radix-8 FFTs
7/17/2015
DFT - Hafiz Malik
36
Direct DFT & FFT Comparison
• Time (1 multiplication per microsec)
26
29
212
215
218
7/17/2015
N
64
512
4096
32768
262144
Direct DFT
FFT
.02 sec
.002 sec
1
.02 sec
67
.2
1 hr 11 mins 2
3 days 4 hrs 19
DFT - Hafiz Malik
37