ms21_5.m

来自「《数字信号处理实验(MATLAB版)》这」· M 代码 · 共 75 行

%双线性变换法设计数字低通\高通\带通\带阻滤波器
%数字滤波器指标
wp=0.4*pi;           %滤波器的通带截止频率
ws=0.6*pi;           %滤波器的阻带截止频率
w0=sqrt(wp*ws);B=0.4*pi;  %
Rp=1;As=20;          %输入滤波器的通阻带衰减指标
ripple=10^(-Rp/20);  %计算通带衰减对应的幅度值
Attn=10^(-As/20);    %计算阻带衰减对应的幅度值
%转换为模拟原型滤波器指标
Fs=2;T=1/Fs;
Omgp=(2/T)*tan(wp/2); %原型通带频率预修正
Omgs=(2/T)*tan(ws/2); %原型阻带频率预修正
Omg0=(2/T)*tan(w0/2); 
B0=(2/T)*tan(B/2);
%模拟原型滤波器计算
[n,Omgc]=cheb1ord(Omgp,Omgs,Rp,As,'s') %计算阶数n和截止频率
[z0,p0,k0]=cheb1ap(n,Rp); %归一化切比雪夫1型原型设计
ba=k0*real(poly(z0));     %求原型滤波器系数b
aa=real(poly(p0));        %求原型滤波器系数a
%由模拟原型滤波器变换为须求解的模拟滤波器
[ba1,aa1]=lp2lp(ba,aa,Omgc);    %变换为模拟低通滤波器
[ba2,aa2]=lp2hp(ba,aa,Omgs);    %变换为模拟高通滤波器
[ba3,aa3]=lp2bp(ba,aa,Omg0,B0);  %变换为模拟带通滤波器
[ba4,aa4]=lp2bs(ba,aa,Omg0,B0);  %变换为模拟带阻滤波器
%求模拟滤波器的频率特性
[Ha1,Omg1]=freqs(ba1,aa1);
[Ha2,Omg2]=freqs(ba2,aa2);
[Ha3,Omg3]=freqs(ba3,aa3);
[Ha4,Omg4]=freqs(ba4,aa4);    
figure(1);
subplot(2,2,1),plot(Omg1/pi,abs(Ha1));axis([0,Fs,0,1.1]);  
set(gca,'XTickMode','manual','XTick',[0,Omgp/pi,Omgs/pi,Fs]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
title('模拟低通滤波器');
subplot(2,2,2),plot(Omg2/pi,abs(Ha2));axis([0,Fs,0,1.1]);    
set(gca,'XTickMode','manual','XTick',[0,Omgp/pi,Omgs/pi,Fs]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
title('模拟高通滤波器');
subplot(2,2,3),plot(Omg3/pi,abs(Ha3));axis([0,Fs,0,1.1]);   
set(gca,'XTickMode','manual','XTick',[0,(Omg0-B0/2)/pi,Omg0/pi,(Omg0+B0/2)/pi,Fs]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
title('模拟带通滤波器');
subplot(2,2,4),plot(Omg4/pi,abs(Ha4));axis([0,Fs,0,1.1]);  
set(gca,'XTickMode','manual','XTick',[0,(Omg0-B0/2)/pi,Omg0/pi,(Omg0+B0/2)/pi,Fs]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
title('模拟带阻滤波器');
%用双线性变换法将模拟滤波器变换成数字滤波器, 计算数字滤波器系数
[bd1,ad1]=bilinear(ba1,aa1,Fs)   %双线性变换为数字低通滤波器
[bd2,ad2]=bilinear(ba2,aa2,Fs)   %双线性变换为数字高通滤波器
[bd3,ad3]=bilinear(ba3,aa3,Fs)   %双线性变换为数字带通滤波器
[bd4,ad4]=bilinear(ba4,aa4,Fs)   %双线性变换为数字带阻滤波器
%求数字滤波器的频率特性
[H1,w1]=freqz(bd1,ad1);   
[H2,w2]=freqz(bd2,ad2);   
[H3,w3]=freqz(bd3,ad3);   
[H4,w4]=freqz(bd4,ad4);   
%显示数字滤波器的幅频特性
figure(2);
subplot(2,2,1),plot(w1/pi,abs(H1));axis([0,1,0,1.1]);
set(gca,'XTickMode','manual','XTick',[0,wp/pi,ws/pi,1]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
title('数字低通滤波器');
subplot(2,2,2),plot(w2/pi,abs(H2));axis([0,1,0,1.1]);
set(gca,'XTickMode','manual','XTick',[0,wp/pi,ws/pi,1]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
title('数字高通滤波器');
subplot(2,2,3),plot(w3/pi,abs(H3));axis([0,1,0,1.1]);
set(gca,'XTickMode','manual','XTick',[0,(w0-B/2)/pi,w0/pi,(w0+B/2)/pi,1]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
title('数字带通滤波器');
subplot(2,2,4),plot(w4/pi,abs(H4));axis([0,1,0,1.1]);
set(gca,'XTickMode','manual','XTick',[0,(w0-B/2)/pi,w0/pi,(w0+B/2)/pi,1]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
title('数字带阻滤波器');