MATLAB实现数字FIR的高通-和带通等滤波器的源程序.doc

MATLAB实现数字FIR的高通_和带通等滤波器的源程序 作者： 日期：15 个人收集整理 勿做商业用途利用汉宁窗设计型数字高通滤波器clear all;Wp=0。6pi;Ws=0.4pi；tr_width=WpWs; 过渡带宽度N=ceil(6.2pi/tr_width） %滤波器长度n=0：1:N-1;Wc=(Ws+Wp）/2; %理想低通滤波器的截止频率hd=ideal_hp1（Wc,N）; %理想低通滤波器的单位冲激响应w_han=（hanning（N）)； %汉宁窗h=hd.*w_han; 截取得到实际的单位脉冲响应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:Ws/delta_w+1）) 实际阻带纹波subplot(221）stem（n,hd)title(理想单位脉冲响应hd(n)subplot(222）stem（n，w_han）title（汉宁窗w（n)subplot(223)stem(n,h）title（实际单位脉冲响应hd（n)）subplot(224）plot（w/pi，db)title(幅度响应(dB）axis（0,1，100,10）clear all；Wp=0。6*pi;Ws=0。4pi；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_han=(hanning（N)； %汉宁窗h=hd。*w_han； 截取得到实际的单位脉冲响应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：Ws/delta_w+1)）) %实际阻带纹波subplot（221）stem（n,hd）title（理想单位脉冲响应hd（n)subplot(222）stem（n,w_han）title（汉宁窗w(n）)subplot（223)stem（n,h)title(实际单位脉冲响应hd(n）subplot（224）plot(w/pi,db)title（幅度响应（dB)axis（0,1，-100，10)基于切比雪夫一致逼近法设计FIR数字低通滤波器clear all;f=0 0。6 0.7 1; %给定频率轴分点A=1 1 0 0; %给定在这些频率分点上理想的幅频响应weigh=1 10； 给定在这些频率分点上的加权b=remez(32，f,A，weigh)； 设计出切比雪夫最佳一致逼近滤波器h,w=freqz(b,1，256，1)；h=abs（h)；h=20*log10（h）；subplot（211）stem（b,。);grid；title（切比雪夫逼近滤波器的抽样值）subplot(212）plot（w，h);grid；title（滤波器幅频特性(dB)利用汉宁窗设计型数字带阻滤波器clear all；Wpl=0.2*pi;Wph=0。8pi；Wsl=0。4*pi；Wsh=0。6*pi；tr_width=min(WslWpl）,（WphWsh）)； 过渡带宽度N=ceil（6。2pi/tr_width） 滤波器长度n=0：1:N1;Wcl=(Wsl+Wpl)/2; %理想低通滤波器的截止频率Wch=（Wsh+Wph)/2;hd=ideal_bs(Wcl,Wch,N）; 理想低通滤波器的单位冲激响应w_hann=(hanning(N）； 汉宁窗h=hd.*w_hann； 截取得到实际的单位脉冲响应db，mag，pha，w=freqz_m2（h,1)； %计算实际滤波器的幅度响应delta_w=2pi/1000;Ap=（min（db（1:1：Wpl/delta_w+1）) %实际通带纹波As=-round(max（db（Wsl/delta_w+1：1：Wsh/delta_w+1)）) 实际阻带纹波subplot(221)stem(n,hd)title（理想单位脉冲响应hd(n)）subplot（222）stem（n,w_hann）title(汉宁窗w（n)subplot(223)stem(n，h)title（实际单位脉冲响应hd(n）)subplot(224）plot(w/pi,db)title(幅度响应（dB）axis（0，1,100,10）本文为互联网收集，请勿用作商业用途个人收集整理，勿做商业用途利用三角窗设计型数字带通滤波器clear all;Wpl=0.4pi；Wph=0。6pi；Wsl=0.2*pi;Wsh=0。8*pi;tr_width=min（(WplWsl),（Wsh-Wph）; 过渡带宽度N=ceil（6.1*pi/tr_width） 滤波器长度n=0:1:N-1；Wcl=(Wsl+Wpl)/2; 理想低通滤波器的截止频率Wch=(Wsh+Wph）/2；hd=ideal_bp2（Wcl,Wch，N); 理想低通滤波器的单位冲激响应w_tri=(triang（N）)； %三角窗h=hd。w_tri; %截取得到实际的单位脉冲响应db，mag,pha，w=freqz_m2(h,1); 计算实际滤波器的幅度响应delta_w=2*pi/1000;Ap=-(min(db(Wpl/delta_w+1:1：Wph/delta_w+1)）) 实际通带纹波As=-round(max（db（Wsh/delta_w+1:1：501)） %实际阻带纹波subplot（221）stem(n,hd)title(理想单位脉冲响应hd(n）)subplot(222）stem（n,w_tri）title(三角窗w(n)subplot（223)stem(n,h）title（实际单位脉冲响应hd（n）subplot(224）plot(w/pi，db）title(幅度响应（dB）axis（0，1,-100，10)文档为个人收集整理，来源于网络文档为个人收集整理，来源于网络利用布拉克曼窗设计型数字带通滤波器clear all;Wpl=0。4pi；Wph=0。6*pi；Wsl=0。2*pi；Wsh=0。8pi;tr_width=min(（WplWsl）,（WshWph）； %过渡带宽度N=ceil(11*pi/tr_width）+1 %滤波器长度n=0:1：N1；Wcl=(Wsl+Wpl）/2； %理想低通滤波器的截止频率Wch=(Wsh+Wph)/2;hd=ideal_bp1(Wcl,Wch,N）； %理想低通滤波器的单位冲激响应w_bman=(blackman(N）; 布拉克曼窗h=hd。*w_bman; 截取得到实际的单位脉冲响应db,mag，pha，w=freqz_m2(h，1); %计算实际滤波器的幅度响应delta_w=2*pi/1000;Ap=(min（db（Wpl/delta_w+1:1：Wph/delta_w+1) %实际通带纹波As=round（max（db(Wsh/delta_w+1:1:501）) %实际阻带纹波subplot(221）stem（n，hd）title（理想单位脉冲响应hd（n）subplot（222）stem(n，w_bman）title(布拉克曼窗w（n）subplot(223)stem(n，h）title（实际单位脉冲响应hd（n)subplot(224）plot（w/pi,db）title(幅度响应(dB）axis(0,1,-100，10)文档为个人收集整理,来源于网络本文为互联网收集，请勿用作商业用途利用海明窗设计型数字低通滤波器clear all;Wp=0。2*pi；Ws=0.4pi;tr_width=Ws-Wp; 过渡带宽度N=ceil（6。6pi/tr_width)+1 滤波器长度n=0：1：N1；Wc=(Ws+Wp)/2; 理想低通滤波器的截止频率hd=ideal_lp1（Wc,N)； 理想低通滤波器的单位冲激响应w_ham=(hamming(N）); %海明窗h=hd.*w_ham； 截取得到实际的单位脉冲响应db，mag,pha,w=freqz_m2（h，1); %计算实际滤波器的幅度响应delta_w=2pi/1000；Ap=（min（db（1:1:Wp/delta_w+1)）) 实际通带纹波As=-round（max（db（Ws/delta_w+1：1:501）) %实际阻带纹波subplot(221）stem(n，hd）title(理想单位脉冲响应hd(n）)subplot(222)stem（n，w_ham)title(海明窗w（n)subplot(223)stem(n,h)title（实际单位脉冲响应hd（n）)subplot（224)plot（w/pi，db)title(幅度响应（dB）axis(0，1,-100,10）%-functiondb，mag，pha，w=freqz_m2(b,a)滤波器的幅值响应（相对、绝对）、相位响应%db:相对幅值响应mag:绝对幅值响应%pha： 相位响应w 采样频率；b 系统函数H（z)的分子项（对FIR，b=h）%a 系统函数H（z)的分母项(对FIR,a=1)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）； %相位响应本文为互联网收集,请勿用作商业用途文档为个人收集整理，来源于网络利用模拟Butterworth滤波器设计数字低通滤波器 exa4-8_pulseDF for example4-8 using Butterworth analog lowpass filter to design digital lowpass filter%利用模拟Butterworth滤波器设计数字低通滤波器脉冲响应不变法wp=0.2*pi;ws=0。3*pi;Rp=1;As=15;T=1;%性能指标Rip=10（Rp/20）；Atn=10（-As/20)；OmgP=wpT；OmgS=ws*T;N,OmgC=buttord（OmgP,OmgS,Rp,As,s);%选取模拟滤波器的阶数cs,ds=butter(N,OmgC,s）; 设计出所需的模拟低通滤波器b，a=impinvar(cs，ds，T）; %应用脉冲响应不变法进行转换求得相对、绝对频响及相位、群迟延响应db，mag,pha，grd，w=freqz_m（b，a）;%下面绘出各条曲线subplot(2,2,1)；plot（w/pi，mag）；title（幅频特性)；xlabel（w(/pi）);ylabel（|H(jw)；axis(0，1,0,1.1);set(gca,XTickMode,manual,XTick,0 0.2 0.3 0.5 1）;set(gca，YTickMode,manual，YTick，0 Atn Rip 1）;gridsubplot(2，2,2）;plot(w/pi,db);title（幅频特性(dB)）；xlabel(w(/pi);ylabel(dB）;axis（0，1,-40，5)；set(gca,XTickMode,manual，XTick,0 0.2 0。3 0.5 1）;set(gca，YTickMode,manual，YTick,-40 -As Rp 0）；gridsubplot(2,2，3)；plot(w/pi，pha/pi);title（相频特性）;xlabel（w（/pi）;ylabel(pha（/pi）);axis(0，1,1，1）;set(gca，XTickMode，manual，XTick,0 0。2 0。3 0。5 1）;gridsubplot（2,2，4);plot（w/pi，grd）;title（群延迟）；xlabel(w(/pi）);ylabel(Sample）;axis(0，1,0,12）;set（gca,XTickMode,manual,XTick，0 0.2 0。3 0.5 1)；gridfunctiondb，mag，pha，grd，w=freqz_m(b,a）%滤波器幅值响应（绝对、相对)、相位响应及群延迟Usage： db,mag，pha,grd，w=freqz_m(b，a) %500点对应0，pi%db 相对幅值响应；mag 绝对幅值响应；pha 相位响应；grd 群延迟响应w 采样频率；b 系统函数H（z)的分子项（对FIR，bh)%a 系统函数H(z）的分母项(对FIR，a1）H，w=freqz(b，a,500)；%500点的复频响应mag=abs（H);db=20*log10（mag+eps）/max（mag)）；pha=angle（H);grd=grpdelay(b,a,w）;本文为互联网收集，请勿用作商业用途文档为个人收集整理，来源于网络基于频域抽样法的FIR数字带阻滤波器设计clear all；N=41；T1=0.598；alpha= （N-1)/2；l=0：N-1；wl= （2*pi/N)l;Hrs=ones(1,6)，T1，zeros(1，7）,T1，ones（1,11）,T1，zeros（1，7），T1,ones(1，6); %理想振幅采样响应Hdr=1，1,0，0,1,1;wdl=0，0。3,0.3，0。7,0.7,1；k1=0:floor（(N1）/2）;k2=floor（N-1）/2）+1:N1；angH=pi/2-alpha*（2pi）/N*(k1+0。5),-pi/2+alpha*（2*pi)/N*(Nk2-0。5）; %相位约束条件Hdk=Hrs.exp（jangH）； 构成Hd（k）h1=ifft(Hdk，N）；n=0：1:N-1；h=real（h1。*exp(j*pi*n/N）)； %实际单位冲激响应db,mag，pha,w=freqz_m2(h,1);Hr,ww，a，L=hr_type3（h）； 实际振幅响应subplot（221）plot(wl/pi+1/N,Hrs，.,wdl，Hdr）title（频率样本Hd（k) :N=41)axis（0 1 -0.1 1。2)subplot（222）stem(l，h）title(实际单位脉冲响应h(n)subplot(223)plot(ww/pi，Hr，wl/pi+1/N,Hrs,.）title（实际振幅响应H（w）axis（0 1 0.1 1。2）subplot（224）plot(w/pi,db)title（幅度响应（dB）)axis（0 1 -80 10）function db，mag,pha，w = freqz_m(b，a）；滤波器的幅值响应(相对、绝对）、相位响应db：相对幅值响应%mag：绝对幅值响应%pha: 相位响应w 采样频率；b 系统函数H（z)的分子项(对FIR，b=h)a 系统函数H（z）的分母项(对FIR，a=1)H，w = freqz（b,a,1000,whole）;H = (H（1：1：501）)； w = （w(1：1:501)）;mag = abs（H);db = 20log10(（mag+eps）/max(mag);pha = angle（H）;% pha = unwrap(angle(H）；function Hr，w，c，L=hr_type3(h）；计算所设计的3型滤波器的振幅响应%Hr=振幅响应%b=3型滤波器的系数L=Hr的阶次h=3型滤波器的单位冲击响应M=length(h）;L=（M-1）/2；c= 2*h(L+1:-1:1）;n=0:1：L；w=0：1:500*2pi/500；Hr=sin(w*n)*c；本文为互联网收集，请勿用作商业用途本文为互联网收集,请勿用作商业用途基于频域抽样法的FIR数字带通滤波器设计wsl=0.12pi；%低阻带边缘wsh=0。82*pi；%高阻带边缘wpl=0.32pi;低通带边缘wph=0.62*pi;高通带边缘delta=(wpl-wsl）;过度带M=ceil（2*pi3/delta)；抽样点数al=（M1）/2；wl=(2pi/M); %抽样间隔k=0：M1；T1=0。12; T2=0。6;过渡带样本点Hrs=zeros(1，ceil（0。12pi/wl）+1)，T2，T1，ones(1,ceil（0。3*pi/wl）),T1，T2,zeros（1，ceil(0.3734*pi/wl）)，T2，T1,ones（1，ceil(0。3pi/wl)）,T1,T2，zeros(1,ceil（0.12pi/wl)+1)；wdl=0 0。12 0。32 0。62 0。82 1;k1=0:floor(（M1）/2）;k2=floor（M1）/2)+1：M1;angH=al(2pi)/M*k1，al*（2*pi)/M*(M-k2);H=Hrs。exp(j*angH）;h=real（ifft(H);%傅立叶反变换figure(1)；冲击响应图stem（k，h）；title（impulse response）；xlabel（n）；ylabel（h(n）)；grid；figure(2）;幅频曲线图Hf=abs(H）;w=kwl/pi；plot（w，Hf，b-）axis(0 1 -0。1 1。1);title（amplitude response）;xlabel（frequency in pi units);ylabel（Hr(w)；set(gca，xtickmode,manual，xtick,wdl）；set(gca,ytickmode，manual，ytick，0 0.12 0.6 1)；grid;figure（3）;fs=15000；c，f3=freqz(h,1);f3=f3/pi*fs/2；plot（f3,20*log10(abs（c)）；title(频谱特性）；xlabel(频率/HZ）；ylabel(衰减/dB）；grid；t=（0:100）/fs;x=sin(2pit*700)+sin（2pit*3200）+sin（2pit6200）；q=filter(h,1,x）;a,f1=freqz（x);f1=f1/pifs/2；b，f2=freqz(q)；f2=f2/pifs/2；figure（4);subplot（2,1，1)；plot（f1,abs(a)）；title（输入波形频谱图）；xlabel（频率)；ylabel(幅度）subplot（2，1,2）；plot（f2，abs（b）)；title（输出波形频谱图);xlabel(频率）;ylabel(幅度)个人收集整理，勿做商业用途本文为互联网收集，请勿用作商业用途基于汉宁窗的FIR数字高通滤波器设计function s2Fs=15000；t=（0：100）/Fs；x=sin(2*pi500t)+sin（2*pi*3000t)subplot(245)；stem（x)；title（原始信号);axis（0，100，2,2);Ws=7pi/30；Wp=13pi/30；tr_wid=WpWs； 过渡带宽度N=ceil(11*pi/tr_wid） 滤波器长度n=0：1:N1;Wc=（Ws+Wp)/2； 理想高通滤波器的截止频率hd=ideal_hp1(Wc,N); 理想高通滤波器的单位冲激响应w_bla=（blackman(N); 布拉克曼h=hd。*w_bla; %截取得到实际的单位脉冲响应db,mag,pha,grd，w=freqz_m(h,1)； %计算实际滤波器的幅度响应delta_w=2*pi/1000；As=round(max（db（1:1:Ws/delta_w+1) %实际阻带纹波,round是取整函数y=filter(h,1，x)subplot（246）plot（y)title(滤波后的信号）;axis(0,100，-1,1)subplot（241)stem（n,hd)title(理想单位脉冲响应hd（n)）subplot（242)stem（n，w_bla)title(布拉克满窗w（n）)subplot(243）stem（n，h）title(实际单位脉冲响应hd（n）)subplot（244）plot（w/pi,db)title（幅度响应（dB）)axis(0，1,-100,10）function db，mag，pha，grd，w = freqz_m(b,a）;H,w = freqz（b,a，1000,whole)；H = (H(1:1:501)）； w = (w(1：1：501）;mag = abs(H);db = 20log10(（mag+eps)/max（mag）;pha = angle(H);grd = grpdelay(b,a,w)；subplot（247）;plot（pha)title(相频响应)function hd=ideal_hp1(Wc，N）alp=(N1)/2；n=0：1：N1;m=nalp+eps； %eps是一个很小很小的数hd=sin（pim）sin（Wcm）./（pim)；文档为个人收集整理，来源于网络文档为个人收集整理，来源于网络用双线性法设计巴特沃斯高通数字滤波器clear all; clc； close allfs=120; T=1/fs；rp=1; rs=30;Wp=0.35pi; Ws=0。65pi； %数字滤波器指标wp=2*tan(Wp/2)/T; ws=2*tan(Ws/2)/T；预畸变，将数字滤波器的指标变为模拟滤波器的指标N,w=buttord(wp,ws,rp,rs，s）； 求滤波器阶数和3dB截止频率Z，P，K=buttap（N); %设计模拟低通滤波器Md,Nd=zp2tf(Z，P,K); 将零极点形式转换为传输函数形式M,N=lp2hp（Md,Nd,w)； 对低通滤波器进行频率变换h,w=freqs（M,N，512)； 模拟滤波器的幅频响应subplot(2，1，1)；plot（w，abs(h）)； grid;xlabel(Hz)；ylabel(幅度)； title(模拟高通滤波器）；Mh,Nh=bilinear(M，N,1/T); %对模拟滤波器双线性变换h1,w1=freqz（Mh,Nh)； 数字滤波器的幅频响应subplot（2,1，2）;plot（w1/pi,20*log10（abs(h1)）); grid;xlabel(/）;ylabel（幅度(dB）); title(数字高通滤波器); 图5 模拟滤波器与设计的滤波器的单位冲击响应k=0：2000； k2=1：1001;x=10sin（pi/10*k/fs）+5*sin（10*pi*k/fs)+3sin(30*pi*k/fs）;figuresubplot(2,1，1)X=fft（x)*2/2001;y=filter（Mh,Nh，x）;plot（k，y)； ylim(-5 5); title（高通数字滤波器输出);Y=fft(y)2/2001；df=fs/2001； ff=（k2-1）*df；subplot(2,1，2); plot(ff,abs(X（k2）)，r,linewidth,2）； hold onplot（ff，abs（Y（k2）)，b）; title(输入输出频谱比较)； grid;个人收集整理，勿做商业用途个人收集整理，勿做商业用途