没有合适的资源?快使用搜索试试~ 我知道了~
资源推荐
资源详情
资源评论
第 7 章 FIR 滤波器设计
第六章我们介绍了无限冲激响应(IIR)滤波器的设计方法。其中最常用的由模拟滤波
器转换为数字滤波器的方法为双线性变换法,因为这种方法无混叠效应,效果较好。但通
过前面的例子我们看到,IIR 数字滤波器相位特性不好(非线性,如图 6-11、图 6-13、图
6-15 ),也不易控制。然而在现代信号处理中,例如图像处理、数据传输、雷达接收以及
一些要求较高的系统中对相位特性要求较为严格,这种滤波器就无能为力了。改善相位特
性的方法是采用有限冲激响应滤波器。本章首先对 FIR 滤波器原理及其使用函数作基本介
绍,然后重点介绍窗函数法设计 FIR 滤波器,并对最优滤波器设计函数进行介绍。
7.1 FIR 滤波器原理概述及滤波函数
7.1.1 FIR 滤波器原理及设计方法分类
根据第 6 章对数字滤波器的介绍,我们知道 FIR 滤波器的传递函数为:
1
0
)(
)(
)(
)(
N
n
n
znh
zX
zY
zH
(7-1)
可得 FIR 滤波器的系统差分方程为:
1
0
)()()()(
)1()1()1()1()()0()(
N
m
nxnbmnxmb
NnxNbnxbnxbny
因此,FIR 滤波器又称为卷积滤波器。根据第 4 章中所描述的系统频率响应,FIR 滤波
器的频率响应表达式为:
1
0
)(
N
n
jnj
enbeH
(7-2)
信号通过 FIR 滤波器不失真条件与(6-6)式所描述的相同,即滤波器在通带内具有恒
定的幅频特性和线性相位特性。理论上可以证明(这里从略):当 FIR 滤波器的系数满足
下列中心对称条件:
)1()()1()( nNbnbnNbnb 或
(7-3)
时,滤波器设计在逼近平直幅频特性的同时,还能获得严格的线性相位特性。线性相位
FIR 滤波器的相位滞后和群延迟在整个频带上是相等且不变的。对于一个 N 阶的线性相位
FIR 滤波器,群延迟为常数,即滤波后的信号简单地延迟常数个时间步长。这一特性使通
带频率内信号通过滤波器后仍保持原有波形形状而无相位失真。
本章主要介绍的 FIR 数字滤波器设计方法及 MATLAB 信号处理工具箱提供的 FIR 数
字滤波器设计函数,见表 7-1。由于篇幅所限,本章我们主要介绍窗函数法和最优化设计方
法。
表 7-1 FIR 滤波器设计的主要方法
函数设计方法
说明 工具函数
窗函数法 理想滤波器加窗处理 fir1(单频带) , fir2(多频
带) , kaiserord
最优化设计 平方误差最小化逼近理
想 幅 频 响 应 或 Park-
McClellan 算法产生等波纹滤
波器
firls , remez,remezord
约束最小二乘逼近 在满足最大误差限制条
件下使整个频带平方误差最
小化
fircls,fircls1
升余弦函数 具有光滑、正弦过渡带
的低通滤波器设计
Fircos
7.1.2 FIR 数字滤波器滤波函数
相对于 IIR 滤波器的滤波函数,FIR 数字滤波器滤波函数除了 dimpulse 和 dstep 仅
适用于 IIR 滤波器外,其他各种函数可直接应用于 FIR 滤波器,只是输入的分母多项式向
量 a=1。另外,MATLAB 还提供了一个函数 tlt,该函数利用效率高的基于 FFT 算法实
现对数据的滤波,该函数只适用于 FIR 滤波器,调用形式为:
y=tlt(b,x[,n])
式中,b 为 FIR 滤波器的系数向量;x 为输入数据;n 为 FFT 长度,缺省时,函数选用最
佳的 FFT 长度,y 为滤波器的输出。该函数执行下面的操作:
n=length(x);
y=it(t(x).*t(b,n)./t(a,n));
应注意,y=fftfilt(b,x)等价于 y=filter(b,a,x)。
7.2 FIR 滤波器的窗函数设计
7.2.1 窗函数的基本原理
FIR 滤波器设计的主要任务是根据给定的性能指标确定滤波器的系数 b,即系统单位
脉冲序列 h(n),它是一个有限长序列。
FIR 滤波器的理想频率响应,可写成复数形式的 Fourier 级数形式:
n
nj
d
j
d
enheH
(7-4)
式中,h
d
(n)是对应的单位脉冲响应序列。这说明滤波器的频率响应和单位脉冲响应互为
Fourier 变换对。因此其单位脉冲响应可由下式求得,
deeHnh
njj
dd
2
1
(7-5)
求得序列
nh
d
后,通过 z 变换,可得到
zH
d
n
n
dd
znhzH )()(
(7-6)
注意,这里
nh
d
为无限长序列,因此
zH
d
是物理上不可实现的。如何变成物理上
可实现呢?一个自然的想法是只取其中的某些项,即只截取
nh
d
中的一部分,比如 n=0,
…,N-1,N 为正整数。这种处理相当于将
nh
d
,n=-∞~∞与函数 w(n)相乘,w(n)具有下
列形式:
Nn
Nnn
nw
0,1
,0,0
)(
w(n)相当于一个矩形,我们称之为矩形窗。即我们可采用矩形窗函数 w(n)将无限脉冲响
应
nh
d
截取一段 h(n)来近似为
nh
d
,这种截取在数学上表示为:
h(n)=
nh
d
w(n) (7-7)
这里应该强调的是,加窗函数不是可有可无的,而是将设计变为物理可实现所必须的。
截取之后的滤波器传递函数变为:
1
0
)()(
N
n
n
znhzH
(7-8)
式中,N 为窗口宽度,H(z)是物理可实现系统。
为了获得线性相位,FIR 滤波器 h(n)必须满足中心对称条件(即 7-3 式),序列 h(n)
的延迟为
2/1 N
。
这种方法的基本原理是用一定宽度的矩形窗函数截取无限脉冲响应序列获得有限长的
脉冲响应序列,从而得到 FIR 滤波器的脉冲响应,故称为 FIR 滤波器的窗函数设计法。
经过加矩形窗后所得的滤波器实际频率响应能否很好地逼近理想频率响应呢?图 7-1
示意给出了理想滤波器加矩形窗后的情况。理想低通滤波器的频率响应如图中左上角图,
矩形窗的频率响应为左下角图。时间域内的乘积(7-7)式要求实际频率响应为这两个频
率响应函数在频域内的卷积(卷积定理),即得到图形为图 7-1(右图)。
图 7-1 FIR 滤波器理想与实际频率响应
由图可看出,加矩形窗后使实际频率响应偏离理想频率响应,主要影响有三个方面:
(1)理想幅频特性陡直边缘处形成过渡带,过渡带宽取决于矩形窗函数频率响应的
主瓣宽度。
(2)过渡带两侧形成肩峰和波纹,这是矩形窗函数频率响应的旁瓣引起的,旁瓣相
对值越大,旁瓣越多,波纹越多。
(3)随窗函数宽度 N 的增大,矩形窗函数频率响应的主瓣宽度减小,但不改变旁瓣
的相对值。
为了改善 FIR 滤波器性能,要求窗函数的主瓣宽度尽可能窄,以获得较窄的过渡带;
旁瓣相对值尽可能小,数量尽可能少,以获得通带波纹小,阻带衰减大,在通带和阻带内均
平稳的特点,这样可使滤波器实际频率响应更好地逼近理想频率响应。
这里我们明确两个概念:截断和频谱泄漏。信号是无限长的,而在进行信号处理时只
能采取有限长信号,所以需要将信号“截断”。在信号处理中, “截断”被看成是用一个有限
长的“窗口”看无限长的信号,或者从分析的角度是无限长的信号 x(t)乘以有限长的窗函数
w(t)。由傅立叶变换性质可知,时间域内的乘积对应于频率域的卷积,即
(7-9)
这里,x(t)是频宽有限信号,而 w(t)是频宽无限信号,
表示互为 Fourier 变换对。截断
后的信号也必须是频宽无限信号,这样就是有限频带的信号分散到无限频带中去,这样就
产生了所谓频谱泄漏。从能量的角度来看,频谱泄漏也是能量的泄漏,因为加窗后使原来
信号集中的窄频带内的能量分散到无限的频带宽度范围内。频谱泄漏是不可避免的,但要
尽量减小。
上边只考虑了矩形窗,如果我们使窗的主瓣宽度尽可能地窄,旁瓣尽可能地小,可以
获得性能更好的滤波器,能否改变窗的形状而达到这个目的呢?回答是肯定的。其实数字
信号处理的前驱者们设计了不同于矩形窗的很多窗函数,这些窗函数在主瓣和旁瓣特性方
面各有特点,可满足不同的要求。为此,用窗函数法设计 FIR 数字滤波器时,要根据给定
的滤波器性能指标选择窗口宽度 N 和窗函数 w(n)。下面我们介绍窗函数。
7.2.2 MATLAB 信号处理中提供的窗函数
(1)矩形窗:前面分析中所用的矩形窗可用下面函数来实现 w=boxcar (N),N 为
窗的长度(以下函数与此同),w 为返回的窗函数序列。
(2)汉宁窗:w=hanning(N)
汉宁窗的表达式为:
(7-10)
(3)哈明窗:w=hamming(N)
哈明窗的表达式为:
(7-11)
(4)Bartlett 窗:w=bartlett(N)
Bartlett 窗的表达式为:
当 N 为奇数时,
(7-12)
当 N 为偶数时,
(7-13)
(5) Blackman 窗:w= blackman(N)
Blackman 窗的表达式为:
(7-14)
Blackman 窗比其他相同尺寸窗 (哈明窗,汉宁窗) 具有主瓣较宽和旁瓣泄漏较小的特
点。
(6)三角窗:w=triang(N)
三角窗的表达式为:
当 N 为奇数时,
(7-15)
当 N 为偶数时,
(7-16)
三角窗和 Bartlett 窗十分类似。三角窗的两端值不为零,而 Bartlett 窗则为零,这一
点可从例 7-1 中看出。
(7)Kaiser 窗:w=kaiser(n,beta)
其中,beta 是 Kaiser 窗参数,影响窗旁瓣幅值的衰减率。
Kaiser 窗表达式:
(7-17)
式中, I
0
[.]是修正过的零阶 Bessel 函数。
Kaiser 窗用于滤波器设计时,若旁瓣幅值为
,则
( 7-18 )
(8) Chebyshev 窗:w=chebwin(n,r)
式中, r 是窗口的旁瓣幅值在主瓣以下的分贝数。
Chebyshev 窗具有主瓣宽度最小,而旁瓣等高,高度可调整的特点。
下面我们在 MATLAB 观看各种窗函数的形状和频率域图象来验证上述所讲特点。
【例 7-1】 用 MATLAB 编程绘制各种窗函数的形状。窗函数的长度为 21。
%Samp7_1
clf
Nwin=21;n=0:Nwin-1; %数据总数和序列序号
gure(1)
for ii=1:4
switch ii
case 1
w=boxcar(Nwin); %矩形窗
stext='矩形窗';
case 2
w=hanning(Nwin); %汉宁窗
stext='汉宁窗';
case 3
w=hamming(Nwin); %哈明窗
stext='哈明窗';
case 4
w=bartlett(Nwin); %Bartlett 窗
stext='Bartelett 窗';
end
posplot=['2,2,' int2str(ii)]; %指定绘制窗函数的图形位置
subplot(posplot);
stem(n,w); %绘出窗函数
hold on
plot (n ,w,'r'); %绘制包络线
xlabel('n'); ylabel('w(n)'); title(stext);
hold o; grid on
end
gure(2)
for ii=1:4
switch ii
case 1
w=blackman(Nwin); %Blackman 窗
stext='Blackman 窗';
case 2
w=triang(Nwin); %三角窗
stext='三角窗';
case 3
w=kaiser(Nwin,4); %Kaiser 窗
stext='Kaiser 窗(Beta=4)';
case 4
w=chebwin(Nwin,40); %Chebyshev 窗
stext='Chebyshev 窗(r=40)';
end
posplot=['2,2,' int2str(ii)];
subplot(posplot);
stem(n,w); %绘出窗函数
hold on
plot (n,w,'r'); %绘出包络线
xlabel('n');ylabel('w(n)');title(stext);
hold o;grid on;
end
程序运行结果见图 7-2 。可以看到各种窗函数的形状。
剩余25页未读,继续阅读
资源评论
- bilei1177391952013-03-12方法和原理介绍得很清楚,谢谢!
hanshuiyiheng
- 粉丝: 5
- 资源: 32
上传资源 快速赚钱
- 我的内容管理 展开
- 我的资源 快来上传第一个资源
- 我的收益 登录查看自己的收益
- 我的积分 登录查看自己的积分
- 我的C币 登录后查看C币余额
- 我的收藏
- 我的下载
- 下载帮助
安全验证
文档复制为VIP权益,开通VIP直接复制
信息提交成功