第七讲 常用DSP算法
合集下载
相关主题
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
失，相反，由于FFT的计算步骤少，硬件位 数有限引起的积累误差小。
FFT算法原理
FFT算法节省运算时间的关键在于，它把数 据组x(n)分解为奇数下标和偶数下标，然后 利用指数函数（即旋转因子）的周期性消
去多余操作。
一个N点DFT可分解成两个N/2点DFT，同理， 逐步分解，经过log2N步之后可分解为N/2 个两点的计算（基-2碟式运算）。
取x(k)的共轭，得x*(k)； 求x*(k)的FFT； 再求共轭，乘以1/N，即得。
不同基的FFT算法
基-2算法：序列长度N为2的整数幂，即N=2L， 其精确复乘次数为（考虑省去非乘法的次数）：
N
3
2 log2 N 2 N 2
基-4算法：序列长度N为4的整数幂，即N=22L， 其精确复乘次数为（考虑省去非乘法的次数）：
时不变系统可以实现输入信号与系统冲激 响应之间的线性卷积。
两个序列卷积的傅里叶变换等于它们的傅 里叶变换相乘，即可在频域计算卷积。
频域取样导致信号时域周期重复，理论上 利用DFT只能计算圆周（循环）卷积，而不 能实现线性卷积。
圆周卷积等价于线性卷积的条件
N L K 1
其中，N为DFT的长度，序列x(n)的 长度为L，序列h(n)的长度为K。 为了利用DFT计算线性卷积，必须 选择满足上式的DFT长度N，并对 x(n)至少补K-1个零，对h(n)至少补 L-1个零。
FIR滤波器
数字滤波是数字信号处理中的基本处 理算法。
FIR（finite impulse response）滤波 器是最常用的一种数字滤波器。
FIR滤波器的差分方程与Z变换
N 1
y(n) ak x(n k) k 0
其中，x(n)为输入，y(n)为输出，ak为滤波器系数。
H (Z )
1
N 1
X (k)e j2nk / N (n 0,1,
,N
1)
N k0
定义旋转因子: WN e j 2nk / N
DFT的一些重要结论
有限长离散时间信号的频域离散表述可对 傅里叶变换取样得到；
只有当傅里叶变换一个周期内[0，2π）的 取样数N大于等于原信号长度L时，表述才 有实用意义。
当信号x(n)长度L小于N时，则对x(n)补零 构成N点序列，再进行DFT运算。
数字信号处理方法与实现
贺知明 副教授
电子科技大学 四川•成都
常用DSP算法
DFT(离散付里叶变换) FFT算法 FIR滤波器 IIR滤波器 自适应滤波器
DFT(离散付里叶变换)
DFT的定义为:
N 1
X (k) x(n)e j2nk / N (k 0,1, , N 1) n0
x(n)
窗函数法——基本方法
FIR滤波器的窗函数设计法
设所需滤波器的理想响应Hd(ejw) H d (e jw ) hd (n)e jnw
hd
(n)
1
2
2
0
Hd
(e
jw
)e
jnw dw
对hd(n)作截尾处理
Q
Ha (e jw ) hd (n)e jnw nQ
Q为FIR滤波器 的阶数，Q越 大，近似程度 越高
1965年，Cooley和Tukey提出快速付里叶 变换（FFT），极大地提高了运算速度。
运算量对比
N2
N 2
log 2
N
基-2 FFT运算近 似的复乘次数
例：1024个采样点 直接DFT：需1048576次乘法； FFT：仅需5120次乘法。
注意
FFT是DFT的一种快速算法； FFT没有对DFT作任何近似，精度无任何损
采用循环寻址功能，使长度为N的缓冲区足 以存储N阶滤波器所需的输入数据。
运算量为N/2 ·log2N
抽取（decimation）
两种形式的蝶式运算
DIT（按时间抽取）
C A WNk B D A WNk B
DIF （按频率抽取）
C AB D ( A B)WNk
位反序排列（以8点FFT为例）
0 1 234 56 7 000 001 010 011 100 101 110 111
如果时间窗长度正好为整数倍个信号周期， 除一点有值外，其余全为0，变换结果没有 产生谱泄漏。
通常情况，时间窗不等于信号周期的整数 倍，频谱中，除中心值外，还有以1/f（倍 频程）缓慢衰减的旁瓣。
时间窗加权
为消除谱泄漏，加权，压低旁瓣。 时间窗本身的作用相当于宽度与之相等的
一个矩形窗函数的加权。
对h(n)补零，使长度变为N； 用长度为N的DFT计算每个序列块的圆周卷
积。 拼接构成输出信号。
两种方法的区别
在重叠相加法中，计算长度N+K-1的DFT， 并将结果相加；
在重叠保留法中，只需计算长度N的DFT， 并通过舍弃、保留手段拼接构成最后输出 结果。
DFT
FFT
FFT算法
DFT运算归结为一系列的乘加运算，谱的 每一个频率点都要对N个乘法项求和，对N 个谱点总共作N2次乘法，计算量大，不利 于实时处理。
FIR滤波器的设计过程
通过计算机编程设计FIR滤波器，常用 Matlab仿真软件进行。
根据任务要求，首先确定所需FIR滤波器阶 数，然后利用仿真软件计算滤波器系数， 即完成FIR滤波器的设计。
FIR滤波器的计算过程
输入数据存在RAM缓冲区中，N个系数存 在系数存储器中，用程序循环执行
N 1
y(n) h(k)x(n k) k 0
000 100 010 110 001 101 011 111 0 4 26 15 37
反向FFT（IFFT）
x(n)
1 N
N 1
X (k )WNnk
k 0
1 N
N 1
[X
k 0
* (k)WNnk ]*
1 N
{FFT [x *(k)]}*
可以使用同一计算程序（硬件结构）作 正反两个方向的计算。三个步骤：
将序列x(n)分解成长度为N的序列块xm(n)； 分别对h(n)和xm(n)补零，使长度变为N+K-1； 用N+K-1点DFT计算每块的圆周卷积； 将结果相加，其中序列块卷积结果ym(n)的后K-
1个取样和ym+1(n)的前K-1个取样重叠，故称重 叠相加法。
重叠保留法
将序列x(n)分解成长度为N的序列块xm(n)， 其中有K-1个取样与相邻序列重叠，第一个 序列块在开头补K-1个零；
FIR滤波器的窗函数设计法
对hd(n)截尾，实际上是对hd(n)乘上一个矩 形窗，这是FIR滤波器设计的最直接和简便 的方法，但采用矩形窗存在较大的Gibbis现 象（有限项近似所产生的振荡），这是由 于矩形窗函数的边缘截止得太突然所引起， 加合适的窗可减小。
通常可采用前面介绍的其它经典窗函数 （Hamming窗、Hanning窗、Blackman窗 等）。
3 8
N
log
2
N
N
1
基-4算法比基-2算法更高效。
谱泄漏
DFT是对有限时间间隔（时间窗）里的采样 数据的变换，导致谱泄漏。
单位矩形序列r(n)乘以无限正弦波f(n)
f '(n) f (n)r(n) 单位矩形序列的
F '(k) F (k) R(k) 谱，SinC函数
无限正弦波谱，δ函数
说明
频域系数通常为固定系数，可提前算好存 起来，实时处理时不必再算，可省一次FFT 运算。
FFT在DSP实现过程中应注意的问题
DSP作FFT运算时，应尽量使用片内RAM， 速度快，并利用位反转（位反序）寻址方 式——专门为FFT运算提供的DSP功能。
定点DSP实现FFT，应注意防止中间结果溢 出，其方法是对中间数值归一化（以因子2 进行归一化，可运用DSP芯片的移位功能， 不会增加运算量），如作2m点FFT，可以2m 作归一化。
其它窗函数形式
矩形窗 三角窗 布莱克曼窗 凯塞窗 切比雪夫窗
FFT用于快速卷积
时域直接卷积 频域处理方法
实用中，这两种方法究竟哪一种速度快， 取决于输入序列的性质，通常采用FFT的频 域方法较快，但存在固有的处理延时。
FFT在DSP实现过程中应注意的问题
应调整信号的采样点数，使之正好或接近 为2（或4，由基数定）的幂次，尽量避免 大量补零的冗余。
对N个点的一次循环计算产生一个输出 y(n)，也就是说，FIR滤波计算就是循环卷积 计算，即移动的加权平均处理。
用前面方法计算得到的一组系数h(n) （即滤波器的设计），即可用于FIR计算。
FIR滤波器的DSP实现
注意输入数据x(n)与滤波器系数h(n)的合理 存放。
数据存储区
程序存储区
均为连续的存储区
选择窗函数的原则
使信号在窗的边缘为0，以减小截断所产生 的不连续效应。
信号经窗函数加权处理后，不应丢失太多 的信息。
选择什么样的窗函数，总是在旁瓣的衰减 程度、中心峰的宽度和能量损失等因素间 取折衷（一组矛盾）。
常用加权函数的一般形式
w(n) (1 ) cos ( 2n)
M
β= 2，α=0.5时，称为汉宁（Hanning）窗； β= 2，α=0.54时，称为汉明（Hamming）窗。
长序列与短序列卷积的DFT实现方法
在大多数实际情况中，长序列对应于系统 输入，短序列对应于系统冲激响应。
将长序列分割成长度为N的序列块，计算每 一块与短序列的卷积。
每一个序列块的卷积必须适当地合并，以 得到长序列和短序列卷积的最后结果。
有两种合并的方法，分别是重叠相加法和 重叠保留ห้องสมุดไป่ตู้。
重叠相加法
输入
h(N-1) x h(N-2) x
…...
h(0) x
+
Z-1
+
…...
Z-1
+
转置型横向滤波器
输出
FIR滤波器的特点
精确严格的线性相位 FIR滤波器是稳定系统 非递归结构 可采用FFT算法，提高运算效率
FIR滤波器的设计方法
频率抽样法——FIR滤波器的逼近方法
设计一个长度为N的FIR滤波器，可以对要求 达到的频率响应在N个点取样，把这些取样值进 行反傅氏变换，即可得到所要设计的滤波器的冲 激响应。
Y (Z )
N 1
h(k)Z k
X (Z ) k0
展开，得
H (Z ) h(0) Z 1 h(1) Z 1(h(2) h(N 1))
标准FIR滤波器结构
输入
Z-1
Z-1 …...
Z-1
h(0) x h(1) x …... h(N-1) x
+
输出
标准横向滤波器结构
转置FIR滤波器结构
DFT的一些重要结论
N越大，则它的DFT与傅里叶变换越近似， 因为在区间[0，2π）的取样数增加了。
N的取值通常需要根据实际应用中允许的运 算复杂度决定，因为补零个数越多，则DFT 计算所需要的运算和存储器件越多。
DFT实现数字滤波器
线性和圆周卷积 重叠相加法 重叠保留法
线性和圆周（循环）卷积
FFT算法原理
FFT算法节省运算时间的关键在于，它把数 据组x(n)分解为奇数下标和偶数下标，然后 利用指数函数（即旋转因子）的周期性消
去多余操作。
一个N点DFT可分解成两个N/2点DFT，同理， 逐步分解，经过log2N步之后可分解为N/2 个两点的计算（基-2碟式运算）。
取x(k)的共轭，得x*(k)； 求x*(k)的FFT； 再求共轭，乘以1/N，即得。
不同基的FFT算法
基-2算法：序列长度N为2的整数幂，即N=2L， 其精确复乘次数为（考虑省去非乘法的次数）：
N
3
2 log2 N 2 N 2
基-4算法：序列长度N为4的整数幂，即N=22L， 其精确复乘次数为（考虑省去非乘法的次数）：
时不变系统可以实现输入信号与系统冲激 响应之间的线性卷积。
两个序列卷积的傅里叶变换等于它们的傅 里叶变换相乘，即可在频域计算卷积。
频域取样导致信号时域周期重复，理论上 利用DFT只能计算圆周（循环）卷积，而不 能实现线性卷积。
圆周卷积等价于线性卷积的条件
N L K 1
其中，N为DFT的长度，序列x(n)的 长度为L，序列h(n)的长度为K。 为了利用DFT计算线性卷积，必须 选择满足上式的DFT长度N，并对 x(n)至少补K-1个零，对h(n)至少补 L-1个零。
FIR滤波器
数字滤波是数字信号处理中的基本处 理算法。
FIR（finite impulse response）滤波 器是最常用的一种数字滤波器。
FIR滤波器的差分方程与Z变换
N 1
y(n) ak x(n k) k 0
其中，x(n)为输入，y(n)为输出，ak为滤波器系数。
H (Z )
1
N 1
X (k)e j2nk / N (n 0,1,
,N
1)
N k0
定义旋转因子: WN e j 2nk / N
DFT的一些重要结论
有限长离散时间信号的频域离散表述可对 傅里叶变换取样得到；
只有当傅里叶变换一个周期内[0，2π）的 取样数N大于等于原信号长度L时，表述才 有实用意义。
当信号x(n)长度L小于N时，则对x(n)补零 构成N点序列，再进行DFT运算。
数字信号处理方法与实现
贺知明 副教授
电子科技大学 四川•成都
常用DSP算法
DFT(离散付里叶变换) FFT算法 FIR滤波器 IIR滤波器 自适应滤波器
DFT(离散付里叶变换)
DFT的定义为:
N 1
X (k) x(n)e j2nk / N (k 0,1, , N 1) n0
x(n)
窗函数法——基本方法
FIR滤波器的窗函数设计法
设所需滤波器的理想响应Hd(ejw) H d (e jw ) hd (n)e jnw
hd
(n)
1
2
2
0
Hd
(e
jw
)e
jnw dw
对hd(n)作截尾处理
Q
Ha (e jw ) hd (n)e jnw nQ
Q为FIR滤波器 的阶数，Q越 大，近似程度 越高
1965年，Cooley和Tukey提出快速付里叶 变换（FFT），极大地提高了运算速度。
运算量对比
N2
N 2
log 2
N
基-2 FFT运算近 似的复乘次数
例：1024个采样点 直接DFT：需1048576次乘法； FFT：仅需5120次乘法。
注意
FFT是DFT的一种快速算法； FFT没有对DFT作任何近似，精度无任何损
采用循环寻址功能，使长度为N的缓冲区足 以存储N阶滤波器所需的输入数据。
运算量为N/2 ·log2N
抽取（decimation）
两种形式的蝶式运算
DIT（按时间抽取）
C A WNk B D A WNk B
DIF （按频率抽取）
C AB D ( A B)WNk
位反序排列（以8点FFT为例）
0 1 234 56 7 000 001 010 011 100 101 110 111
如果时间窗长度正好为整数倍个信号周期， 除一点有值外，其余全为0，变换结果没有 产生谱泄漏。
通常情况，时间窗不等于信号周期的整数 倍，频谱中，除中心值外，还有以1/f（倍 频程）缓慢衰减的旁瓣。
时间窗加权
为消除谱泄漏，加权，压低旁瓣。 时间窗本身的作用相当于宽度与之相等的
一个矩形窗函数的加权。
对h(n)补零，使长度变为N； 用长度为N的DFT计算每个序列块的圆周卷
积。 拼接构成输出信号。
两种方法的区别
在重叠相加法中，计算长度N+K-1的DFT， 并将结果相加；
在重叠保留法中，只需计算长度N的DFT， 并通过舍弃、保留手段拼接构成最后输出 结果。
DFT
FFT
FFT算法
DFT运算归结为一系列的乘加运算，谱的 每一个频率点都要对N个乘法项求和，对N 个谱点总共作N2次乘法，计算量大，不利 于实时处理。
FIR滤波器的设计过程
通过计算机编程设计FIR滤波器，常用 Matlab仿真软件进行。
根据任务要求，首先确定所需FIR滤波器阶 数，然后利用仿真软件计算滤波器系数， 即完成FIR滤波器的设计。
FIR滤波器的计算过程
输入数据存在RAM缓冲区中，N个系数存 在系数存储器中，用程序循环执行
N 1
y(n) h(k)x(n k) k 0
000 100 010 110 001 101 011 111 0 4 26 15 37
反向FFT（IFFT）
x(n)
1 N
N 1
X (k )WNnk
k 0
1 N
N 1
[X
k 0
* (k)WNnk ]*
1 N
{FFT [x *(k)]}*
可以使用同一计算程序（硬件结构）作 正反两个方向的计算。三个步骤：
将序列x(n)分解成长度为N的序列块xm(n)； 分别对h(n)和xm(n)补零，使长度变为N+K-1； 用N+K-1点DFT计算每块的圆周卷积； 将结果相加，其中序列块卷积结果ym(n)的后K-
1个取样和ym+1(n)的前K-1个取样重叠，故称重 叠相加法。
重叠保留法
将序列x(n)分解成长度为N的序列块xm(n)， 其中有K-1个取样与相邻序列重叠，第一个 序列块在开头补K-1个零；
FIR滤波器的窗函数设计法
对hd(n)截尾，实际上是对hd(n)乘上一个矩 形窗，这是FIR滤波器设计的最直接和简便 的方法，但采用矩形窗存在较大的Gibbis现 象（有限项近似所产生的振荡），这是由 于矩形窗函数的边缘截止得太突然所引起， 加合适的窗可减小。
通常可采用前面介绍的其它经典窗函数 （Hamming窗、Hanning窗、Blackman窗 等）。
3 8
N
log
2
N
N
1
基-4算法比基-2算法更高效。
谱泄漏
DFT是对有限时间间隔（时间窗）里的采样 数据的变换，导致谱泄漏。
单位矩形序列r(n)乘以无限正弦波f(n)
f '(n) f (n)r(n) 单位矩形序列的
F '(k) F (k) R(k) 谱，SinC函数
无限正弦波谱，δ函数
说明
频域系数通常为固定系数，可提前算好存 起来，实时处理时不必再算，可省一次FFT 运算。
FFT在DSP实现过程中应注意的问题
DSP作FFT运算时，应尽量使用片内RAM， 速度快，并利用位反转（位反序）寻址方 式——专门为FFT运算提供的DSP功能。
定点DSP实现FFT，应注意防止中间结果溢 出，其方法是对中间数值归一化（以因子2 进行归一化，可运用DSP芯片的移位功能， 不会增加运算量），如作2m点FFT，可以2m 作归一化。
其它窗函数形式
矩形窗 三角窗 布莱克曼窗 凯塞窗 切比雪夫窗
FFT用于快速卷积
时域直接卷积 频域处理方法
实用中，这两种方法究竟哪一种速度快， 取决于输入序列的性质，通常采用FFT的频 域方法较快，但存在固有的处理延时。
FFT在DSP实现过程中应注意的问题
应调整信号的采样点数，使之正好或接近 为2（或4，由基数定）的幂次，尽量避免 大量补零的冗余。
对N个点的一次循环计算产生一个输出 y(n)，也就是说，FIR滤波计算就是循环卷积 计算，即移动的加权平均处理。
用前面方法计算得到的一组系数h(n) （即滤波器的设计），即可用于FIR计算。
FIR滤波器的DSP实现
注意输入数据x(n)与滤波器系数h(n)的合理 存放。
数据存储区
程序存储区
均为连续的存储区
选择窗函数的原则
使信号在窗的边缘为0，以减小截断所产生 的不连续效应。
信号经窗函数加权处理后，不应丢失太多 的信息。
选择什么样的窗函数，总是在旁瓣的衰减 程度、中心峰的宽度和能量损失等因素间 取折衷（一组矛盾）。
常用加权函数的一般形式
w(n) (1 ) cos ( 2n)
M
β= 2，α=0.5时，称为汉宁（Hanning）窗； β= 2，α=0.54时，称为汉明（Hamming）窗。
长序列与短序列卷积的DFT实现方法
在大多数实际情况中，长序列对应于系统 输入，短序列对应于系统冲激响应。
将长序列分割成长度为N的序列块，计算每 一块与短序列的卷积。
每一个序列块的卷积必须适当地合并，以 得到长序列和短序列卷积的最后结果。
有两种合并的方法，分别是重叠相加法和 重叠保留ห้องสมุดไป่ตู้。
重叠相加法
输入
h(N-1) x h(N-2) x
…...
h(0) x
+
Z-1
+
…...
Z-1
+
转置型横向滤波器
输出
FIR滤波器的特点
精确严格的线性相位 FIR滤波器是稳定系统 非递归结构 可采用FFT算法，提高运算效率
FIR滤波器的设计方法
频率抽样法——FIR滤波器的逼近方法
设计一个长度为N的FIR滤波器，可以对要求 达到的频率响应在N个点取样，把这些取样值进 行反傅氏变换，即可得到所要设计的滤波器的冲 激响应。
Y (Z )
N 1
h(k)Z k
X (Z ) k0
展开，得
H (Z ) h(0) Z 1 h(1) Z 1(h(2) h(N 1))
标准FIR滤波器结构
输入
Z-1
Z-1 …...
Z-1
h(0) x h(1) x …... h(N-1) x
+
输出
标准横向滤波器结构
转置FIR滤波器结构
DFT的一些重要结论
N越大，则它的DFT与傅里叶变换越近似， 因为在区间[0，2π）的取样数增加了。
N的取值通常需要根据实际应用中允许的运 算复杂度决定，因为补零个数越多，则DFT 计算所需要的运算和存储器件越多。
DFT实现数字滤波器
线性和圆周卷积 重叠相加法 重叠保留法
线性和圆周（循环）卷积