基于频率采样法FIR带通滤波器设计.docx

基于频率采样法FIR带通滤波器设计.docx

《基于频率采样法FIR带通滤波器设计.docx》由会员分享，可在线阅读，更多相关《基于频率采样法FIR带通滤波器设计.docx（11页珍藏版）》请在冰点文库上搜索。

基于频率采样法FIR带通滤波器设计

**大学

电子信息工程学院

DSP课程设计报告

题目：

基于频率采样法的FIR带通滤波器的设计

专业班级：

通信工程

组员：

***

指导老师：

***

2011年6月13日

一、设计目的………………………………………………3

二、设计要求………………………………………………3

三、设计原理………………………………………………4

四、程序代码………………………………………………6

五、调试分析………………………………………………12

六、设计总结………………………………………………13

七、设计心得………………………………………………13

八、参考文献………………………………………………

一、设计目的：

1、掌握用频率取样法设计FIRDF的方法，并掌握该方法的MATLAB编程。

2、熟悉频率取样理论，熟悉内插函数及其应用。

3、了解FIRDF的频率特性和相位特性，观察过渡带取样点对滤波器幅频特性的影响。

二、设计内容:

基于MATLAB结合FFT和IFFT，利用频率采样法设计FIR数字带通滤波器，然后用自己设计的滤波器对采集的加噪后的语音信号进行滤波，并将滤波前后的信号进行比较，回放语音信号。

三、设计原理：

FIR滤波器的单位脉冲响应

是有限长的（

），其z变换为

的（N-1）阶多项式：

可得FIR滤波器的系统差分方程：

因此，FIR滤波器又称为卷积滤波器。

FIR滤波器的频率响应表达式为：

信号通过FIR滤波器不失真条件是在通带内具有恒定的幅频特性和相位线性特性。

理论上可以证明：

当FIR滤波器的系数满足下列中心对称条件：

时，滤波器设计在逼近平直幅频特性的同时，还能获得严格的线性相位特性。

线性相位FIR滤波器的相位滞后和群延迟在整个频带上是相等且不变的。

对于一个N阶的线性相位FIR滤波器，群延迟为常数，即滤波后的信号简单地延迟常数个时间步长。

这一特性使通带频域内信号通过滤波器后仍保持原有波形形状而无相位失真。

FIR滤波器设计的方法很多，如窗函数法、频率采样法以及其他的各种优化设计方法，本次课程设计使用频率采样法设计FIR带通滤波器。

频率采样法是从频域出发。

因为有限长序列h（n）又可用其离散傅立叶变换H（k）来唯一表示，H（k）与所要求的FIR滤波器系统函数Hd（z）之间存在着频率取样关系。

即Hd（z）在Z平面单位圆上按角度等分的取样值等于Hd（k）的各相应值，就以此Hd（k）值作为实际FIR数字滤波器频率特性的取样值H（k）,或者说H（k）正是所要求的频率响应H（ejw）的N各等间隔的取样值。

频率采样法就是根据频域采样理论，由滤波特性指标构造希望逼近的滤波器频响函数Hd（ejω），对其在[0,2π]上采样得到:

然后，就可求出单位脉冲响应h（n），或是系统函数H（z）。

这样，h（n）或是H（z）就是滤波器的设计结果。

频率取样法设计的基本思想：

把给出的理想频率响应进行取样，通过IDFT从频谱样点直接求得有限脉冲响应。

其设计过程如下

频率取样法的关键是正确确定数字频域系统函数H（k）在Ω∈[0，2π]内的N个样点，其约束条件为

H（k）=H（N-k）

ϕ（m）=-ϕ（N-m）0≤k≤N-1

频率采样法的优点是可以在频域直接设计，并且适合最优化设计；缺点是采样频率只能等于

的整数倍，因而不能确保截止频率

的自由取值，要想实现自由地选择截止频率，必须增加采样点数N，但是这又使计算量加大。

四、程序代码:

语音程序：

filename='111';

[s,fs,nbits]=wavread（filename）;

sound（s,fs,nbits）;%回放语音信号

n=length（s）;%求出语音信号的长度

t=0:

1/fs:

（n-1）/fs;

Y=fft（s）;%傅里叶变换

figure

（1）

subplot（2,1,1）;plot（s）;title（'原始信号波形'）;gridon

subplot（2,1,2）;plot（abs（Y））;title（'原始信号频谱'）

sound（s,fs）;

gridon

仿真结果：

图一

带通滤波器程序：

N=40;

alfa=（40-1）/2;

k=0:

N-1;

w1=（2*pi/N）*k;

T1=0.109021;T2=0.59417456;

hrs=[zeros（1,5）,T1,T2,ones（1,7）,T2,T1,zeros（1,9）,T1,T2,ones（1,7）,T2,T1,zeros（1,4）];

hdr=[0,0,1,1,0,0];wd1=[0,0.2,0.35,0.65,0.8,1];

k1=0:

floor（（N-1）/2）;k2=floor（（N-1）/2）+1:

N-1;

angH=[-alfa*（2*pi）/N*k1,alfa*（2*pi/N*（N-k2））];

H=hrs.*exp（j*angH）;

h=real（ifft（H））;

[db,mag,pha,grd,w]=freqz_m（h,1）;

[Hr,ww,a,L]=Hr_Type2（h）;

subplot（2,2,1）

plot（w1（1:

21）/pi,hrs（1:

21）,'o',wd1,hdr）

axis（[0,1,-0.1,1.1]）;

title（'带通：

N=40，T1=0.109021，T2=0.59417456'）

ylabel（'Hr（k）'）;

set（gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1]）

set（gca,'YTickMode','manual','YTick',[0,0.059,0.109,1]）;

grid%绘制带网格的图像

subplot（2,2,2）;

stem（k,h）;

axis（[-1,N,-0.4,0.4]）

title（'脉冲响应'）;ylabel（'h（n）'）;text（N+1,-0.4,'n'）

subplot（2,2,3）;plot（ww/pi,Hr,w1（1:

21）/pi,hrs（1:

21）,'o'）;

axis（[0,1,-0.1,1.1]）;title（'振幅响应'）

xlabel（'频率（单位：

pi）'）;ylabel（'Hr（w）'）

set（gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1]）;

set（gca,'YTickMode','manual','YTick',[0,0.059,0.109,1]）;

grid

subplot（2,2,4）;plot（w/pi,db）;axis（[0,1,-100,10]）;

grid

title（'幅度响应'）;xlabel（'频率（单位：

pi）'）;ylabel（'分贝'）

set（gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1]）

set（gca,'YTickMode','manual','YTick',[-60;0]）;

set（gca,'YTickLabelMode','manual','YTickLabels',[60;0]）;

仿真结果：

图二

总程序：

%设计一个通带为0.35π到0.65π，阻带下截止为0.2π，上截止为0.8π

Clear;

N=40;

alfa=（40-1）/2;

k=0:

N-1;

w1=（2*pi/N）*k;

T1=0.109021;T2=0.59417456;

hrs=[zeros（1,5）,T1,T2,ones（1,7）,T2,T1,zeros（1,9）,T1,T2,ones（1,7）,T2,T1,zeros（1,4）];

hdr=[0,0,1,1,0,0];wd1=[0,0.2,0.35,0.65,0.8,1];

k1=0:

floor（（N-1）/2）;k2=floor（（N-1）/2）+1:

N-1;

angH=[-alfa*（2*pi）/N*k1,alfa*（2*pi/N*（N-k2））];

H=hrs.*exp（j*angH）;

h=real（ifft（H））;

[db,mag,pha,grd,w]=freqz_m（h,1）;

[Hr,ww,a,L]=Hr_Type2（h）;

figure

（1）

subplot（2,2,1）

plot（w1（1:

21）/pi,hrs（1:

21）,'o',wd1,hdr）

axis（[0,1,-0.1,1.1]）;

title（'带通：

N=40，T1=0.109021，T2=0.59417456'）

ylabel（'Hr（k）'）;

set（gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1]）

set（gca,'YTickMode','manual','YTick',[0,0.059,0.109,1]）;

grid%绘制带网格的图像

subplot（2,2,2）;

stem（k,h）;

axis（[-1,N,-0.4,0.4]）

title（'脉冲响应'）;ylabel（'h（n）'）;text（N+1,-0.4,'n'）

subplot（2,2,3）;plot（ww/pi,Hr,w1（1:

21）/pi,hrs（1:

21）,'o'）;

axis（[0,1,-0.1,1.1]）;title（'振幅响应'）

xlabel（'频率（单位：

pi）'）;ylabel（'Hr（w）'）

set（gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1]）;

set（gca,'YTickMode','manual','YTick',[0,0.059,0.109,1]）;

grid

subplot（2,2,4）;plot（w/pi,db）;axis（[0,1,-100,10]）;

grid

title（'幅度响应'）;xlabel（'频率（单位：

pi）'）;ylabel（'分贝'）

set（gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1]）

set（gca,'YTickMode','manual','YTick',[-60;0]）;

set（gca,'YTickLabelMode','manual','YTickLabels',[60;0]）;

[s,fs,nbits]=wavread（'sj.wav'）;%信号de取样频率为44100HZ

x=s（:

1）;

sound（x,fs）;

L=length（x）;

f=fs*（0:

L-1）/L;

t=0:

1/fs:

（L-1）/fs;%将所加噪声信号的点数调整到与原始信号相同

d=0.03*abs（max（x））*[cos（2*pi*22000*t）]';%噪声为500和3300Hz的余弦信号

dz=[cos（2*pi*11025*t）]';

xd=x.*dz;

xz=xd+d;

sound（xz,fs）;%播放加噪声后的语音信号

X=fft（x）;%求信号的频谱

XD=fft（xd）;%信号调制后的频谱

XZ=fft（xz）;

figure

（2）

subplot（3,1,1）;plot（t,x）

title（'未加噪的信号'）;xlabel（'times'）;ylabel（'幅度'）;

subplot（3,1,2）;plot（t,xd）

title（'调制后的信号'）;xlabel（'times'）;ylabel（'幅度'）;

subplot（3,1,3）;plot（t,xz）

title（'调制加噪后的信号'）;xlabel（'timen'）;ylabel（'fuzhin'）;

figure（3）

subplot（3,1,1）;

plot（f,abs（X））;

title（'原始语音信号频谱'）;

xlabel（'频率（单位：

Hz）'）;

ylabel（'幅度'）;

subplot（3,1,2）;

plot（f,abs（XD））;

title（'调制后的信号频谱'）;

xlabel（'频率（单位：

Hz）'）;

ylabel（'幅度'）;

subplot（3,1,3）;

plot（f,abs（XZ））;

title（'加噪后的信号频谱'）;

xlabel（'频率（单位：

Hz）'）;

ylabel（'幅度'）;

y=fftfilt（h,xd）;

Y=fft（y）;

sound（3*y,fs）;

figure（4）

subplot（2,1,1）;plot（t,3*y）

title（'滤波后的信号'）;xlabel（'times'）;ylabel（'幅度'）;

subplot（2,1,2）;plot（f,abs（Y））;

title（'滤波后的信号频谱'）;xlabel（'频率（单位：

Hz）'）;ylabel（'幅度'）;

仿真结果：

图三

图四

图五

图一为原始语音的时域和频域图：

语音信号采样频率为44100HZ；（注：

为了让信号不失真，采样频率大于信号频率的2倍。

）

图二为带通滤波器图：

通带为0.35π到0.65π，阻带下截止为0.2π，上截止为0.8π；

图三、图四分别为未加噪声和调制加入噪声后的语音信号时域图和频域图：

d=0.03*abs（max（x））*[cos（2*pi*22000*t）]';噪声为22000Hz的余弦信号

图五为经过滤波器滤除噪声之后的信号时域和频域波形。

五、调试分析

1、语音部分：

我们先利用电脑自带的录音筒录了一段语音，语音的采样频率是44100HZ，然后用MATLAB编写程序进行测试，对其进行时域和频域的分析，得到图一。

2、带通滤波器：

我们所要设计的基于MATLAB的频率采样法的FIR带通滤波器的参数是通带为0.35π到0.65π，阻带下截止为0.2π，上截止为0.8π。

在正式设计滤波器之前，我们组员认真、仔细地查阅了很多资料，并参考了滤波器设计等应用图书，最终经过不断调试、修改得到了符合要求的仿真程序，在这个过程中，出现了不少困难和问题：

首先，在对话音信号进行时域和频谱分析时没有对话音信号取单通道，由于矩阵大小不一致致使在运算中始终得不到正确结果，调试无法通过。

起初对滤波器的通带和阻带的上下截止频率的设置和加噪后语音信号的中心频率出现了偏差，不能很完美地滤掉噪声。

此外，由于滤波器是根据相对频率设计的，再根据ω=Ω*Τ转换到模拟频率，当语音信号的采样频率变化时，滤波器的通带相应的发生改变，而我们在更换语音信号时忘记了这一点，以至于对有些语音信号不能达到消噪的目的而苦恼不解。

当然，在整个过程中，由于不细心出现了各种小问题，最终在队友的齐心努力下一一纠正。

3、加噪声过程：

根据语音信号的中心频率和频谱分布情况，将话音信号调制到中心频率为11025HZ的载波上，再加入高频噪声，噪声为22000Hz的余弦信号。

原先所加的噪声频率在语音信号的的频带内，不利于带通滤波器的滤波处理，这样就会导致滤波器无法滤除和语音信号相混杂的噪声，达不到滤波的效果，经过小组内的讨论之后，将噪声改为高频信号。

六、设计总结：

频率采样法设计带通滤波器相对来说比较简单，但是阻带失真衰减很大，增加采样点数N不能改善阻带最小衰耗。

改善阻带衰耗的唯一办法是加宽过渡带。

频率采样法特别适用于设计窄带选频滤波器。

因为这时只有少数几个非零值的

，计算量大为降低。

但由于频率抽样点的分布必须符合一定规律，在规定通、阻带截止频率方面不够灵活。

比如当截止频率

不是

整数倍时会产生较大逼近误差，因而该方法的应用不及窗口法普遍。

七、设计心得：

课程设计是培养我们综合运用所学知识,发现、提出、分析和解决实际问题,锻炼实践能力的重要环节。

本次课程设计的题目是基于频域抽样法的FIR数字带通滤波器设计，通过仔细阅读课本相关章节和借阅MATLAB教程书籍，为实际设计打好了理论基础，在此基础上，通过自己动手设计完成了课程设计要求。

通过这次课设，我们更进一步理解数字滤波器设计原理，学会了数字滤波器设计的方法和一般步骤，能够独立设计一个数字滤波器，实现了把理论知识转化为实际动手能力的过程。

我还从本次课程设计中体会到了MATLAB软件的强大功能，了解到它在各种工程计算中的重要作用，为我以后进一步学习打下了良好的基础。

当然这次课程设计也暴露了我的一些问题，比如学习程序设计教程不够快，数字信号处理基础不够扎实，虽然MATLAB使用的语言和语法都继承于C语言，但还是花了不少时间学习其中的函数，最后才能把课程设计顺利完成。

九、参考文献：

●《应用matlab语言处理数字信号与数字图像》科学出版社

●《Matlab7.x数字信号处理》人民邮电出版社

●《数字信号处理》北京理工大学出版社

●matlab中文论坛

●

●