数字带通巴特沃思滤波器的设计-实验报告-中国地质大学(北京)-数字信号处理-陈玉东

中国地质大学(北京)实验报告题目：数字低通巴特沃思滤波器的设计学生姓名学号院（系）专业指导教师陈玉东提交日期201？年1月9日1实验四数字带通巴特沃思滤波器的设计1．实验目的基于MATLAB环境，掌握数字滤波器的设计方法。2.实验环境硬件环境:计算机,软件环境：MATLAB平台。3.实验原理利用原型设计法，对满足指标要求的模拟低通原型滤波器，通过双线形变换和数字频率变换，得到满足指标要求的数字带通滤波器。4.设计流程图5.实验内容和要求1．编写设计数字带通巴特沃思滤波器的相应函数文件。2．编写.m程序文件，通过调用相应函数，计算满足下列指标要求的数字带通滤波器：通带允许起伏Rp=1dB，p1=0.2，p2=0.3。阻带衰减As=25dB，s1=0.1，s2=0.4。预畸变求模拟带通滤波器指标输入给定的数字滤波器指标转换成模拟低通原型指标确定模拟低通原型阶数和截止频率[N,c]=buttord(p,s,Rp,As)确定模拟低通原型的H(s)[bs,as]=butter(N,Wn,'s')用[b,a]=lp2bp(bs,as,W0,Bw)在模拟域进行频带转换变成模拟带通滤波器，然后用双线性变换[b,a]=bilinear(b,a,T)求得数字带通滤波器。调用impz.m求h(n)，并求出H(z)的mag,db,phase,grpdelay。保存所写.m文件打印所有结果显示N，Ωc显示b，a，并整理成文档，以备写报告时用26.程序实现%%由巴特沃斯模拟低通滤波器映射成数字带通滤波器的程序%相应参数及指标见下面注释%%%%%%%%%%%%%%%%%%%%%%%%%%日期班级姓名学号%wp1=0.2*pi;%通带的下边缘截止频率%wp2=0.3*pi;%通带的上边缘截止频率%ws1=0.1*pi;%阻带的下边缘截止频率%ws2=0.4*pi;%阻带的上边缘截止频率%Rp=1;%通带的波纹指标%As=25;%阻带的衰减指标%设计得到的数字带通滤波器的系统函数分子系数为%%b=[0.00110-0.004400.00660-0.004400.0011];%设计得到的数字带通滤波器的系统函数分母系数为%%a=[1-4.961012.2137-18.816319.7331-14.35117.1024-2.19880.3386]；clc,clearall;wp1=0.2*pi;wp2=0.3*pi;ws1=0.1*pi;ws2=0.4*pi;T=1;%约定%预畸变Omgp1=(2/T)*tan(wp1/2);Omgp2=(2/T)*tan(wp2/2);Omgs1=(2/T)*tan(ws1/2);Omgs2=(2/T)*tan(ws2/2);%确定模拟带通滤波器的指标if(Omgp1*Omgp2-Omgs1*Omgs2)epsOmgs1=Omgp1*Omgp2/Omgs2elseif(Omgp1*Omgp2-Omgs1*Omgs2)epsOmgs2=Omgp1*Omgp2/Omgs1endOmgpWidth=Omgp2-Omgp1;%通带的宽度Omg0=sqrt(Omgp2*Omgp1);%通带的中心频率3%模拟低通滤波器原型确定Wplp=1;Wslp=(Omg0^2-Omgs1^2)*Wplp/(OmgpWidth*Omgs1)fprintf('\nWp=%fWs=%f\n',Wplp,Wslp);Rp=1;%通带的波纹指标As=25;%阻带的衰减指标Ripple=10^(-Rp/20)%通带的波纹Attn=10^(-As/20)%阻带的衰减[N,Wn]=buttord(Wplp,Wslp,Rp,As,'s');[bs,as]=butter(N,Wn,'s')WRange=3.5;w=[0:500]*WRange/500;H=freqs(bs,as,w);mag0=abs(H);w_edges=[0Wplp,WslpWRange];%由模拟低通映射生成模拟带通滤波器[num_s,den_s]=lp2bp(bs,as,Omg0,OmgpWidth)bs=num_s;as=den_s;H=freqs(bs,as,w);mag=abs(H);db=20*log10((mag+eps)/max(mag));phase=angle(H);%plotsOmega=[Omgs1Omgp1Omgp2Omgs2];omega=[ws1wp1wp2ws2]/pi;figure(1);clf;XTickLimit=4;dBLimit=40;subplot(221);plot(w,mag0);title('模拟低通滤波器原型频率响应');ylabel('|Ha(s)|');xlabel('模拟角频率（弧度/秒）');axis([0,XTickLimit,-0.1,1.1]);set(gca,'XTickMode','manual','XTick',[0,Wplp,Wslp4,XTickLimit]);set(gca,'YTickMode','manual','Ytick',[0,Attn,Ripple,1]);grid;subplot(222);plot(w,mag);title('模拟带通滤波器原型频率响应');ylabel('|Ha(s)|');xlabel('模拟角频率（弧度/秒）');axis([0,XTickLimit,-0.1,1.1]);set(gca,'XTickMode','manual','XTick',[0,Omega,XTickLimit]);set(gca,'YtickMode','manual','YTick',[0,Attn,Ripple,1]);grid;subplot(223);plot(w,db);title('模拟带通滤波器原型频率响应（dB）');xlabel('模拟角频率（弧度/秒）');ylabel('分贝');axis([0,XTickLimit,-dBLimit,1])set(gca,'XTickMode','manual','XTick',[0,Omega,XTickLimit]);set(gca,'YTickMode','manual','Ytick',[-dBLimit,-As,-Rp,0]);grid;subplot(224);plot(w,phase/pi);title('模拟带通滤波器相位响应');ylabel('\pi弧度');xlabel('模拟角频率（弧度/秒）');axis([0,XTickLimit,-1.1,1.1]);set(gca,'XTickMode','manual','XTick',[0,Omega,XTickLimit]);grid;h=findobj(gcf,'Color','b');set(h,'Color','r','linewidth',1.5)%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%[b,a]=bilinear(bs,as,1)%%plots%[H,w]=freqz(b,a,1000,'whole');5H=(H(1:501))';w=(w(1:501))';mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);grd=grpdelay(b,a,w);figure(2);clf;subplot(221);plot(w/pi,mag);title('数字带通滤波器频率响应');ylabel('|H(z)|');xlabel('数字频率\pi（弧度）');axis([0,1,-0.1,1.1]);set(gca,'XTickMode','manual','XTick',[0,omega,1]);set(gca,'YTickMode','manual','YTick',[0,Attn,Ripple,1]);grid;subplot(223);plot(w/pi,db)title('数字带通滤波器频率响应（dB）');xlabel('数字频率\pi（弧度）');ylabel('分贝');axis([0,1,-dBLimit,1.1]);set(gca,'XTickMode','manual','XTick',[0,omega,1]);set(gca,'YTickMode','manual','YTick',[-dBLimit,-As,-Rp,0]);grid;subplot(222);plot(w/pi,phase/pi);title('数字带通滤波器相位响应');xlabel('数字频率\pi(弧度)')ylabel('\pi(弧度)');axis([0,1,-1.1,1.1]);set(gca,'XTickMode','manual','XTick',[0,omega,1]);set(gca,'YTickMode','manual','YTick',[-1,0,1]);grid;subplot(224);plot(w/pi,grd);title('群延时');6xlabel('数字频率\pi（弧度）');ylabel('样本数');delayLimit=max(grd);axis([0,1,0,delayLimit]);set(gca,'XTickMode','manual','XTick',[0,omega,1]);set(gca,'YTickMode','manual','YTick',[0:4:delayLimit]);grid;h=findobj(gcf,'Color','b');set(h,'Color','m','linewidth',2)%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%figure(3)clf;subplot(211);zplane(b,a);title('零极点分布图');subplot(212);[h,t]=impz(b,a,70);stem(h);title('单位冲击响应');h=findobj(gcf,'Color','b');set(h,'Color','m','linewidth',1.5)%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%[N,Wn]=buttord([Omgp1,Omgp2],[Omgs1,Omgs2],Rp,As,'s')[bs0,as0]=butter(N,Wn,'bandpass','s');w=[0:500]*WRange/500;H=freqs(bs,as,w);mag0=abs(H);fileID=fopen('Test_4_H,txt','w');fprintf(fileID,'%s','b=');fprintf(fileID,'%12.8f',b);fprintf(fileID,'\n\n%s','a=');fprintf(fileID,'%12.8f',a);fclose(fileID);77.图形绘制012.7013400.05620.89131模拟低通滤波器原型频率响应|Ha(s)|模拟角频率（弧度/秒）00.45570.64981.01911.4531400.05620.89131模拟带通滤波器原型频率响应|Ha(s)|模拟角频率（弧度/秒）00.45570.64981.01911.45314-40-25-10模拟带通滤波器原型频率响应（dB）模拟角频率（弧度/秒）分贝00.45570.64981.01911.45314-1-0.500.51模拟带通滤波器相位响应弧度模拟角频率（弧度/秒）800.10.20.30.4100.05620.89131数字带通滤波器频率响应|H(z)|数字频率（弧度）00.10.20.30.41-40-25-10数字带通滤波器频率响应（dB）数字频率（弧度）分贝00.10.20.30.41-101数字带