实验5 频率采样法设计FIR数字滤波器(预习报告)

实验5 频率采样法设计FIR 数字滤波器

1．实验程序及运行结果

实验内容1：直接频率采样设计法

自定义一个能够产生已绝对值的幅度响应、相对dB 标尺的幅度响应、相位响应和群时延响应。函数名称定义为freqz_m 。

函数代码：freqz_m

%功能扩展的频响特性函数freqz_m.m function [db,mag,pha,grd,w]=freqz_m(b,a) %H ：频响特性的样值向量 %w ：频响特性的位置向量 %db ：幅度衰减（dB ） %mag ：H 的模 %pha ：H 的相位 %grd ：H 的群时延

%[H,w]=freqz(b,a,1000,'whole'); %H=(H(1:1:501))';w=(w(1:1:501))'; [H,w]=freqz(b,a,1024,'whole'); H=(H(1:1:512))';w=(w(1:1:512))'; mag=abs(H);

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

grd=grpdelay(b,a,w);% grpdelay 群时延函数 设计分析：

c k 取2c c k N ωπ≤

的最大整数，因此=22c c k N ω

π

=，故频率采样值()g H k 为 15()1,1,1,0,,0,1,1g H k ??

=??????

个零

源程序：shzxhchlshiyan8_1.m

%shzxhchlshiyan8_1

%频率采样法的直接设计方法（FIR 低通滤波器） clear all;close all;clc;clf;

M=20;tao=(M-1)/2;%M ：采样点数（20点），tao ：群时延（9.5） k=0:M-1;wl=(2*pi/M)*k;%w1描述各采样点频率值

Hrs=[1,1,1,zeros(1,15),1,1];%频域采样值（频率采样样本向量） %Hg(0)=Hg(1)=Hg(2)=Hg(18)=Hg(19)=1 %kc=(wc/2π)N,floor(kc)=2

Hdr=[1,1,0,0];wdl=[0,0.25,0.25,1];%Hdr 与wdl 共同确定1和0出现的频率位置边界 %当设置wdl=[0,0.25,0.25,1]时，是锐截止，无过渡带

%当设置wdl=[0,0.2,0.3,1]时，是设置频率在0.2π到0.3π之间的过渡带

k1=0:floor((M-1)/2);k2=floor((M-1)/2)+1:M-1;%（频率采样位置向量，0-9、10-19） %k2=ceil((M-1)/2):M-1

%floor:向小的整数方向取整，ceil ：向大的整数方向取整

%fix：取小数的整数部分，round：四舍五入取整

angH=[-tao*(2*pi)/M*k1,tao*(2*pi)/M*(M-k2)];% 计算-T(2π/N)kc（相位计算）

H=Hrs.*exp(j*angH);%求传输函数

h=real(ifft(H,M)); %求脉冲响应（单位冲激响应），即系数值

%h=abs(ifft(H,M));%%%不能使用绝对值函数

[db,mag,pha,grd,w]=freqz_m(h,1);%分子系数即h，分母系数恒等于1（FIR滤波器）%w：频率的位置向量

%db：幅度衰减（dB）

%mag：H的模

%pha：H的相位

%grd：H的群时延

[Hr,ww,a,L]=Hr_type2(h);

%Hr：振幅响应的样值向量

%ww：在[0,6π]区间内计算Hr的512个频率点

%a：第Ⅱ类线性相位FIR低通滤波器的系数

%L：Hr的阶次

%-----------------------------------

%subplot(2,2,1);plot(wl(1:11)/pi,Hrs(1:11),'o',wdl,Hdr);%画频率样本（半个周期）

%axis([0,1,-0.1,1.1]);grid;set(gca,'Xtick',[0,0.2,0.3,1]);

%title('频率样本:M=20');xlabel('w(单位:π)');ylabel('Hg(k)');

%title('Frequency Samples:M=20');xlabel('frequency in pi units');ylabel('Hr(k)');

%-----------------------------------

%----------------------------

%%subplot(2,2,1);plot(wl(1:11)/pi,Hrs(1:11));

%%subplot(2,2,1);stem(wl(1:11)/pi,Hrs(1:11));

%%subplot(2,2,1);plot(wl(1:11)/pi,Hrs(1:11),'o');

%%subplot(2,2,1);plot(wl(1:11)/pi,Hrs(1:11),wdl,Hdr);

%%subplot(2,2,1);plot('o',wdl,Hdr);%错

%%subplot(2,2,1);plot(wdl,Hdr,'o');

%%subplot(2,2,1);plot(wl/pi,Hrs,'o',wdl,Hdr);

Hdr=[1,1,0,0,1,1];wdl=[0,0.25,0.25,1.75,1.75,2];%画频率样本（1个周期）

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

axis([0,2,-0.1,1.1]);grid;set(gca,'Xtick',[0,0.2,0.3,1.8,1.9,2]);

title('频率样本:M=20');xlabel('w(单位:π)');ylabel('Hg(k)');

%%%%title('Frequency Samples:M=20');xlabel('frequency in pi units');ylabel('Hr(k)'); %----------------------------

subplot(2,2,2);stem(k,h,'.');%画脉冲响应

axis([-1,M,-0.1,0.3]);

title('脉冲响应');xlabel('n');ylabel('h(n)');

%title('Impulse Response');xlabel('n');ylabel('h(n)');

subplot(2,2,3);plot(ww/pi,Hr,wl(1:11)/pi,Hrs(1:11),'o');%%%画振幅响应

axis([0,1,-0.2,1.2]);grid;set(gca,'Xtick',[0,0.2,0.3,1]);

title('振幅响应');xlabel('w(单位:π)');ylabel('H_d(w)');

%title('Amplitude Response');xlabel('frequency in pi units');ylabel('Hr(w)');

%--------------------------------------------

%%subplot(2,2,3);plot(ww/pi,Hr,wl/pi,Hrs,'o');%%%画振幅响应 %%%axis([0,2,-0.2,1.2]);grid;set(gca,'Xtick',[0,0.2,0.3,1.7,1.8,2]); %%%axis([0,6,-0.2,1.2]);

%%title('振幅响应');xlabel('w(单位:π)');ylabel('H_d(w)'); %--------------------------------------------

subplot(2,2,4);plot(w/pi,db);%画幅度响应(dB) axis([0,1,-60,10]);grid;set(gca,'Xtick',[0,0.2,0.3,1]); title('幅度响应(dB)');xlabel('w(单位:π)');ylabel('dB');

%title('Magnitude Response');xlabel('frequency in pi units');ylabel('Decibels'); 程序运行结果：

频率样本:M=20

w(单位:π)

H g (k )

-0.1

0.10.2

0.3脉冲响应

n

h (n

)

0.20.3

1

振幅响应

w(单位:π)

H d (w

)

0.20.3

1

-60-40

-200幅度响应(dB)

w(单位:π)

d B

图8-1 直接频率采样设计法

程序运行结果分析：

阻带最小衰减大约是16dB ，通带、阻带衰减均不符合设计要求。

实验内容2：频率采样设计法（过渡带设置一个采样点） 设计分析：

设T1=0.5。当采样点数40N =时，频率采样值()g H k 为

1129()1,1,1,1,1,T 0,,0,T ,1,1,1,1g H k ??

=??????

个零，

源程序（T 1=0.5

）：shzxhchlshiyan8_2 %shzxhchlshiyan8_2

%频率采样法的最优设计方法（FIR 低通滤波器） %T1=0.5(非最优) T1为过渡带内的样本值 clear all;close all;clc;clf;

M=40;tao=(M-1)/2;%M ：采样点数（40点），tao ：群时延（19.5） k=0:M-1;wl=(2*pi/M)*k;

Hrs=[ones(1,5),0.5,zeros(1,29),0.5,ones(1,4)];%

%前面五个1加一个过渡点0.5，后面一个过渡点0.5加四个1 Hdr=[1,1,0,0];wdl=[0,0.25,0.25,1];

k1=0:floor((M-1)/2);k2=floor((M-1)/2)+1:M-1; angH=[-tao*(2*pi)/M*k1,tao*(2*pi)/M*(M-k2)]; H=Hrs.*exp(j*angH); h=real(ifft(H,M));

[db,mag,pha,grd,w]=freqz_m(h,1); [Hr,ww,a,L]=Hr_type2(h);

subplot(2,2,1);plot(wl(1:21)/pi,Hrs(1:21),'o',wdl,Hdr); axis([0,1,-0.1,1.1]);grid;set(gca,'Xtick',[0,0.2,0.3,1]);

title('频率样本:M=40');xlabel('w(单位:π)');ylabel('Hg(k)');

subplot(2,2,2);stem(k,h,'.');axis([-1,M,-0.1,0.3]);title('脉冲响应');xlabel('n');ylabel('h(n)'); subplot(2,2,3);plot(ww/pi,Hr,wl(1:21)/pi,Hrs(1:21),'o'); axis([0,1,-0.2,1.2]);grid;set(gca,'Xtick',[0,0.2,0.3,1]); title('振幅响应');xlabel('w(单位:π)');ylabel('H_d(w)');

subplot(2,2,4);plot(w/pi,db);axis([0,1,-80,10]);grid;set(gca,'Xtick',[0,0.2,0.3,1],'Ytick',[-30,0]); title('幅度响应(dB)');xlabel('w(单位:π)');ylabel('dB'); 程序运行结果：

0.20.3

1

频率样本:M=40

w(单位:π)H g (k )

010

203040

-0.1

0.10.2

0.3脉冲响应

n

h (n )

振幅响应

w(单位:π)

H d (w

)

-30

0幅度响应(dB)

w(单位:π)

d B

图8-2 最优频率采样设计法（T1=0.5）

程序运行结果分析：

通带衰减符合设计要求，阻带最小衰减提高到约30dB ，阻带衰减仍不符合设计要求。 源程序（T1=0.39）：shzxhchlshiyan8_3 %shzxhchlshiyan8_3

%频率采样法的最优设计方法（FIR 低通滤波器） %T1=0.39(最优) T1为过渡带内的样本值 clear all;close all;clc;clf;

M=40;tao=(M-1)/2;%M ：采样点数（40点），tao ：群时延（19.5） k=0:M-1;wl=(2*pi/M)*k;

Hrs=[ones(1,5),0.39,zeros(1,29),0.39,ones(1,4)];%

%前面五个1加一个过渡点0.39，后面一个过渡点0.39加四个1 Hdr=[1,1,0,0];wdl=[0,0.25,0.25,1];

k1=0:floor((M-1)/2);k2=floor((M-1)/2)+1:M-1; angH=[-tao*(2*pi)/M*k1,tao*(2*pi)/M*(M-k2)]; H=Hrs.*exp(j*angH); h=real(ifft(H,M));

[db,mag,pha,grd,w]=freqz_m(h,1); [Hr,ww,a,L]=Hr_type2(h);

subplot(2,2,1);plot(wl(1:21)/pi,Hrs(1:21),'o',wdl,Hdr); axis([0,1,-0.1,1.1]);grid;set(gca,'Xtick',[0,0.2,0.3,1]);

title('频率样本:M=40');xlabel('w(单位:π)');ylabel('Hg(k)'); subplot(2,2,2);stem(k,h,'.');axis([-1,M,-0.1,0.3]); title('脉冲响应');xlabel('n');ylabel('h(n)');

subplot(2,2,3);plot(ww/pi,Hr,wl(1:21)/pi,Hrs(1:21),'o'); axis([0,1,-0.2,1.2]);grid;set(gca,'Xtick',[0,0.2,0.3,1]); title('振幅响应');xlabel('w(单位:π)');ylabel('H_d(w)'); subplot(2,2,4);plot(w/pi,db);

axis([0,1,-60,10]);grid;set(gca,'Xtick',[0,0.2,0.3,1],'Ytick',[-60,-43,0]); title('幅度响应(dB)');xlabel('w(单位:π)');ylabel('dB'); 程序运行结果：

0.20.3

1

频率样本:M=40

w(单位:π)H g (k )

010

203040

-0.1

0.10.2

0.3脉冲响应

n

h (n )

振幅响应

w(单位:π)

H d (w

)

-60-43

幅度响应(dB)

w(单位:π)

d B

图8-3 最优频率采样设计法（T1=0.39）

程序运行结果分析：

通带衰减符合设计要求，阻带最小衰减提高到约43dB ，阻带衰减仍不符合设计要求。 实验内容3：频率采样设计法（过渡带设置两个采样点） 设计分析：

设T1=0.5925、T2=0.1099。当采样点数60N =时，频率采样值()g H k 为

12217"1"6"1"43"0"()1,,1,T T 0,,0,T ,T ,1,,1g H k ??

=??????

个个个,,

源程序：shzxhchlshiyan8_4

%shzxhchlshiyan8_4

%频率采样法的最优设计方法（FIR 低通滤波器）

%T1=0.5925,T2=0.1099(最优：可查表得到) T1、T2为过渡带内的两个过渡点样本值 clear all;close all;clc;clf;

M=60;tao=(M-1)/2;%M ：采样点数（60点），tao ：群时延（29.5） k=0:M-1;wl=(2*pi/M)*k;

Hrs=[ones(1,7),0.5925,0.1099,zeros(1,43),0.1099,0.5925,ones(1,6)];%

%前面七个1加二个过渡点0.5925、0.1099，后面二个过渡点0.1099、0.5925加六个1 Hdr=[1,1,0,0];wdl=[0,0.2,0.3,1];%设置频率在0.2π到0.3π之间的过渡带 k1=0:floor((M-1)/2);k2=floor((M-1)/2)+1:M-1; angH=[-tao*(2*pi)/M*k1,tao*(2*pi)/M*(M-k2)]; H=Hrs.*exp(j*angH); h=real(ifft(H,M));

[db,mag,pha,grd,w]=freqz_m(h,1); [Hr,ww,a,L]=Hr_type2(h);

subplot(2,2,1);plot(wl(1:31)/pi,Hrs(1:31),'o',wdl,Hdr); axis([0,1,-0.1,1.1]);grid; set(gca,'Xtick',[0,0.2,0.3,1]);

title('频率样本:M=60');xlabel('w(单位:π)');ylabel('Hg(k)'); subplot(2,2,2);stem(k,h,'.'); axis([-1,M,-0.1,0.3]); title('脉冲响应');

xlabel('n');ylabel('h(n)');

subplot(2,2,3);plot(ww/pi,Hr,wl(1:31)/pi,Hrs(1:31),'o'); axis([0,1,-0.2,1.2]);grid; set(gca,'Xtick',[0,0.2,0.3,1]); title('振幅响应');

xlabel('w(单位:π)');ylabel('H_d(w)'); subplot(2,2,4);plot(w/pi,db); axis([0,1,-100,10]);grid;

set(gca,'Xtick',[0,0.2,0.3,1],'Ytick',[-100,-63,0]); title('幅度响应(dB)');

xlabel('w(单位:π)');ylabel('dB'); 程序运行结果：

0.20.3

1

频率样本:M=60

w(单位:π)H g (k )

0204060

-0.1

0.10.2

0.3脉冲响应

n

h (n )

0.20.3

1

振幅响应

w(单位:π)

H d (w

)

0.20.3

1

-100

-63

幅度响应(dB)

w(单位:π)

d B

图8-4 最优频率采样设计法（T1=0.5925、T2=0.1099）

程序运行结果分析：

阻带最小衰减提高到约63dB ，通带、阻带衰减均符合设计要求。

2．思考题简答

（1）简要说明用直接频率采样设计法设计FIR 数字滤波器的基本方法。 答：略（见指导书）

（2）简要说明用最优频率采样设计法设计FIR 数字滤波器的基本方法。 答：略（见指导书）

（3）简要说明频率采样法设计步骤。 答：略（见指导书）