.专业整理.
实验三:FIR数字滤波器的设计
实验目的
1) 掌握用窗函数法,频率采样法及优化设计法设计FIR 滤波器的原理及方法。 2) 熟悉线性相位FIR 滤波器的幅频特性和相频特性。 3) 了解各种不同窗函数对滤波器性能的影响。
一、 实验内容
1. N=45,计算并画出矩形窗、汉明窗、布莱克曼窗的归一化的幅度谱,并比较各自的主要
特点。
clear all; N=45;
wn1=kaiser(N,0); wn2=hamming(N); wn3=blackman(N); [h1,w1] = freqz(wn1,N); [h2,w2] = freqz(wn2,N); [h3,w3] = freqz(wn3,N);
plot(w1/pi,20*log10(abs(h1)),'r-',w2/pi,20*log10(abs(h2)),'b-',w3/pi,20*log10(abs(h3)),'g-');
axis([0,1,-120,10]);grid;
xlabel('归一化频率/\\pi'); ylabel('幅度/dB'); title('三种窗口函数');
legend('矩形窗','汉明窗','布莱克曼窗',3);
三种窗口函数 0矩形窗汉明窗布莱克曼窗-20-40幅度/dB-60-80-100-120 00.10.20.30.40.5归一化频率/0.60.70.80.91
分析:阻带衰减和过渡带带宽是相互矛盾的,矩形窗过渡带带宽窄,但是阻带衰减比较少;布莱克曼窗过渡带带宽宽,但是阻带衰减比较大
.学习帮手.
.专业整理.
2. N=15,带通滤波器的两个通带边界分别是ω1=0.3π,ω2=0.5π。用汉宁窗设计此线性
相位带通滤波器,观察它的实际3dB 和20dB 带宽。N=45,重复这一设计,观察幅频和相位特性的变化,注意长度N 变化的影响。
N=15;
h= fir1(N-1,[0.3 0.5],'bandpass',hanning(N)); figure(1)
freqz(h,1);axis([0,1,-60,10]); title('N=15,汉宁窗'); N=45;
h= fir1(N-1,[0.3 0.5],'bandpass',hanning(N)); figure(2)
freqz(h,1);axis([0,1,-60,10]); title('N=45,汉宁窗');
N=15,汉宁窗0-20-40-60Magnitude (dB)00.10.20.30.40.50.60.70.8Normalized Frequency  ( rad/sample)0.91500Phase (degrees)0-500-100000.10.20.30.40.50.60.70.8Normalized Frequency  ( rad/sample)0.91
.学习帮手.
.专业整理.
N=45,汉宁窗0-20-40-60Magnitude (dB)00.10.20.80.70.60.50.40.3Normalized Frequency  ( rad/sample)0.91500Phase (degrees)0-500-1000-150000.10.20.80.70.60.50.40.3Normalized Frequency  ( rad/sample)0.91
结论:增加窗口函数的长度能够在幅度频谱和相位频谱上获得较好的特性。但代价是增加了计算量和系统的阶数。
3. 分别改用矩形窗和布莱克曼窗,设计(2)中的带通滤波器,观察并记录窗函数对滤波器
幅频特性的影响,比较三种窗的特点。
clear all; %矩形窗
N=15;h= fir1(N-1,[0.3 0.5],'bandpass',kaiser(N,0)); [h1,w1]=freqz(h,1);
subplot(2,1,1); plot(w1/pi,20*log10(abs(h1)));axis([0,1,-60,10]); grid;
xlabel('归一化频率/\\pi'); ylabel('幅度/dB');title('N=15,矩形窗'); N=45;h= fir1(N-1,[0.3 0.5],'bandpass',kaiser(N,0)); [h1,w1]=freqz(h,1);
subplot(2,1,2); plot(w1/pi,20*log10(abs(h1)));axis([0,1,-60,10]); grid;
xlabel('归一化频率/\\pi'); ylabel('幅度/dB');title('N=45,矩形窗'); clear all; %布莱克曼窗
N=15;h= fir1(N-1,[0.3 0.5],'bandpass',blackman(N)); [h1,w1]=freqz(h,1);
subplot(2,1,1);plot(w1/pi,20*log10(abs(h1))); axis([0,1,-60,10]);grid;
xlabel('归一化频率/\\pi'); ylabel('幅度/dB');title('N=15,布莱克曼窗'); N=45;h= fir1(N-1,[0.3 0.5],'bandpass',blackman(N)); [h1,w1]=freqz(h,1);
subplot(2,1,2);plot(w1/pi,20*log10(abs(h1)));
.学习帮手.
.专业整理.
axis([0,1,-60,10]);grid;
xlabel('归一化频率/\\pi'); ylabel('幅度/dB');title('N=45,布莱克曼窗');
N=15,矩形窗0幅度/dB-20-40-6000.10.20.30.40.50.6归一化频率/N=45,矩形窗0.70.80.910幅度/dB-20-40-6000.10.20.30.40.50.6归一化频率/N=15,布莱克曼窗0.70.80.91
0幅度/dB-20-40-6000.10.20.30.40.50.6归一化频率/N=45,布莱克曼窗0.70.80.910幅度/dB-20-40-6000.10.20.30.40.50.6归一化频率/0.70.80.91
总结:同实验内容1的结论,除此之外,利用各窗口函数构造的带通滤波器的特性也有
.学习帮手.
.专业整理.
些差别。汉宁窗在这方面具有较好的特性。
4.  用Kaiser 窗设计一专用线性相位滤波器,N=40,当β=4、6、10 时,分别设计、比较
它们的幅频和相频特性,注意β取不同值时的影响。
clear all; N=40;
f = [0 0.2 0.2 0.4 0.4 0.6 0.6 0.8 0.8 1]; a = [0 0 1 1 0 0 1 1 0 0];
beta=4;h = fir2(N-1,f,a,kaiser(N,beta)); [h1,w1]=freqz(h,1);
figure;plot(w1/pi,20*log10(abs(h1)));axis([0,1,-60,10]);grid; xlabel('归一化频率/\\pi'); ylabel('幅度/dB'); title('beta=4 时凯塞窗专用线性相位滤波器'); beta=6;h = fir2(N-1,f,a,kaiser(N,beta)); [h1,w1]=freqz(h,1);
figure;plot(w1/pi,20*log10(abs(h1)));axis([0,1,-60,10]);grid; xlabel('归一化频率/\\pi'); ylabel('幅度/dB'); title('beta=6 时凯塞窗专用线性相位滤波器'); beta=10;h = fir2(N-1,f,a,kaiser(N,beta)); [h1,w1]=freqz(h,1);
figure;plot(w1/pi,20*log10(abs(h1)));axis([0,1,-60,10]);grid; xlabel('归一化频率/\\pi'); ylabel('幅度/dB'); title('beta=10 时凯塞窗专用线性相位滤波器');
.学习帮手.
.专业整理.
beta=4 时凯塞窗专用线性相位滤波器100-10幅度/dB-20-30-40-50-6000.10.20.30.60.50.4归一化频率/0.70.80.91
beta=6 时凯塞窗专用线性相位滤波器100-10幅度/dB-20-30-40-50-6000.10.20.30.40.50.6归一化频率/0.70.80.91
.学习帮手.
.专业整理.
beta=10 时凯塞窗专用线性相位滤波器100-10幅度/dB-20-30-40-50-6000.10.20.30.60.50.4归一化频率/0.70.80.91
总结:Kaiser 窗的beta值越大,过渡带的带宽越宽,阻带的衰减越是厉害;beta值越小,过渡带的带宽越窄,阻带的衰减越是却有所降低。所以在实际应用的时候,要权衡过渡带带宽和阻带衰减,以获得最优性能!
5. 用频率采样法设计(4)中的滤波器,过渡带分别设一个过渡点,令H(k)=0.5。比较两种
不同方法的结果。
clear all; N=40;
Hk=[zeros(1,3) 0.5 ones(1,5) 0.5 zeros(1,1) 0.5 ones(1,5) 0.5 zeros(1,5) -0.5 -ones(1,5) -0.5 zeros(1,1) -ones(1,5) -0.5 zeros(1,3)]; k=0:N-1;
hn=real(ifft(Hk.*exp(-j*pi*(N-1)*k/N))); [H w]=freqz(hn, 1);
plot(w/pi, 20*log10(abs(H))); axis([0 1 -80 10]);grid; xlabel('归一化频率/\\pi') ylabel('幅度/dB')
title('频率采样法设计专用线性相位滤波器');
.学习帮手.
.专业整理.
频率采样法设计专用线性相位滤波器100-10-20幅度/dB-30-40-50-60-70-8000.10.20.30.40.50.6归一化频率/0.70.80.91
总结:利用频率采样,获得的滤波器通带内波动较好,但衰减有所降低,过渡带也比较宽!
6. 用雷米兹(Remez)交替算法设计(4)中的滤波器,并比较(4)、(5)、(6)三种不同方法的
结果。
clear all; N=40;
f=[0 0.15 0.2 0.4 0.45 0.55 0.6 0.8 0.85 1]; a=[0 0 1 1 0 0 1 1 0 0]; wt=[2 1 2 1 2]; b=remez(N-1,f,a,wt); [h,w]=freqz(b,1);
plot(w/pi,20*log10(abs(h))); axis([0 1 -70 10]);grid; xlabel('归一化频率/\\pi') ylabel('幅度/dB')
title('雷米兹交替算法设计专用线性相位滤波器');
.学习帮手.
.专业整理.
雷米兹交替算法设计专用线性相位滤波器100-10幅度/dB-20-30-40-50-6000.10.20.30.60.50.4归一化频率/0.70.80.91
总结:利用雷米兹交替算法可以在通带波动和阻带的衰减上选择一个较好的平衡。
7. 利用雷米兹交替算法,设计一个线性相位高通FIR 数字滤波器,其指标为:通带边界频
率f c=800Hz,阻带边界频率f r=500Hz,通带波动δ=1dB,阻带最小衰减At=40dB,采样频率f s=5000Hz。
clear all; fedge=[500 800]; mval=[0 1]; dev=[0.01 0.109]; fs=5000;
[N,fpts,mag,wt]=remezord(fedge,mval,dev,fs); b=remez(N,fpts,mag,wt); [h,w]=freqz(b,1);
plot(w*2500/pi,20*log10(abs(h))); axis([0 2500 -60 10]);grid; xlabel('频率/Hz') ylabel('幅度/dB')
title('雷米兹交替算法设计线性相位高通FIR 数字滤波器');
.学习帮手.
.专业整理.
雷米兹交替算法设计线性相位高通FIR 数字滤波器100-10幅度/dB-20-30-40-50-6005001000频率/Hz150020002500
结论:雷米兹算法是建立在频域采样的基础上,利用最大值最小化原理,在带内波动和阻带衰减上进行最优化处理,所以得到的滤波器的效果很好。
二、 思考题
1) 定性地说明用本实验程序设计的FIR 滤波器的3dB 截止频率在什么位置?它等于理想
频率响应Hd (ejω)的截止频率吗?
答:如果在时域通过矩形窗截取的方法,3dB截止频率并不等于理想频率响应的截止频率,但如果通过在时域采样的方法,3dB 截止频率和理想的截止频率有较好的吻合
2) 如果没有给定h(n)的长度N,而是给定了通带边缘截止频率ωc 和阻带临界频率ωp,以
及相应的衰减,能根据这些条件用窗函数法设计线性相位FIR 低通滤波器吗? 答:可以的,利用matlab fir1函数,将M取得大一些,就可以获得比较好的线性相位低通滤波器。
.学习帮手.