四、实验的程序代码 MATLAB提供了[z,p,k]=buttap(N)函数来设计归一化的巴特沃思模拟滤波器，它返回零点数组z和极点数组p,以及增益k。我们借助该函数重新定义一个未归一化的巴特沃思滤波器函数。然后编写一个这个函数文件，以文件名u_buttap.m保存在当前路径，具体代码如下： function [b,a] = u_buttap(N,Omegac); [z,p,k] = buttap(N); p = p*Omegac; k = k*Omegac^N; B = real(poly(z)); b = k*B; b0= k; a = real(poly(p)); 我们定义一个函数[b,a] =afd_butt(Omegp,Omegs,Ap,As)，功能为：由参数Omegp,Omegs,Ap,As计算OmegC和阶数N。编写该函数文件，以文件名afd_butt.m保存在当前路径，具体代码如下： function [b,a] =afd_butt(Omegp,Omegs,Ap,As); if Omegp<= 0 error('Passband edge must be larger than 0') end if Omegs<=Omegp error('Stopband edge must be larger than Passband edge') end if (Ap <= 0) | (As < 0) error('PB ripple and/or SB attenuation must be larger than 0') end N = ceil((log10((10^(Ap/10)-1)/(10^(As/10)-1)))/(2*log10(Omegp/Omegs)) ); fprintf('\n*** Butterworth Filter Order = %2.0f \n',N) OmegaC = Omegp/((10^(Ap/10)-1)^(1/(2*N))); [b,a]=u_buttap(N,OmegaC); 根据实验内容，首先编写一个脚本文件，用脉冲响应不变法设计数字低通滤波器。以文件名shiyan4_1.m保存在当前路径，具体代码如下： %----------------- DF指标------------------------------------------ wp=0.4*pi; ws=0.6*pi; Ap=0.5; As=50; %-----------------模拟滤波器原形指标的频率逆映射------------------- T=1;Fs=1/T; OmegaP=wp/T; OmegaS=ws/T; [cs,ds]=afd_butt(OmegaP,OmegaS,Ap,As); %巴特沃思原型滤波器计算----- %------------模拟滤波器到数字滤波器的脉冲响应不变法变换------------ [b,a]=impinvar(cs,ds,Fs); [h,w]=freqz(b,a); subplot(2,2,1);plot(w/pi,abs(h));title('幅度响应');grid; subplot(2,2,2);plot(w/pi,angle(h));title('相位响应');grid; subplot(2,2,3);plot(w/pi,20*log10(abs(h)));title('幅度响应dB');grid; n=[0:1:59];imp=[1;zeros(59,1)];y=filter(b,a,imp); subplot(2,2,4);plot(n,y);title('脉冲响应');grid; 然后编写一个脚本文件，用双线性法设计数字低通滤波器。以文件名shiyan4_2.m保存在当前路径，具体代码如下： clear; close all; wp=0.4*pi; ws=0.6*pi; Ap=0.5; As=50; T=1;Fs=1/T; OmegaP=(2/T)*tan(wp/2); OmegaS=(2/T)*tan(ws/2); [cs,ds]=afd_butt(OmegaP,OmegaS,Ap,As); %------------模拟滤波器到数字滤波器的双线性法变换------------ [b,a]=bilinear(cs,ds,Fs); [h,w]=freqz(b,a); subplot(2,2,1);plot(w/pi,abs(h));title('幅度响应');grid on; subplot(2,2,3);plot(w/pi,20*log10(abs(h)));title('幅度响应dB');axis([0,1,-80,5]);grid; n=[0:1:59];imp=[1;zeros(59,1)];y=filter(b,a,imp); subplot(2,2,4);plot(n,y);title('脉冲响应');grid; 然后编写一个脚本文件，用双线性法设计数字高通滤波器。以文件名shiyan4_3.m保存在当前路径，具体代码如下： clear; close all; f=10;fp=2;fs=1.5; wp=2*pi*fp*(1/f); ws=2*pi*fs*(1/f); Ap=3; As=40; [N,wn]=buttord(wp/pi,ws/pi,Ap,As); [b,a]=butter(N,wn,'high'); [h,w]=freqz(b,a); subplot(2,2,1);plot(w/pi,abs(h));title('幅度响应');grid; subplot(2,2,3);plot(w/pi,20*log10(abs(h)));title('幅度响应dB');axis([0,1,-80,5]);grid; n=[0:1:59];imp=[1;zeros(59,1)];y=filter(b,a,imp); subplot(2,2,4);plot(n,y);title('脉冲响应');grid; n=[0:1000];t=0.0001*n; x=2*sin(2*pi*40*t)+5*sin(2*pi*4000*t); y=filter(b,a,x); figure(2); subplot(211); plot(n,x,'p-'); subplot(212); plot(n,y,'p-'); 本文来源：https://www.wddqw.com/doc/4dfcf883f78a6529647d538a.html