ms21_2.m

《数字信号处理实验(MATLAB版)》这

%双线性变换法设计数字高通滤波器
%数字滤波器指标
fs=150;fp=250;Fs=1000;T=1/Fs;
wp=fp/Fs*2*pi;             %数字滤波器的通带截止频率
ws=fs/Fs*2*pi;             %数字滤波器的阻带截止频率
Rp=1;As=20;                %输入滤波器的通阻带衰减指标
ripple=10^(-Rp/20);
Attn=10^(-As/20);
%转换为模拟滤波器指标
Omgp=(2/T)*tan(wp/2);
Omgs=(2/T)*tan(ws/2);
%模拟原型滤波器计算
[n,Omgc]=ellipord(Omgp,Omgs,Rp,As,'s') %计算阶数n和截止频率
[z0,p0,k0]=ellipap(n,Rp,As);   %归一化椭圆原型设计
ba=k0*real(poly(z0));      %求原型滤波器系数b
aa=real(poly(p0));         %求原型滤波器系数a
[ba1,aa1]=lp2hp(ba,aa,Omgc);%变换为模拟高通滤波器
%注意,以上4行求滤波器系数ba1、aa1的程序，可由下一条程序替代
%[ba1,aa1]=butter(n,Omgc,'s');  %直接求模拟滤波器系数
%用双线性变换法计算数字滤波器系数
[bd,ad]=bilinear(ba1,aa1,Fs) %双线性变换
%求数字系统的频率特性
[H,w]=freqz(bd,ad);   
dbH=20*log10(abs(H)/max(abs(H))); %化为分贝值
%
subplot(2,2,1),plot(w/2/pi*Fs,abs(H),'k'); 
ylabel('|H|');title('幅度响应');axis([0,Fs/2,0,1.1]);
set(gca,'XTickMode','manual','XTick',[0,fs,fp,Fs/2]);
set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid
subplot(2,2,2),plot(w/2/pi*Fs,angle(H)/pi*180,'k');
ylabel('\phi');title('相位响应');axis([0,Fs/2,-180,180]);
set(gca,'XTickMode','manual','XTick',[0,fs,fp,Fs/2]);
set(gca,'YTickMode','manual','YTick',[-180,0,180]);grid
subplot(2,2,3),plot(w/2/pi*Fs,dbH,'k');
title('幅度响应(dB)');axis([0,Fs/2,-40,5]); 
ylabel('dB');xlabel('频率(\pi)');
set(gca,'XTickMode','manual','XTick',[0,fs,fp,Fs/2]);
set(gca,'YTickMode','manual','YTick',[-50,-20,-1,0]);grid
subplot(2,2,4),zplane(bd,ad);
axis([-1.1,1.1,-1.1,1.1]);title('零极图');