FIR數(shù)字濾波器的設(shè)計--等波紋最佳逼近法一、等波最佳逼近的原理簡介等波紋最佳逼近法是一種優(yōu)化設(shè)計法，即最大誤差最小化準則，它克服了窗函數(shù)設(shè)計法和頻率采樣法的缺點，使最大誤差(即波紋的峰值)最小化，并在整個逼近頻段上均勻分布。用等波紋最佳逼近法設(shè)計的FIR數(shù)字濾波器的幅頻響應(yīng)在通帶和阻帶都是等波紋的，而且可以分別控制通帶和阻帶波紋幅度，這就是等波紋的含義。最佳逼近是指在濾波器長度給定的條件下，使加權(quán)誤差波紋幅度最小化。與窗函數(shù)設(shè)計法和頻率采樣法比較，由于這種設(shè)計法使濾波器的最大逼近誤差均勻分布，所以設(shè)計的濾波器性能價格比最高。階數(shù)相同時，這種設(shè)計法使濾波器的最大逼近誤差最小，即通帶最大衰減最小，阻帶最小衰減最大；指標相同時，這種設(shè)計法使濾波器階數(shù)最低。等波紋最佳逼近法的設(shè)計思想。用Hd()表示希望逼近的幅度特性函數(shù)，要求設(shè)計線性相位FIR數(shù)字濾波器時,Hd()必須滿足線性相位約束條件。用H表示實際設(shè)計的濾波器的幅度特性函數(shù)。定義加權(quán)誤差函數(shù)為WHdH式中，W為幅度誤差加權(quán)函數(shù)，用來控制不同頻帶(一般指通帶和阻帶)的幅度逼近精度。等波紋最佳逼近法的設(shè)計在于找到濾波器的系數(shù)向量hn，使得在通帶和阻帶內(nèi)的最大絕對值幅度誤差為最小，這也就是最大誤差最小化問題。二、等波紋逼近法設(shè)計濾波器的步驟和函數(shù)介紹1.根據(jù)濾波器的設(shè)計指標的要求：邊界頻率，通帶最大衰減，阻帶最大衰等估計濾波器階數(shù)n，確定幅度誤差加權(quán)函數(shù)W2.采用Parks-McClellan算法，獲得所設(shè)計濾波器的單位脈沖響應(yīng) hn實現(xiàn)FIR數(shù)字濾波器的等波紋最佳逼近法的MATLAB信號處理工具函數(shù)為firpm和firpmord。firpm函數(shù)采用數(shù)值分析中的多重交換迭代算法求解等波紋最佳逼近問題，求的滿足等波紋最佳逼近準則的FIR數(shù)字濾波器的單位脈沖響應(yīng)hn。firpmord根據(jù)逼近指標,計算采用Parks-McClellan算法等波紋最佳逼近濾波器的最低階數(shù)，誤差加權(quán)向量 w,歸一化邊界頻率向量f。3對firpm和firpmord的說明firpm函數(shù)功能：采用Parks-McClellan算法設(shè)計FIR濾波器函數(shù)格式：hn=firpm(n,f,m,w)n是濾波器的階hn是數(shù)字濾波器的單位脈沖響應(yīng)，其長度為n+1f是希望濾波器的邊界頻率向量，要求f是單調(diào)增向量，并且從0開始，以1結(jié)束，1對于數(shù)字頻率m是與f對應(yīng)的希望濾波器的幅度向量，m和f的長度相等，mk表示希望濾波器在頻率點fk上的幅頻響應(yīng)，m和f給出了希望濾波器的幅度特性。w是誤差加權(quán)向量，其長度為f的一半。w(i)表示對m中第i個頻率段幅度逼近精度的加權(quán)值。w缺省時，函數(shù)默認w全為1，即每個頻率段的逼近誤差加權(quán)值相同。firpmord函數(shù)函數(shù)功能：根據(jù)逼近指標，計算采用Parks-McClellan算法等波紋最佳逼近濾波器的最低階數(shù) n，誤差加權(quán)向量w和歸一化邊界頻率f。其返回參數(shù)作為firpm函數(shù)的調(diào)用參數(shù)。函數(shù)格式：[n,f,m,w]=firpmord(f,m,rip,fs)f可以是歸一化邊界數(shù)字頻率向量，也可以是模擬邊界頻率向量，但必須以0開始，以1結(jié)束或fs/2結(jié)束，并且其中省略了0和fs/2兩個頻率點。fs是時域采樣頻率，單位Hz。fs缺省時,函數(shù)默認fs=2Hz.但這是f的長度(包括省略0和fs/2兩個頻率點)是m的兩倍，即m中的每個元素表示f給定的一個逼近頻段上的希望逼近的幅度值。rip表示f和m描述的各逼近段允許的波紋振幅(幅頻響應(yīng)最大偏差)。4.注意事項：省略fs時，f必須是歸一化的數(shù)字頻率有時計算的階數(shù)n略小，使設(shè)計結(jié)果達不到指標要求，這時要取n+1或n+2三、程序低通濾波器設(shè)計%Lowpassfilter:wp=0.4pi,ws=0.6pi,N=26%peakpassbandrippleis0.01,peakstopbandrippleis0.001clc;clear;clearall;wp_l=0.4*pi;ws_l=0.6*pi;N1=26;deltal_1=0.01;deltal_2=0.001; %初始化參數(shù)f_l=[0wp_lws_lpi]/pi;%設(shè)定歸一化的firpm函數(shù)參數(shù)deltal=[deltal_1,deltal_2];a_l=[1100];%幅值[n,f,m,w_l]=firpmord(f_l(2:3),[10],deltal);%計算權(quán)值fil_l=firpm(N1,f_l,a_l,w_l);%生成濾波器figure(1);stem(fil_l,'.');title('impulseresponseoflowpassfilter');xlabel('timesequence');ylabel('amplitude');fft_l=fft(fil_l,1024);db_l=20*log10(abs(fft_l));%對數(shù)表示figure(2);plot([0:length(fft_l)-1]/length(fft_l),db_l); %將橫坐歸化到0-1,單位為1/1024title('Log-magnituderesponseoflowpassfilter');xlabel('Normalizedfrequency');ylabel('Magnituderesponse(dB)');axis([00.5-10010]); %設(shè)坐標軸的范圍，只取0-2pi的一半afl=abs(fft_l);%取傅變的幅值afl_p=zeros(1,1024/2);afl_s=zeros(1,1024/2);%初始化兩個一行512列的零序列afl_p(1:fix((wp_l/(2*pi))*1024))=afl(1:fix((wp_l/(2*pi))*1024))-1; %將通帶幅值賦給序列一afl_s(fix((ws_l/(2*pi))*1024)+3:end)=afl(fix((ws_l/(2*pi))*1024)+3:end/2); %將阻帶幅值賦給序列二figure(3);plot([0:1024/2-1]/1024,afl_p+afl_s);%將通帶和阻帶的波形顯示到一幅圖像上title('passbandandstopbandapproximationerroroflowpassfilter');xlabel('Normalizedfrequency');ylabel('Passbandandstopbandripple');帶通濾波器設(shè)計%Bandpassfilter:ws1=0.2pi,wp1=0.28pi,wp2=0.72pi,ws2=0.8pi,N=80%peakpassbandrippleis0.01,eachpeakstopbandrippleis0.001ws_b1=0.2*pi;wp_b1=0.28*pi;wp_b2=0.72*pi;ws_b2=0.8*pi;N2=80;deltab_1=0.001;deltab_2=0.01;deltab_3=0.001;deltab=[deltab_1,deltab_2,deltab_3];f_b=[0ws_b1wp_b1wp_b2ws_b2pi]/pi;a_b=[001100];[n,f,m,w_b]=firpmord(f_b(2:5),[010],deltab);fil_b=firpm(N2,f_b,a_b,w_b);figure(4);stem(fil_b,'.');title('impulseresponseofbandpassfilter');xlabel('timesequence');ylabel('amplitude');fft_b=fft(fil_b,1024);db_b=20*log10(abs(fft_b));figure(5);plot([0:length(fft_b)-1]/length(fft_b),db_b);title('Log-magnituderesponseofbandpassfilter');xlabel('Normalizedfrequency');ylabel('Magnituderesponse(dB)');axis([00.5-10010]);afb=abs(fft_b);afb_p=zeros(1,1024/2);afb_s=zeros(1,1024/2);afb_p(fix((wp_b1/(2*pi))*1024)+3:fix((wp_b2/(2*pi))*1024))=afb(fix((wp_b1/(2*pi))*1024)+3:fix((wp_b2/(2*pi))*1024))-1;afb_s(fix((ws_b2/(2*pi))*1024)+2:end)=afb(fix((ws_b2/(2*pi))*1024)+2:end/2);afb_s(1:fix((ws_b1/(2*pi))*1024))=afb(1:fix((ws_b1/(2*pi))*1024));figure(6);plot([0:1024/2-1]/(1024),afb_p+afb_s);title('passbandandstopbandapproximationerrorofbandpassfilter');xlabel('Normalizedfrequency');ylabel('Passbandandstopbandripple');高通濾波器設(shè)計%Highpassfilter:ws=0.4pi,wp=0.55pi,N=34

%peakpassbandrippleis0.01,peakstopbandrippleis0.001ws_h=0.4*pi;wp_h=0.55*pi;N3=34;deltah_1=0.001;deltah_2=0.01;deltah=[deltah_1,deltah_2];f_h=[0ws_hwp_hpi]/pi;a_h=[0011];[n,f,m,w_h]=firpmord(f_h(2:3),[01],deltah)fil_h=firpm(N3,f_h,a_h,w_h);figure(7);stem(fil_h,'.');title('impulseresponseofhighpassfilter');xlabel('timesequence');ylabel('amplitude');fft_h=fft(fil_h,1024);db_h=20*log10(abs(fft_h));figure(8);plot([0:length(fft_h)-1]/length(fft_h),db_h);title('Log-magnituderesponseofhighpassfilter');xlabel('Normalizedfrequency');ylabel('Magnituderesponse(dB)');axis([00.5-10010]);afh=abs(fft_h);afh_p=zeros(1,1024/2);afh_s=zeros(1,1024/2);afh_p(fix((wp_h/(2*pi))*1024)+5:end)=afh(fix((wp_h/(2*pi))*1024)+5:end/2)-1;afh_s(1:fix((ws_h/(2*pi))*1024))=afh(1:fix((ws_h/(2*pi))*1024));figure(9);plot([0:1024/2-1]/(1024),afh_p+afh_s);title('passbandandstopbandapproximationerrorofhighpassfilter');xlabel('Normalizedfrequency');ylabel('Passbandandstopbandripple');仿真圖