采用频率采样法的FIR滤波器要点Word文件下载.docx
《采用频率采样法的FIR滤波器要点Word文件下载.docx》由会员分享,可在线阅读,更多相关《采用频率采样法的FIR滤波器要点Word文件下载.docx(18页珍藏版)》请在冰豆网上搜索。
2.1FIR滤波器的设计
FIR(FiniteImplseResponse[3]滤波器:
有限长单位冲激响应滤波器,又称为非递归型滤波器,是数字信号处理系统中最基本的元件,他可以在保证任意幅频特性的同时具有严格的线性相频特性,同时其单位抽样响应是有限长的,因而滤波器是稳定的系统。
因此,FIR滤波器在通信、图像处理、模式识别等领域都有着广泛的应用。
有限长单位冲激响应(FIR)滤波器有以下特点:
(1)系统的单位冲激响应h(n)在有限个n值处不为0;
(2)系统函数H(z)在|z|>
0处收敛,极点全部在z=0处(因果系统);
(3)结构上主要是非递归结构,没有输出到输入的反馈,但有些结构中(例如频率抽样结构)也包含有反馈的递归部分。
2.2频率采样型结构
把一个有限长序列(长度为N点)的z变换H(z)在单位圆上作N等分抽样,就得
到H(k),其主值序列就等于h(n)的离散傅里叶变换H(k)。
那里也说到用H(k)表示的H(z)的内插公式为
(2-1)
H(k)
1-WnV1
这个公式就为FIR滤波器提供了另外一种结构,这种结构由两部分级联组成
第二部分由N各谐振器组成的谐振柜
它是由N个一阶网络并联而成,而这每一个一阶网络都是一个谐振器
其结构如下图所示:
频率抽样型结构特点:
(1)它的系数H(k)直接就是在滤波器在Wkk处得频率响应。
因此,控制得
N
频率响应是很直接得。
(2)结构有两个主要缺点:
a.所有的相乘系数及H(k)都是复数,应将它们先化成二阶实数,这样乘起来较复杂,增加乘法次数,存储量。
b.所有谐振器的极点都是在单位圆上,由Wn"
决定考虑到系数量化的影响,当系数量
化时,极点会移动,有些极点就不能被梳状滤波器的零点所抵消。
(零点由延时单元决定,不受量化的影响)系统就不稳定了。
(3)将一阶网络合并为二阶网络:
a.第k和第N-k个谐振器合并为一个实系数的二阶网络,因为h(n)是实数,他的DFT
也是圆周共轭对称的。
H(k)二H*(N-k)k=1,2,3,N-1(2-5)
因此,可以将第k和第N-k个谐振器合并为一个二阶网络。
(2-6)
H(k)九邨以-1
1-z」[w3Wf]r『wEw^z』〔—z'
2rcos(2k)r2z
其中:
5二?
Re[H(k)],-1k=-2rRe[H(k)WN]
b.第k和第N-k个谐振器合并为一个二阶网络的极点在单位圆内,而不是在单位圆
上,因而从频率响应的几何解释可知,它相当于一个有限Q的谐振器。
其谐振频率为
2二wkk。
r3
“才]J
▼
L
—r
1■
图2-3二阶网络结构图
除了共轭复根外,还有实根。
当N=偶数时,有一对实根,它们分别是k=0,k=N/2两个点
当N=奇数时,只有一个实根z=r(k=0),即只有H0(z)
c.修正频率抽样结构流图(N=偶数)
d.
2办
Zrcus(-)
V
<
尸22齐叽
2.3频率采样法
2.3.1设计思路:
这种设计方法是从频域进行设计的一种方法,首先给定一个希望逼近的频率响应
知道H(k)后,由IDFT定义,可以用这N个采样值H(k)来唯一确定有限长序列h(n),即
4N」
h(n)-1'
H(k)WN』kn=0,1,2...,N_1(2-11)
Nk_Q
NJ
H(z)='
h(n)z』(2-12)
n=0
N」
H(ej)=h(n)ejn(2-13)
k-Q
内插公式:
:
sin(人/2~2“2(2-15)
'
Nsin(o/2)
四种线性相位的FIR滤波器如下表2-1所示:
表2-1四种线性相位的FIR滤波器
爭㈣
1型
鮒N为奇数
同⑴)关于5=0、JU加偶对称
0
2型
凰"
=加MW),N为偶数
耳关于u>
=0、2x偶对8?
关于(o-ju奇对称
2
3型
烈jg占Ml迪』N炖数
耳⑴)关于g加奇对称
y-i打
4型
如砂g>
M对禹数
预⑴)关于2兀奇对称,关于斫禺对称
22
2.3.2逼近误差及其改进措施:
(1)采样点上滤波器的实际响应是严格地和理想频率响应数值相等的;
(2)在采样点之间的频率响应则是由各采样点的加权内插函数的延伸叠加而成的,因而有一定的逼近误差,误差大小取决于理想频率响应曲线形状。
(3)理想频率响应特性变化越平缓,则内插值越接近理想值,逼近误差越小;
(4)如果采样点之间的理想频率特性变化越陡,则内插值与理想值的误差就越大,因而在理想频率特性的不连续点附近,就会产生尖峰和起伏。
2.3.3滤波器性能的改善:
(1)增加过滤带采样点,它可以大大减少震荡,阻带衰减也可以得到进一步改善。
一般一点到二点的过滤带采样即可得到满意的结果。
(2)增加采样点密度,过渡带的宽度与采样点数N成反比。
但N值意味着或长度的增加,滤波器运算量必然增大[4]。
3.1设计流程图
3设计步骤
图3-1流程图
3.2音频信号获取
在网上下载一段吹管乐,并截取适当片段。
把音频放到电脑录音机里面,把音频属性
设置为8000Hz。
如图3-2所示:
图3-2音频信号设置
然后在MATLAB软件平台下,利用函数wavread对音频信号进行采样,源程序为:
[x,fs,bits]=wavread('
梅花三弄.wav'
),记住采样频率和采样点数,MATLAB实现得:
fs=8000;
bits=8。
3.3音频信号的频谱分析
梅花三弄.wav'
);
%输入参数为文件的全路径和文件
名,输出的第一个参数是每个样本的值,fs是生成该波形文件时的采样率,bits是波形文件每样本的编码位数。
sound(x,fs,bits);
%按指定的采样率和每样本编码位数回放
n=length(x);
%计算信号x的长度
fn=3200;
%单频噪声频率,此参数可改
t=0:
1/fs:
(n-1)/fs;
%计算时间范围,样本数除以采样频率
X=X(:
1)'
;
%将双声道转为单声道
y=x+0.1*sin(fn*2*pi*t);
%加噪声
sound(y,fs,bits);
%应该可以明显听出有尖锐的单频啸叫声
X=abs(fft(x));
Y=abs(fft(y));
%对原始信号和加噪信号进行fft变换,取幅度谱
X=X(1:
n/2);
Y=Y(1:
%截取前半部分
deltaf=fs/n;
%计算频谱的谱线间隔
f=0:
deltaf:
fs/2-deltaf;
%计算频谱频率范围
运行结果如下图:
图3-3加噪前后时域频域对比图
由图3-3可以看出,在频域为3200Hz处加入一个单频噪声,而加入噪声之后,时域
的波形图出现了明显失真,通过听取原音频信号x和加噪声音频信号y,可以明显听到y
音频信号中有一明显尖锐噪声。
3.4滤波器设计
设计一个低通滤波器,将单频信号滤出去,源程序如下所示:
fcd=2000;
fcu=2500;
%定义设计指标,通带截止频率fed和阻带
截止频率feu
T1=0.6025,T2=0.127;
Ns=fix((fcu-fcd)*2*pi/(2*pi*N));
Np=ceil(N-Ns);
计算通带内样本数
Ak=[ones(1,Np),zeros(1,2*Ns+1),ones(1,Np-1)];
%得到理想滤波器幅度特性采样值
Ak(Np+1)=T;
Ak(2*N-Np+1)=T;
%设置过渡带样本为最优值
M=length(Ak);
alpha=(M-1)/2;
k1=0:
floor((M-1)/2);
k2=floor((M-1)/2)+1:
M-1;
angH=[-alpha*(2*pi)/M*k1,alpha*(2*pi)/M*(M-k2)];
%计算分段线性相位
Hk=Ak.*exp(j*angH);
%幅度和相位相乘得实际滤波器频率特性Hk
hn=real(ifft(Hk));
%对HK作IFFT得到脉冲响应hn
[db,mag,pha,grd,w]=freq乙m(hn,1);
%计算hn对应系统的频率响应
%作图检查幅度特性是否满足预定指标
程序运行结果如下图3-4所示:
图3-4滤波器参数图
由图3-4可以看出,滤波器的衰减大于设定值As=40,满足性能指标,滤波器的衰减
可由增加过度带宽来得到。
3.5信号滤波处理
源程序如下所示:
y_fil=fftfilt(hn,y);
Y_fi匸fft(y_fil);
Y_fi匸abs(Y_fil);
Y_fil=Y_fil(1:
程序运行结果如下图3-5所示:
图3-5滤波前后时域、频域对比图
频谱取前一半由上图3-5可以看出,在滤波之后时域图能得到恢复,频域图中的单频噪声信号也得到滤除,说明了设计的滤波器能滤除加入的噪声信号,因此说达到了设计要求。
3.6结果分析
开始通过分析决定设计一个带阻的滤波器来滤除加入的单频噪声,根据噪声的频率来设计阻带的范围。
在收集音频信号后,按照步骤用频率采样法设计频率采样型滤波器。
由图3-4可知,
设计的滤波器达到要求。
我们观察到图3-5滤波前后音频型号的波形对比图,发现时域波
形中加干扰噪声后有明显变化,不过经过滤波后几乎没有变化,说明设计的滤波器达到要求。
再通过听取原始语音信号,加噪信号y和滤波之后的信号y_fil。
对比之后,发现滤波器确实滤除了噪声。
从理想的角度考虑,该带阻滤波器的阻带带宽应该可以变得更窄,让滤波效果更好,但是这样的采样值会变得非常大。
考虑到实际的情况,权衡之后,决定牺牲带宽来使得滤波器的阶数降低,因此在上图3-5中我们可以看出,在噪声频谱左右两边的信号也被滤除了。
4出现的问题及解决办法
在这次课程设计当中,由于基础不扎实,出现了很多问题,即MATLAB软件操作不
当,也有知识掌握程度不够出现的各种问题。
1.在一开始找音频后,没有修改参数值,导致频率抽样过高,在老师提醒下,改为8000Hz;
2.在调用音乐文件时,没有将文件放在MATLAB的工作文件夹下面,到时文件找不
到。
还有有程序中用到的各种函数没有放入正确的位置;
3.在绘制加噪前后频率对比图时,留白过多,对比不明显,通过axis函数对横纵坐标
进行限定;
4.在设计滤波器的过程中,没有准确理解每一条指令代表的含义,导致程序前后不对应,出现更多错误,通过MATLAB中的错误提示,准确找到错误的那一行代码,进行修改;
5.在设计滤波器的过程中,滤波器的衰减小于开始所设置的值,通过牺牲过渡带和调节过渡带的采样值,即T1和T2来使得衰减大于所设定值AS;
6.最后听取滤波后声音过程当中,没有将其保存,请教同学之后,学会了如何保存滤波后的声音文件。
5结束语
这是第二次课程设计,在前面的课程设计当中我们学习到如何使用MATLAB,所以对于MATLAB软件的使用并没有那么陌生了,尽管如此,在使用MATLAB的过程中还是出现了很多错误,比如说忘记添加函数文件,参数前后不对应之类的低级错误。
此次课程设计,让我更深入的了解到频率抽样法以及频率抽样型的滤波器,由开始的无从下手,再翻阅书上的例子。
首先决定使用低通的滤波器将噪声滤除,但由于是加入的是一个单频的信号噪声,就决定使用带阻滤波器的例子,所以就通过书籍《信号与系统》上高通的例题来进行修改,当然这过程不是一蹴而就的。
经历过一次又一次的错误,猜得出来最后的模型。
还有这也是我和其他同学一起讨论来的。
这次课设让我对滤波器的类型有了一个更加完整了解,在设计中也使我对一些概念有了更深刻的认识。
比如:
在滤波器分类方面,我深刻的了解了低通,高通滤波器与带通,带阻滤波器的特性区别。
还有在课程设计中每一次的数据输入都有其重要意义,用MATLAB编译程序时,可以根据滤波器指标的要求实时知道对滤波器的影响。
通过一次次
的调试和权衡使滤波器的性能达到最佳。
课程设计不仅要求对滤波器理论的研究,更重要的是培养一种遇到问题解决问题的思维。
因为有了这次课程设计,我懂得了书本知识只是实际应用的理论指导。
如果仅仅只学习书本知识,不去在实践中运用,那只是停留在知识表面,不知其因的层面。
比如在数学计算上,可以将噪声完全滤除。
而在这次课设中,若要完全滤除噪音,滤波器的阶数就会增高,在现实生活中是很难实现的,所以噪声是不能完全滤除的。
课程设计结束了,我相信这次课程设计对今后的学习是很有帮助的,它让我将理论更好的和实践相结合,提高了动手的能力,也填补了自己学习上的一些不足。
这次课设的成功,不仅仅是我一个人的努力的结果,更离不开指导老师与同学的帮助,在此向老师和同学们表示衷心的感谢。
参考文献
[1]吴镇扬,数字信号处理[M],高等教育出版社,2004
[2]张圣勤,MATLAB7.C实用教程[M],北京:
机械工程出版社,2006
[3]程佩青,数字信号处理教程[M],北京:
清华大学出版社,2002
1994
[4]高西全,丁玉美,数字信号处理[M],第三版,西安:
西安科大出版社,
附录:
源程序
%输入参数为文件的全路径和文件名,输出的第一个参数是每个样本的值,fs是生成该波形文件时的采样率,bits是波形文件每样本的编码位数。
sound(x,fs,bits);
(n-1)/fs;
x=x(:
;
%将双声道转为单声道
sound(y,fs,bits);
%对原始信号和加噪信号进行fft变换,取幅度谱X=X(1:
%计算频谱的谱线间隔
%计算频谱频率范围
figure
(1)subplot(2,2,1);
plot(t,x);
axis([0,3,-1.3,1.3]);
gridonxlabel('
时间(单位:
s)'
ylabel('
幅度'
title('
原始音乐信号'
subplot(2,2,2);
plot(f,X(1:
n/2));
axis([0,5000,0,4000]);
频率(单位:
Hz)'
幅度谱'
音乐信号幅度谱图'
subplot(2,2,3);
plot(t,y);
gridon
xlabel('
加入单频干扰的音乐信号'
subplot(2,2,4);
plot(f,Y);
频率(单位:
Hz)'
加入干扰后的音乐信号幅度谱图’);
fcd=2000;
%定义设计指标,通带截止频率fcd和阻带截止频率fcuRp=0.5;
As=40;
wp=fcd/fs*2*pi;
ws=fcu/fs*2*pi;
%将模拟指标转换为数字指标
delta_w=2*pi/1000;
Rp=-(min(db(1:
1:
wp/delta_w+1)));
As=-round(max(db(ws/delta_w+1:
501)));
%设置自由样本数m
%计算通带内样本数
Ak=[ones(1,Np),zeros(1,2*Ns+1),ones(1,Np-1)];
k1=0:
angH=[-alpha*(2*pi)/M*k1,alpha*(2*pi)/M*(M-k2)];
Hk=Ak.*exp(j*angH);
%幅度和相位相乘得实际滤波器频率特性Hk
%对HK作IFFT得到脉冲响应hn
figure
(2)
subplot(2,2,1);
plot(w/pi,db)
gridon;
w/pi'
db'
滤波器db'
X_l=[0,0,wp/pi,ws/pi;
1,1,wp/pi,ws/pi];
Y_l=[-As,-Rp,-100,-100;
-As,-Rp,0,0];
%在wp,ws,Rp,As处画线以更直观判断设计是否达标,每列参数是每个线条的端点坐
标
line(X」,Y_l,'
Color'
'
r'
LineWidth'
2,'
LineStyle'
--'
)%添加线条,红色,线宽为2
axis([01-805])
plot(w/pi,mag);
mag'
滤波器幅度响应'
axis([0101.2])
plot(w/pi,pha);
pha'
滤波器相位响应'
axis([01-44])
subplot(2,2,4);
stem(l,hn);
n'
y(n)'
滤波器单位冲击响应'
axis([040-11])
Y_fil二fft(y_fil);
Y_fil二abs(Y_fil);
figure(3)
subplot(3,2,1);
时间t'
原始语音信号'
axistight;
subplot(3,2,2);
plot(f,X);
频率f);
原始语音信号幅度谱'
axis([0500004000]);
subplot(3,2,3);
加干扰后的语音信号'
subplot(3,2,4);
加干扰后的语音信号幅度谱
subplot(3,2,5);
plot(t,y_fil);
滤波后语音信号时间'
axis([04-1.51.5]);
subplot(3,2,6);
plot(f(1:
n/2),Y_fil);
滤波后语音信号的信号频谱图'