數(shù)字信號處理實驗報告.doc

《數(shù)字信號處理實驗報告.doc》由會員分享，可在線閱讀，更多相關(guān)《數(shù)字信號處理實驗報告.doc（25頁珍藏版）》請在裝配圖網(wǎng)上搜索。

《數(shù)字信號處理》 實驗報告 課程名稱：《數(shù)字信號處理》 學(xué) 院：信息科學(xué)與工程學(xué)院 專業(yè)班級：通信1502班 學(xué)生姓名：侯子強 學(xué) 號：0905140322 指導(dǎo)教師：李宏 2017年5月28日 實驗一 離散時間信號和系統(tǒng)響應(yīng) 一. 實驗?zāi)康? 1. 熟悉連續(xù)信號經(jīng)理想采樣前后的頻譜變化關(guān)系，加深對時域采樣定理的理解 2. 掌握時域離散系統(tǒng)的時域特性 3. 利用卷積方法觀察分析系統(tǒng)的時域特性 4. 掌握序列傅里葉變換的計算機實現(xiàn)方法，利用序列的傅里葉變換對離散信號及系統(tǒng)響應(yīng)進(jìn)行頻域分析 二、實驗原理 1. 采樣是連續(xù)信號數(shù)字化處理的第一個關(guān)鍵環(huán)節(jié)。對采樣過程的研究不僅可以了解采樣前后信號時域和頻域特性的變化以及信號信息不丟失的條件，而且可以加深對離散傅里葉變換、Z變換和序列傅里葉變換之間關(guān)系式的理解。 對連續(xù)信號以T為采樣間隔進(jìn)行時域等間隔理想采樣，形成采樣信號： 式中為周期沖激脈沖，為的理想采樣。 的傅里葉變換為： 上式表明將連續(xù)信號采樣后其頻譜將變?yōu)橹芷诘?，周期為Ωs=2π/T。也即采樣信號的頻譜是原連續(xù)信號xa(t)的頻譜Xa(jΩ)在頻率軸上以Ωs為周期，周期延拓而成的。因此，若對連續(xù)信號進(jìn)行采樣，要保證采樣頻率fs≥2fm，fm為信號的最高頻率，才可能由采樣信號無失真地恢復(fù)出原模擬信號 計算機實現(xiàn)時，利用計算機計算上式并不方便，因此我們利用采樣序列的傅里葉變換來實現(xiàn)，即 而為采樣序列的傅里葉變換 2. 時域中，描述系統(tǒng)特性的方法是差分方程和單位脈沖響應(yīng)，頻域中可用系統(tǒng)函數(shù)描述系統(tǒng)特性。已知輸入信號，可以由差分方程、單位脈沖響應(yīng)或系統(tǒng)函數(shù)求出系統(tǒng)對于該輸入信號的響應(yīng)。本實驗僅在時域求解，對于差分方程可用Matlab中的工具箱函數(shù)filter()函數(shù)求解 一個時域離散線性時不變系統(tǒng)的輸出與輸入間的關(guān)系為： 可用Matlab中的工具箱函數(shù)conv()函數(shù)求解 三、實驗內(nèi)容及步驟 1. 時域采樣定理的驗證 給定模擬信號： 式中 。 其幅頻特性如圖所示： 選擇三種采樣頻率Fs=1kHz, 300Hz, 200Hz, 生成采樣序列 分別用序列表示。編寫程序計算三個序列的幅頻特性曲線，并繪圖顯示。觀察在折疊頻率附近與連續(xù)信號頻譜有無明顯差別，分析頻譜混疊現(xiàn)象。 實驗程序如下 %時域采樣定理的驗證 %Fs=1KHz Tp=64/1000; %Tp=64ms Fs=1000;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5; xnt=A*exp(-alph*n*T).*sin(omega*n*T); Xk=T*fft(xnt,M); %MFFT yn=xa(nT);subplot(3,2,1); stem(xnt); % box on;title((a) Fs=1000Hz); k=0:M-1;fk=k/Tp; subplot(3,2,2);plot(fk,abs(Xk));title((a) T*FT[xa(nT)],Fs=1000Hz); xlabel(f(Hz));ylabel(幅度); axis([0,Fs,0,1.2*max(abs(Xk))]) %Fs=300Hz Tp=64/1000; Fs=300;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5; xnt=A*exp(-alph*n*T).*sin(omega*n*T); Xk=T*fft(xnt,M); yn=xa(nT);subplot(3,2,1); stem(xnt); box on;title((a) Fs=300Hz); k=0:M-1;fk=k/Tp; subplot(3,2,2);plot(fk,abs(Xk),r);title((a)T*FT[xa(nT)],Fs=300Hz’); xlabel(f(Hz));ylabel(‘幅度’); axis([0,Fs,0,1.2*max(abs(Xk))]) %Fs=200Hz Tp=64/1000; %64ms Fs=300;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5; xnt=A*exp(-alph*n*T).*sin(omega*n*T); Xk=T*fft(xnt,M); yn=xa(nT);subplot(3,2,1); stem(xnt,.); box on;title((a) Fs=200Hz); k=0:M-1;fk=k/Tp; subplot(3,2,2);plot(fk,abs(Xk));title((a) T*FT[xa(nT)],Fs=200Hz); xlabel(f(Hz));ylabel(幅度); axis([0,Fs,0,1.2*max(abs(Xk))]); 2. 給定一個低通濾波器的差分方程為： 輸入序列 （1）分別求出和的系統(tǒng)響應(yīng)，并畫出其波形 （2） 求出系統(tǒng)的單位脈沖響應(yīng)，畫出其波形 A=[1,-0.9]; B=[0.05,0.05]; x1n=[ones(1,8),zeros(1,50)] x2n=ones(1,200); hn=impz(B,A,50); subplot(3,1,1);stem(hn); title(‘（1）系統(tǒng)單位脈沖響應(yīng)h(n)’); y1n=filter(B,A,x1n); subplot(3,1,2);stem(y1n); title((2)系統(tǒng)對R（8）的響應(yīng) y1(n)); y2n=filter(B,A,x2n); subplot(3,1,3);stem(y2n); title(‘(3)系統(tǒng)對u(n)的響應(yīng))y2(n)’); 3. 給定系統(tǒng)的單位脈沖響應(yīng)為 用線性卷積法求分別對系統(tǒng)和的輸出響應(yīng)，并畫出波形 x1n=ones(1,8); h1n=[ones(1,10) zeros(1,20)]; h2n=[1,2.5,2.5,1,zeros(1,10)]; y11n=conv(h1n,x1n); y22n=conv(h2n,x1n); subplot(2,2,1);stem(h1n,’.b’); title(‘(4)系統(tǒng)單位脈沖響應(yīng)h1(n)’); subplot(2,2,2);stem(y11n,.b); title(‘(5)h1(n)與R8(n)的卷積y11(n)’); subplot(2,2,3);stem(h2n,.b); title(‘(6)系統(tǒng)單位脈沖響應(yīng)h2(n)’); subplot(2,2,4);stem(y22n,.b); title(‘(7)h2(n)與R8（n）的卷積y22(n)’); 四、實驗思考 1. 在分析理想采樣序列特性的實驗中，采樣頻率不同時，相應(yīng)理想采樣序列的傅里葉變換頻譜的數(shù)字頻率度量是否都相同？它們所對應(yīng)的模擬頻率是否相同？為什么？ 答：當(dāng)采樣頻率不同時，數(shù)字度量不同，但是模擬頻率相同。 因為數(shù)字頻率W是模擬角頻率Ω用采樣頻率FS歸一化頻率。數(shù)字頻率和模擬角頻率之間的關(guān)系是W=ΩT,模擬信號的模擬角頻率Ω不變，當(dāng)采樣頻率不同時，T不同，所以數(shù)字頻率Ω不同。因此，采樣頻率不同時，相應(yīng)理想采樣序列的傅里葉變換頻譜的數(shù)字頻率度量不相同，但是它們所對應(yīng)的模擬頻率相同。 2. 如果輸入信號為無線長序列，系統(tǒng)的單位脈沖響應(yīng)是有限長序列，可否用線性卷積法求系統(tǒng)的響應(yīng)？如何求？ 答：（1）對輸入信號序列分段； （2）求單位脈沖響應(yīng)與各段的卷積； （3）將各段卷積結(jié)果相加。 3. 如果信號經(jīng)過低通濾波器，把信號的高頻分量濾掉，時域信號會有何變化？用前面第二個實驗結(jié)果進(jìn)行分析說明 答：把信號經(jīng)過低通濾波器，把信號的高頻成分濾掉，時域信號的劇烈將變得平滑。 五、 實驗心得及體會 通過本次實驗我重新溫習(xí)了MATLAB這個軟件的使用方法，運行環(huán)境。通過這款軟件使我們的學(xué)習(xí)更加便利。 實驗二 用FFT對信號作頻譜分析 一、實驗?zāi)康? 1. 進(jìn)一步加深DFT算法原理和基本性質(zhì)的理解 2. 掌握用FFT對連續(xù)信號和時域離散信號進(jìn)行頻譜分析的方法 3. 了解用FFT進(jìn)行頻譜分析時可能出現(xiàn)的分析誤差及其原因，以便在實際中正確應(yīng)用FFT 二、實驗原理 用FFT對信號作頻譜分析是學(xué)習(xí)數(shù)字信號處理的重要內(nèi)容，經(jīng)常需要進(jìn)行譜分析的信號是模擬信號和時域離散信號。對信號進(jìn)行譜分析的重要問題是頻譜分辨率F和分析誤差。頻譜分辨率直接和FFT的變換區(qū)間N有關(guān)，F(xiàn)FT能夠?qū)崿F(xiàn)的頻率分辨率是2p/N，因此要求2p/NF?？梢愿鶕?jù)此式選擇FFT的變換區(qū)間N。誤差主要來自于用FFT作頻譜分析時，得到的是離散譜，而信號（周期信號除外）是連續(xù)譜，只有當(dāng)N較大時，離散譜的包絡(luò)才能逼近于連續(xù)譜，因此N要適當(dāng)選擇大一些。 周期信號的頻譜是離散譜，只有用整數(shù)倍周期的長度作FFT，得到的離散譜才能代表周期信號的頻譜。如果不知道信號周期，可以盡量選擇信號的觀察時間長一些。 對模擬信號進(jìn)行譜分析時，首先要按照采樣定理將其變?yōu)闀r域離散信號。如果是模擬周期信號，也應(yīng)該選取整數(shù)倍周期的長度，經(jīng)過采樣后形成周期序列，按照周期序列的譜分析進(jìn)行。 三、實驗步驟及內(nèi)容 1. 對以下給出的各序列進(jìn)行譜分析： 選擇FFT的變換區(qū)間N為8和16兩種情況進(jìn)行頻譜分析。分別打印其幅頻特性曲線，并進(jìn)行對比、分析、討論。 x1n=[ones(1,4)]; %產(chǎn)生R4(n)序列向量 X1k8=fft(x1n,8); %計算x1n的8點DFT X1k16=fft(x1n,16); %計算x1n的16點DFT N=8; f=2/N*(0:N-1); figure(1); subplot(1,2,1);stem(f,abs(X1k8),.); %繪制8點DFT的幅頻特性圖 title((1a) 8點DFT[x_1(n)]);xlabel(ω/π);ylabel(幅度); N=16; f=2/N*(0:N-1); subplot(1,2,2);stem(f,abs(X1k16),.); %繪制8點DFT的幅頻特性圖 title((1a) 16點DFT[x_1(n)]);xlabel(ω/π);ylabel(幅度); %x2n 和 x3n M=8;xa=1:(M/2); xb=(M/2):-1:1; x2n=[xa,xb]; %產(chǎn)生長度為8的三角波序列x2(n) x3n=[xb,xa]; X2k8=fft(x2n,8); X2k16=fft(x2n,16); X3k8=fft(x3n,8); X3k16=fft(x3n,16); figure(2); N=8; f=2/N*(0:N-1); subplot(2,2,1);stem(f,abs(X2k8),.); %繪制8點DFT的幅頻特性圖 title((2a) 8點DFT[x_2(n)]);xlabel(ω/π);ylabel(幅度); subplot(2,2,3);stem(f,abs(X3k8),.); %繪制8點DFT的幅頻特性圖 title((3a) 8點DFT[x_3(n)]);xlabel(ω/π);ylabel(幅度); N=16; f=2/N*(0:N-1); subplot(2,2,2);stem(f,abs(X2k16),.); %繪制8點DFT的幅頻特性圖 title((2a) 16點DFT[x_2(n)]);xlabel(ω/π);ylabel(幅度); subplot(2,2,4);stem(f,abs(X3k16),.); %繪制8點DFT的幅頻特性圖 title((3a) 16點DFT[x_3(n)]);xlabel(ω/π);ylabel(幅度); 2. 對以下各周期序列進(jìn)行頻譜分析 選FFT的變換區(qū)間N為8和16兩種情況分別對以上序列進(jìn)行頻譜分析。分別打印其幅頻特性曲線，并進(jìn)行對比、分析、討論。 n=0:8; xn4=cos((pi.*n)/4); subplot(2,3,1);stem(n,xn4,.); X8k4 = fft(xn4,8); n21 = 0:length(X8k4)-1; subplot(2,3,2);stem(n21,X8k4,.); X16k4 = fft(xn4,16); n22 = 0:length(X16k4)-1; subplot(2,3,3);stem(n22,X16k4,.);%end n=0:16; xn5=cos((pi.*n)/4)+cos((pi.*n)/8); subplot(2,3,4);stem(n,xn5,.); X8k5 = fft(xn5,8); n21 = 0:length(X8k5)-1; subplot(2,3,5);stem(n21,X8k5,.); X16k5 = fft(xn5,16); n22 = 0:length(X16k5)-1; subplot(2,3,6);stem(n22,X16k5,.); 3. 對模擬周期信號進(jìn)行頻譜分析 選擇樣頻率Fs=64Hz，對變換區(qū)間N=16，32，64三種情況進(jìn)行譜分析。分別打印其幅頻特性曲線，并進(jìn)行對比、分析、討論。 程序如下： Fs=64;T=1/Fs; N=16;n=0:N-1; %FFT的變換區(qū)間N=16 x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %對x6(t)16點采樣 X6k16=fft(x6nT); %計算x6nT的16點DFT X6k16=fftshift(X6k16); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %頻率分辨率F k=-N/2:N/2-1;fk=k*F; %產(chǎn)生16點DFT對應(yīng)的采樣點頻率（以零頻率為中心） subplot(3,1,1);stem(fk,abs(X6k16),.); box on %繪制8點DFT的幅頻特性圖 title((6a) 16點|DFT[x_6(nT)]|);xlabel(f(Hz));ylabel(幅度); axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k16))]) N=32;n=0:N-1; %FFT的變換區(qū)間N=16 x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %對x6(t)32點采樣 X6k32=fft(x6nT); %計算x6nT的32點DFT X6k32=fftshift(X6k32); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %頻率分辨率F k=-N/2:N/2-1;fk=k*F; %產(chǎn)生16點DFT對應(yīng)的采樣點頻率（以零頻率為中心） subplot(3,1,2);stem(fk,abs(X6k32),.); box on %繪制8點DFT的幅頻特性圖 title((6b) 32點|DFT[x_6(nT)]|);xlabel(f(Hz));ylabel(幅度); axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k32))]) N=64;n=0:N-1; %FFT的變換區(qū)間N=16 x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %對x6(t)64點采樣 X6k64=fft(x6nT); %計算x6nT的64點DFT X6k64=fftshift(X6k64); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %頻率分辨率F k=-N/2:N/2-1;fk=k*F; %產(chǎn)生16點DFT對應(yīng)的采樣點頻率（以零頻率為中心） subplot(3,1,3);stem(fk,abs(X6k64),.); box on%繪制8點DFT的幅頻特性圖 title((6a) 64點|DFT[x_6(nT)]|);xlabel(f(Hz));ylabel(幅度); axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k64))]) 四、實驗思考 1. 在N=8時， 和的幅頻特性會相同嗎？為什么？N=16時呢？ 答：在N=8時， 和 的幅頻特性相同，而N=16時不相同。 因為 = ,所以 和 的8點DFT的模相等。但當(dāng)N=16時， 和 不滿足循環(huán)移位關(guān)系，所以兩者幅頻特性不相同。 2. 對于周期序列，如果周期不知道，如何用FFT進(jìn)行譜分析？ 答：周期信號的周期預(yù)先不知道時,可先截取M點進(jìn)行DFT,再將截取長度擴(kuò)大1倍截取,比較結(jié)果,如果二者的差別滿足分析誤差要求,則可以近似表示該信號的頻譜,如果不滿足誤差要求就繼續(xù)將截取長度加倍,重復(fù)比較,直到結(jié)果滿足要求 。 五、實驗總結(jié)及心得體會 通過實驗進(jìn)一步加深DFT算法原理和基本性質(zhì)的理解，掌握用FFT對連續(xù)信號和時域離散信號進(jìn)行頻譜分析的方法，了解用FFT進(jìn)行頻譜分析時可能出現(xiàn)的分析誤差及其原因，以便在實際中正確應(yīng)用FFT。 實驗三 用雙線性變換法設(shè)計IIR數(shù)字濾波器 一、實驗?zāi)康? 1. 熟悉用雙線性變換法設(shè)計IIR數(shù)字濾波器的原理和方法 2. 掌握IIR數(shù)字濾波器的Matlab實現(xiàn)方法 3. 通過觀察對實際心電圖信號的濾波作用，獲得數(shù)字濾波的感性認(rèn)識 二、實驗原理 設(shè)計IIR數(shù)字濾波器一般采用間接設(shè)計法——脈沖響應(yīng)不變法和雙線性變換法，應(yīng)用最廣泛的是雙線性變換法。 脈沖響應(yīng)不變法的基本思想是：使數(shù)字濾波器的單位脈沖響應(yīng)h(n)近似于模擬濾波器的單位脈沖響應(yīng)ha(t)，即使 其S平面和Z平面的映射關(guān)系為： 雙線性變換法的基本思想是：使描述數(shù)字濾波器的差分方程近似描述模擬濾波器的微分方程 S平面和Z平面的映射關(guān)系為： 雙線性變換法中的頻率變換是一種非線性變換，這種非線性引起的幅頻特性畸變可通過預(yù)變形矯正法而得到校正。 設(shè)計IIR 數(shù)字濾波器的一般步驟: （1）確定所需類型數(shù)字濾波器的技術(shù)指標(biāo)：通帶截止頻率ωp、通帶衰減αp、阻帶截止頻率ωs、阻帶衰減αs。 （2）將所需類型數(shù)字濾波器的技術(shù)指標(biāo)轉(zhuǎn)換成相應(yīng)類型模擬濾波器的技術(shù)指標(biāo)。 （3）設(shè)計該類型模擬濾波器 （4）通過復(fù)頻率變換將模擬濾波器轉(zhuǎn)換成所需類型的數(shù)字濾波器。 三、實驗內(nèi)容 1. 分別用脈沖響應(yīng)不變法和雙線性變換法設(shè)計一個巴特沃斯低通IIR數(shù)字濾波器，設(shè)計指標(biāo)參數(shù)為：在通帶內(nèi)頻率低于0.2p時，最大衰減小于1dB，在阻帶內(nèi)[0.3p，p]頻率區(qū)間上，最小衰減大于15dB。觀察并畫出所設(shè)計數(shù)字濾波器的幅頻特性曲線和相頻特性曲線，記錄帶寬和衰減量，檢查是否滿足要求。比較這兩種方法的優(yōu)缺點。 Matlab程序為： %脈沖響應(yīng)法 T=1;wp=0.2*pi;ws=0.3*pi;rp=1;as=15; %輸入低通濾波器要求 [n,wpo]=buttord(wp,ws,rp,as,s); %計算階數(shù) [B,A]=butter(n,wpo,s); %計算表達(dá)式分子分母的系數(shù)矩陣 [B1,A1]=impinvar(B,A); [Hk,w]=freqz(B1,A1); subplot(2,1,1); plot(w/pi,20*log10(abs(Hk))); %畫出濾波器損耗函數(shù)曲線 grid on; title((1)脈沖響應(yīng)不變法衰減曲線); xlabel(頻率（w/pi）); ylabel(幅度（dB）); %雙線性法 wp=2*tan(0.2*pi/2); ws=2*tan(0.3*pi/2); rp=1;as=15; [n,wpo]=buttord(wp,ws,rp,as,s); %計算階數(shù) [B,A]=butter(n,wpo,s); %計算表達(dá)式分子分母的系數(shù)矩陣 [B1,A1]=bilinear(B,A,1); [Hk,w]=freqz(B1,A1); subplot(2,1,2); plot(w/pi,20*log10(abs(Hk))); %畫出濾波器損耗函數(shù)曲線 title((1)雙線性法衰減曲線); grid on; xlabel(頻率（w/pi）);ylabel(幅度（dB）); 優(yōu)缺點比較： （1）脈沖響應(yīng)不變法會產(chǎn)生頻譜混疊，但具有很好的線性特性，其單位脈沖響應(yīng)完全模仿模擬濾波器的單位沖激響應(yīng)波形，時域逼近性好。適合于帶通、低通濾波器的設(shè)計。 （2）雙線性法很好地消除了頻譜混疊，但是其數(shù)字頻率與模擬頻率之間不具有線性關(guān)系。 2. 用雙線性變換法設(shè)計一個切比雪夫高通IIR數(shù)字濾波器，設(shè)計指標(biāo)參數(shù)為：在通帶內(nèi)頻率高于0.3KHz, 最大衰減小于1dB，在阻帶內(nèi)頻率低于0.2 KHz，最小衰減大于20dB，T=1ms。 畫出所設(shè)計數(shù)字濾波器的幅頻特性曲線和相頻特性曲線，觀察其通帶損耗和阻帶衰減是否滿足要求。 Matlab程序如下： %切比雪夫高通濾波器的設(shè)計 fp=2*pi*300*0.001; fs=2*pi*200*0.001; wp=2000*tan(fp/2); ws=2000*tan(fs/2); %進(jìn)行頻率變換 rp=1;as=20; %輸入高通濾波器要求 [n,wpo]=buttord(wp,ws,rp,as,s); %計算階數(shù) [B,A]=butter(n,wpo,high,s); %計算表達(dá)式分子分母的系數(shù)矩陣 [B1,A1]=impinvar(B,A); [Hk,w]=freqz(B,A); plot(w/pi,20*log10(abs(Hk))); %畫濾波器損耗函數(shù)曲線 grid on; title((3)切比雪夫高通IIR數(shù)字濾波器); xlabel(頻率（HZ）);ylabel(幅度（dB）); 3. 人體心電圖信號在測量過程中往往受到工業(yè)高頻干擾，所以必須經(jīng)過低通濾波處理后，才能作為判斷心臟功能的有用信息。下面給出一實際心電圖信號采樣序列樣本x(n)，其中存在高頻干擾。 用1所計的濾波器對心電圖信號采樣序列x(n)進(jìn)行仿真濾波處理，畫出處理前后的信號波形。 xn=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0]; subplot(2,1,1); plot(xn); title(xn); y1n=filter(B1,A1,xn); subplot(2,1,2); plot(y1n); title(xn); 四、實驗思考 1. 用雙線性變換法設(shè)計數(shù)字濾波器過程中，變換公式 中T取值，對設(shè)計結(jié)果有無影響？為什么？ 答：沒有。因為在第一步和第二步中從數(shù)字角頻率W變到模擬角頻率Ω，再從模擬角頻率Ω變到數(shù)字角頻率W，兩次變換是對稱的，只要兩次變換過程的T是相同的即可，T的取值是無關(guān)緊要的 2. 雙線性變換法中Ω和ω之間的關(guān)系是非線性的，在實驗中你注意到這種非線性關(guān)系了嗎？從哪幾種數(shù)字濾波器的幅頻特性曲線中可以觀察到這種非線性關(guān)系？ 答：注意到了。雙線性變換是從S平面映射到S1平面，再從S1平面映射到Z平面，一個線性相位的模擬濾波器經(jīng)過雙線性法變換后，就變成了非線性的了。切比雪夫的幅頻特性是非線性的。 五、實驗總結(jié)及心得體會 通過實驗學(xué)會了用雙線性變換法設(shè)計IIR數(shù)字濾波器的原理和方法，并且掌握IIR數(shù)字濾波器的Matlab實現(xiàn)方法，通過觀察對實際心電圖信號的濾波作用，獲得數(shù)字濾波的感性認(rèn)識。 實驗四 用窗函數(shù)法設(shè)計FIR數(shù)字濾波器 一、實驗?zāi)康? 1. 掌握用窗函數(shù)法設(shè)計FIR數(shù)字濾波器的原理和方法 2. 熟悉線性相位FIR數(shù)字濾波器特性 3. 了解各種窗函數(shù)對濾波特性的影響 二、實驗原理 窗函數(shù)法設(shè)計 FIR 濾波器的步驟為： （1）構(gòu)建希望逼近的理想頻率響應(yīng)函數(shù)及技術(shù)指標(biāo) （2）求濾波器的單位脈沖響應(yīng) 如果復(fù)雜，可對從采樣M個點，采樣值為，則： （3）根據(jù)對過渡帶及阻帶衰減的要求，選擇窗函數(shù)的形式，并估計窗口寬度N，設(shè)要求的過渡帶寬為，則 （4）計算濾波器的單位脈沖響應(yīng)： （5）求H(ejω)，分析其幅頻特性，若不滿足要求，可適當(dāng)改變窗函數(shù)形式或長度N，重復(fù)上述設(shè)計過程，以得到滿意的結(jié)果。 窗函數(shù)傅里葉變換W(ejω)的主瓣決定了H(ejω)過渡帶寬，W(ejω)的旁瓣大小和多少決定了H(ejω)在通帶和阻帶范圍內(nèi)波動幅度，常用的幾種窗函數(shù)有： 矩形窗；Hanning窗；Hamming窗；Blackmen窗；Kaiser窗 三、實驗內(nèi)容及步驟 1. 用升余弦窗設(shè)計一線性相位低通FIR 數(shù)字濾波器，截止頻率 。窗口長度N=15，33。要求在兩種窗口長度情況下，分別求出h(n)，打印出相應(yīng)的幅頻特性和相頻特性曲線，觀察3dB帶寬和20dB 帶寬，總結(jié)窗口長度N對濾波特性的影響 %用升余弦窗設(shè)計一線性相位低通FIR 數(shù)字濾波器 %N=15 hn1=fir1(14,pi/4,hanning(15)); figure(1); [Hn1,w]=freqz(hn1,1); subplot(2,1,1);plot(w/pi,abs(Hn1)); title(N=15的h(n)的幅頻曲線) xlabel(w);ylabel(幅度); w1=angle(Hn1); subplot(2,1,2);plot(w/pi,w1); title(N=15的h(n)的相頻曲線) xlabel(w);ylabel(angle(Hn1)) %N為33 hn2=fir1(32,pi/4,hanning(33)); figure(2); [Hn2,w]=freqz(hn2,1); subplot(2,1,1);plot(w/pi,abs(Hn2)); title(N=33的h(n)的幅頻曲線) title(N=33的h(n)的相頻曲線) xlabel(w);ylabel(幅度); w2=angle(Hn2); subplot(2,1,2);plot(w/pi,w2); xlabel(w);ylabel(angle(Hn2)) 對窗口長度N對濾波特性的影響的總結(jié)： 調(diào)整窗口長度N只能有效地控制過渡帶的寬度，當(dāng)N增大，主瓣幅度加高，同時旁瓣也加高，WRg（ω）的主瓣和旁瓣幅度變窄。導(dǎo)致波動頻率加快。因此，加大N，并不是解決吉布斯效應(yīng)的有效方法。 對窗口長度N對濾波特性的影響的總結(jié)： 調(diào)整窗口長度N只能有效地控制過渡帶的寬度，當(dāng)N增大，主瓣幅度加高，同時旁瓣也加高，WRg（ω）的主瓣和旁瓣幅度變窄。導(dǎo)致波動頻率加快。因此，加大N，并不是解決吉布斯效應(yīng)的有效方法。 2. N=33, , 用四種窗函數(shù)設(shè)計線性相位低通濾波器，繪制相應(yīng)的幅頻特性曲線，觀察3dB帶寬和20dB 帶寬以及阻帶最小衰減，比較四種窗函數(shù)對濾波器特性的影響 %矩形窗 hn1=fir1(32,pi/4,boxcar(33)); figure(1); [Hn1,w]=freqz(hn1,1); subplot(4,1,1);plot(w/pi,abs(Hn1)); title(（1）（2）矩形窗的幅頻曲線) xlabel(w);ylabel(幅度); w1=angle(Hn1); subplot(4,1,2); plot(w/pi,w1); xlabel(w);ylabel(angle(Hn)) %三角窗 hn2=fir1(32,pi/4,bartlett(33)); [Hn2,w]=freqz(hn2,1); subplot(4,1,3);plot(w/pi,abs(Hn2)); title(（3）(4）三角窗的幅頻曲線) xlabel(w);ylabel(幅度); w2=angle(Hn2); subplot(4,1,4);plot(w/pi,w2); xlabel(w);ylabel(angle(Hn2)) %漢寧窗 hn3=fir1(32,pi/4,hanning(33)); figure(2); [Hn3,w]=freqz(hn3,1); subplot(4,1,1);plot(w/pi,abs(Hn3)); title(（1）（2）漢寧窗的幅頻曲線) xlabel(w);ylabel(h(n）幅頻曲線); w3=angle(Hn3); subplot(4,1,2);plot(w/pi,w3); xlabel(w);ylabel(angle(Hn3)) %哈明窗 hn4=fir1(32,pi/4,hamming(33)); [Hn4,w]=freqz(hn4,1); subplot(4,1,3);plot(w/pi,abs(Hn4)); title(（3）(4）哈明窗的幅頻曲線) xlabel(w);ylabel(h(n）幅頻曲線); w4=angle(Hn4); subplot(4,1,4);plot(w/pi,w4); xlabel(w);ylabel(angle(Hn4)) 3. 調(diào)用信號產(chǎn)生函數(shù)xtg（見教材P295）產(chǎn)生具有加性噪聲的信號xt，顯示xt波形及其頻譜。設(shè)計一FIR 數(shù)字低通濾波器，從高頻噪聲中提取xt中的單頻調(diào)幅信號，要求信號幅頻失真小于0.1dB,噪聲頻譜衰減不小于60dB。 （1）觀察xt的頻譜，確定濾波器指標(biāo)參數(shù)。 （2）根據(jù)濾波器指標(biāo)選擇合適的窗函數(shù)，計算窗函數(shù)的長度N，設(shè)計一個FIR低通濾波器，繪圖顯示濾波器的頻響特性曲線 （3）用設(shè)計的FIR低通濾波器對xt進(jìn)行濾波，繪圖顯示濾波器輸出信號的時域和頻域波形圖 Matlab程序為： function xt=xtg %信號x(t)產(chǎn)生函數(shù)，并顯示信號的時域波形和幅頻特性曲線 %xt=xtg產(chǎn)生一個長度為N，有加性高頻噪聲的單調(diào)調(diào)幅信號xt，n=1000， %采樣頻率fs=1000HZ %載波頻率fc=fs/10=100HZ,調(diào)制正弦波頻率f0=fc/10=10HZ N=1000;Fs=1000;T=1/Fs;Tp=N*T;t=0:T:(N-1)*T; fc=Fs/10;f0=fc/10; mt=cos(2*pi*f0*t); ct=cos(2*pi*fc*t); xt=mt.*ct;nt=2*rand(1,N)-1; %設(shè)計高通濾波器hn,用于濾除噪聲nt中的低頻成分，生成高通噪聲 fp=150;fs=200; rp=0.1;as=70; %濾波器指標(biāo) fb=[fp,fs];m=[0,1]; dev=[10^(-as/20),(10^(rp/20)-1)/(10^(rp/20)+1)]; [n,fo,mo,W]=remezord(fb,m,dev,Fs); %確定remez函數(shù)所需參數(shù) hn=remez(n,fo,mo,W); %調(diào)用remez函數(shù)進(jìn)行設(shè)計，用于濾除噪聲nt中的低頻成分 yt=filter(hn,1,10*nt); %濾除隨機噪聲中的低頻成分，生成高通噪聲yt %以下為繪圖成分 xt=xt+yt;%噪聲加信號 fst=fft(xt,N);k=0:N-1;f=k/Tp; subplot(2,1,1);plot(t,xt); grid; xlabel(t/s);ylabel(x(t)); axis([0,Tp/5,min(xt),max(xt)]);title((a)xt信號加噪聲波形); subplot(2,1,2);plot(f,abs(fst)/max(abs(fst))); grid; title((b)xt信號加噪聲頻譜); axis([0,Fs/2,0,1.2]);xlabel(f/HZ); ylabel(幅度); %阻帶衰減as不低于60，所以選用布萊克曼窗 hn=fir1(119,pi/10,blackman(120)); figure(2); [Hn,w]=freqz(hn,1); subplot(2,1,1);plot(w/pi,abs(Hn)); title(濾波器的頻響特性); xlabel(w);ylabel(h(n）頻響特性曲線); %濾波過程 yn=conv(xt,hn); figure(3); subplot(2,1,1);plot(yn); title(xt經(jīng)過低通后的時域波形); [Hn,w]=freqz(yn,1); subplot(2,1,2);plot(w/pi,abs(Hn)/max(abs(Hn))); xlabel(w);title(xt經(jīng)過低通后的頻域波形); 四、實驗思考 1. 如果給定通帶截止頻率和阻帶截止頻率以及阻帶最小衰減，如何用窗函數(shù)法設(shè)計線性相位低通濾波器？ 答：1）根據(jù)阻帶衰減確定窗函數(shù)，然后根據(jù)通帶與阻帶截止頻率確定窗的長度N 2）寫出低通濾波器的頻域函數(shù) 3）求低通濾波器的時域函數(shù) 4）求= 2. 如果要求用窗函數(shù)法設(shè)計帶通濾波器，且給定上、下邊帶截止頻率為和，試求理想帶通的單位脈沖響應(yīng) 答：1）根據(jù)阻帶衰減確定窗函數(shù)，然后根據(jù)通帶與阻帶截止頻率確定窗的長度N 2）寫出帶通濾波器的頻域函數(shù) 3）求帶通濾波器的時域函數(shù) 4）求= 五. 實驗總結(jié) IIR濾波器的優(yōu)點是可利用模擬濾波器設(shè)計的結(jié)果，缺點是相位是非線性的。而FIR濾波器具有良好的線性相位。FIR濾波器線性相位的特點： 如果FIR濾波器的單位抽樣響應(yīng)H(N)為實數(shù)，而且滿足以下任一條件： 偶對稱H(N)＝H(N-1-N) 奇對稱H(N)＝-H(N-1-N) 且其對稱中心在N＝(N-1)/2處 則濾波器具有準(zhǔn)確的線性相位。 以下對FIR及IIR進(jìn)行對比： FIR濾波器 IIR濾波器 設(shè)計方法 一般無解析的設(shè)計公式，要借助計算機程序完成 利用AF的成果，可簡單、有效地完成設(shè)計 設(shè)計結(jié)果 線性相位(最大優(yōu)點) 只能得到幅頻特性，相頻特性未知， 如需要線性相位，須用全通網(wǎng)絡(luò)校 準(zhǔn)，但增加濾波器階數(shù)和復(fù)雜性 穩(wěn)定性 極點全部在原點(永遠(yuǎn)穩(wěn)定) 有穩(wěn)定性的問題 心得體會： 1.更加熟悉掌握用窗函數(shù)法設(shè)計FIR數(shù)字濾波器的原理和方法； 2. 熟悉了線性相位FIR數(shù)字濾波器特性； 3. 了解了各種窗函數(shù)對濾波特性的影響； 4.對FIR和IIR 有了更加深刻的理解. 5.實驗過程中也遇到了很多困難，在同學(xué)幫助下得到解決。