第五次DSP實(shí)驗(yàn)

1、第五次DSP實(shí)驗(yàn)FIR濾波一、 實(shí)驗(yàn)?zāi)康?、熟悉 FIR 濾波器設(shè)計(jì)的基本方法。2、掌握用窗函數(shù)設(shè)計(jì) FIR 數(shù)字濾波器的原理及方法，熟悉相應(yīng)的計(jì)算機(jī)高級(jí)語(yǔ)言編程。3、熟悉線性相位 FIR 濾波器的幅頻特性和相位特性。4、了解各種不同窗函數(shù)對(duì)濾波器性能的響應(yīng)。二、實(shí)驗(yàn)原理與方法1.FIR 濾波器的設(shè)計(jì)在前面的實(shí)驗(yàn)中， 我們介紹了 IIR 濾波器的設(shè)計(jì)方法并實(shí)踐了其中的雙線性變換法， IIR具有許多誘人的特性； 但與此同時(shí)， 也具有一些缺點(diǎn)。 例如： 若想利用快速傅立葉變換技術(shù)進(jìn)行快速卷積實(shí)現(xiàn)濾波器，則要求單位脈沖響應(yīng)是有限列長(zhǎng)的。此外， IIR 濾波器的優(yōu)異幅度響應(yīng)， 一般是以相位的非線性為代

2、價(jià)的， 非線性相位會(huì)引起頻率色散。 FIR 濾波器具有嚴(yán)格的相位特性， 這對(duì)于語(yǔ)音信號(hào)處理和數(shù)據(jù)傳輸是和重要的。 目 前 FIR 濾波器的設(shè)計(jì)方法主要有三種： 窗函數(shù)法、 頻率取樣法和切比雪夫等波紋逼近的最優(yōu)化設(shè)計(jì)方法。 常用的是窗函數(shù)法和切比雪夫等波紋逼近的最優(yōu)化設(shè)計(jì)方法。 本實(shí)驗(yàn)項(xiàng)目 中的窗函數(shù)法比較簡(jiǎn)單， 可應(yīng)用現(xiàn)成的窗函數(shù)公式， 在技術(shù)指標(biāo)要求不高的時(shí)候是比較靈活方便的。 它是從時(shí)域出發(fā)， 用一個(gè)窗函數(shù)截取理想的 hd(n) 得到h(n)，以有限長(zhǎng)序列 h(n)近似理想的hd(n) ；如果從頻域出發(fā)， 用理想的h(ejw) 在單位圓上等角度取樣得到 H(k)，根據(jù) h(k)得到 H(z

3、)將逼近理想的 ，這就是頻率取樣法。2.窗函數(shù)設(shè)計(jì)法同其他的數(shù)字濾波器的設(shè)計(jì)方法一樣， 用窗函數(shù)設(shè)計(jì)濾波器也是首先要對(duì)濾波提出性能指標(biāo)。一般是給定一個(gè)理想的頻率響應(yīng) ， 使所設(shè)計(jì)的 FIR 濾波器的頻率響應(yīng)去逼近所要求的理性的濾波器的響應(yīng) 。 窗函數(shù)法設(shè)計(jì)的任務(wù)在于尋找一個(gè)可實(shí)現(xiàn)（有限長(zhǎng)單位脈沖響應(yīng)）的傳遞函數(shù) 三、實(shí)驗(yàn)內(nèi)容1. 傅立葉反變換數(shù)值計(jì)算子程序：代碼如下：function xa = myifft(x,m,n) na = 0:n-1; xa=x * ( exp( 1i * pi * 2 /length(m) ) ) . (m' * na)/length(m);end2.窗函數(shù)

4、產(chǎn)生子程序function r = window(n) r = ones(1,n-1);endfunction r = hamming(n) r = ones(1,n-1); nr = 0:n-2; r = r.*(0.54-0.46*cos(2*pi*nr/(n-1); %stem(nr,r);endfunction r = hanning(n) r = ones(1,n-1); nr = 0:n-2; r = r.*(0.5-0.5*cos(2*pi*nr/(n-1); %stem(nr,r);endfunction r = blackman(n) r = ones(1,n-1); nr

5、= 0:n-2; r = r.*(0.42-0.5*cos(2*pi*nr/(n-1)+0.08*cos(4*pi*nr/(n-1); %stem(nr,r);endfunction r = kaiser(n,b) r = ones(1,n-1); nr = 0:n-2; tr = b*sqrt(1-(2*nr/(n-1)-1).2); r = r.*besseli(0,tr)/besseli(0,b); %stem(nr,r);end3.主函數(shù)function xn = main(hd,n,m) nhd = 0:m-1; hd1 = zeros(1,m); hd1(1:m/2) = hd(1

6、:m/2); hd1(end:-1:m/2+2) = hd1(2:m/2); if mod(n,2) = 0 hd1(1:m/2) = hd1(1:m/2) .* exp(-1i*pi*nhd(1:m/2)/m*(n-1); hd1(m/2+1) = 0; hd1(m/2+2:m) = hd1(m/2+2:m) .* (-exp(-1i*pi*nhd(m/2+2:m)/m*(n-1); xn = real(myifft(hd1,nhd,n); else hd1 = hd1 .* exp(-1i*pi*nhd/m*(n-1); xn = real(myifft(hd1,nhd,n); enden

7、d3上機(jī)部分1. 用窗函數(shù)法設(shè)計(jì)一個(gè)長(zhǎng)度 N 等于 8 的線性相位 FIR 濾波器。其理想的幅頻特性為 分別用矩形窗、Hanning 窗、Hamming 窗、Blackman 窗、Kaiser 窗（ =8.5）設(shè)計(jì)該濾波器，比較設(shè)計(jì)結(jié)果代碼如下：function test() n =8; m = 10*n; nhd = 0:m-1; hd = zeros(1,m); hd(1:floor(m*0.2)+1) = 1; subplot(2, 2, 1); stem(nhd(1:m/2)/m*2,abs(hd(1:m/2); xlabel( ' frequency in pi units&

8、#39; ); title( '濾波器頻域設(shè)計(jì)' ); ylabel( 'Hd' ); xn = main(hd,n,m); nx = 0:n-2; xn = xn(1:n-1).*window(n); subplot(2, 2, 2); stem(nx,xn); xlabel( ' n' ); title( '濾波器時(shí)域設(shè)計(jì)' ); ylabel( 'x' ); k = 0:1:1000; %置初值 X, w = dtft ( xn, nx, k ); %調(diào)用dtft函數(shù)，完成變換 magX = abs( X )

9、; angX = angle( X ); %取模和相位 subplot( 2, 2, 3 ); plot( w(1:501) / pi, magX(1:501) ); grid %分別作出幅度，相位，實(shí)部，虛部對(duì)的歸一化w的圖 %這里只做出w取值為0,pi的圖像，并對(duì)w做歸一化 xlabel( ' frequency in pi units' ); title( 'Magnitude part' ); ylabel( 'Magnitude' ); subplot( 2, 2, 4 ); plot( w(1:501) / pi, angX(1:50

10、1) ); grid xlabel( ' frequency in pi units' ); title( 'Angle part' ); ylabel( 'Radians' ); end實(shí)驗(yàn)結(jié)果：窗函數(shù)類型h(0)h(1)h(2)h(3)h(4)h(5)h(6)矩形窗-0.0898 -0.01250.19770.38430.38430.1977-0.0125Hanning窗0 -0.00240.12090.36520.36520.1209-0.0024Hamming窗-0.0072 -0.00320.12700.36680.36680.1270

11、-0.0032Blackman窗0.0000 -0.00110.09080.35370.35370.0908-0.0011Kaiser(=8.5)-0.0001 -0.00120.09180.35410.35410.0918-0.0012可以看到，實(shí)驗(yàn)結(jié)果和參考答案是基本一樣的，所以認(rèn)為實(shí)驗(yàn)程序調(diào)試成功2.用 Hanning 窗設(shè)計(jì)一個(gè)線性相位帶通濾波器，其長(zhǎng)度 N=15上下邊帶截至頻率分別 ， ，求 h(n)，繪制它的幅頻和相位特性曲線，觀察它的實(shí)際3dB和20dB帶寬。如果N=45，重復(fù)這一設(shè)計(jì)，觀察幅頻和相位特性的變化，注意長(zhǎng)度 N 變化對(duì)結(jié)果的影響。代碼如下：function ex1(

12、) n =15; m = 20*n; nhd = 0:m-1; hd = zeros(1,m); hd(floor(m*0.15)+1:floor(m*0.25)+1) = 1; subplot(2, 2, 1); stem(nhd(1:m/2)/m*2,abs(hd(1:m/2); xlabel( ' frequency in pi units' ); title( '濾波器頻域設(shè)計(jì)' ); ylabel( 'Hd' ); xn = main(hd,n,m); nx = 0:n-2; xn = xn(1:n-1).*hanning(n); su

13、bplot(2, 2, 2); stem(nx,xn); xlabel( ' n' ); title( '濾波器時(shí)域設(shè)計(jì)' ); ylabel( 'x' ); k = 0:1:1000; %置初值 X, w = dtft ( xn, nx, k ); %調(diào)用dtft函數(shù)，完成變換 magX = abs( X ); angX = angle( X ); %取模和相位 subplot( 2, 2, 3 ); plot( w(1:501) / pi, magX(1:501) ); grid %分別作出幅度，相位，實(shí)部，虛部對(duì)的歸一化w的圖 %這里只做出

14、w取值為0,pi的圖像，并對(duì)w做歸一化 xlabel( ' frequency in pi units' ); title( 'Magnitude part' ); ylabel( 'Magnitude' ); subplot( 2, 2, 4 ); plot( w(1:501) / pi, angX(1:501) ); grid xlabel( ' frequency in pi units' ); title( 'Angle part' ); ylabel( 'Radians' ); end運(yùn)行

15、結(jié)果如下： 這里n=15，如圖可以得到一個(gè)帶通的濾波器，可以測(cè)得其通帶為0.288pi0.512pi，這要大于0.3pi0.5pi的設(shè)計(jì)要求，而且可以看到，其過(guò)渡帶較寬，通帶增益衰減快，因此性能比較一般。相位有時(shí)會(huì)有+-pi的跳變。 這里n=45，如圖可以得到一個(gè)帶通的濾波器，可以測(cè)得其通帶為0.318pi0.48pi，這要小于0.3pi0.5pi的設(shè)計(jì)要求，但是其過(guò)渡帶較窄，通帶增益平穩(wěn)，因此性能比較好。這是由于函數(shù)階次高，長(zhǎng)度長(zhǎng)，頻譜泄露小，但是可以看到為了達(dá)到要求的帶寬，設(shè)計(jì)師應(yīng)放松對(duì)截止頻率的要求。相位會(huì)有+-pi的跳變。而去變化頻繁。3.改用矩形窗和 Blackman窗，設(shè)計(jì)步驟（1

16、）中的帶通濾波器，觀察并記錄窗函數(shù)對(duì)濾波器幅頻和相位特性的影響，比較這三種窗函數(shù)的特點(diǎn)實(shí)驗(yàn)結(jié)果如下：矩形窗，n=15矩形窗，n=45Blackman窗，n=15Blackman窗，n=45綜合比較這幾個(gè)函數(shù)，可以發(fā)現(xiàn)，隨著n的增大，將導(dǎo)致系統(tǒng)函數(shù)幅頻響應(yīng)變窄，過(guò)渡帶變陡，阻帶變平，更趨于0，而通帶性狀變好，變平，更趨于1，通帶變窄，已經(jīng)不能滿足設(shè)計(jì)要求。而相頻響應(yīng)的變化將變得更加劇烈。對(duì)比兩個(gè)不同的窗函數(shù)，可以發(fā)現(xiàn)blackman窗產(chǎn)生的系統(tǒng)函數(shù)具有更好的性能。而矩形窗函數(shù)產(chǎn)生的系統(tǒng)函數(shù)的幅頻響應(yīng)在通帶和阻帶內(nèi)均有波紋，不管n=15還是n=45，通帶均過(guò)窄，不滿足設(shè)計(jì)要求。相頻響應(yīng)也不是線性的

17、。3.用 Kaiser 窗設(shè)計(jì)一個(gè)專用的線性相位濾波器。N=40，理想的幅頻特性。當(dāng) 值分別 4,6,8 時(shí)，設(shè)計(jì)相應(yīng)的濾波器，比較它們的幅頻和相位特性，觀察并分析值不同的時(shí)候?qū)Y(jié)果有什么影響。代碼如下：function ex3() n =40; m = 20*n; nhd = 0:m-1; hd = zeros(1,m); hd(floor(m*0.1)+1:floor(m*0.2)+1) = 1; hd(floor(m*0.3)+1:floor(m*0.4)+1) = 1; subplot(2, 2, 1); stem(nhd(1:m/2)/m*2,abs(hd(1:m/2); xlabe

18、l( ' frequency in pi units' ); title( '濾波器頻域設(shè)計(jì)' ); ylabel( 'Hd' ); xn = main(hd,n,m); nx = 0:n-2; xn = xn(1:n-1).*kaiser(n,4); subplot(2, 2, 2); stem(nx,xn); xlabel( ' n' ); title( '濾波器時(shí)域設(shè)計(jì)' ); ylabel( 'x' ); k = 0:1:1000; %置初值 X, w = dtft ( xn, nx, k

19、); %調(diào)用dtft函數(shù)，完成變換 magX = abs( X ); angX = angle( X ); %取模和相位 subplot( 2, 2, 3 ); plot( w(1:501) / pi, magX(1:501) ); grid %分別作出幅度，相位，實(shí)部，虛部對(duì)的歸一化w的圖 %這里只做出w取值為0,pi的圖像，并對(duì)w做歸一化 xlabel( ' frequency in pi units' ); title( 'Magnitude part' ); ylabel( 'Magnitude' ); subplot( 2, 2, 4

20、); plot( w(1:501) / pi, angX(1:501) ); grid xlabel( ' frequency in pi units' ); title( 'Angle part' ); ylabel( 'Radians' ); end結(jié)果分析如下：Kaiser窗，b=4Kaiser窗，b=6Kaiser窗，b=8Kaiser窗，b=4，6，8對(duì)比可以發(fā)現(xiàn)，隨著b的增大，過(guò)渡帶的范圍也將增大，導(dǎo)致阻帶和通帶的性能變差。綜合來(lái)看，b=4時(shí)有最好的性能。，但其實(shí)它們都沒(méi)辦法保證應(yīng)有的帶寬。思考題1. 對(duì)不同的窗函數(shù)和不同的n，3db截止頻率都有不同，根據(jù)上面的實(shí)驗(yàn)，可以發(fā)現(xiàn)，隨著n的增大，將導(dǎo)致3db截止頻率向通帶中心頻率靠近，也就是說(shuō)通帶變窄。它顯然不等于理想頻率響應(yīng)的截止頻率。2.由1中的討論可知，n的大小將直接影響設(shè)計(jì)通帶邊緣截至頻率和阻帶臨界頻率，所以可以通過(guò)選定合適的n已經(jīng)窗函數(shù)以達(dá)到設(shè)計(jì)要求。但這些關(guān)系實(shí)