實(shí)驗(yàn)四窗函數(shù)法設(shè)計(jì)FIR數(shù)字濾波器.doc

.實(shí)驗(yàn)四 窗函數(shù)法設(shè)計(jì)FIR數(shù)字濾波器一、實(shí)驗(yàn)?zāi)康?、掌握窗函數(shù)法設(shè)計(jì)FIR數(shù)字濾波器的原理及具體方法。 2、掌握頻率取樣法設(shè)計(jì)FIR數(shù)字濾波器的原理和基本方法。 3、學(xué)習(xí)利用窗函數(shù)法和頻率取樣法設(shè)計(jì)低通、帶通、高通、帶阻數(shù)字濾波器。二、實(shí)驗(yàn)環(huán)境計(jì)算機(jī)、MATLAB軟件三、實(shí)驗(yàn)基礎(chǔ)理論窗函數(shù)設(shè)計(jì)FIR濾波器1.基本原理窗函數(shù)設(shè)計(jì)法的基本思想為，首先選擇一個(gè)適當(dāng)?shù)睦硐氲臑V波器，然后用窗函數(shù)截取它的單位脈沖響應(yīng)，得到線性相位和因果的FIR濾波器。這種方法的重點(diǎn)是選擇一個(gè)合適的窗函數(shù)和理想濾波器，使設(shè)計(jì)的濾波器的單位脈沖響應(yīng)逼近理想濾波器的單位脈沖響應(yīng)。2.設(shè)計(jì)步驟（1）給定理想濾波器的頻率響應(yīng)，在通帶上具有單位增益和線性相位，在阻帶上具有零響應(yīng)。一個(gè)帶寬為的低通濾波器由下式給定： 其中為采樣延遲，其作用是為了得到一個(gè)因果系統(tǒng)。（2）確定這個(gè)濾波器的單位脈沖響應(yīng)為了得到一個(gè)長(zhǎng)度為N的因果的線性相位FIR濾波器，我們令 （3）用窗函數(shù)截取得到所設(shè)計(jì)FIR數(shù)字濾波器： 3.窗函數(shù)的選擇常用的窗函數(shù)有矩形（Rectangular）窗，漢寧（Hanning）窗，海明（Hamming）窗、布萊克曼（Blackman）窗、凱瑟（Kaiser）窗等表4-1 MATLAB中產(chǎn)生窗函數(shù)的命令MATLAB函數(shù)窗函數(shù)MATLAB函數(shù)窗函數(shù)Boxcar矩形窗函數(shù)Blackman布萊克曼窗Hanning漢寧窗函數(shù)Kaiser凱瑟窗函數(shù)Hamming海明窗表4-2 常用窗函數(shù)的特性窗函數(shù)窗函數(shù)頻率特性加窗后濾波器指標(biāo)旁瓣峰值dB主瓣寬度過(guò)渡帶寬最小阻帶衰減dB矩形窗-134/N1.8/N-21漢寧窗-318/N6.2/N-44海明窗-418/N6.6/N-53布萊克曼窗-5712/N11/N-74凱瑟窗是一種廣泛在實(shí)際中廣泛應(yīng)用的窗函數(shù)，它由下式給定：其中是修正的零階貝塞爾函數(shù)，參數(shù)控制最小阻帶衰減，這種窗函數(shù)對(duì)于相同的N可以提供不同的過(guò)渡帶寬。由于貝塞爾函數(shù)比較復(fù)雜，這種窗函數(shù)的設(shè)計(jì)方程很難推導(dǎo)，然而幸運(yùn)的是，有一些經(jīng)驗(yàn)設(shè)計(jì)方程可以直接使用。已知給定的指標(biāo) ，濾波器長(zhǎng)度N和凱瑟窗參數(shù)可以按如下凱瑟窗方程給出過(guò)渡帶帶寬：頻率取樣設(shè)計(jì)FIR濾波器1.基本原理頻率取樣法從頻域出發(fā)，把理想的濾波器等間隔采樣得到,將作為實(shí)際設(shè)計(jì)濾波器的：得到以后可以由來(lái)確定唯一確定濾波器的單位脈沖響應(yīng),可以由求得：其中為內(nèi)插函數(shù)：有求得的頻率響應(yīng)將逼近。如果我們?cè)O(shè)計(jì)的是線性相位FIR濾波器，則的幅度和相位滿足線性相位濾波器的約束條件。我們將表示為如下形式當(dāng)為實(shí)數(shù)，則由此得到即為中心偶對(duì)稱。在利用線性相位條件可知，對(duì)于1型和2型線性相位濾波器：對(duì)于3型和4型線性相位濾波器2.設(shè)計(jì)步驟（1）由給定的理想濾波器給出和。（2）由求得（3）根據(jù)求得或四、實(shí)驗(yàn)內(nèi)容1、設(shè)計(jì)一個(gè)數(shù)字低通FIR濾波器，其技術(shù)指標(biāo)如下：分別采用矩形窗、漢寧窗、海明窗、布萊克曼窗、凱瑟窗設(shè)計(jì)該濾波器。結(jié)合實(shí)驗(yàn)結(jié)果，分別討論采用上述方法設(shè)計(jì)的數(shù)字濾波器是否都能滿足給定指標(biāo)要求。（1） 矩形窗程序代碼：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp;N=ceil(1.8*pi/tr_width)n=0:N-1;wc=(wst+wp)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc(wc/pi)*(n-alpha);w_boxcar=boxcar(N);h=hd.*w_boxcar;subplot(221);stem(n,hd,filled);axis tight;xlabel(n);ylabel(hd(n);Hr,w1=zerophase(h);subplot(222);plot(w1/pi,Hr);axis;xlabel(omega/pi);ylabel(H(omega);subplot(223);stem(n,h,filled);axis tight;xlabel(n);ylabel(h(n);H,w=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(H);axis tight;xlabel(omega/pi);ylabel(dB);grid on;MATLAB圖形：（2）漢寧窗程序代碼：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp;N=ceil(6.2*pi/tr_width)n=0:N-1;wc=(wst+wp)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc(wc/pi)*(n-alpha);w_boxcar=hanning(N);h=hd.*w_boxcar;subplot(221);stem(n,hd,filled);axis tight;xlabel(n);ylabel(hd(n);Hr,w1=zerophase(h);subplot(222);plot(w1/pi,Hr);axis;xlabel(omega/pi);ylabel(H(omega);subplot(223);stem(n,h,filled);axis tight;xlabel(n);ylabel(h(n);H,w=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(H);axis tight;xlabel(omega/pi);ylabel(dB);grid on;MATLAB圖形：（3）海明窗程序代碼：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp;N=ceil(6.6*pi/tr_width)n=0:N-1;wc=(wst+wp)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc(wc/pi)*(n-alpha);w_boxcar=hamming(N);h=hd.*w_boxcar;subplot(221);stem(n,hd,filled);axis tight;xlabel(n);ylabel(hd(n);Hr,w1=zerophase(h);subplot(222);plot(w1/pi,Hr);axis;xlabel(omega/pi);ylabel(H(omega);subplot(223);stem(n,h,filled);axis tight;xlabel(n);ylabel(h(n);H,w=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(H);axis tight;xlabel(omega/pi);ylabel(dB);grid on;MATLAB圖形：（4）布萊克曼窗程序代碼：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp;N=ceil(11*pi/tr_width)n=0:N-1;wc=(wst+wp)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc(wc/pi)*(n-alpha);w_boxcar=blackman(N);h=hd.*w_boxcar;subplot(221);stem(n,hd,filled);axis tight;xlabel(n);ylabel(hd(n);Hr,w1=zerophase(h);subplot(222);plot(w1/pi,Hr);axis;xlabel(omega/pi);ylabel(H(omega);subplot(223);stem(n,h,filled);axis tight;xlabel(n);ylabel(h(n);H,w=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(H);axis tight;xlabel(omega/pi);ylabel(dB);grid on;MATLAB圖形為：（5）凱瑟窗程序代碼：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp;As=50;N=ceil(As-7.95)/(2.285*tr_width)+1;beta=0.1102*(As-8.7);n=0:N-1;wc=(wst+wp)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc(wc/pi)*(n-alpha);w_boxcar=kaiser(N,beta);h=hd.*w_boxcar;subplot(221);stem(n,hd,filled);axis tight;xlabel(n);ylabel(hd(n);Hr,w1=zerophase(h);subplot(222);plot(w1/pi,Hr);axis;xlabel(omega/pi);ylabel(H(omega);subplot(223);stem(n,h,filled);axis tight;xlabel(n);ylabel(h(n);H,w=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(H);axis tight;xlabel(omega/pi);ylabel(dB);grid on;MATLAB圖形：2、設(shè)計(jì)一個(gè)數(shù)字帶通FIR濾波器，其技術(shù)指標(biāo)如下：下阻帶邊緣：下通帶邊緣：上通帶邊緣：上阻帶邊緣：程序代碼：wp1=0.2*pi;Rp1=1;wst1=0.35*pi;A1=60;width1=wst1-wp1;N1=ceil(11*pi/width1)+1;n1=0:(N1-1);wc1=(wp1+wst1)/2;alpha=(N1-1)/2;wp2=0.65*pi;Rp2=1;wst2=0.8*pi;A2=60;width2=wst2-wp2;N2=ceil(11*pi/width2)+1;n2=0:(N2-1);wc2=(wp2+wst2)/2;alpha=(N2-1)/2;hd=(wc2/pi)*sinc(wc2/pi)*(n2-alpha)-(wc1/pi)*sinc(wc1/pi)*(n1-alpha);w_w=blackman(N1); h=hd.*w_w; subplot(221);stem(n1,h,filled);subplot(222);H,w=freqz(h,1);plot(w/pi,20*log10(abs(H)/max(abs(H);subplot(223);Hr,w1=zerophase(h);plot(w1/pi,Hr);subplot(224);stem(n1,hd,filled);Hr,wl=zerophase(h); grid on;MATLAB圖形： 3.采用頻率取樣法設(shè)計(jì)FIR數(shù)字低通濾波器，滿足以下指標(biāo)（1）取N=20，過(guò)渡帶沒(méi)有樣本。（2）取N=40，過(guò)渡帶有一個(gè)樣本，T=0.39。（3）取N=60，過(guò)渡帶有兩個(gè)樣本，T1=0.5925，T2=0.1099。（4）分別討論采用上述方法設(shè)計(jì)的數(shù)字濾波器是否都能滿足給定的指標(biāo)要求。（1）程序代碼：N=20;alpha=(N-1)/2;L=0:N-1;wL=(2*pi/N)*L; Hrs=1,1,1,zeros(1,15),1,1; Hdr=1,1,0,0; wdL=0,0.25,0.25,1;k1=0:floor(N-1)/2);k2=floor(N-1)/2)+1:N-1;angH=-alpha*(2*pi)/N*k1,alpha*(2*pi)/N*(N-k2); H=Hrs.*exp(j*angH); h=ifft(H,N); w=0:500*pi/500;Hr,wr=zerophase(h); subplot(221); plot(wdL,Hdr,wL(1:N/2+1)/pi,Hrs(N/2+1);axis(0,1,-0.1,1.1);xlabel(omega (pi);ylabel(Hr(k);subplot(222); stem(L,h,filled);axis(0,N-1,-0.1,0.3);xlabel(n);ylabel(h(n); subplot(223);plot(wr/pi,Hr,wL(1:N/2+1)/pi,Hrs(1:N/2+1);axis(0,1,-0.2,1.2);grid on; xlabel(omega (pi);ylabel(Hr(omega);subplot(224); plot(wr/pi,20*log10(abs(Hr)/max(abs(Hr);axis(0,1,-50,5);grid on;xlabel(omega (pi);ylabel(dB)MATLAB圖形如下：（2）程序代碼：N=40;alpha=(N-1)/2;L=0:N-1;wL=(2*pi/N)*L; Hrs=1,1,1,1,1,0.39,zeros(1,29),0.39,1,1,1,1;Hdr=1,1,0.39,0,0;wdL=0,0.2,0.25,0.3,1;k1=0:floor(N-1)/2);k2=floor(N-1)/2)+1:N-1;angH=-alpha*(2*pi)/N*k1,alpha*(2*pi)/N*(N-k2); H=Hrs.*exp(j*angH); h=ifft(H,N); w=0:500*pi/500;Hr,wr=zerophase(h); subplot(221); plot(wdL,Hdr,wL(1:N/2+1)/pi,Hrs(N/2+1);axis(0,1,-0.1,1.1);xlabel(omega (pi);ylabel(Hr(k);subplot(222); stem(L,h,filled);axis(0,N-1,-0.1,0.3);xlabel(n);ylabel(h(n);subplot(223); plot(wr/pi,Hr,wL(1:N/2+1)/pi,Hrs(1:N/2+1);axis(0,1,-0.2,1.2);grid on;xlabel(omega (pi);ylabel(Hr(omega);subplot(224); plot(wr/pi,20*log10(abs(Hr)/max(abs(Hr);axis(0,1,-50,5);grid on;xlabel(omega (pi);ylabel(dB)grid on;MATLAB圖形如下：3、程序代碼：N=60;alpha=(N-1)/2;L=0:N-1;wL=(2*pi/N)*L; Hrs=1,1,1,1,1,1,1,0.5925,0.1099,zeros(1,43),0.1099,0.5925,1,1,1,1,1,1;Hdr=1,1,0.5925,0.1099,0,0;wdL=0,0.2,7/30,8/30,0.3,1;k1=0:floor(N-1)/2);k2=floor(N-1)/2)+1:N-1;angH=-alpha*(2*pi)/N*k1,alpha*(2*pi)/N*(N-k2); H=Hrs.*exp(j*angH); h=ifft(H,N);w=0:500*pi/500;Hr,wr=zerophase(h); subplot(221); plot(wdL,Hdr,wL(1:N/2+1)/pi,Hrs(N/2+1);axis(0,1,-0.1,1.1);xlabel(omega (pi);ylabel(Hr(k);subplot(222); stem(L,h,filled);axis(0,N-1,-0.1,0.3);xlabel(n);ylabel(h(n);subplot(223); plot(wr/pi,Hr,wL(1:N/2+1)/pi,Hrs(1:N/2+1);axis(0,1,-0.2,1.2);grid on; xlabel(omega (pi);ylabel(Hr(omega);subplot(224); plot(wr/pi,20*log10(abs(Hr)/max(abs(Hr);axis(0,1,-50,5);grid on;xlabel(omega (pi);ylabel(dB)grid on;MATLAB圖形:各階數(shù)的通阻帶指標(biāo)如下表：階數(shù)通帶波動(dòng)0.25dB阻帶衰減50dB是否滿足指標(biāo)N=20不滿足不滿足否N=40不滿足滿足否N=60滿足滿足是4.采用頻率取樣技術(shù)設(shè)計(jì)下面的高通濾波器對(duì)于高通濾波器，N必須為奇數(shù)（或1型濾波器）。選擇N=33，過(guò)渡帶有兩個(gè)樣本，過(guò)渡帶的最優(yōu)值為T(mén)1=0.1095，T2=0.598.程序代碼：N=33;alpha=(N-1)/2;L=0:N-1;wL=(2*pi/N)*L;Hrs=0,0,0,0,0,0,0,0,0,0,0,0.1095,0.598,1,1,1,1,1,1,1,0.598,0.1095,0,0,0,0,0,0,0,0,0,0,