基于频率抽样设计法线性相位型FIR数字低通滤波器设计报告

目录

任务书 .................................... 2 摘要： .................................... 3 关键字： .................................. 3 FIR数字滤波器设计的原理 ................... 3 设计过程： ................................ 5 实验结果 .................................. 7 1.程序： ................................ 7 2.截图 ................................. 11 结束语

任务书

题目8 主要 内容 基于频率抽样设计法线性相位型FIR数字低通滤波器设计 1、根据设计指标，利用频率抽样法设计线性相位型FIR数字低通滤波器； 2、对比分析改变抽样点数及过渡点数对数字滤波器性能的影响。 设计 要求 1、根据设计指标，构建线性相位型FIR滤波器之频响抽样值序列； 2、确定数字滤波器的冲击响应及系统函数； 3、对比分析抽样点数及过渡点数的增减对系统频响的影响； 4、生成一个时间序列，利用已设计的FIR滤波器进行滤波，分析滤波效果。 主要 仪器 设备 主要 参考 文献 1、计算机1台，安装MATLAB软件 [美]维纳.K.恩格尔，约翰.G.普罗科斯著，刘树棠译.数字信号处理——使用MATLAB[M].西安：西安交通大学出版社，2002. 飞思科技产品研发中心编著.MATLAB7辅助信号处理技术与应用[M].北京：电子工业出版社，2005. 课程设计进度安排（起止时间、工作内容） 课程设计共设16个设计题目，每班3至4人为1组，1人1套设备，每组选作不同的题目，4个班共分4批。完整课程设计共20学时，为期1周，具体进度如下： 5学时 学习题目相关知识，掌握实现原理； 5学时 用MATLAB语言实现题目要求； 5学时 进一步完善功能，现场检查、答辩； 5学时 完成并提交课程设计报告。 课程设计开始日期 2012.12.31 课程设计完成日期 地 点 2013.1.6 实验楼501 课程设计实验室名称 信号处理实验室 资料下载地址

各班公共邮箱 2

摘要：

论文从基于频率抽样设计法线性相位型FIR数字低通滤波器的设计原理和设计过程方面来阐述我们的设计思路，并用matlab软件运行后得到了论文中的实验结果，我们还对课程设计中遇到的问题及解决方法进行了阐述，最后描述了自己的心得体会。

关键字：

频率抽样法 线性相位型FIR数字低通滤波器 matlab 有限长冲激响应

一、FIR数字滤波器设计的原理 频率抽样设计法

FIR低通滤波器的设计一般方法有两种，即频率抽样法和窗函数法，频率抽样法设计不同于窗函数法，窗函数是从时域出发，把理想的hd(n)用一定形状得窗函数截取成有限长的h(n)，以此h(n)来近似理想的hd(n)，这样得到的频率响应H(ejw)逼近于所要求的理想的频率响应Hd(ejw)。

频率抽样法则是从频域出发，把给定的理想频率响应Hd(ejw)加以等间隔抽样，即Hd(ejw)|w?2?kN?Hd(k)然后以此Hd(k)作为实际FIR数字滤波器的频率特性

w?2?kN的抽样值H(k)，即令H(k)?Hd(k)?Hd(ejw)|k?0,1,,N?1，知道H(k)后，由DFT定义，可以用频域的这N个抽样值H(k)来唯一确定有限长序列h(n)，而由X(z)的内插公式知道，利用这N个频域抽样值H(k)同样可求得FIR滤波器的系统函数H(z)及频率响应H(ejw)。这个H(z)或H(ejw)将逼近Hd(z)或

Hd(ejw)，H(z)和H(ejw)的内插公式为

1?z?NH(z)?NjwN?1k?0H(k) ??k?11?Wzk?0N2?k) （2－2） NN?1H(e)??H(k)?(w?

3

其中?(w)是内插函数

wN)?jw(N?1)122 （2－3） ?(w)?ewNsin()2将式（2－3）代入（2－2）式，化简后可得

wNsin()?kN?1N?1?j1?j(2)wjw2 H(e)?eH(k)eN?w?kNk?0sin(?)2Nw?ksin[N(?)]?kN?1N?1?j()w1jN(N?1)jw2N 2即 H(e)?eH(k)?e??w?kNk?0sin(?)2N

sin(H(k)??h(n)en?0N?1?j2?nkN

从内插公式（2－2）看到，在各频率抽样点上，滤波器的实际频率响应是严格地和理想频率响应数值相等，即H(ej2?kN)?H(k)?Hd(k)?Hd(ej2?kN)。但是在

抽样点之间的频率响应则是由各抽样点的加权内插函数的延伸叠加而形成的，因而有一定的逼近误差，误差大小取决于理想频率响应曲线形状，理想频率响应特性变化越平缓，则内插值越接近理想值，逼近误差越小，如下图梯形理想频率特性所示。

反之，如果抽样点之间的理想频率特性变化越陡，则内插值与理想值之误差就越大，因而在理想频率特性的不连续点附近，就会产生肩峰和波纹，

4

当h(n)为实数时，满足H(k)?H*((N?k))NRN(k)?H*(N?k)，由此得出

H(k)?H(N?k，)Q(k)??Q(N?k)，也就是说，H(k)的模H(k)以k?N2为对称中心呈偶对称，H(k)的相角Q(k)以k?N2为对称中心呈奇对称。再利用线性相位的条件Q(ejw)??N?1w，即可得到（N为偶数）： 2?2??N?1?k?0,??Nk?2?????0k? Q(k)????2??N?1??N(N?k)k?????N2???2??N?,??1??2?N

2??1?,,N?1?2?N?1?jk()??N?N2H(k)ek?0,,?1???2????N0k? H(k)??

2?2?N?1?j(N?k)()?N?N2k???1?,,N?1?H(N?k)e?2??

???N?H(0)sinN?1????j??2??j?2H?e??e ?????Nsin????2??N?12k?1??????k??????k?????sinN????sin?N??????H?k???2N2N?????????????

?k??k?N???????sin?sin???????2N2N????????

二、设计过程：

（一）设计思想：

根据给定DLPF幅频特性要求（通带截止频率ωp=0.5π，通带最大衰减αp=0.5 dB，阻带截止频率ωs=0.6π，阻带最小衰减αs=50 dB）取得DLPF的H(k)。设冲激响应h(n)为偶对称，N为偶数。根据线性相位约束条件设计FIR线性相位数字滤波器的H(k)。根据H(k)生成DLPF的h(n)。利用频率抽样法设计加一个

5

过渡带，加两个过渡带的低通滤波器，看其对滤波器性能的影响。设计与之相对应的DLPF，给出窗函数及所设计滤波器的幅度特性，对比分析DLPF幅频特性是否符合要求最后利用x=rand(1,sizex)函数随随机生成一个序列来验证设计的滤波器是否具有低通滤波器的特性。 （二）功能结构：

FIR数字滤波器具有严格的线性相位，低通滤波器只能让低频的通过而把高频的部分滤掉。 设计步骤： （三）设计步骤

① 根据给定DLPF幅频特性要求（通带截止频率ωp=0.5π，通带最大衰减αp=0.5 dB，阻带截止

频率ωs=0.6π，阻带最小衰减αs=50 dB）取得DLPF的X(K)； ② 根据线性相位型数字滤波器条件，构建线性相位型DLPF的X(K)； ③ 根据X(K)生成DLPF的h(n)；

④ 设计与之相对应的DLPF，给出窗函数及所设计滤波器的幅度特性，对比分析DLPF幅频特 性是否符合要求；

⑤ 试说明过渡点对所设计数字滤波器性能的影响；

⑥ 产生一个有干扰频率的时域序列（借助FFT分析说明其有干扰），使之通过所设计的DLPF，

对滤波输出结果作出分析，说明输出结果。

⑦ 扩展部分：自拟指标，设计一个DBPF，追求最佳性能，并检验设计效果。

6

三、实验结果

（一）程序：

% wp=0.5*pi;Rp=0.5dB; ws=0.6*pi;As=50dB;

频率抽样法

wp=0.5*pi;ws=0.6*pi; tr_width=ws-wp;

% 用频率抽样法设计FIR滤波器,过渡带内一个样本T1，N=40。 alpha=(N-1)/2; %N为偶数 l=[0:1:N-1]; wl=(2*pi/N)*l;

Hrs=[ones(1,11),T1,zeros(1,17),T1,ones(1,10)]; %偶对称

Hdr=[1 1 0 0 ];wdl=[0 0.5 0.6 1]; k1=0:(N/2-1); k2=(N/2+1):N-1;%依据公式7-107 angH=[-alpha*(2*pi)/N*k1,0,alpha*(2*pi)/N*(N-k2)]; H=Hrs.*exp(1i*angH);

h=real(ifft(H,N));[H,w]=freqz(h,1,1000,'whole'); %[db mag pha grd w]=freqz_m(h,1);

db=20*log10((abs(H)+eps)/max(abs(H)));%求FIR滤波器频响的dB值 delta_w=2*pi/1000; %将2pi等分1000份 %plot figure(1);

subplot(221);plot(wl(1:21)/pi,Hrs(1:21),'o',wdl,Hdr,'linewidth',2); title('理想滤波器频域波形');

axis([0 1 -0.1 ,1.2]); ylabel('Hr(k)');

set(gca,'XTickMode','manual','XTick',[0 0.5 0.6 1]); set(gca,'YTickMode','manual','YTick',[0 T1 1]); grid; subplot(222); stem(l,h,'m');

title('单位脉冲响应');axis([-1,N,-0.15,0.5]);ylabel('h(n)');

subplot(223); plot(w/pi,abs(H),wl(1:31)/pi,Hrs(1:31),'o','linewidth',2);axis([0 1 -0.2,1.2]); title('频域抽样');xlabel('频率');ylabel('Hr(w)');grid; set(gca,'XTickMode','manual','XTick',[0 0.5 0.6 1]); set(gca,'YTickMode','manual','YTick',[0 T1 1]);

subplot(224); plot(w/pi,db,'r');axis([0,1,-100,10]);grid; ylabel('分贝'); title('幅频响应');xlabel('频率');

set(gca,'XTickMode','manual','XTick',[0,0.5,0.6,1]); set(gca,'YTickMode','manual','YTick',[-50,0]);

% 用频率抽样法设计FIR滤波器,过渡带内无样本，N=20。 tr_width=ws-wp;

alpha=(N-1)/2; %N为偶数 l=[0:1:N-1]; wl=(2*pi/N)*l;

Hrs=[ones(1,6),zeros(1,9),ones(1,5)]; %偶对称

7

Hdr=[1 1 0 0 ];wdl=[0 0.5 0.6 1]; k1=0:(N/2-1); k2=(N/2+1):N-1;

angH=[-alpha*(2*pi)/N*k1,0,alpha*(2*pi)/N*(N-k2)]; H=Hrs.*exp(j*angH); h=real(ifft(H,N));[H,w]=freqz(h,1,1000,'whole'); %[db mag pha grd w]=freqz_m(h,1);

db=20*log10((abs(H)+eps)/max(abs(H)));%求FIR滤波器频响的dB值 delta_w=2*pi/1000; %将2pi等分1000份 %plot

figure(2);clf;

subplot(221);plot(wl(1:11)/pi,Hrs(1:11),'o',wdl,Hdr,'linewidth',2);

title('理想滤波器频域波形');

axis([0 1 -0.1 ,1.2]); ylabel('Hr(k)');

set(gca,'XTickMode','manual','XTick',[0 0.5 0.6 1]); %set(gca,'YTickMode','manual','YTick',[0 T1 1]); grid;

subplot(222); stem(l,h,'m');

title('单位脉冲响应');axis([-1,N,-0.15,0.55]);ylabel('h(n)');

subplot(223); plot(w/pi,abs(H),wl(1:11)/pi,Hrs(1:11),'o','linewidth',2);axis([0 1 -0.2,1.2]); title('频域抽样');xlabel('频率');ylabel('Hr(w)');grid; %set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 1]); %set(gca,'YTickMode','manual','YTick',[0 T1 1]);

subplot(224); plot(w/pi,db,'r');axis([0,1,-100,10]);grid; ylabel('分贝'); title('幅频响应');xlabel('频率');

set(gca,'XTickMode','manual','XTick',[0,0.5,0.6,1]); set(gca,'YTickMode','manual','YTick',[-50,0]);

% 用频率抽样法设计FIR滤波器,过渡带内两个样本T1、T2，N=60。 tr_width=ws-wp;

alpha=(N-1)/2; %N为偶数 l=[0:1:N-1]; wl=(2*pi/N)*l;

Hrs=[ones(1,16),T1,T2,zeros(1,25),T2, T1,ones(1,15)]; %偶对称 Hdr=[1 1 0 0 ];wdl=[0 0.5 0.6 1]; k1=0:(N/2-1); k2=(N/2+1):N-1;

angH=[-alpha*(2*pi)/N*k1,0,alpha*(2*pi)/N*(N-k2)]; H=Hrs.*exp(j*angH); h=real(ifft(H,N));[H,w]=freqz(h,1,1000,'whole'); %[db mag pha grd w]=freqz_m(h,1);

db=20*log10((abs(H)+eps)/max(abs(H)));%求FIR滤波器频响的dB值 delta_w=2*pi/1000; %将2pi等分1000份 %plot

figure(3);clf;

subplot(221);plot(wl(1:31)/pi,Hrs(1:31),'o',wdl,Hdr,'linewidth',2); title('理想滤波器频域波形');

axis([0 1 -0.1 ,1.2]); ylabel('Hr(k)');

set(gca,'XTickMode','manual','XTick',[0 0.5 0.6 1]); set(gca,'YTickMode','manual','YTick',[0 T2 T1 1]); grid; subplot(222); stem(l,h,'m');

title('单位脉冲响应');axis([-1,N,-0.15,0.55]);ylabel('h(n)');

8

subplot(223); plot(w/pi,abs(H),wl(1:31)/pi,Hrs(1:31),'o','linewidth',2);axis([0 1 -0.2,1.2]); title('频域抽样');xlabel('frequency in pi units');ylabel('Hr(w)');grid; set(gca,'XTickMode','manual','XTick',[0 0.5 0.6 1]); set(gca,'YTickMode','manual','YTick',[0 T2 T1 1]);

subplot(224); plot(w/pi,db,'r');axis([0,1,-100,10]);grid; ylabel('Decibels'); title('幅频响应');xlabel('频率');

set(gca,'XTickMode','manual','XTick',[0,0.5,0.6,1]); set(gca,'YTickMode','manual','YTick',[-50,0]); figure(4);

subplot(211); plot(w/pi,db,'r');axis([0,1,-100,10]);grid; ylabel('分贝'); title('幅频响应');xlabel('频率');

set(gca,'XTickMode','manual','XTick',[0,0.5,0.6,1]); set(gca,'YTickMode','manual','YTick',[-50,0]);

subplot(212); plot(w/pi,angle(H),'b');axis([0,1,-3.5,3.5]);grid; ylabel('rad'); title('相频响应');

%生成一个时间序列，利用已设计的FIR滤波器进行滤波，分析滤波效果 alpha=(N-1)/2; %N为偶数 l=[0:1:N-1]; wl=(2*pi/N)*l;

Hrs=[ones(1,16),T1,T2,zeros(1,25),T2, T1,ones(1,15)]; %偶对称 Hdr=[1 1 0 0 ];wdl=[0 0.5 0.6 1]; k1=0:(N/2-1); k2=(N/2+1):N-1;

angH=[-alpha*(2*pi)/N*k1,0,alpha*(2*pi)/N*(N-k2)]; H=Hrs.*exp(j*angH); h=real(ifft(H,N)); x=rand(1,60); LL=0:29; X=fft(x,60); XX=[X(1:30)];

y=conv(x,h); %x(n)和h(n)做卷积 yL=length(y);

disp('y(n)序列长度为N=');disp(yL);%求卷积后的序列长度 Y=fft(y,60); YY=[Y(1:30)]; %plot figure(5)

m=0:1:length(x)-1; subplot(2,2,1);

title('序列时域波形'); stem(m,x,'.');

xlabel('时间序列 n');ylabel('x(n)'); subplot(2,2,2); title('序列频域'); stem(LL,abs(XX),'.');

xlabel('k=0,1,2,...,29');ylabel('|X(k)|');axis([0,30,0,30]); m=0:1:length(y)-1; subplot(2,2,3);

9

title('卷积后波形'); stem(m,y,'.');

xlabel('时间序列 n');ylabel('y(n)'); subplot(2,2,4);

title('滤波后的频域抽样'); stem(LL,abs(YY),'.');

xlabel('k=0,1,2,...,29');ylabel('|Y(k)|');axis([0,30,0,30]);

窗函数法

% wp=0.5*pi;Rp=0.5dB; ws=0.6*pi;As=50dB;

% 查第237页表7.1可知用Hamming窗、Blackman窗均可(最小阻带衰减>=50dB)， % 但Hamming窗具有较小的过渡带(6.6pi/M)，故选择Hamming窗。 wp=0.5*pi;ws=0.6*pi; tr_width=ws-wp;

M=ceil(6.6*pi/tr_width)+1; n=[0:1:M-1]; wc=(ws+wp)/2; hd=ideal_lp(wc,M);

w_ham=(hamming(M))'; h=hd.*w_ham;

[H,w]=freqz(h,1,1000,'whole'); H=(H(1:1:501))';w=(w(1:1:501))'; mag=abs(H);

db=20*log10((mag+eps)/max(mag)); pha=angle(H);

grd=grpdelay(h,1,w); delta_w=2*pi/1000; Rp=0.5; As=50; figure(1);

subplot(2,2,1);stem(n,hd);title('理想脉冲响应'); axis([0 M-1 -0.1 0.56]);xlabel('n');ylabel('hd(n)'); subplot(2,2,2);plot(w/pi,pha);

axis([0,1,-3.5,3.5]); ylabel('角度');xlabel('频率'); title('相位频率响应');

set(findobj(gcf,'Type','line','Color',[0 0 1]), 'Color','b', 'LineWidth',2); grid;

subplot(2,2,3);stem(n,h);title('实际脉冲响应'); axis([0 M-1 -0.1 0.56]);xlabel('n');ylabel('h(n)'); subplot(2,2,4); plot(w/pi,db);

axis([0,1,-100,10]); ylabel('幅度');title('幅度响应'); set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 1]); set(gca,'YTickMode','manual','YTick',[-50 0]);grid; pause

10

实验中用到的.m文件 hd=ideal_lp(wc,M)

function hd=ideal_lp(wc,M) alpha=(M-1)/2; n=[0:1:(M-1)]; m=n-alpha+eps; hd=sin(wc*m)./(pi*m);

end

（二）截图

图1.添加一个过渡点所设计的滤波器、冲激响应h(n)、幅度特性（N=40，过渡点T1=0.37）

11

图2.不添加过渡点设计的低通滤波器、冲激响应h(n)、幅度特性（N=20，无过渡点）

12

图3.添加两个过渡点设计的低通滤波器、冲激响应h(n)、幅度特性（N=60，过渡点T1=0.59 T2=0.11）

13

图4.添加两个过渡点设计的低通滤波器幅频和相频响应

14

图5.添加两个过渡点设计的低通滤波器的冲激响应h(n)

15

本文来源：https://www.bwwdw.com/article/uo62.html