数字信号处理 第四章资料
数字信号处理DSP第4章
k 0,1, , N 1
2
13
4.2 按时间抽取(DIT)的基2–FFT算法
将系数统一为 WNk 2 WN2k ,则可得
x[0]
N 4点
x[4]
DFT
G[0]
X [0]
G[1]
X [1]
x[2]
N 4点
WN0
x[6]
DFT
WN2
G[2]
1 G[3]
1
X [2] X [3]
x[1]
N 4点
X m1[i] WNr X m1[ j] , X m1[i] WNr X m1[ j]
m 1, 2 ,
每一个蝶形需要一次复数乘法和两次复数加法。
17
4.2 按时间抽取(DIT)的基2–FFT算法
N点的DIT-FFT计算量为
复数乘法:
1
N 2
log2
N
N 2
复数加法:
2
N 2
log2
N
N
例: 如果每次复数乘法需要100us,每次复数加法需要20us,来 计算N=1024点DFT,则需要
12
4.2 按时间抽取(DIT)的基2–FFT算法
同理
( N 4)1
( N 4)1
G[k] DFT[g[r]]
g[2l]WN2lk2
g[2l 1]WN(22l1)k
l 0
l 0
( N 4)1
( N 4)1
g[2l]WNlk 4 WNk 2
g[2l 1]WNlk 4 ,
l 0
l 0
k 0,1,
(3) WN0 WN4 WN8 WN12 WN16 WN20 WN24 WN28
或 WN4i i 0,1, 2, 3, 4, 5, 6, 7 (dm 1)
数字信号处理 第四章
4.2利用模拟滤波器理论的设计方法 4.2利用模拟滤波器理论的设计方法: 利用模拟滤波器理论的设计方法:
IIR滤波器的设计是基于模拟滤波器的成熟技术完成 IIR滤波器的设计是基于模拟滤波器的成熟技术完成 滤波器的设计是基于模拟滤波器的成熟技术 最先设计的模拟低通滤波器称为原型滤波器 最先设计的模拟低通滤波器称为原型滤波器 常用的模拟滤波器有: 常用的模拟滤波器有: 巴特沃思(Butterworth)滤波器 巴特沃思(Butterworth)滤波器 切比雪夫(Chebyshev)滤波器 切比雪夫(Chebyshev)滤波器 椭圆(Ellipse) (Ellipse)滤波器 椭圆(Ellipse)滤波器 贝塞尔(Bessel) (Bessel)滤波器 贝塞尔(Bessel)滤波器
8
3、数字滤波器的设计方法 、
IIR滤波器的设计方法 IIR滤波器的设计方法
借助模拟滤波器设计 IIR滤波器 IIR滤波器 脉冲响应不变法 双线性变换法 最优化设计法 直接设计数字IIR滤波器 直接设计数字IIR滤波器 IIR
零极点累试法
Notice: IIR滤波器与 : 滤波器与FIR滤波器的设计方法大不 滤波器与 滤波器的设计方法大不 相同,IIR滤波器设计结果是系统函数 相同, 滤波器设计结果是系统函数H(z), FIR , 滤波器设计结果是系统函数 滤波器设计结果是单位脉冲响应h(n)。 滤波器设计结果是单位脉冲响应 。
20
1、巴特沃思Butterworth低通滤波器 巴特沃思Butterworth低通滤波器
----用巴特沃思函数近似滤波器的频响函数 用
Butterworth低通模平方函数
1 | Ha ( jΩ ) | = 1 + (Ω / Ω c ) 2 N
数字信号处理 第4章 FFT基本思想和2种基本的FFT
= −W
W的对称性
W的可约性
2 rk WN rk = WN / 2
长序列变成短序列 若N → 2个N / 2
2 则N 2次复述乘法 →(N / 2)= N 2 / 2次复数乘法 2
从信号的特殊性上考虑
– 如奇、偶、虚、实性
W 0 X (0) X (1) W 0 = X (2) W 0 0 X (3) W
对 N = 2M , 共可分 M 次,即 m = 0,1,L , M − 1,
8点FFT时间抽取算法信号流图
每一级有 N/2 个如下的“蝶形”单元:
xm ( p )
xm +1 ( p )
W
r N
xm (q)
−1
xm +1 (q )
算法讨论( “级”的概念、碟形单元、 “组” 的概念、旋转因子的分布、码位倒置)
r =2l ,r =2l +1
A(k ), B(k )
C(k) = D(k) =
N / 4−1 l =0
∑x(4l)W
l =0
lk N/4
, k = 0,1,..., N / 4 −1
N / 4−1
lk x(4l + 2)WN / 4 , k = 0,1,..., N / 4 −1 ∑
k A(k) = C(k) +WN / 2 D(k), k = 0,1,..., N / 4 −1 k A(k + N / 4) = C(k) −WN / 2 D(k), k = 0,1,..., N / 4 −1
x(6)
n N
N n = 0,1,L , 2
由此得到基本 运算单元
g (0) g (1) g (2) g (3)
数字信号处理课件第四章资料
5、时间抽取蝶形运算流图符号
X1(k)
X1(k) WNk X 2 (k)
X 2 (k )
WNk
1 X1(k) WNk X 2 (k)
返回DIF 返回例题
设 N 23 8
X1(k)
X 2 (k )
WNk
k 0
W80
1
W81
2
W82
3
W83
X (k)
k 0,1,,7
l0
l 0
X1(k) X 3(k) WNk X 4 (k)
2
X1(k
N 4
)
X 3 (k )
W Nk
2
X
4
(k)
k 0,1,..., N 1 4
x2(r)也进行同样的分解:
x5 (l) x2 (2l)
x6 (l) x2 (2l 1)
l 0,1,..., N 1 4
)
N
/ 21
x1(r)WNrk/ 2
X1(k)
r 0
r 0
X2(k N / 2) X2(k) X (k) X1(k) WNk X 2 (k)
W (kN N
/
2)
WNkWNN
/
2
WNk
N点X(k)可以表示成前 N点和后 点N 两部分:
2
2
前半部分X(k):
X (k) X1(k) WNk X 2 (k)
N 1
X (k) x(n)WNnk k = 0, 1, …, N-1
n0
x(n)
1 N
N 1
X (k )WNnk
k 0
n = 0, 1, …, N-1
二者的差别只在于WN 的指数符号不同,以及差一 个常数因子1/N,所以IDFT与DFT具有相同的运算量。
数字信号处理 第四章
x5(l)WNkl/ 4 x6 (l)WNkl/ 4
= =
DFT[x5(l)] = DFT[x6 (l)] =
DFT [ x2 (2l )] DFT[x2(2l + 1)]
l = 0,1,..., N / 4 − 1 l = 0,1,..., N / 4 − 1
统一系数:
Wk N /2
→ WN2k
第四章 快速傅里叶变换 jingqilu@
N/2仍为整数,进一步分解:N/2→ N/4
将x1
(
r
)按奇偶分解:⎧ ⎨ ⎩
x1 x1
(2l (2l
) = x3(l) + 1) = x4
(l
)
l = 0,1,
, N −1 4
⎧⎪ X1(k ) = X 3(k ) + WNk /2 X 4 (k )
⎨ 则有:⎪⎩ X1(k
+
N 4
)
=
X 3 (k )
X
2
(k
)
= =
N / 2−1
x1(r)WNrk/ 2
r=0
N / 2−1
x2 (r)WNrk/ 2
r =0
= =
DFT [x1(r)] DFT [x2 (r)]
第四章 快速傅里叶变换 jingqilu@
4.2.2 时域抽取法基2FFT基本原理
再利用周期性求X(k)的后半部分
∵ X1 (k ), X 2 (k )是以N / 2为周期的
两个N/2点DFT
N 2/ 2
N (N / 2 –1)
总计
N2/2+N/2≈N2/2 N(N/2-1)+N ≈N2/2
分解前的计算量
精品课件-数字信号处理(第四版)(高西全)-第4章
点DFT和(4.2.10)式或(4.2.11)式所示的N/4个蝶形运算,
如图4.2.3所示。依次类推,经过M次分解,最后将N点DFT
分解成N个1点DFT和M级蝶形运算,而1点DFT就是时域序列
本身。一个完整的8点DIT-FFT运算流图如图4.2.4所示。
图中用到关系式
。W图N中k / m输入W序Nmk列不是顺序排
In Time FFT,简称DIT-FFT ); 频域抽取法FFT (Decimation In Frequency FFT,简称DIF-FFT)。本节介 绍DIT-FFT
设序列x(n)的长度为N,且满足N=2M,M为自然数。按n 的奇偶把x(n)分解为两个N/2点的子序列
x1(r) x(2r), x2 (r) x(2r 1),
x1
(2l
1)WNk
( /
2l 2
1)
l 0
l 0
N / 41
N / 41
x3 (l)WNkl/ 4 WNk / 2
x4
(l
)WNk
l /
4
l 0
l 0
X 3 (k ) WNk/ 2 X 4 (k )
k 0, 1, , N 1 2
(4.2.9)
第4章 快速傅里叶变换(FFT)
式中
N / 41
r0
2
(4.2.6)
由于X1(k)和X2(k)均以N/2为周期,
kN
WN 2
WNk
且
,因此X(k)又可表示为
第4章 快速傅里叶变换(FFT)
X (k) X1(k) WNk X 2 (k),
X
(k
N 2
)
X1(k)
WNk
X
数字信号处理第四章 fft
第四章 快速傅立叶变换(FFT)
基2FFT算法
FFT的基本思想 长为N的序列x(n)的) N
N 1
x ( n )W N
nk
0 k N 1
nk
n0 N 1
X ( k )W N
0 n N 1
x ( n )
n0
N 2
x(
N 2
rn
N 2
1
x1 ( n )W N
2
rn
n0
基2频域抽取FFT(Sande-Tukey算法 ,DIF-FFT)
X ( 2 r 1)
x ( n ) x ( N 2
n0 1
N 2
N 2
1
n ) W N
( 2 r 1) n
nk
将k分成偶和奇数,即将X(k)分解成奇偶两组,在偶数组中k=2r, e 则 ; jk 1 在奇数组中k=2r+1,则 e 1;有:
jk
X (2r )
x ( n )
n0 1
N 2
1
x ( N n ) W N 2 n) N W
2
2 rn
实序列的FFT算法
Z(k)与H(k)、G(k)的关系 一般而言,H(k),G(k)均是复函数,因此,关键是怎样从Z(k)中分 离出H(k)和G(k)。Z(k)可写作:
WN W
N 2
N 2
k N
e
j
2 N N 2
W N W N
k
k
) X 1 (k
) WN
k
(k
N 2
)
X 2 (k
数字信号处理-第四章(加绪论共八章)精品PPT课件
21
4.1.2 小结
• 四种FIR数字滤波器的相位特性只取决于 h(n)的对称性,而与h(n)的值无关。 •幅度特性取决于h(n)。 •设计FIR数字滤波器时,在保证h(n)对称的 条件下,只要完成幅度特性的逼近即可。
注意:当H(ω)用│H(ω)│表示时,当H(ω)为 奇对称时,其相频特性中还应加一个固定 相移π。
22
4.1.3 线性相位FIR滤波器的零点特性
h(n) h(N 1 n)
H z zN1H z1
N 1
N 1
H z hnzn hN 1 nzn
N 1 n0
n0
N 1
h(m)z N 1m z N 1 h m z m
m0
m0
zN 1H z1
1)若 z = zi 是H(z)的零点,则 z = zi-1 也是零点
4.1.1 线性相位的条件 线性相位意味着一个系统的相频特性是 频率的线性函数,即 ()
式中为常数,此时通过这一系统的各频率分
量的时延为一相同的常数,系统的群时延为
g
d () d
5
线性相位FIR滤波器的DTFT为
N1
H e j h n e jn H e j () H ()e j
26
截短并移位的脉冲响应
过渡带带宽=阻带边缘频率-通带边缘频率 设计中用的通带边缘频率=所要求的通带边缘频率+(过渡带带宽/2) 27
20
例1 N=5, h (0) = h (1) = h (3) = h (4) = -1/2, h (2) = 2,求 幅度函数H (ω)。
解 N为奇数并且h(n)满足偶
对称关系 a (0) = h (2) = 2 a (1) = 2 h (3) = -1 a (2) = 2 h (4) = -1 H (ω) = 2 - cosω- cos2ω
数字信号处理第四章..
第四章线性时不变离散时间系统的频域分析一、传输函数和频率响应例4.1传输函数分析Q4.1clear;M = input('Enter the filter length M: ');w = 0:2*pi/1023:2*pi;num = (1/M)*ones(1,M);den = [1];h = freqz(num, den, w);subplot(2,1,1)plot(w/pi,abs(h));gridtitle('Magnitude Spectrum |H(e^{j\omega})|')xlabel('\omega /\pi');ylabel('Amplitude');subplot(2,1,2)plot(w/pi,angle(h));gridtitle('Phase Spectrum arg[H(e^{j\omega})]')xlabel('\omega /\pi');ylabel('Phase in radians');M=2 M=10 M=15幅度谱为偶对称,相位谱为奇对称,这是一个低通滤波器。
M越大,通带越窄且过渡带越陡峭。
Q4.2使用修改后的程序P3.1,计算并画出当w=[0,pi]时传输函数的因果线性时不变离散时间系统的频率响应。
它表示哪种类型的滤波器?w = 0:pi/511:pi;num = [0.15 0 -0.15];den = [1 -0.5 0.7];如下图1这是一个带通滤波器。
图1 图2Q4.3对下面的传输函数重做习题Q4.2:,式(4.36)和式(4.37)给出的两个滤波器之间的区别是什么?你将选择哪一个滤波器来滤波,为什么?w = 0:pi/511:pi;num = [0.15 0 -0.15];den = [0.7 -0.5 1];如上图2也是一个带通滤波器,这两个滤波器的幅度谱是一样的,相位谱不太一样,我会选择第一个带通滤波器,因为它的相位谱更加平滑,相位失真小。
精品课件-数字信号处理-第4章
第4章 快速傅里叶变换(FFT)
可见,1次复数乘法相当于4次实数乘法和2次实数加法,1次 复数加法相当于2次实数加法。因此,计算一个X(k)需要4N次 实数乘法和2N+2(N-1)=4N-2次复数加法,整个DFT运算需要 4N2次实数乘法和N(4N-2)次实数加法。
当N >>1时,N(N-1)≈N2,N(4N-2)≈4N2,因此,不管 以复数还是实数进行统计,直接计算N点DFT的乘法或加法次 数都与N2成正比,随着N的增加,运算量增加越来越快,特别 是N很大时,运算量将非常可观。例如,N=8时,次数为N2=64; N=1024时,次数为N2=1 048 576,超过100万次。对于各种实 时性很强的信号处理应用来说,要求计算速度特别快,因此必 须改进DFT的计算方法,减少运算量。
第4章 快速傅里叶变换(FFT) 第4章 快速傅里叶变换(FFT)
4.1 直接计算DFT的运算量及改进途径 4.2 时间抽取法(DIT)基-2 FFT算法 4.3 频域抽取法(DIF)基-2 FFT算法 4.4 快速傅里叶逆变换(IFFT)算法 4.5 FFT算法的工程实现考虑 习题
第4章 快速傅里叶变换(FFT)
DIF-FFT
第4章 快速傅里叶变换(FFT)
4.2 时间抽取法(DIT)基-2 FFT算法 4.2.1 DIT-FFT算法原理
设序列x(n)长度N=2M,N点DIT-FFT算法对应序列奇偶分解, 令x1(r)、x2(r)
x1(r) x(2r), r 0,1, 2,
, N 1 2
x2 (r) x(2r 1), r 0,1, 2,
4.1 直接计算DFT的运算量及改进途径 4.1.1 直接计算DFT的运算量
设x(n)是长度为N的有限长序列,其N点DFT为
数字信号处理课件第4章
N −1
2 ( = ∑ x(2r )WN rk + ∑ x(2r + 1)WN2 r +1) k r =0
N 2
−1
N 2
−1
r =0
2 k 2 = ∑ x1 (r )(WN ) rk + WN ∑ x2 (r )(WN ) rk r =0 r =0
−1
N 2
−1
根据可约性,W = e
2 N
N 2
X1(k + N ) = X3 (k) −WNk X4 (k), k = 0,1,L, N −1 4 4
2
(二) N/4点DFT
同样对n为奇数时 , 点分为两个N/4点的 同样对 为奇数时, N/2点分为两个 为奇数时 点分为两个 点的 序列: 序列 x5 (l) = x2 (2l), l = 0,1,L, N −1 4
3
k 则有:X ( k ) = X 1 (k ) + WN X 2 (k ) k X (k + 4) = X 1 (k ) − WN X 2 (k ), k = 0,1,2,3
(一) N/2点DFT 一 点
整个过程如下图所示: 整个过程如下图所示
x1(0)=x(0) )= ( x1(1)= (2) )=x( x1(2)= (4) )=x( x1(3)= (6) )=x( x2(0)= (1) )=x( x2(1)= (3) )=x( x2(2)= (5) )=x( x2(3)= (7) )=x( X1(0) N/2点 N/2点 X1(1) X1(2) DFT X1(3) X2(0) X2(1) X2(2) X2(3)
2
N 2 N2 ( ) = 2 4 N N ( − 1) 2 2
精品文档-数字信号处理(吴瑛)-第4章
第4章 快速傅里叶变换(FFT)
由于1点序列的DFT值为其序列本身,因此在最后一次分
解后,流图中已经没有直接计算DFT的环节。第三次分解旋转
因子为
W20 ,WN运0 算流图如图4.3.4
由以上例子我们可以看出,由于每一次分解都是按输入序
列在时域上的次序是偶数还是奇数来抽取的,最终分解成N个1
点DFT,因此称为基2时分FFT。
第4章 快速傅里叶变换(FFT)
4.1 引言 4.2 提高DFT运算效率的基本途径 4.3 基2时分FFT算法 4.4 基2频分FFT算法 4.5 IDFT的快速算法 4.6 实序列DFT的有效计算方法 4.7 线性调频Z(Chirp-Z)变换算法 习题与上机题
第4章 快速傅里叶变换(FFT)
4.1 引 言
利用旋转因子的可约性,即
W nkm Nm
WNnk
,
X(k)
X
(k)
N 2
1
x1
(r
)WNrk/
2
WNk
N 2
1
x2
(r
)WNrk/
2
,
0 k N 1
r 0
r 0
第4章 快速傅里叶变换(FFT)
由于X1(k)和X2(k)都隐含周期性,周期为 N ,因此上式 2
X (k) X~1(k) WNk X~2 (k), 0 k N 1
N 2
)
X1(k
)
WNk
X
2
(k
)
第4章 快速傅里叶变换(FFT)
X (k) X1(k) WNk X 2 (k),
0 k N 1 2
X
(k
N 2
)
X1(k
数字信号处理 第4章 信号与系统的复频域分析
极点的分布反映了系统的各种特征。
系统函数往往用零点和极点在S平面上的分 布图来表示,以”○”表示零点,以”×” 表示极点,以“⊙”表示重零点,以”*” 表示重极点。
jω
×
1
○
*
-2
-1
○
01
○
2
σ
×
-1
H
(s)
s(s (s2 2s
求上式的拉氏反变换,就可以得到系统的
冲激响应为:
n
h(t) bm kie pit i 1
每一极点对应一分量 epit ,(有r重极点时对 应 t e r1 pit ),极点位置就决定了该分量 的时域性质。
在H(s)的系数都为实数时,如果有一极点
为复数,必有另一极点是该极点的共轭复 数,同时系数k也将为共轭复数,一对共轭 极点组成的响应分量仍然为实数。
系统稳定性:对于任何一个有界的激励, 稳定系统产生的响应在任何时候都是有界 的。也就是要求系统的冲激响应有界(随 着t→∞,|h(t)|将逐渐衰减到零)。系统的 冲激响应的时域性质可由系统函数的极点 位置确定,因此,系统的稳定性可由系统 函数的极点位置来判断。
1、系统函数的极点全部位于左半S平面时, 随着t→∞将逐渐衰减到零,系统稳定。因
1
F (s)estds F (s)estds
2 j C0 Ci
Ci
0
k
Re
s(sk
)
1
2
j
Ci
F
(s)e st ds
F (s)estds 0 t 0
C1
F (s)estds 0 t 0
C2
数字信号处理 第四章
线性相关的FFT算法
1. 2. 3. 4.
计算步骤: X 求N点FFT, ( k ) DFT x ( n ) ; Y 求N点FFT, ( k ) DFT y ( n ) ; 求乘积,R ( k ) X ( k )Y ( k ) ; r 求N点IFFT, ( n ) IDFT R ( k ) 。 同样,可以利用已有的FFT程序计算IFFT, 求 1 1
mF 3 2 N log
2
3 N N N 1 log 2
2
N
线性卷积的FFT算法
[结论]:用线性相位FIR滤波器来比较直接计算 线性卷积和FFT法计算线性卷积这两种方法 的乘法次数,得
Km md mF ML 3 2 N 1 log 2 N 2 ML 3 2M L 11 log 2 M L 1 2
算法原理 运算量 按时间抽选的FFT算法的特点
按时间抽选的FFT算法的特点
原位运算(同址运算)
运算规律:每级(每列)计算都是由N/2个蝶形运 算构成,每一个蝶形结构完成下述基本迭代运算:
X m ( k ) X m 1 ( k ) X m 1 ( j )W Nr r X m ( j ) X m 1 ( k ) X m 1 ( j )W N
第4章 快速Fourier变换(FFT)
4.1 引言 4.2 直接计算DFT的问题及改进的途径 4.3 按时间抽选(DIT)的基-2FFT算法(库利- 图基算法) 4.4 按频率抽选(DIF)的基-2 FFT算法(桑德 -图基算法) 4.5 离散Fourier反变换(IDFT)的快速计算方法 4.10 线性卷积与线性相关的FFT算法
数字信号处理第四章
第四章线性时不变离散时间系统的频域分析一、传输函数和频率响应例4。
1传输函数分析Q4.1clear;M = input('Enter the filter length M: ');w = 0:2*pi/1023:2*pi;num = (1/M)*ones(1,M);den = [1];h = freqz(num, den, w);subplot(2,1,1)plot(w/pi,abs(h));gridtitle(’Magnitude Spectrum |H(e^{j\omega})|')xlabel(’\omega /\pi');ylabel('Amplitude’);subplot(2,1,2)plot(w/pi,angle(h));gridtitle('Phase Spectrum arg[H(e^{j\omega})]')xlabel(’\omega /\pi’);ylabel('Phase in radians’);M=2 M=10 M=15幅度谱为偶对称,相位谱为奇对称,这是一个低通滤波器。
M越大,通带越窄且过渡带越陡峭。
Q4。
2使用修改后的程序P3。
1,计算并画出当w=[0,pi]时传输函数的因果线性时不变离散时间系统的频率响应.它表示哪种类型的滤波器?w = 0:pi/511:pi;num = [0.15 0 -0.15];den = [1 —0。
5 0.7];如下图1这是一个带通滤波器.图1 图2Q4。
3对下面的传输函数重做习题Q4。
2:,式(4。
36)和式(4.37)给出的两个滤波器之间的区别是什么?你将选择哪一个滤波器来滤波,为什么?w = 0:pi/511:pi;num = [0.15 0 —0。
15];den = [0。
7 -0。
5 1];如上图2也是一个带通滤波器,这两个滤波器的幅度谱是一样的,相位谱不太一样,我会选择第一个带通滤波器,因为它的相位谱更加平滑,相位失真小。
《数字信号处理》 第4章
右图为描述倒位序的树状图(N=8)
5 倒位序的实现
对照表
变址功能
产生倒序数的十进制运算规律 N=2M,用M位二进制数表示,则从左至右的十进制权值为:
N 1 4
x1(2l)WNk22l
N 1 4
x1(2l
1)WNk22l1
r0
l0
l0
N1
N1
4
4
x3(l)WN kl4WN k2 x4(l)WN kl4
l0
l0
X 3(k) W N k2X 4(k),k0 ,1 ,
,N 1 2
式中
N1 4
N1 4
X3(k)DFTx3(l) x3(l)WN kl4 X4(k)DFTx4(l) x4(l)WN kl4
47线性调频变换chirp变换算法471算法原理已知序列xn0nn1是有限长序列其z变换为为适应z可沿z平面更一般的路径取值就沿z平面上的一段螺线作等分角的采样z的这些采样点zk为因此有其中a决定起始采样点z0的位置a0表示z0的矢量半径长度通常取a010表示z0的相角0表示两相邻采样点之间的角度差w0一般为正值表示螺线的伸展率图471线性调频变换在平面的螺线采样当mn即时各采样点zk就均匀等间隔地分布在单位圆上这就是求序列的dft
N
W N k(N n)W N (N k)nW N kn,
W
2 N
1
N
k
WN 2
WNk
利用这些特性,使DFT运算中有些项可以合并,并且可以 将长序列的DFT分解为几个短序列的DFT,以减少DFT的运算 次数。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
例如 10点 DFT 100次 复数乘;
1024点 DFT 1,048,576次 复数乘,即100万次的复数乘运算!
二、DFT的高效计算
j 2π
W
N
W
k N
e
e
N j 2πk
N
N
WN2
•
• WNk
0
2 k N
•WN0
WN•k
周期性 WWN0NkrNWNrNWNk1
对称性:
N WWN2N(k N2)W
n0
n0
k 0,1,....N 1
N次 复数乘 计算一个X(k)工作量: 或4N次 实数乘
( N-1)次 复数加 2N+2(N-1)=4N-2次 实数加
全部计算N个X(k) N²次 复数乘
N( N-1)次 复数加
或 4N ²次 实数乘
N(4N-2)次 实数加
结论:直接计算DFT的计算量和N的平方成正比
用下面的蝶形图也可清楚地说明这种运算。
A
A+BC X1(k)☉ •
C
B
A-BC
X2(k)☉
W
k
N•
蝶形运算符号
运算量:几次乘?几次加? 一个蝶形运算:一次复数乘、两次复数加
经过一次时域抽取计算量:
X(k)
•☉
N
-1
X(k ) •☉ 2
复数乘: 2 ( N / 2)2 N / 2 N ( N 1) / 2
(2l
1)W
k N
(2l /2
1)
l0
l0
N / 41
N / 41
x3
(l
)W
kl N/
4
WNk
/
2
x4
(l
)W
kl N/
4
l0
l0
于是
X3(k) WNk / 2 X4(k) ,
k 0,1,... N 1 2
X1(k) X 3(k) WNk /2 X4(k)
X1(k
N
/
4)
X3(k)
碟形运算
复数加: 2N / 2( N / 2 1) 2N / 2 N 2 / 2
运算量减少近一半
以N=8为例
x(0) ☉ x(1) ☉ x(2) ☉
.☉ .☉ .☉
☉
DFT
(N=8)
x(7) ☉
☉ X(0) ☉ X(1) ☉ X(2) ☉. ☉. ☉.
☉
☉ X(7)
直接计算需要8×8次复数乘、8×7次复数加
r0
r0
其中X1 (k)和X2 (k)分别为 x(2r)和x(2r+1) 的N/2点DFT:
N 1 2
X1(k) x(2r)WNkr/2 DFT[x(2r)] r0
, k 0,1,... N 1 2
N 1 2
X2(k) x(2r 1)WNkr/2 DFT[x(2r 1)] r0
0 N
1
WNk
可约性:
WWNnNnkk
WmmNnk
W
nk N/
/m m
利用WN 因子的周期性和对称性,可导出一个高效的快速算法
1965年 Cooley & Tukey 奠定FFT,把长序列短分解,
使得乘法计算量由N 2 次降为
N 2 log2 N
次。
以1024点为例,计算量降为5120 次,仅为原来的4.88%。
x(0) ☉
X1(0)
x(2) ☉ DFT X1(1)
x(4) ☉ (N=4) X1(2)
•
x(6) ☉ x(1) ☉ x(3) ☉
DFT
X1(3) W80
X2(0) X2(1)
W81
•
x(5) ☉ (N=4) X2(2)
x(7) ☉
X2(3)
☉ X(0) ☉ X(1) ☉ X(2) ☉ X(3) ☉ X(4) ☉ X(5) ☉ X(6) ☉ X(7)
36次复数乘 32次复数加
3、第二次抽取
将 x1(r) 按奇偶分解成两个N/4长的子序列x3(l ) 和 x4(l )
x3 (l )
x1(2l )
,
x4 (l ) x1(2l 1)
l 0,1,... N 1 4
N / 41
N / 41
X1(k)
x1 (2l
)W
2kl N /2
x1
6
(k
)
X2(k
N
/
4)
X5(k)
W
k N
/
2
X
6
(
k
)
k 0,1,... N 1 4
经过第二次分解,将N/2点的DFT分解成两个N/4点的DFT和N/4个蝶形运算。
数字信号处理
第四章 快速傅立叶变换(FFT)
§4.1 引言 § 4.2 按时间抽取(DIT)的FFT算法 § 4.3 按频率抽取(DIF)的FFT算法 § 4.4 离散傅立叶反变换(IDFT)的快速计算方法 § 4.5 进一步减少运算量的措施
§4.1 引言
一、DFT直接计算工作量很大
N 1
N 1
X (k) x(n)WNkn {(Re x(n) ReW Im x(n) Im W ) j(Re x(n) Im W Im x(n) ReW )}
r 0,1,... N 1 2
x(n)的DFT为:
N 1
N 1
2
2
X (k)
x(n)WNkn
x(n)WNkn= x(2r )WN2kr x(2r 1)WNk(2r1)
n偶 数
n奇 数
r0
r0
N 1
N 1
2
2
x1
(r
)W
kr N/
2
W
k N
x2(r)WNkr/ 2 X1(k) WNk X 2(k) , k 0,1,...N 1
, k 0,1,... N 1 2
由于X1 (k)和X2 (k)都是N/2点DFT,而X (k)有N点,所以得计算后N/2点.
由于它们均以N/2为周期,且
X1(k
N) 2
N 1 2
(k N )r
x1(r )WN / 2 2
r0
X1
(k
WN
(k)
N 2
)
W 同理
k N
,因此
N X2(k 2
)
本章以基2 的FFT算法为重点
§4.2 按时间抽取(DIT)的FFT算法(库利-图基算法)
一、算法原理(时域奇偶分,频域前后分)
设x(n)长度N, N=2M, M为自然数
1、第一次抽取:
将x(n)按偶、奇分成两组,可得两各自长度为N/2的奇偶序列
x1(r) = x(2r)
x2(r) = x(2r+1)
பைடு நூலகம்
W
k N
/
2
X4(k)
k 0,1,... N 1 4
同样,将 x2 (r ) 按奇偶分解成两个N/4长的子序列 x5 (l ) 和 x6 (l )
x5 (l )
x2 (2l )
,
x6 (l ) x2 (2l 1)
l 0,1,... N 1 4
X 2 (k )
X
5
(
k
)
W
k N
/
2
X
X2(k)
X (k) X1(k) WNk X2(k),
k 0,1,... N 1 2
X1(k)☉ •
X(k
N 2
)
X1(k)
WNk
X2(k),
k 0,1,... N 1 2
X2(k)☉
W
k
N•
-1
X(k)
•☉
N X(k ) •☉ 2
这样,将一个N点DFT分解成两个N/2点 DFT。
2、蝶形运算