基于MATLAB布莱克曼窗FIR数字低通滤波器设计程序.docVIP

程序： clear all; Wp=1.34*pi/3;%input通带频率; Ws=7.4*pi/30;%input阻带频率; Ar=75; tr_width =Wp-Ws;%过渡带宽度 N=ceil(6.2*pi/tr_width);%滤波器长度，按布莱克曼窗计算所需的滤波器长度n=0:1:N-1; Wc=(Ws+Wp)/2;%理想低通滤波器的截止频率，截止频率取为两边缘频率的平均值hd=ideal_hp1(Wc,N);%理想低通率的单位冲击响应 w_ham=( blackman(N));string=[blackman,N=,num2str(N)];%调用布莱克曼窗函数 h=hd.*w_ham;%截取得到实际的单位冲击响应，为理想脉冲响应与窗函数乘积[db,mag,pha,w]=freqz_m2(h,[1]);%计算实际滤波器的幅度响应 %delta_w=2*pi/1000;% 频率分辨率Ap=-(min(db(Wp/delta_w+1:1:501)));%实际通带波纹 %As=-round(max(db(1:1:501)));%实际阻带波纹 subplot(4,2,1);stem(n,hd);title(理想脉冲响应hd(n)); axis([0,N-1,-0.5,1]);xlabel(n);ylabel(hd(n)) subplot(4,2,2);stem(n,w_ham);axis([0,N-1,0,1.1]); xlabel(n);ylabel(w(n));text(1.5,1.3,string); subplot(4,2,3);stem(n,h);title(实际脉冲响应h(b)); axis([0,N,-1.4,1.4]);xlabel(n);ylabel(h(n)); subplot(4,2,4);plot(w,pha);title(相频特性); axis([0,3.15,-4,4]);xlabel(频率(pi));ylabel(相位(Φ)); subplot(4,2,5);plot(w/pi,db);title(幅度特性(dB)); axis([0,1,-100,10]);xlabel(频率(pi));ylabel((分贝数)); subplot(4,2,6);plot(w,mag);title(频率特性); axis([0,3.15,0,1.5]);xlabel(频率(rad));ylabel(幅值); 调用子函数1：ideal_hp1 function hd=ideal_hp1(Wc,N) alpha=(N-1)/2; n=0:1:N-1; m=n-alpha+eps; hd=[sin(pi*m)-sin(Wc*m)]./(pi*m);%全通-低通； 调用子函数2：ideal_hp1%计算z域频率响应； function[db,mag,pha,w]=freqz_m2(b,a)% 频率的取值范围在0—1之间，标准化频率1对应 %数字频率为π，对应的模拟频率为采样频率的一半； [H,w]=freqz(b,a,1000,whole);% b直接型滤波器分子系数;a直接型滤波器分母系数 %db为[0:pi]之间的相对振幅;mag为[0 pi]区间的绝对值;pha为[0 pi]区间的相位响应H=(H(1:1:501)); w=(w(1:1:501)); mag=abs(H);%幅度 db=20*log10((mag+eps)/max(mag)); pha=angle(H);%相位 %用含多种频率成份离散时间信号作输入信号，得到通过滤波器后与输出信号的对比结果 fs=15000; t=(0:100)/fs; x=sin(2*pi*t*750)+sin(2*pi*t*1850)+sin(2*pi*t*2500)+sin(2*pi*t*3350)+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; subplot(4,2,7); plot(f1,abs(a)); title(输入波形频谱图); xlabel(频率/Hz);ylabel(幅度) subplot(4,2,8); plot(f2,abs(b)); title(输出波形频谱图); xlabel(频率/Hz);ylabel(幅度) % computation of s-domain frequency response:modified version %-------------------------------------------------------- % [dB,mag,pha,grd,w]=fr