第二章 离散傅里叶变换 及其快速算法(8学时 )

Presentation on theme: "第二章 离散傅里叶变换 及其快速算法(8学时 )"— Presentation transcript:

1 第二章 离散傅里叶变换 及其快速算法(8学时 )
数字信号处理 第二章 离散傅里叶变换 及其快速算法(8学时 ) 引言 2.1 离散傅里叶变换（DFT） 2.2 快速傅里叶变换(FFT) 2.3 FFT应用

2 数字信号处理 学习目标 了解周期序列的傅里叶级数及性质，掌握周期卷积过程 理解离散傅里叶变换及性质，掌握圆周移位、共轭对称性，掌握圆周卷积、线性卷积及两者之间的关系 掌握FFT算法的计算量分析，按时间抽选的基-2FFT算法

3 数字信号处理 引言： 离散傅里叶变换不仅具有明确的物理意义，更便于用计算机处理。但是，直至上个世纪六十年代，由于数字计算机的处理速度较低以及离散傅里叶变换的计算量较大，离散傅里叶变换长期得不到真正的应用，快速离散傅里叶变换算法的提出，才得以显现出离散傅里叶变换的强大功能，并被广泛地应用于各种数字信号处理系统中。近年来，计算机的处理速率有了惊人的发展，同时在数字信号处理领域出现了许多新的方法，但在许多应用中始终无法替代离散傅里叶变换及其快速算法。

4 2.1离散傅里叶变换（DFT） 数字信号处理 1. DFT （Discrete Fourier Transform） 是重要的变换
1）分析有限长序列的有用工具。 2）在信号处理的理论上有重要意义。 3）在运算方法上起核心作用，谱分析、 卷积、相关都可以通DFT在计算机上 实现。

5 2. DFT是现代信号处理桥梁 数字信号处理 DFT要解决两个问题： 一是离散与量化， 二是快速运算。 傅氏变换 离散量化 信号处理
DFT(FFT)

6 2.2.1周期序列及其 离散傅里叶级数（DFS） 数字信号处理
周期序列不能进行Z变换，因为其在 n=-到+ 都周而复始永不衰减，即 z 平面上没有收敛域。但是，正象连续时间周期信号可用傅氏级数表达，周期序列也可用离散的傅氏级数来表示，也即用周期为N的正弦序列来表示。

7 数字信号处理 1.离散傅里叶级数 其中：

8 数字信号处理

9 数字信号处理 补充证明

10 数字信号处理

11 数字信号处理

12 数字信号处理 说明 DFS变换对公式表明，一个周期序列虽然是无穷长序列，但是只要知道它一个周期的内容（一个周期内信号的变化情况），其它的内容也就都知道了，所以这种无穷长序列实际上只有N个序列值的信息是有用的，因此周期序列与有限长序列有着本质的联系。

13 数字信号处理 2.DFS的性质 1）线性

14 数字信号处理 2）周期性 时域上周期序列的离散傅里叶级数在频域上仍是一个周期序列。

15 数字信号处理 3）序列的移位

16 数字信号处理

17 数字信号处理 4)共轭对称性

18 数字信号处理 5）周期卷积

19 数字信号处理

20 周期卷积与线性卷积的区别： (1) 线性卷积在无穷区间求和;周期卷积在一个主值周期内求和
数字信号处理 周期卷积与线性卷积的区别： (1) 线性卷积在无穷区间求和;周期卷积在一个主值周期内求和 (2) 两个不同长度的序列可以进行线性卷积;只有同周期的两个序列才能进行周期卷积，且周期不变

21 数字信号处理 例 m 求两个序列的周期卷积 m 计算区

22 数字信号处理 m 计算区 m

23 数字信号处理 m m

24 数字信号处理

25 数字信号处理 n 1 3 4 计算区

26 周 期 卷 积

27 周期卷积计算方法 数字信号处理 1 1 1 1 0 0 y(n) 0 1 2 1 0 0 0 0 0 1 2 1 Y(0)=1
y(n) Y(0)=1 Y(1)=1 Y(2)=3 Y(3)=4 Y(4)=4 Y(5)=3

28 数字信号处理 例 m 求两个序列的线性卷积 m 计算区

29 线性卷积计算方法 数字信号处理 1 1 1 1 0 0 y(n) 0 1 2 1 0 0 0 0 1 2 1 Y(0)=0 0 0 1 2
y(n) Y(0)=0 1 0 Y(1)=1 Y(2)=3 0 0 Y(3)=4 Y(4)=4 Y(5)=3 Y(6)=1 1 Y(7)=0

30 2.1.2离散傅里叶变换DFT 数字信号处理 一.预备知识 1.余数运算表达式 如果 ， m为整数；则有：
如果 ， m为整数；则有： 此运算符表示n被N除，商为m，余数为 。 是 的解,或称作取余数,或说作n对N取模值, 或简称为取模值,n模N。 数字信号处理

31 数字信号处理 例

32 数字信号处理 2. 先取模值，后进行函数运作; 而 将 视作周期延 拓。

33 数字信号处理

34 定义从n=0 到（N-1）的第一个周期为主值序列或区间。
数字信号处理 N-1 n x(n) ... ... n N-1 定义从n=0 到（N-1）的第一个周期为主值序列或区间。

35 周期延拓Matlab程序 数字信号处理 subplot(1,1,1) % 画x((n+4))11的图
n = 0:10; x = 10*(0.8) .^ n; n1 = -11:21; x1 = [zeros(1,11), x, zeros(1,11)]; subplot(2,2,1); stem(n1,x1); title('初始序列 x(n)') axis([-10,17,-1,11]);text(18,-1,'n') x2 = [x, x, x]; subplot(2,2,3); stem(n1,x2); title('周期延伸') axis([-11,21,-1,11]);text(18,-1,'n')

36 数字信号处理

37 例周期延拓 数字信号处理 (1)周期延拓：N=5时 2 1 3 0.5 n x(n) 2 1 3 x(n) 0.5 n

38 数字信号处理

39 数字信号处理 有限长序列隐含着周期性。

40 数字信号处理

41 数字信号处理

42 数字信号处理 例

43 数字信号处理 6、DFT的性质

44 数字信号处理

45 n N-1 1 2 3 4 5 n=0 N=6 数字信号处理 n 周期延拓 左移2

46 循环移位程序 数字信号处理 subplot(1,1,1)% a) 画x((n-6))15的图
n = 0:10; x = 10*(0.8) .^ n; y = cirshftt(x,6,15); %右移 n = 0:14; x = [x, zeros(1,4)]; subplot(2,1,1); stem(n,x); title('初始序列') ylabel('x(n)'); axis([-1,15,-1,11]);text(15.5,-1,'n') subplot(2,1,2); stem(n,y); title('循环移位序列, N=15') ylabel('x((n-6) mod 15)'); axis([-1,15,-1,11]);text(15.5,-1,'n')

47 function y = cirshftt(x,m,N)
% 长度为 N 的x序列: (时域)作m采样点圆周移位 % [y] = cirshftt(x,m,N)% y = 包含圆周移位的输出序列 % x = 长度 <= N的输入序列% m = 移位采样数 % N = 圆周缓冲器长度 % 方法: y(n) = x((n-m) mod N) % Check for length of x if length(x) > N error('N 必须 >= x的长度') end x = [x zeros(1,N-length(x))]; n = [0:1:N-1]; n = mod(n-m,N); y = x(n+1);

48 数字信号处理

49 数字信号处理 有限长序列的循环移位导致频谱线性相移 而对频谱幅度无影响。

50 数字信号处理

51 数字信号处理

52 数字信号处理 b)时域循环卷积过程： 1）补零 2）周期延拓 3）翻褶，取主值序列 4）循环移位 5）相乘相加

53 数字信号处理 N-1 n

54 数字信号处理 y(0) y(1) y (2) y (3)

55 数字信号处理

56 数字信号处理 2 3 1 N-1 n *

57 线性卷积与循环卷积步骤比较 线性卷积 翻转、移位、相乘求和 2 3 1 x(n) 5 4 n N1=5 2 1 3 h(n) n N2=3
N1=5 2 1 3 h(n) n N2=3 线性卷积 翻转、移位、相乘求和 X(m) y(n) h(m) h(-m) 3 2 1 Y(0)=5 h(1-m) 3 2 1 Y(1)=14 h(2-m) Y(2)=26 h(3-m) Y(3)=20 h(4-m) Y(4)=14 h(5-m) 3 2 Y(5)=8 h(6-m) 2 1 Y(6)=3 h(7-m) Y(7)=0

58 数字信号处理 得到线性卷积结果的示意图 14 26 5 y(n) 20 8 3 N=7 n

59 循环卷积 数字信号处理 2）其中一个序列周期延拓 3）翻褶，取主值序列 4）循环移位 5）相乘相加 2 3 1 x(n) 5 4 n 2 1
2 1 3 h(n) n N2=3 N1=5 循环卷积 1）循环卷积：(N=7)不足的，补零加长 2）其中一个序列周期延拓 3）翻褶，取主值序列 4）循环移位 5）相乘相加

60 2 1 3 h(k) 2 3 1 5 4 N=7 k 2 3 1 h((k))N k k X(m) 5 4 3 2 1 0 0 y(n)
2 3 1 5 4 N=7 k X(m) y(n) h(m) h((m))NRN h((-m))NRN Y(0)=5 h((1-m))NRN Y(1)=14 h((2-m))NRN Y(2)=26 h((3-m))NRN Y(3)=20 h((4-m))NRN Y(4)=14 h((5-m))NRN Y(5)=8 h((6-m))NRN Y(6)=3 h((7-m))NRN Y(7)=5 2 3 1 h((k))N k k

61 可见,线性卷积与循环卷积相同(当N≥[N1(5)+N2(3)-1]=7时)
数字信号处理 得到循环卷积的示意图 14 26 5 n y(n) 20 8 3 可见,线性卷积与循环卷积相同(当N≥[N1(5)+N2(3)-1]=7时)

62 数字信号处理 N=5 X(m) 5 4 3 2 1 y(n) h(m) 1 2 3 0 0 h((m))NRN h((-m))NRN
y(n) h(m) h((m))NRN h((-m))NRN Y(0)=13 h((1-m))NRN Y(1)=17 h((2-m))NRN Y(2)=26 h((3-m))NRN Y(3)=20 h((4-m))NRN Y(4)=14 h((5-m))NRN Y(5)=13 h((6-m))NRN Y(6)=17 h((7-m))NRN Y(7)=26 N=5

63 可见,线性卷积与循环卷积不同(当N<[N1(5)+N2(3)-1]=7时)
数字信号处理 得到循环卷积的示意图 17 26 13 n y(n) 20 14 可见,线性卷积与循环卷积不同(当N<[N1(5)+N2(3)-1]=7时)

64 数字信号处理 总结

65 数字信号处理 循环卷积

66 循环卷积程序 x1 = [1,2,2]; x2 = [1,2,3,4]; y = circonvt(x1,x2,4) Stem(y)

67 N=4

68 N=7

69 function y = circonvt(x1,x2,N)
% x1 = 长度 N1 <= N 的输入序列 % x2 = 长度 N2 <= N 的输入序列 % N = 循环缓冲器的大小 % 方法 y(n) = sum (x1(m)*x2((n-m) mod N)) % Check for length of x1 if length(x1) > N error('N 必须 >= x1的长度') end % Check for length of x2 if length(x2) > N error('N 必须 >= x2的长度') x1=[x1 zeros(1,N-length(x1))]; x2=[x2 zeros(1,N-length(x2))]; m = [0:1:N-1]; x2 = x2(mod(-m,N)+1); H = zeros(N,N); for n = 1:1:N H(n,:) = cirshftt(x2,n-1,N); y = x1*H';

70 线性卷积 function [y,ny]=conv_improve(x,nx,h,nh)
%conv(x,h)可以实现两个有限长度序列的卷积 ny1=nx(1)+nh(1); ny2=nx(length(x))+nh(length(h)); ny=[ny1:ny2]; y=conv(x,h); 在命令窗口调用卷积函数。x=[ ]; nx=[-3:3]; h=[ ]; nh=[N:N+5];

71 例

72 X1(n) n X2(n) n X2((m))10 n X1((m))10 n n

73 周期卷积 X1((m))10 X2((8-m))10 X2((9-m))10 X2((10-m))10 周期性 1111100000
X2((-m))10 Y(0)=-3 X2((1-m))10 Y(1)=-1 X2((2-m))10 Y(2)=1 X2((3-m))10 Y(3)=3 X2((4-m))10 Y(4)=5 X2((5-m))10 Y(5)=3 X2((6-m))10 Y(6)=1 X2((7-m))10 Y(7)=-1 X2((8-m))10 Y(8)=-3 X2((9-m))10 Y(9)=-5 X2((10-m))10 Y(10)=-3 周期性

74 N=10周期卷积结果 Y(n) n 9

75 X1(n) n n X2(n) X2((m))10 n X2((-m))10 n

76 循环卷积 X2((8-m))10 X2((9-m))10 X2((10-m))10 X1(m) 1111100000 X2(m)
X2((m))10 X2((-m))10 Y(0)=-3 X2((1-m))10 Y(1)=-1 X2((2-m))10 Y(2)=1 X2((3-m))10 Y(3)=3 X2((4-m))10 Y(4)=5 X2((5-m))10 Y(5)=3 X2((6-m))10 Y(6)=1 X2((7-m))10 Y(7)=-1 X2((8-m))10 Y(8)=-3 X2((9-m))10 Y(9)=-5 X2((10-m))10 Y(10)=-3

77 N=10点的循环卷积结果 n 9 Y(n)

78 线性卷积 X2(8-m) X2(9-m) X2(10-m) X2(11-m) X2(12-m) X2(13-m) X1(m)
X2(m) X2(-m) 1 Y(0)=1 X2(1-m) 11 Y(1)=2 X2(2-m) 111 Y(2)=3 X2(3-m) 1111 Y(3)=5 X2(4-m) 11111 Y(4)=5 X2(5-m) Y(5)=3 X2(6-m) -1-1-1 Y(6)=1 X2(7-m) -1-1 Y(7)=-1 X2(8-m) -1 Y(8)=-3 X2(9-m) Y(9)=-5 X2(10-m) Y(10)=-4 X2(11-m) Y(11)=-3 X2(12-m) Y(12)=-2 X2(13-m) Y(13)=-1

79 Y(n) n Y(n) n 9 循环卷积 Y(0)=-3 Y(1)=-1 Y(2)=1 Y(3)=3 Y(4)=5 Y(5)=3
线性卷积 Y(0)=1 Y(1)=2 Y(2)=3 Y(3)=5 Y(4)=5 Y(5)=3 Y(6)=1 Y(7)=-1 Y(8)=-3 Y(9)=-5 Y(10)=-4 Y(11)=-3 Y(12)=-2 Y(13)=-1 9 n Y(n)

80 数字信号处理 2.2 利用DFT进行连续信号的频谱分析 1、混叠 2、泄漏 3、栏栅效应 4、DFT的分辨率 5、周期信号的谱分析

81 数字信号处理

82 数字信号处理 离散傅里叶变换不仅具有明确的物理意义，更便于用计算机处理。但是，直至上个世纪六十年代，由于数字计算机的处理速度较低以及离散傅里叶变换的计算量较大，离散傅里叶变换长期得不到真正的应用，快速离散傅里叶变换算法的提出，才得以显现出离散傅里叶变换的强大功能，并被广泛地应用于各种数字信号处理系统中。近年来，计算机的处理速率有了惊人的发展，同时在数字信号处理领域出现了许多新的方法，但在许多应用中始终无法替代离散傅里叶变换及其快速算法

83 数字信号处理 2.3 快速傅里叶变换 有限长序列通过离散傅里叶变换 (DFT)将其频 域离散化成有限长序列.但其计算量太大(与N的平方成正比）, 很难 实时地处理问题 , 因 此 引 出 了 快 速 傅 里 叶 变 换(FFT) . FFT 并 不 是 一 种 新 的 变 换 形 式 ,它 只 是 DFT 的 一 种 快 速 算 法 . 并 且 根 据 对 序 列 分 解 与 选 取 方 法 的 不 同 而 产 生 了 FFT 的 多 种 算 法 . FFT 在 离 散 傅 里 叶 反 变 换 、 线 性 卷 积 和 线 性 相 关 等 方 面 也 有 重 要 应 用。

84 FFT产生故事 虽然频谱分析和DFT运算很重要，但在很长一段时间里，由于DFT运算复杂，并没有得到真正的运用，而频谱分析仍大多采用模拟信号滤波的方法解决。 当时加文(Garwin)在自已的研究中极需要一个计算傅里叶变换的快速方法。他注意到图基(J.W.Turkey)正在写有关傅里叶变换的文章，因此详细询问了图基关于计算傅里叶变换的技术知识。图基概括地对加文介绍了一种方法，它实质上就是后来的著名的库利(Cooley J.W)图基算法。在加文的迫切要求下，库利很快设计出一个计算机程序。1965年库利--图基在<计算数学>、Mathematic of Computation 杂志上发表了著名的“机器计算傅里级数的一种算法”文章，提出一种快速计算DFT的方法和计算机程序--揭开了FFT发展史上的第一页。 促使FFT算法产生原因还有1967年至1968年间FFT的数字硬件制成，电子数字计算机的条件， 使DFT的运算大简化了。 1984年，法国的杜梅尔(P.Dohamel)和霍尔曼(H.Hollmann)将基2分解和基4分解糅合在一起，提出了分裂基FFT算法。其运算量比前几种算法都有所减少，运算流图却与基2FFT很接近，运算程序也很短。它是目前一种实用的高效新算法。

85 直接计算DFT算法存在的问题及改进途径 问题提出： 设有限长序列x(n),非零值长度为N, 计算对x(n)进行一次DFT运算，共
数字信号处理 直接计算DFT算法存在的问题及改进途径 问题提出： 设有限长序列x(n),非零值长度为N, 计算对x(n)进行一次DFT运算，共 需多大的运算工作量?

86 数字信号处理 1.比较DFT与IDFT的运算量 其中x(n)为复数， 也为复数 所以DFT与IDFT二者计算量相同。

87 2.以DFT为例复数运算量 N*N次复数相乘和N*(N-1)次复数加法 数字信号处理 计算一个X(k)(一个频率成分)值，运算量为

88 3.一次复数乘法换算成实数运算量 数字信号处理 一个复数乘法包括4个实数乘法和2个实数相法。
(a+jb)(c+jd)=(ac-bd)+j(bc+ad) 所以所有X(k)就要4N2次实数乘法运算,2N2＋2N(N-1)＝ N(4N-2)次实数加法运算.当N很大时,运算量将是惊人的, 这样,难以做到实时处理。 4次实数乘法 2次实数加法

89 例子 数字信号处理 例1:当N=1024点时，直接计算DFT需要： N2=220=1048576次，即一百多万次的复乘运算
这对实时性很强的信号处理(如雷达信号处理)来 讲，它对计算速度有十分苛刻的要求-->迫切需 要改进DFT的计算方法，以减少总的运算次数。 例2：石油勘探，24道记录，每道波形记录长度5秒，若每秒抽样500点/秒， 每道总抽样点数=500*5=2500点 24道总抽样点数=24*2500=6万点 DFT运算时间=N2=(60000)2=36*108次

90 数字信号处理 5、 FFT的计算工作量 FFT算法对于N点DFT,仅需(N/2)log2N 次复数乘法运算和Nlog2N 次复数加法。

91 例 数字信号处理 如果计算机的速度为平均每次复数乘需要5×10-6秒 ，每次复加需要1×10-6秒，用来计算N=1024点DFT,
问1)直接计算需要多少时间？2）用FFT算法计算需要多少时间？

92 6. FFT算法分类： 数字信号处理 1.按抽取方法分： 时间抽取法(DIT Decimation-In-Time);
频率抽取法(DIF Decimation-In-Frequency) 2.按“基数”分： 基-2FFT算法；基-4FFT算法；混合基FFT算法；分裂基FFT算法 3.其它方法： 线性调频Z变换(Chrip-z法)

93 数字信号处理 2.3.1 按时间抽取的DFT 1、 的特性

94 数字信号处理 例子 利用以上特性,得到改进DFT直接算法的方法.

95 2、DFT的基本思想 数字信号处理 快速傅里叶变换( FFT ) 就是在此特性基础上发展起来的：
（1）利用DFT系数的对称性和周期性,合并DFT运算中的某些项； （2）将长序列分解为短序列,从而减少其运算量。    因合并与分解方法的不同产生了多种DFT的快速算法。

96 2.3.1 时域抽取法基2FFT基本原理Decimation-in-Time(DIT)
数字信号处理 2.3.1 时域抽取法基2FFT基本原理Decimation-in-Time(DIT) 1、时域抽取算法原理 设输入序列长度为N=2M(M为正整数，将该序列按时间顺序的奇偶分解为越来越短的子序列，称为基2按时间抽取的FFT算法。也称为Coolkey-Tukey算法。 其中基数2----N=2M，M为整数.若不满足这个条件，可以人为地加上若干零值（加零补长）使其达到 N=2M

97 例子 数字信号处理 设一序列x(n)的长度为L=9,应加零补长为 N=24=16 应补7个零值。 x(n)

98 数字信号处理 2.算法步骤 设序列点数 N = 2M，M为整数。若不满足，则补零

99

100 数字信号处理 注

101 数字信号处理 4.蝶形运算 1 -1 1 1 -1

102 *但是,N点DFT的复乘为N2 ;复加N(N-1);与 分解后相比可知,计算工作点差不多减少 一半。
数字信号处理 (1)N/2点的DFT运算量:复乘次数: 复加次数: (2)两个N/2点的DFT运算量:复乘次数: (3)N/2个蝶形运算的运算量:复乘次数: 复乘: 总共运算量: 复加: *但是,N点DFT的复乘为N2 ;复加N(N-1);与 分解后相比可知,计算工作点差不多减少 一半。

103 例 8点FFT的算法 首先可以分解为两个N/2=4点的DFT.具体方法如下:

104 （1）N=8点分解成2个4点的DFT的信号流图 X(0) X(1) X(2) X(3) X(4) X(5) X(6) X(7) x(0)
x(2) N/2点 x(4) DFT x(6) x(1) x(3) N/2点 x(5) DFT x(7) X1(1) X(1) X1(2) X(2) X1(3) X(3) X2(0) N W X(4) -1 X2(1) N W 1 X(5) -1 X2(2) 2 X(6) W N -1 X2(3) X(7) 3 W N X(4)~X(7) 同学们自已写

105 (2) N/2(4点)-->N/4(2点)FFT
由于N=2 L ,所以 N/2仍为偶数,可以进一步把每个N/2点的序列再按其奇偶部分分解为两个N/4的子序列。

106

107 X(0) X(0) N/4 DFT X1(0) X1(1) X1(2) X1(3) X2(0) X2(1) X2(2) X2(3)
-1 W N X(6) X(3) 2 -1 N X(1) X(4) N/4 DFT WN0 WN1 WN2 WN3 -1 X(5) X(5) -1 X(3) N/4 DFT X(6) -1 W N N X(7) X(7) 2 -1

108 (3)N/4（2点）->2个1点DFT 最后剩下两点DFT,它可分解成两个一点DFT，但一点DFT就等于输入信号本身，所以两点DFT可用一个蝶形结表示。取x(0)、x(4)为例。

109 2个1点的DFT蝶形流图 X3(0) X3(1) 1点DFT x(0) 1点DFT x(4) -1

110 一个完整的按时间抽取的8点FFT

111 6、对 N = 2L点FFT流程图画法 数字信号处理 1）输入倒位序，输出自然序
自然顺序n 二进制n n n 倒位序二进制n n n 倒位顺序n

112 数字信号处理 2）蝶形运算 对N = 2L点FFT，共需L级蝶形运算，每级有N/2个蝶形运算组成，蝶形运算两节点的距离:2L-1（L表示级数） 每个蝶形运算有一次复乘和再次复加。如 1 -1

113 数字信号处理 3） WNr 的确定(仅给出方法)

114 数字信号处理 X(0) X(0) X(4) X(1) X(2) X(2) -1 X(6) X(3) X(1) X(4) -1 X(5)
000，000 001，100 010，010 011，110 100，001 101，101 110，011 111，111 X3(0) X3(1) X4(0) X4(1) X5(0) X5(1) X6(0) X6(1) X(1) WN0 X(2) WN0 WN2 -1 X(3) X(4) WN0 WN1 WN2 WN3 -1 X(5) X(6) WN0 WN2 X(7)

115 数字信号处理 N=16 x(0) x(8) x(4) x(12) x(2) x(10) x(6) x(14) x(1) x(9) x(5)
0000，0000 0001，1000 0010，0100 0011，1100 0100，0010 0101，1010 0110，0110 0111，1110 1000，0001 1001，1001 1010，0101 1011，1101 1100，0011 1101，1011 1110，0111 1111，1111 WN0 WN0 WN4 WN0 WN0 WN2 WN4 WN6 WN0 WN0 WN4 WN0 WN1 WN2 WN3 WN4 WN5 WN6 WN7 WN0 WN0 WN0 WN4 WN0 WN0 WN2 WN4 WN6 WN0 WN0 WN4 WN0

116 6、FFT运算量 数字信号处理 由于计算机的乘法运算比加法运算所 需的时间多得多，故以乘法作为比较基准.
由上述分析可知,N=8需三级蝶形运算 N=23 =8,由此可知,N=2L 共需L级蝶形运算,而且每级都由N/2个蝶形运算 组成,每个蝶形运算有一次复乘,两次复加。 因此,N点的FFT的运算量为 复乘: mF =（N/2）L=（N/2） log2 N 复加: aF =N L=N log2 N 由于计算机的乘法运算比加法运算所 需的时间多得多，故以乘法作为比较基准.

117 数字信号处理 比较DFT

118 2.3.2 按频率抽取的FFT 数字信号处理 1、算法原理
设序列点数N=2L，L为整数。将X(k)按k的奇偶分组前，先将输入x(n)按n的顺序分成前后两半：

119 数字信号处理

120 数字信号处理 2.N点DFT按k的奇偶分组可分为两个N/2的DFT 当k为偶数,即 k=2r时,(-1)k =1； 当k为奇数,即 k=2r+1 时 (-1)k =-1 。 这时 X(k) 可分为两部分：

121 数字信号处理 可见,上面两式均为N/2的DFT。

122 数字信号处理 3.蝶形运算 -1

123 4.N=8时,按k的奇偶分解过程先蝶形运算，后DFT:再将N/2点DFT按k的奇偶分解为两个N/4点的DFT,如此进行下去,直至分解为 2点DFT。
W N 1 -1 W N -1 2 W N -1 3 W N -1

124 例如 N=8时频域抽取的FFT流图如下： 数字信号处理 x(1) X(4) x(2) X(2) x(3) X(6) x(4) X(1)
-1 W N 1 2 3

125 5. 时域抽取和频域抽取FFT的区别 数字信号处理 1）基本蝶形不同 时域: 先复乘后加减 频域: 先减后复乘 2）运算量相同
3）都可原位运算 4）两者基本蝶形互为转置

126 2.4 关于FFT应用中的几个问题 数字信号处理 2.4.1 用FFT计算IDFT 比较两式可知,只要DFT的每个系数 换
成 ,最后再乘以常数1/N就可以得到IDFT的快速算法--IFFT。

127 数字信号处理

128 直接调用FFT子程序计算IFFT的方法：
数字信号处理 直接调用FFT子程序计算IFFT的方法： 这就是说,先将X(k)取共轭,即将X(k)的虚部乘-1，直接利用FFT程序计算DFT；然后再取一次共轭；最后再乘1/N,即得 x(n)。所以FFT,IFFT可用一个子程序。 共轭 FFT 共轭 乘1/ N

129 数字信号处理 2.4.2 线性卷积的FFT算法

130

131 数字信号处理

132 数字信号处理

133 利用FFT计算线性卷积 MATLAB程序 数字信号处理 N=1024 x=[2 3 1 4 5 ones(1,N)];
h=[ ones(1,N)]; Lenx=length(x); %求序列x的长度 Lenh=length(h); %求序列h的长度 N=Lenx+Lenh-1; t=cputime(%或tic) Xk=fft(x,N); %计算x序列的DFT Hk=fft(h,N); %计算h序列的DFT Yk=Xk.*Hk; y=ifft(Yk) %求IDFT t1=cputime-t（或toc） stem(y); xlabel('n');

134 数字信号处理

135 直接线性卷积 MATLAB程序 数字信号处理 %[x,nx]为第一个信号 %[h,nh]为第二个信号
%conv(x,h)可以实现两个有限长度序列的卷积 x=[ ]; h=[ ]; ny2=length(x)+length(h)-1; t=cputime ny=[0:ny2]; y=conv(x,h); T1=cputime-t stem(y)

136 数字信号处理

137 数字信号处理 重叠相加法

138 x(n) L-1 x0(n) L+ N2 -1 x1(n) x2(n) N1-1 L-1 h(n) 3 N2 N2 2 N2 N2-1
x(n) 3 N2 N2 2 N2 N2-1 L-1 x0(n) 2 N2 -1 L+ N2 -1 N2 x1(n) 3 N2 -1 2 N2 2 N2+L -1 x2(n)

139 数字信号处理 * L-1 N2 L＋N2-1 L＋2N2-1 2N2

140 数字信号处理 重叠保留法 分段序列中补零的部分不是补零，而是保留原来的输入序列，于是重叠了输入信号段，省略掉输出段的重叠相加。

141 数字信号处理 舍弃yi(n)的前M-1个点，再将yi(n)顺次连接， 即得y(n)。 分段 右移序列 卷积 N *

142

143

144 说明 由于保留法补的不是零值点，而是前段保留下来的（M-1)点有x(n)值，因而每段循环卷积的结果的前（M-1）个点的值不等于线性卷积的值必须舍去。把相临输出段的序列衔接起来就构成了最后正确输出。

145 数字信号处理 2.4.4 用FFT计算相关函数 若L点x(n)，M点y(n)，计算线性相关： x(n)及y(n)的卷积公式

146 相关和卷积的时域关系

147 数字信号处理 FFT计算相关

148 FFT计算相关程序 x=[ ]; y=[ ]; k=length(x); xk=fft(x,2*k); yk=fft(y,2*k); rm=real(ifft(conj(xk).*yk)); rm=[rm(k+2:2*k) rm(1:k)]; m=(-k+1):(k-1); stem(m,rm) xlabel('m'); ylabel('幅度');

149 -8 -6 -4 -2 2 4 6 8 10 12 m 幅度 两个序列的自相关函数

150 通过FFT来分析其信号频率成分。 数字信号处理 t=0:0.001:0.6;
x=sin(2*pi*50*t)+sin(2*pi*120*t); %两个正弦信号合成一个输入 y=x+randn(1,length(t)); %加载噪声信号 Y=fft(y,512); %进行FFT变换 subplot(211);plot(y);title('受噪声污染的信号'); n=0:511;f=1000*n/512; %将信号横向放大，由512放大为1000 subplot(212);plot(f,abs(Y));title('FFT');

151 数字信号处理

152 重点习题 2.2，2.4，2.6（2），2.9，2.19，2.20，2.21（时间抽取）