窗函数法设计数字滤波器

1.基于矩形窗的FIR数字低通滤波器设计参数如下：主程序：% 基于矩形窗的FIR数字低通滤波器设计wp=2*pi*1.75/15;ws=2*pi*3.25/15w=(wp+ws)/2;N=ceil(4*pi/(ws-wp));n=0:1:(N-1);alpha=(N-1)/2;m=n-alpha+eps;hd=sin(w*m)./(pi*m); %得到理想低通滤波器B=boxcar(N);string=['Boxcar,N=',num2str(N)];h=hd.*(B)'; %得到FIR数字滤波器[H,w]=freqz(h,[1],1024); %求其频率响应db=20*log10(abs(H)+eps); %得到幅值pha=angle(H); %得到相位figure(1);subplot(2,2,1);stem(n,h,'r');axis([0,N-1,-0.1,0.35]);xlabel('n');ylabel('h(n)');title('实际低通滤波器的h(n)');text((0.3*N),0.275,string);subplot(2,2,2);plot(w/pi,db,'r');axis([0,1,-100,5]);xlabel('w/pi');ylabel('dB');title('衰减特性(dB)');grid;subplot(2,2,3);plot(w,pha);hold on;plot(0:4,zeros(5),'k');title('相频特性');xlabel('频率(rad)');ylabel('相位(rad)');axis([0,3.2,-4,4]);subplot(2,2,4);plot(w,abs(H));title('频率特性');xlabel('频率W(rad)');ylabel('幅值');axis([0,3.15,0,1.5]);text(0.9,1.3,string);2.基于矩形窗的带通滤波器的设计参数如下：1/8子程序：function hd=ideal_bs(Wcl,Wch,N)alpha=(N-1)/2;n=0:1:N-1;m=n-alpha+eps;hd=[sin(Wch*m)-sin(Wcl*m)]./(pi*m);function hd=ideal_bs(Wcl,Wch,N)alpha=(N-1)/2;n=0:1:N-1;m=n-alpha+eps;hd=[sin(Wch*m)-sin(Wcl*m)]./(pi*m);主程序：clear all;Wph=2*pi*4.75/15;Wpl=2*pi*2.5/15;Wsl=2*pi/15;Wsh=2*pi*6.25/15;tr_width=min((Wpl-Wsl),(Wsh-Wph)); %过渡带宽度N=ceil(4*pi/tr_width); %滤波器长度n=0:1:N-1;Wcl=(Wsl+Wpl)/2; %理想滤波器的截止频率Wch=(Wsh+Wph)/2;hd=ideal_bs(Wcl,Wch,N); %理想滤波器的单位冲击响应w_ham=(boxcar(N))';string=['矩形窗','N=',num2str(N)];h=hd.*w_ham; %截取取得实际的单位脉冲响应[db,mag,pha,w]=freqz_m2(h,[1]); %计算实际滤波器的幅度响应subplot(3,2,1);stem(n,hd);title('理想脉冲响应hd(n)')axis([-1,N,-0.5,0.8]);xlabel('n');ylabel('hd(n)');subplot(3,2,2);stem(n,w_ham);axis([-1,N,0,1.1]);xlabel('n');ylabel('w(n)');text(1.5,1.3,string);subplot(3,2,3);stem(n,h);title('实际脉冲响应h(n)');axis([0,N,-1.4,1.4]);xlabel('n');ylabel('h(n)');subplot(3,2,4);plot(w,pha);title('相频特性');axis([0,3.15,-4,4]);xlabel('频率（rad）');ylabel('相位（?）'); subplot(3,2,5);plot(w/pi,db);title('幅度特性（dB）');axis([0,1,-80,10]);xlabel('频率（pi）');ylabel('分贝数');subplot(3,2,6);plot(w,mag);title('频率特性')axis([0,3.15,0,1.5]);xlabel('频率（rad）');ylabel('幅值');2/8fs=15000;t=(0:100)/fs;x=sin(2*pi*t*750)+sin(2*pi*t*5000)+sin(2*pi*t*6100); q=filter(h,1,x); [a,f1]=freqz(x);f1=f1/pi*fs/2;[b,f2]=freqz(q);f2=f2/pi*fs/2;figure(2);subplot(2,1,1);plot(f1,abs(a));title('输入波形频谱图');xlabel('频率');ylabel('幅度')subplot(2,1,2);plot(f2,abs(b));title('输出波形频谱图');xlabel('频率');ylabel('幅度')3.基于矩形窗的FIR数字带阻滤波器设计参数：图片加载失败点击图片重新加载主程序：clear all;Wph=2*pi*6.25/15;Wpl=2*pi/15;Wsl=2*pi*2.5/15;Wsh=2*pi*4.75/15;tr_width=min((Wsl-Wpl),(Wph-Wsh));%过渡带宽度N=ceil(4*pi/tr_width);%滤波器长度n=0:1:N-1;Wcl=(Wsl+Wpl)/2;%理想滤波器的截止频率Wch=(Wsh+Wph)/2;hd=ideal_bs(Wcl,Wch,N);%理想滤波器的单位冲击响应w_ham=(boxcar(N))';string=['矩形窗','N=',num2str(N)];h=hd.*w_ham;%截取取得实际的单位脉冲响应[db,mag,pha,w]=freqz_m2(h,[1]);%计算实际滤波器的幅度响应delta_w=2*pi/1000;subplot(3,2,1);stem(n,hd);title('理想脉冲响应hd(n)')axis([-1,N,-0.5,0.8]);xlabel('n');ylabel('hd(n)');subplot(3,2,2);stem(n,w_ham);axis([-1,N,0,1.1]);3/8xlabel('n');ylabel('w(n)');text(1.5,1.3,string);subplot(3,2,3);stem(n,h);title('实际脉冲响应h(n)');axis([0,N,-1.4,1.4]);xlabel('n');ylabel('h(n)');subplot(3,2,4);plot(w,pha);title('相频特性');axis([0,3.15,-4,4]);xlabel('频率（rad）');ylabel('相位（Φ）'); subplot(3,2,5);plot(w/pi,db);title('幅度特性（dB）');axis([0,1,-80,10]);xlabel('频率（pi）');ylabel('分贝数'); subplot(3,2,6);plot(w,mag);title('频率特性')axis([0,3,0,2]);xlabel('频率（rad）');ylabel('幅值');fs=15000;t=(0:100)/fs;x=sin(2*pi*t*750)+sin(2*pi*t*3000)+sin(2*pi*t*6100);q=filter(h,1,x);[a,f1]=freqz(x);f1=f1/pi*fs/2;[b,f2]=freqz(q);f2=f2/pi*fs/2;figure(2);subplot(2,1,1);plot(f1,abs(a));title('输入波形频谱图');xlabel('频率');ylabel('幅度')subplot(2,1,2);plot(f2,abs(b));title('输出波形频谱图');xlabel('频率');ylabel('幅度')；调用程序1：function hd=ideal_bs(Wcl,Wch,m);alpha=(m-1)/2;n=[0:1:(m-1)];m=n-alpha+eps;hd=[sin(m*pi)+sin(Wcl*m)-sin(Wch*m)]./(pi*m)调用程序2：function[db,mag,pha,w]=freqz_m2(b,a)[H,w]=freqz(b,a,1000,'whole');H=(H(1:1:501))'; w=(w(1:1:501))';mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);。