數(shù)字信號處理實(shí)驗(yàn)報(bào)告.doc

數(shù)字信號處理實(shí)驗(yàn)報(bào)告課程名稱：數(shù)字信號處理學(xué) 院：信息科學(xué)與工程學(xué)院專業(yè)班級：通信1502班學(xué)生姓名：侯子強(qiáng)學(xué) 號：0905140322指導(dǎo)教師：李宏2017年5月28日實(shí)驗(yàn)一 離散時間信號和系統(tǒng)響應(yīng)一. 實(shí)驗(yàn)?zāi)康?. 熟悉連續(xù)信號經(jīng)理想采樣前后的頻譜變化關(guān)系，加深對時域采樣定理的理解2. 掌握時域離散系統(tǒng)的時域特性3. 利用卷積方法觀察分析系統(tǒng)的時域特性4. 掌握序列傅里葉變換的計(jì)算機(jī)實(shí)現(xiàn)方法，利用序列的傅里葉變換對離散信號及系統(tǒng)響應(yīng)進(jìn)行頻域分析二、實(shí)驗(yà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)行采樣，要保證采樣頻率fs2fm，fm為信號的最高頻率，才可能由采樣信號無失真地恢復(fù)出原模擬信號計(jì)算機(jī)實(shí)現(xiàn)時，利用計(jì)算機(jī)計(jì)算上式并不方便，因此我們利用采樣序列的傅里葉變換來實(shí)現(xiàn)，即而為采樣序列的傅里葉變換2. 時域中，描述系統(tǒng)特性的方法是差分方程和單位脈沖響應(yīng)，頻域中可用系統(tǒng)函數(shù)描述系統(tǒng)特性。已知輸入信號，可以由差分方程、單位脈沖響應(yīng)或系統(tǒng)函數(shù)求出系統(tǒng)對于該輸入信號的響應(yīng)。本實(shí)驗(yàn)僅在時域求解，對于差分方程可用Matlab中的工具箱函數(shù)filter()函數(shù)求解 一個時域離散線性時不變系統(tǒng)的輸出與輸入間的關(guān)系為：可用Matlab中的工具箱函數(shù)conv()函數(shù)求解三、實(shí)驗(yàn)內(nèi)容及步驟1. 時域采樣定理的驗(yàn)證給定模擬信號：式中 。其幅頻特性如圖所示：選擇三種采樣頻率Fs=1kHz, 300Hz, 200Hz, 生成采樣序列分別用序列表示。編寫程序計(jì)算三個序列的幅頻特性曲線，并繪圖顯示。觀察在折疊頻率附近與連續(xù)信號頻譜有無明顯差別，分析頻譜混疊現(xiàn)象。實(shí)驗(yàn)程序如下%時域采樣定理的驗(yàn)證%Fs=1KHzTp=64/1000; %Tp=64ms Fs=1000;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*20.5;omega=pi*50*20.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*FTxa(nT),Fs=1000Hz); xlabel(f(Hz);ylabel(幅度);axis(0,Fs,0,1.2*max(abs(Xk)%Fs=300HzTp=64/1000; Fs=300;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*20.5;omega=pi*50*20.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*FTxa(nT),Fs=300Hz); xlabel(f(Hz);ylabel(幅度);axis(0,Fs,0,1.2*max(abs(Xk)%Fs=200HzTp=64/1000; %64ms Fs=300;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*20.5;omega=pi*50*20.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*FTxa(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);四、實(shí)驗(yàn)思考1. 在分析理想采樣序列特性的實(shí)驗(yàn)中，采樣頻率不同時，相應(yīng)理想采樣序列的傅里葉變換頻譜的數(shù)字頻率度量是否都相同？它們所對應(yīng)的模擬頻率是否相同？為什么？答：當(dāng)采樣頻率不同時，數(shù)字度量不同，但是模擬頻率相同。因?yàn)閿?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)過低通濾波器，把信號的高頻分量濾掉，時域信號會有何變化？用前面第二個實(shí)驗(yàn)結(jié)果進(jìn)行分析說明答：把信號經(jīng)過低通濾波器，把信號的高頻成分濾掉，時域信號的劇烈將變得平滑。五、 實(shí)驗(yàn)心得及體會通過本次實(shí)驗(yàn)我重新溫習(xí)了MATLAB這個軟件的使用方法，運(yùn)行環(huán)境。通過這款軟件使我們的學(xué)習(xí)更加便利。實(shí)驗(yàn)二 用FFT對信號作頻譜分析一、實(shí)驗(yàn)?zāi)康?. 進(jìn)一步加深DFT算法原理和基本性質(zhì)的理解2. 掌握用FFT對連續(xù)信號和時域離散信號進(jìn)行頻譜分析的方法3. 了解用FFT進(jìn)行頻譜分析時可能出現(xiàn)的分析誤差及其原因，以便在實(shí)際中正確應(yīng)用FFT二、實(shí)驗(yàn)原理用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)行。三、實(shí)驗(yà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); %計(jì)算x1n的8點(diǎn)DFTX1k16=fft(x1n,16); %計(jì)算x1n的16點(diǎn)DFTN=8;f=2/N*(0:N-1);figure(1);subplot(1,2,1);stem(f,abs(X1k8),.); %繪制8點(diǎn)DFT的幅頻特性圖title(1a) 8點(diǎn)DFTx_1(n);xlabel(/);ylabel(幅度);N=16;f=2/N*(0:N-1);subplot(1,2,2);stem(f,abs(X1k16),.); %繪制8點(diǎn)DFT的幅頻特性圖title(1a) 16點(diǎn)DFTx_1(n);xlabel(/);ylabel(幅度);%x2n 和 x3nM=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點(diǎn)DFT的幅頻特性圖title(2a) 8點(diǎn)DFTx_2(n);xlabel(/);ylabel(幅度);subplot(2,2,3);stem(f,abs(X3k8),.); %繪制8點(diǎn)DFT的幅頻特性圖title(3a) 8點(diǎn)DFTx_3(n);xlabel(/);ylabel(幅度);N=16;f=2/N*(0:N-1);subplot(2,2,2);stem(f,abs(X2k16),.); %繪制8點(diǎn)DFT的幅頻特性圖title(2a) 16點(diǎn)DFTx_2(n);xlabel(/);ylabel(幅度);subplot(2,2,4);stem(f,abs(X3k16),.); %繪制8點(diǎn)DFT的幅頻特性圖title(3a) 16點(diǎn)DFTx_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,.);%endn=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點(diǎn)采樣 X6k16=fft(x6nT); %計(jì)算x6nT的16點(diǎn)DFT X6k16=fftshift(X6k16); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %頻率分辨率Fk=-N/2:N/2-1;fk=k*F; %產(chǎn)生16點(diǎn)DFT對應(yīng)的采樣點(diǎn)頻率（以零頻率為中心） subplot(3,1,1);stem(fk,abs(X6k16),.);box on %繪制8點(diǎn)DFT的幅頻特性圖 title(6a) 16點(diǎn)|DFTx_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點(diǎn)采樣 X6k32=fft(x6nT); %計(jì)算x6nT的32點(diǎn)DFT X6k32=fftshift(X6k32); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %頻率分辨率F k=-N/2:N/2-1;fk=k*F; %產(chǎn)生16點(diǎn)DFT對應(yīng)的采樣點(diǎn)頻率（以零頻率為中心） subplot(3,1,2);stem(fk,abs(X6k32),.);box on %繪制8點(diǎn)DFT的幅頻特性圖 title(6b) 32點(diǎn)|DFTx_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點(diǎn)采樣 X6k64=fft(x6nT); %計(jì)算x6nT的64點(diǎn)DFT X6k64=fftshift(X6k64); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %頻率分辨率F k=-N/2:N/2-1;fk=k*F; %產(chǎn)生16點(diǎn)DFT對應(yīng)的采樣點(diǎn)頻率（以零頻率為中心） subplot(3,1,3);stem(fk,abs(X6k64),.); box on%繪制8點(diǎn)DFT的幅頻特性圖 title(6a) 64點(diǎn)|DFTx_6(nT)|);xlabel(f(Hz);ylabel(幅度); axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k64)四、實(shí)驗(yàn)思考1. 在N=8時， 和的幅頻特性會相同嗎？為什么？N=16時呢？答：在N=8時， 和 的幅頻特性相同，而N=16時不相同。因?yàn)?= ,所以 和 的8點(diǎn)DFT的模相等。但當(dāng)N=16時， 和 不滿足循環(huán)移位關(guān)系，所以兩者幅頻特性不相同。2. 對于周期序列，如果周期不知道，如何用FFT進(jìn)行譜分析？答：周期信號的周期預(yù)先不知道時,可先截取M點(diǎn)進(jìn)行DFT,再將截取長度擴(kuò)大1倍截取,比較結(jié)果,如果二者的差別滿足分析誤差要求,則可以近似表示該信號的頻譜,如果不滿足誤差要求就繼續(xù)將截取長度加倍,重復(fù)比較,直到結(jié)果滿足要求 。五、實(shí)驗(yàn)總結(jié)及心得體會通過實(shí)驗(yàn)進(jìn)一步加深DFT算法原理和基本性質(zhì)的理解，掌握用FFT對連續(xù)信號和時域離散信號進(jìn)行頻譜分析的方法，了解用FFT進(jìn)行頻譜分析時可能出現(xiàn)的分析誤差及其原因，以便在實(shí)際中正確應(yīng)用FFT。實(shí)驗(yàn)三 用雙線性變換法設(shè)計(jì)IIR數(shù)字濾波器一、實(shí)驗(yàn)?zāi)康?. 熟悉用雙線性變換法設(shè)計(jì)IIR數(shù)字濾波器的原理和方法2. 掌握IIR數(shù)字濾波器的Matlab實(shí)現(xiàn)方法3. 通過觀察對實(shí)際心電圖信號的濾波作用，獲得數(shù)字濾波的感性認(rèn)識二、實(shí)驗(yàn)原理設(shè)計(jì)IIR數(shù)字濾波器一般采用間接設(shè)計(jì)法脈沖響應(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è)計(jì)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è)計(jì)該類型模擬濾波器（4）通過復(fù)頻率變換將模擬濾波器轉(zhuǎn)換成所需類型的數(shù)字濾波器。三、實(shí)驗(yàn)內(nèi)容1. 分別用脈沖響應(yīng)不變法和雙線性變換法設(shè)計(jì)一個巴特沃斯低通IIR數(shù)字濾波器，設(shè)計(jì)指標(biāo)參數(shù)為：在通帶內(nèi)頻率低于0.2p時，最大衰減小于1dB，在阻帶內(nèi)0.3p，p頻率區(qū)間上，最小衰減大于15dB。觀察并畫出所設(shè)計(jì)數(shù)字濾波器的幅頻特性曲線和相頻特性曲線，記錄帶寬和衰減量，檢查是否滿足要求。比較這兩種方法的優(yōu)缺點(diǎn)。Matlab程序?yàn)椋?脈沖響應(yīng)法T=1;wp=0.2*pi;ws=0.3*pi;rp=1;as=15; %輸入低通濾波器要求n,wpo=buttord(wp,ws,rp,as,s); %計(jì)算階數(shù)B,A=butter(n,wpo,s); %計(jì)算表達(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); %計(jì)算階數(shù)B,A=butter(n,wpo,s); %計(jì)算表達(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)缺點(diǎn)比較：（1）脈沖響應(yīng)不變法會產(chǎn)生頻譜混疊，但具有很好的線性特性，其單位脈沖響應(yīng)完全模仿模擬濾波器的單位沖激響應(yīng)波形，時域逼近性好。適合于帶通、低通濾波器的設(shè)計(jì)。（2）雙線性法很好地消除了頻譜混疊，但是其數(shù)字頻率與模擬頻率之間不具有線性關(guān)系。2. 用雙線性變換法設(shè)計(jì)一個切比雪夫高通IIR數(shù)字濾波器，設(shè)計(jì)指標(biāo)參數(shù)為：在通帶內(nèi)頻率高于0.3KHz, 最大衰減小于1dB，在阻帶內(nèi)頻率低于0.2 KHz，最小衰減大于20dB，T=1ms。 畫出所設(shè)計(jì)數(shù)字濾波器的幅頻特性曲線和相頻特性曲線，觀察其通帶損耗和阻帶衰減是否滿足要求。Matlab程序如下：%切比雪夫高通濾波器的設(shè)計(jì)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); %計(jì)算階數(shù)B,A=butter(n,wpo,high,s); %計(jì)算表達(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)過低通濾波處理后，才能作為判斷心臟功能的有用信息。下面給出一實(shí)際心電圖信號采樣序列樣本x(n)，其中存在高頻干擾。 用1所計(jì)的濾波器對心電圖信號采樣序列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);四、實(shí)驗(yàn)思考1. 用雙線性變換法設(shè)計(jì)數(shù)字濾波器過程中，變換公式 中T取值，對設(shè)計(jì)結(jié)果有無影響？為什么？答：沒有。因?yàn)樵诘谝徊胶偷诙街袕臄?shù)字角頻率W變到模擬角頻率，再從模擬角頻率變到數(shù)字角頻率W，兩次變換是對稱的，只要兩次變換過程的T是相同的即可，T的取值是無關(guān)緊要的2. 雙線性變換法中和之間的關(guān)系是非線性的，在實(shí)驗(yàn)中你注意到這種非線性關(guān)系了嗎？從哪幾種數(shù)字濾波器的幅頻特性曲線中可以觀察到這種非線性關(guān)系？答：注意到了。雙線性變換是從S平面映射到S1平面，再從S1平面映射到Z平面，一個線性相位的模擬濾波器經(jīng)過雙線性法變換后，就變成了非線性的了。切比雪夫的幅頻特性是非線性的。五、實(shí)驗(yàn)總結(jié)及心得體會 通過實(shí)驗(yàn)學(xué)會了用雙線性變換法設(shè)計(jì)IIR數(shù)字濾波器的原理和方法，并且掌握IIR數(shù)字濾波器的Matlab實(shí)現(xiàn)方法，通過觀察對實(shí)際心電圖信號的濾波作用，獲得數(shù)字濾波的感性認(rèn)識。實(shí)驗(yàn)四 用窗函數(shù)法設(shè)計(jì)FIR數(shù)字濾波器一、實(shí)驗(yàn)?zāi)康?. 掌握用窗函數(shù)法設(shè)計(jì)FIR數(shù)字濾波器的原理和方法2. 熟悉線性相位FIR數(shù)字濾波器特性3. 了解各種窗函數(shù)對濾波特性的影響二、實(shí)驗(yàn)原理 窗函數(shù)法設(shè)計(jì) FIR 濾波器的步驟為： （1）構(gòu)建希望逼近的理想頻率響應(yīng)函數(shù)及技術(shù)指標(biāo)（2）求濾波器的單位脈沖響應(yīng) 如果復(fù)雜，可對從采樣M個點(diǎn)，采樣值為，則：（3）根據(jù)對過渡帶及阻帶衰減的要求，選擇窗函數(shù)的形式，并估計(jì)窗口寬度N，設(shè)要求的過渡帶寬為，則（4）計(jì)算濾波器的單位脈沖響應(yīng)： （5）求H(ej)，分析其幅頻特性，若不滿足要求，可適當(dāng)改變窗函數(shù)形式或長度N，重復(fù)上述設(shè)計(jì)過程，以得到滿意的結(jié)果。窗函數(shù)傅里葉變換W(ej)的主瓣決定了H(ej)過渡帶寬，W(ej)的旁瓣大小和多少決定了H(ej)在通帶和阻帶范圍內(nèi)波動幅度，常用的幾種窗函數(shù)有：矩形窗；Hanning窗；Hamming窗；Blackmen窗；Kaiser窗三、實(shí)驗(yàn)內(nèi)容及步驟1. 用升余弦窗設(shè)計(jì)一線性相位低通FIR 數(shù)字濾波器，截止頻率 。窗口長度N=15，33。要求在兩種窗口長度情況下，分別求出h(n)，打印出相應(yīng)的幅頻特性和相頻特性曲線，觀察3dB帶寬和20dB 帶寬，總結(jié)窗口長度N對濾波特性的影響%用升余弦窗設(shè)計(jì)一線性相位低通FIR 數(shù)字濾波器%N=15hn1=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為33hn2=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è)計(jì)線性相位低通濾波器，繪制相應(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è)計(jì)一FIR 數(shù)字低通濾波器，從高頻噪聲中提取xt中的單頻調(diào)幅信號，要求信號幅頻失真小于0.1dB,噪聲頻譜衰減不小于60dB。（1）觀察xt的頻譜，確定濾波器指標(biāo)參數(shù)。（2）根據(jù)濾波器指標(biāo)選擇合適的窗函數(shù)，計(jì)算窗函數(shù)的長度N，設(shè)計(jì)一個FIR低通濾波器，繪圖顯示濾波器的頻響特性曲線（3）用設(shè)計(jì)的FIR低通濾波器對xt進(jìn)行濾波，繪圖顯示濾波器輸出信號的時域和頻域波形圖Matlab程序?yàn)椋篺unction 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=10HZN=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è)計(jì)高通濾波器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è)計(jì)，用于濾除噪聲nt中的低頻成分yt=filter(hn,1,10*nt); %濾除隨機(jī)噪聲中的低頻成分，生成高通噪聲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)過低通后的頻域波形);四、實(shí)驗(yàn)思考1. 如果給定通帶截止頻率和阻帶截止頻率以及阻帶最小衰減，如何用窗函數(shù)法設(shè)計(jì)線性相位低通濾波器？答：1）根據(jù)阻帶衰減確定窗函數(shù)，然后根據(jù)通帶與阻帶截止頻率確定窗的長度N2）寫出低通濾波器的頻域函數(shù)3）求低通濾波器的時域函數(shù)4）求=2. 如果要求用窗函數(shù)法設(shè)計(jì)帶通濾波器，且給定上、下邊帶截止頻率為和，試求理想帶通的單位脈沖響應(yīng)答：1）根據(jù)阻帶衰減確定窗函數(shù)，然后根據(jù)通帶與阻帶截止頻率確定窗的長度N2）寫出帶通濾波器的頻域函數(shù)3）求帶通濾波器的時域函數(shù)4）求=五. 實(shí)驗(yàn)總結(jié)IIR濾波器的優(yōu)點(diǎn)是可利用模擬濾波器設(shè)計(jì)的結(jié)果，缺點(diǎn)是相位是非線性的。而FIR濾波器具有良好的線性相位。FIR濾波器線性相位的特點(diǎn)： 如果FIR濾波器的單位抽樣響應(yīng)H(N)為實(shí)數(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è)計(jì)方法一般無解析的設(shè)計(jì)公式，要借助計(jì)算機(jī)程序完成利用AF的成果，可簡單、有效地完成設(shè)計(jì)設(shè)計(jì)結(jié)果線性相位(最大優(yōu)點(diǎn))只能得到幅頻特性，相頻特性未知， 如需要線性相位，須用全通網(wǎng)絡(luò)校 準(zhǔn)，但增加濾波器階數(shù)和復(fù)雜性穩(wěn)定性極點(diǎn)全部在原點(diǎn)(永遠(yuǎn)穩(wěn)定)有穩(wěn)定性的問題 心得體會：1.更加熟悉掌握用窗函數(shù)法設(shè)計(jì)FIR數(shù)字濾波器的原理和方法；2. 熟悉了線性相位FIR數(shù)字濾波器特性；3. 了解了各種窗函數(shù)對濾波特性的影響；4.對FIR和IIR 有了更加深刻的理解.5.實(shí)驗(yàn)過程中也遇到了很多困難，在同學(xué)幫助下得到解決。