醫(yī)學(xué)信號(hào)處理報(bào)告三

西安交通大學(xué)實(shí)驗(yàn)報(bào)告成績(jī)實(shí)驗(yàn)日期：交報(bào)告日期：報(bào)告退發(fā)：教師審批簽字：第頁(yè)（共頁(yè)）實(shí)驗(yàn)日期：交報(bào)告日期：報(bào)告退發(fā)：教師審批簽字：第頁(yè)（共頁(yè)）2013年5月11日2013年5月11日（訂正、重做）課程： 醫(yī)學(xué)信號(hào)處理專業(yè)班號(hào)組別.姓名學(xué)號(hào).同組者 實(shí)驗(yàn)名稱： 陵機(jī)佶號(hào)分析一.實(shí)驗(yàn)要求1）產(chǎn)生1000點(diǎn)的白噪聲信號(hào)，并計(jì)算它的均值、均方值、均方根值、方差。2）計(jì)算該1000點(diǎn)白噪聲的自相關(guān)函數(shù)并作圖顯示白噪聲與它的自相關(guān)函數(shù)。3）計(jì)算腦電信號(hào)的均值、均方值、均方根值、方差，計(jì)算腦電信號(hào)的自相關(guān)函數(shù)并作圖顯示腦電信號(hào)與它的自相關(guān)函數(shù)。（zuoye3_leeg.mat）4）計(jì)算含有噪聲的心電信號(hào)的自相關(guān)函數(shù)并作圖顯示含有噪聲的心電信號(hào)與它的自相關(guān)函數(shù)。（zuoye3_lzshecg.mat）二.實(shí)驗(yàn)原理自相關(guān)函數(shù)是信號(hào)在時(shí)域中特性的平均度量，它用來(lái)描述信號(hào)在一個(gè)時(shí)刻的取值與另一時(shí)刻取值的依賴關(guān)系，其定義式為：三.實(shí)驗(yàn)步驟與結(jié)果分析1）產(chǎn)生1000點(diǎn)的白噪聲信號(hào)，并計(jì)算它的均值、均方值、均方根值、方差。根據(jù)題目的要求，利用randn產(chǎn)生均值為0,方差為1的隨機(jī)白噪聲。程序：clearcic%產(chǎn)生白噪聲信號(hào)fs=1000;t=O:l/fs:l;x=randn(l,length(t));%計(jì)算均值meanx=mean(x)%計(jì)算均方值junfangx=sum(x(l,:).A2)/fs%計(jì)算均方根值junfanggx=junfangxA(l/2)%計(jì)算方差sum((x(l,:)-meanx).A2)/fs程序運(yùn)行完成后，輸出結(jié)果：均值：meanx=0.0305均方值：junfangx=0.9346均方根值：junfanggx=0.9667方差：fangchax=0.9337分析：和設(shè)計(jì)時(shí)均值為0,方差為1略有出入，但是從統(tǒng)計(jì)的角度看，是正確的。2)計(jì)算該1000點(diǎn)白噪聲的自相關(guān)函數(shù)并作圖顯示白噪聲與它的自相關(guān)函數(shù)根據(jù)題目的要求，求出白噪聲的自相關(guān)函數(shù)，延遲為200。程序：%計(jì)算自相關(guān)函數(shù)lag=200;[c,lags]=xcorr(x,lag,'biased');%lag為M,lags為下標(biāo)m%作圖顯示白噪聲subplot(2,l,l)plot(t,x)xlabel('t');ylabel('x(t)');title。白噪聲x(t)')%作圖顯示白噪聲的自相關(guān)subplot(2,l,2)plot(lags/fs,c)xlabel('t');ylabel('Rxx(t)');title('白噪聲的自相關(guān)')程序運(yùn)行完成后，繪制出了白噪聲的時(shí)域圖像和他的自相關(guān)函數(shù)的圖像，如圖(1)：-E圖(1)1000個(gè)點(diǎn)的白噪聲和在maxlags=200時(shí)其自相關(guān)函數(shù)分析：從圖中可以看出白噪聲的自相關(guān)函數(shù)為一個(gè)在0時(shí)刻的沖擊。3)計(jì)算腦電信號(hào)的均值、均方值、均方根值、方差，計(jì)算腦電信號(hào)的自相關(guān)函數(shù)并作圖顯示腦電信號(hào)與它的自相關(guān)函數(shù)加載eeg數(shù)據(jù)，loadzuoye3_leeg.mato其中ecgl變量存放的是原始腦電信號(hào)。根據(jù)要求計(jì)算各個(gè)參數(shù)的值。程序：clearcic%讀入eeg數(shù)據(jù)load('E:\課程'醫(yī)學(xué)信 乍業(yè)3\zuoye3_leeg.mat')fs=250;t=0:l/fs:10-l/fs;%計(jì)算均值meanx=mean(eegl)%計(jì)算均方值junfangx=sum(eegl(l,:).A2)/length(eegl)%計(jì)算均方根值junfanggx=junfangxA(l/2)%計(jì)算方差fangchax=sum((eegl(l,:)-meanx).A2)/length(eegl)%計(jì)算自相關(guān)函數(shù)[c,lags]=xcorr(eegl,'biased');%lag為M,lags為下標(biāo)m%作圖顯示白噪聲subplot(2,l,l)plot(t,eegl)xlabel('t');ylabel('eeg幅值')；title('腦電信號(hào)')%作圖顯示白噪聲的自相關(guān)subplot(2,l,2)plot(lags/fs,c)xlabel('t');ylabel('eeg自相關(guān)')；title('白噪聲的自相關(guān)')程序運(yùn)行完成后，輸出結(jié)果：均值：meanx=-0.0198均方值：junfangx=5.5061e-04均方根值：junfanggx=0.0235均方差：fangchax=1.5968e-04

程序運(yùn)行完成后，繪制出了腦電信號(hào)的時(shí)域圖像和它的自相關(guān)函數(shù)的圖像，如圖(3)：圖（3）腦電信號(hào)的原始波形和其自相關(guān)函數(shù)圖象4）計(jì)算含有噪聲的心電信號(hào)的自相關(guān)函數(shù)并作圖顯示含有噪聲的心電信號(hào)與它的自相關(guān)函數(shù)加載eeg數(shù)據(jù)，loadzuoye3_leeg.mato其中ecgl變量存放的是原始腦電信號(hào)。程序：clearcic%讀入eeg數(shù)據(jù)load('E:\課程、醫(yī)學(xué)信號(hào)處理\作業(yè)3\zuoye3_lzshecg.mat')；fs=250;t=0:l/fs:10-l/fs;%計(jì)算自相關(guān)函數(shù)[c,lags]=xcorr(zshecgl,'biased')；%lag為M,lags為下標(biāo)m%作圖顯示白噪聲subplot(2,l,l)plot(t,zshecgl)xlabel('t')；ylabel('ecg幅值')；title。含有噪聲的心電信號(hào)')%作圖顯示白噪聲的自相關(guān)subplot(2,l,2)plot(lags/fs,c)xlabel('t');ylabel('eeg自相關(guān)')；title('含有噪聲的心電信號(hào)的自相關(guān)')程序運(yùn)行完成后，繪制出了有噪聲的心電信號(hào)和其自相關(guān)函數(shù)的圖像，如圖(4)：含H吟聲的心電偌：；8禽含H吟聲的心電偌：；8禽W號(hào)聲的4屯■:;力■相關(guān)圖(4)含有噪聲的心電信號(hào)和其自相關(guān)函數(shù)分析:可以從圖中看出由于心電信號(hào)由于噪聲的原因使得自相關(guān)函數(shù)在0處的沖擊變得很大。四.實(shí)驗(yàn)總結(jié)通過(guò)本次試驗(yàn)，使我對(duì)隨機(jī)信號(hào)的分析以及自相關(guān)函數(shù)等概念有了更加深入的認(rèn)識(shí)。熟練地掌握了用Matlab求取隨機(jī)信號(hào)的各個(gè)參數(shù)的方法。第頁(yè)(共頁(yè))課程：_醫(yī)學(xué)信號(hào)處理實(shí)驗(yàn)日期：2013年5月11日專業(yè)班號(hào)碩2102組別交報(bào)告日期：2013年5月11日姓名—陳星學(xué)號(hào)3112322010報(bào)告退發(fā)：(訂正、重做)同組者教師審批簽字:西安交通大學(xué)實(shí)驗(yàn)報(bào)告成績(jī)西安交通大學(xué)實(shí)驗(yàn)報(bào)告成績(jī)實(shí)驗(yàn)名稱： 生理佶號(hào)的功率譜估計(jì)一.實(shí)驗(yàn)要求1)已知隨機(jī)信號(hào)x(n)=sin(2*pi*fl*n)+sin(2*pi*f*2n)+w(n).n=l,,,,,N,fl=0.05,12=0.12,w(n)為白噪聲。①N=1024時(shí)，當(dāng)M=1024,M=256時(shí)，求自相關(guān)法功率譜密度并作圖。②當(dāng)N=256,N=1024時(shí)，求周期圖法功率譜密度并作圖。③用welch法計(jì)算功率譜密度并作圖。2)信號(hào)為作業(yè)1獲取的心電信號(hào)，取一個(gè)周期，求周期圖法功率譜密度并作圖。3)信號(hào)為作業(yè)1獲取的心電信號(hào)，取10個(gè)周期，每段一個(gè)周期，用不重疊分段法求平均周期圖法功率譜密度并作圖。4)信號(hào)為作業(yè)1獲取的腦電信號(hào)，當(dāng)N=256,N=1024時(shí)，求周期圖法功率譜密度并作圖。5)信號(hào)為作業(yè)1獲取的腦電信號(hào)，每段長(zhǎng)256時(shí)，用welch法計(jì)算功率譜密度并作圖。生理信號(hào)抽樣頻率250Hz。二.實(shí)驗(yàn)原理D直接法直接法功率譜估計(jì)是間接法功率譜估計(jì)的一個(gè)特例，又稱為周期圖法，它是把隨機(jī)信號(hào)的N個(gè)觀察值直接進(jìn)行傅里葉變換，得到，然后取其幅值的平方，再除以N,作為對(duì)功率譜的估計(jì)。2)改進(jìn)的周期圖法將N點(diǎn)的觀察值分成L個(gè)數(shù)據(jù)段，每段的數(shù)據(jù)為M,然后計(jì)算L個(gè)數(shù)據(jù)段的周期圖的平均，作為功率譜的估計(jì)，以此來(lái)改善用N點(diǎn)觀察數(shù)據(jù)直接計(jì)算的周期圖的方差特性。根據(jù)分段方法的不同，又可以分為Welch法和Bartlett法。3)Welch法所分的數(shù)據(jù)段可以互相重疊，選用的數(shù)據(jù)窗可以是任意窗。三.實(shí)驗(yàn)步驟與結(jié)果分析1)分別利用自相關(guān)法，周期圖法，和welch法計(jì)算隨機(jī)信號(hào)的功率譜。程序：clearde%信號(hào)產(chǎn)生N=1024;n=l:l:N;fl=0.05;f2=0.12；x=sin(2*pi*fl*n)+sin(2*pi*f2*n)+randn(l,N);%相關(guān)函數(shù)估計(jì)c=xcorr(x,'biased');%譜估計(jì)M=Nlenth=length(x);px2047=abs(fft(c,2047));%譜估計(jì)M=256c2=[c(lenth-255:lenth-l),c(lenth:lenth+255)];px511=abs(fft(c2));%繪圖p=fftshift(px2047);p2=fftshift(px511);ml=-1023:1023;m2=-255:255;figure(l)subplot(2,l,l)plot(ml*N/2047,20*logl0(p))title('M=1024')subplot(2,l,2)plot(m2*N/511,20*logl0(p2))title('M=512')

mm圖（1）M=1024和M=256的自相關(guān)法功率譜程序：clearcic%信號(hào)產(chǎn)生N=1024;n=l:l:N;fl=0.05;f2=0.12;x=sin(2*pi*fl*n)+2*sin(2*pi*f2*n)+randn(l,N);%譜估計(jì)N=1024pxl024=abs(fft(x,1024)),A2/1024;%N=256px256=abs(fft(x,256)).A2/256;%繪圖p=fftshift(pxl024);p2=fftshift(px256);ml=-511:512;m2=-127:128;figure(l)subplot(2,l/l)plot(ml,20*logl0(p))title('N=1024')subplot。,1,2)plot(m2,20*logl0(p2))title('N=256')程序運(yùn)行結(jié)果如圖(2)：圖(2)N=1024和N=256的周期圖法功率譜程序：clearde%信號(hào)產(chǎn)生N=1024;n=l:l:N;fl=0.05;f2=0.12;xn=sin(2*pi*fl*n)+sin(2*pi*f2*n)+randn(l,N);%welch法海寧窗256個(gè)數(shù)據(jù)為一段重合128個(gè)數(shù)據(jù)L=7w=hanning(256),;pxx=(abs(fft(w.*xn(l:256))).A2+abs(fft(w.*xn(129:384))).A2+abs(fft(w.*xn(257:512))).A2+abs(fft(w.*xn(385:640))).A2+abs(fft(w.*xn(513:768))).A2+abs(fft(w.*xn(641:896))).A2+abs(fft(w.*xn(769:1024))).A2)/(256*7);p=fftshift(pxx);plot(-127:128,10*logl0(p));title('welch法')；程序運(yùn)行結(jié)果如圖(3)：法法圖(3)welch法功率譜2)對(duì)心電信號(hào)取…個(gè)周期，用周期圖法求功率譜密度并作圖程序：clearcic%讀入ecg數(shù)據(jù)load('E:\課程'醫(yī)學(xué)信號(hào)處理、作業(yè)3\zuoye3_lzshecg.mat')；N=2500fs=250;m=-1249:1250;px2500=abs(fft(zshecgl,2500)).A2/N;p=fftshift(px2500);plot(m*fs/N,10*logl0(p));title('心電信號(hào)周期圖法功率譜(一個(gè)周期)')；程序運(yùn)行結(jié)果如圖(4)：心電俱；青明陰法正章謂＜7WW!）圖（4）心電信號(hào)一個(gè)周期周期圖法功率譜3)對(duì)心電信號(hào)取10個(gè)周期，每段…個(gè)周期，用不重疊分段法求平均周期圖法功率譜密度并作圖。程序：clearcic%讀入ecg數(shù)據(jù)load('E:\課程'醫(yī)學(xué)信號(hào)處理、作業(yè)3\zuoye3_lzshecg.mat')；N=2500;fs=250;m=-124:125%t=0:l/fs:10-l/fs;figure(2)pxx=(abs(fft(zshecgl(l:250))).A2+abs(fft(zshecgl(251:500))).A2+abs(fft(zshecgl(501:750))).A2+abs(fft(zshecgl(751:1000))).A2+abs(fft(zshecgl(1001:1250))).A2+abs(fft(zshecgl(1251:1500))).A2+abs(fft(zshecgl(1501:1750))).A2+abs(fft(zshecgl(1751:2000))).A2+abs(fft(zshecgl(2001:2250))).A2+abs(fft(zshecgl(2251:2500))).A2)/250*10;p=fftshift(pxx);plot(m*fs/250,10*logl0(p));title('心電信號(hào)周期圖法功率譜(10個(gè)周期)')；

程序運(yùn)行完成后，結(jié)果如圖（5）：圖（5）心電信號(hào)一個(gè)周期周期圖法功率譜分析：可以看出10個(gè)周期的功率譜密度圖像要比一個(gè)周期的平滑許多。4）求腦電信號(hào)當(dāng)N=256,N=1024時(shí)的周期圖法功率譜密度并作圖。程序：clearcic%讀入eeg數(shù)據(jù)load('E:\課程'醫(yī)學(xué)信號(hào)處理、作業(yè)3\zuoye3_leeg.mat')fs=250%譜估計(jì)N=1024pxl024=abs(fft(eegl,1024)).A2/1024;%N=256px256=abs什ft(eegl,256)).△2/256;%繪圖p=fftshift(pxl024);p2=fftshift(px256);ml=-5U:512;m2=-127:128；figure(l)subplot(2,l,l)plot(ml*fs/1024,20*logl0(p))title('N=1024周期圖法eeg功率譜')subplot(2,l,2)plot(m2*fs/256,20*logl0(p2))title('N=256周期圖法eeg功率譜')程序運(yùn)行完成后，結(jié)果如圖(6)：圖（6）N=1024和N=256時(shí)的腦電信號(hào)周期圖法功率譜5）求腦電信號(hào)每段長(zhǎng)256時(shí)welch法的功率譜密度并作圖。程序:clearcic%讀入eeg數(shù)據(jù)load('E:\課程'醫(yī)學(xué)信號(hào)處理\作業(yè)3\zuoye3_leeg.mat')fs=250;xn=eegl;%welch法海寧窗256個(gè)數(shù)據(jù)為一段重合128個(gè)數(shù)據(jù)L=18w=hanning(256),;pxx=(abs(fft(w.*xn(l:256))).A2+abs(fft(w.*xn(129:384))),A2+abs(fft(w.*xn(257:512))).A2+abs(fft(w.*xn(385:640))).A2+abs(fft(w.*xn(513:768))).A2+abs(fft(w.*xn(641:896))).A2+abs(fft(w.*xn(769:1024))).A2+abs(fft(w.*xn(897:1152))).A2+abs(fft(w.*xn(l025:1280))).A2+abs(fft(w.*xn(1153:1408))).A2+abs(fft(w.*xn(1281:1536))).A2+abs(fft(w.*xn(1409:1664))).A2+abs(fft(w.*xn(1537:1792))).A2+abs(fft(w.*xn(1665:1920))).A2+abs(fft(w.*xn(1793:2048))).A2+abs(fft(w.*xn(1921:2176))).A2+abs(fft(w.*xn(2049:2304))).A2+abs(fft(w.*xn(2177:2432))).A2)/(256*7);p=fftshift(pxx);plot((-127:128)*fs/256,10*logl0(p));title('welch法eeg功率譜')；程序運(yùn)行完成后，結(jié)果如圖(7)：圖（7）N=1024和N=256時(shí)的腦電信號(hào)周期圖法功率譜四.實(shí)驗(yàn)總結(jié)通過(guò)本次的實(shí)驗(yàn)，使我對(duì)自相關(guān)法、周期圖法和welch法求功率譜有了較為深刻的理解，同時(shí)使我對(duì)利用matlab進(jìn)行功率譜估計(jì)方法有了比較深入地認(rèn)識(shí)，對(duì)以后的學(xué)習(xí)具有非常大的幫助。實(shí)驗(yàn)日期：交報(bào)告日期：報(bào)告退發(fā)：教師審批簽字：第頁(yè)（共頁(yè)）2013年實(shí)驗(yàn)日期：交報(bào)告日期：報(bào)告退發(fā)：教師審批簽字：第頁(yè)（共頁(yè)）2013年5月11日2013年5月H日（訂正、重做）實(shí)驗(yàn)名稱:腦電佶號(hào)的AR參數(shù)模型譜佶計(jì)西安交通大學(xué)實(shí)驗(yàn)報(bào)告課程：醫(yī)學(xué)信號(hào)處理 專業(yè)班號(hào)碩2102組別姓名陳星學(xué)號(hào)3112322010同組者 一.實(shí)驗(yàn)要求1）已知32點(diǎn)觀測(cè)值x（n）=[0.4282 1.14541.55971.89941.6854 2.30752.4679 1.97901.60631.2804 -0.20830.0577 0.02060.3572 1.65720.74881.6666 1.98302.6914 1.25211.8691 1.68550.62420.1763 1.34900.6955 1.2941-0.6177]A①計(jì)算自相關(guān)序列氏-團(tuán)）1.0475 0.4319并作圖顯示。0.03120.5802②用求出的自相關(guān)序列Rx（⑼來(lái)估計(jì)3階AR模型參數(shù)a3⑴，a3（2）,a3（3）及A2方差0③求該觀測(cè)序列的3階AR模型功率譜估計(jì)并作圖顯示。④若要求Ep<=0.3783,求AR模型的階數(shù)并計(jì)算功率譜，要求作圖顯示。2）信號(hào)為作業(yè)1獲取的腦電信號(hào)，求3階AR模型功率譜估計(jì)并作圖顯示。二.實(shí)驗(yàn)原理AR模型譜估計(jì)原理：該模型除bo=l外，勺（'=1八1%）的取值均為零。其系統(tǒng)函數(shù)為“（z）=-—k=\x(n)=w(n)-2\akx(n-k)其差分方程為 *=i ,可知其功率譜密度為' A⑵A—)所以可推出其功率譜為2P(°)=噓——AR模型又稱為p階自回歸模型，簡(jiǎn)稱AR模型。其系統(tǒng)函數(shù)只有極點(diǎn)，沒(méi)有零點(diǎn),因此也稱為全極點(diǎn)模型。三.實(shí)驗(yàn)步驟與結(jié)果分析A1)計(jì)算32點(diǎn)觀測(cè)值的自相關(guān)序列，然后利用公式計(jì)算3階AR模型參數(shù)a3⑴，A A Aa3(2),a3(3)及方差繼而可以求出3階AR模型功率譜估計(jì)并作圖顯示，最后求解在限制條件下Ep<=0.3783,求出AR模型的階數(shù)并計(jì)算功率譜，作圖顯示。將AR模型計(jì)算功率譜寫(xiě)成函數(shù)puguji.m。程序：functionP=puguji(a,p/E);%a為AR模型參數(shù)，p為階數(shù),E為方差w="pi:0.01:pi;forn=l:l:length(w)sum=l;fork=l:l:psum=sum+a(k+l)*exp(-j*w(n)*k);endP(n)=E/abs(sum)A2;endplot(w/2/pi,P)clearcicxn=[0.42821.14541.55971.89941.68542.30752.46791.97901.60631.2804-0.20830.05770.02060.35721.65720.74881.66661.98302.69141.25211.86911.68550.62420.17631.34900.69551.29411.04750.43190.03120.5802-0.6177];%求自相關(guān)函數(shù)并作圖[rx,lags]=xcorr(xn,'biased');subplot(3,l,l)plot(rx)title('自相關(guān)函數(shù)')；%m=0e0=rx(32);a0=l;%m=lall=-rx(33)/e0;el=eO*(l-allA2);%m=2a22=-(rx(34)+all*rx(33))/el;a21=all*(l+a22);e2=el*(l-a22A2);%m=3a33=-(rx(35)+a21*rx(34)+a22*rx(33))/e2;a31=a21+a33*a22;a32=a22+a33*a21;e3=e2*(l-a33A2);%顯示所要求的參數(shù)a31=a31a32=a32a33=a33e3=e3%求3階功率譜并作圖A=[l,a31,a32,a33];subplot(3,l,2)puguji(A,3,e3);title('3階AR模型譜估計(jì)')；

%利用matlab的函數(shù)來(lái)求解階數(shù)nE=l;i=4;while(E>03783)[aE]=aryule(xn,i);i=i+l;endi=i-l%求i階Ar功率subplot(3,l,3)puguji(a,i,E);title('12階AR模型譜估計(jì)')程序運(yùn)行結(jié)果：Aa3(l)：a31=-0.6984Aa3(2)：a32=-0.2748Aa3(3)：a33=0.0915A方差6？：e3=0.4678在Ep〈=0.3783最小階數(shù)為:i=12圖(1)中做出了序列的自相關(guān)函數(shù)，3階和12階的AR模型功率譜:圖(1)AR模型功率譜估計(jì)2）信號(hào)為腦電信號(hào)，求3階AR模型功率譜估計(jì)并作圖顯示。程序：cleardeload（'E:\課程'醫(yī)學(xué)信號(hào)處理'作業(yè)3\zuoye3_leeg.mat'）xn=eegl;[a,Ep]=aryule（xn,3）%求3階功率譜并作圖puguji（a,3,Ep）;title（'3階AR模型譜估計(jì)，）；程序運(yùn)行結(jié)果如圖（2）：圖（1）腦電信號(hào)的3階AR模型功率譜估計(jì)四.實(shí)驗(yàn)討論經(jīng)典譜估計(jì)法可以利用FFT計(jì)算，因而有計(jì)算效率高的優(yōu)點(diǎn)，在譜分辨力要求不是太高的地方常用這種方法。但頻率分辨率地是經(jīng)典譜估計(jì)的一個(gè)無(wú)法回避的缺點(diǎn)。如周期圖法在計(jì)算中把觀測(cè)到的有限長(zhǎng)的N個(gè)數(shù)據(jù)以外的數(shù)據(jù)認(rèn)為是零,而B(niǎo)T法僅利用N個(gè)有限的觀測(cè)數(shù)據(jù)作自相關(guān)函數(shù)估計(jì)，實(shí)質(zhì)上也就是假設(shè)除已知數(shù)據(jù)外的自相關(guān)函數(shù)全為零，這些顯然都是與事實(shí)不符的。為了克服以上缺點(diǎn),人們提出了平均，加窗平滑等方法，在一定程度上改善了經(jīng)典譜估計(jì)的性能。但是，經(jīng)典譜估計(jì)，始終無(wú)法解決，頻率分辨率與譜估計(jì)穩(wěn)定性之間的矛盾，特別是在數(shù)據(jù)記錄長(zhǎng)度比較短時(shí)，這一矛盾尤其突出?，F(xiàn)代譜估計(jì)理論也就是在這種背景下產(chǎn)生的，以1967年Burg提出的最大燧譜分析法為代表的現(xiàn)代譜估計(jì)法，不認(rèn)為在觀察到的N個(gè)數(shù)據(jù)以外的數(shù)據(jù)全為零。因此克服了經(jīng)典法的這個(gè)缺點(diǎn)，提高了譜估計(jì)的分辨率。后來(lái)發(fā)現(xiàn)線性預(yù)測(cè)自回歸模型法（簡(jiǎn)稱AR模型法）與Burg的最大烯譜分析法是等價(jià)的，它們都可歸結(jié)為通過(guò)Yule-Walker方程求解自回歸模型的系數(shù)問(wèn)題。五.實(shí)驗(yàn)總結(jié)通過(guò)本次試驗(yàn)，使我對(duì)AR模型有了更加深入的認(rèn)識(shí)和了解，讓我對(duì)現(xiàn)代譜估計(jì)法的優(yōu)勢(shì)有了更深層次的認(rèn)識(shí)，并掌握了利用matlab來(lái)計(jì)算AR模型譜估計(jì)的參數(shù)的方法對(duì)今后的學(xué)習(xí)大有裨益。西安交通大學(xué)實(shí)驗(yàn)報(bào)告第頁(yè)(共頁(yè))課程：_醫(yī)學(xué)信號(hào)處理實(shí)驗(yàn)日期：2013年5月11日專業(yè)班號(hào)碩2102組別交報(bào)告日期：2013年5月11日姓名—陳星學(xué)號(hào)3112322010報(bào)告退發(fā)：(訂正、重做)同組者_(dá)教師審批簽字：實(shí)驗(yàn)名稱： 自適應(yīng)濾波春設(shè)計(jì)與仿真一.實(shí)驗(yàn)要求1)噪聲抵消器圖⑴2)按照以上噪聲抵消器框圖設(shè)計(jì)自適應(yīng)濾波器。其中s(k)為心電信號(hào)，如ecg0410,x(k)=sin(2n50nT)為50Hz干擾，T=l/250o3)使用以下迭代公式求濾波器系數(shù)h(n)及e(k)0ep)ad(n)~戈，巾He-H"+/xe(i)寫(xiě)成分量形式：hn+,(k)=h?(k)+ux(n-k)e(n) k=0,…M可取M=5其中步長(zhǎng)入max,入max為矩陣R的最大特征值，口可取0.1左右。初始值可設(shè)h°(k)=Ok=0,-Mo當(dāng)n足夠大時(shí)，有e(n)=s(n)。要求顯示觀測(cè)數(shù)據(jù)d(k)與濾波器輸出e(k)的波形?？蓞⒖糾atlab幫助―>Demos—ToolboxesfFilterDesignfAdaptiveFiIter中的有關(guān)例子。二.實(shí)驗(yàn)原理D維納濾波從連續(xù)的(或離散的)輸入數(shù)據(jù)中濾除噪聲和干擾以提取有用信息的過(guò)程稱為濾波，而相應(yīng)的裝置稱為濾波器。根據(jù)濾波器的輸出是否為輸入的線性函數(shù)，可將它分為線性濾波器和非線性濾波器兩種。濾波器研究的一個(gè)基本課題就是：如何設(shè)計(jì)和制造最佳的或最優(yōu)的濾波器。所謂最佳濾波器是指能夠根據(jù)某?最佳準(zhǔn)則進(jìn)行濾波的濾波器。20世紀(jì)40年代，維納奠定了關(guān)于最佳濾波器研究的基礎(chǔ)。即假定線性濾波器的輸入為有用信號(hào)和噪聲之和，兩者均為廣義平穩(wěn)過(guò)程且知它們的二階統(tǒng)計(jì)特性，維納根據(jù)最小均方誤差準(zhǔn)則(濾波器的輸出信號(hào)與需要信號(hào)之差的均方值最小)，求得了最佳線性濾波器的參數(shù)，這種濾波器被稱為維納濾波器。在維納研究的基礎(chǔ)上，人們還根據(jù)最大輸出信噪比準(zhǔn)則、統(tǒng)計(jì)檢測(cè)準(zhǔn)則以及其他最佳準(zhǔn)則求得的最佳線性濾波器。實(shí)際上，在一定條件下，這些最佳濾波器與維納濾波器是等價(jià)的。因而，討論線性濾波器時(shí)，一般均以維納濾波器作為參考。維納濾波理論用于解決最小均方誤差下的線性濾波問(wèn)題。設(shè)接收到(或觀測(cè)到)的信號(hào)為隨機(jī)信號(hào)x@)■項(xiàng))+4) ⑴其中s(t)是未知的實(shí)隨機(jī)信號(hào)，n(t)是噪聲。要設(shè)計(jì)的線性濾波器，其沖擊響應(yīng)為h(t,t),輸入為x(t),輸出為W),即喇?04小》 ⑵令為估計(jì)誤差。沖擊響應(yīng)h(t,C按最小均方誤差準(zhǔn)則確定，即h(t,T)必須滿足使

達(dá)到最小。根據(jù)最小均方誤差估計(jì)的正交條件，有以下關(guān)系成立TOC\o"1-5"\h\z啡(0-匚坤加網(wǎng)則｝.0 -84-48 ⑷令凡(5)見(jiàn)1)?小仲崎 ⑹則有口k款&。< -8$- ⑦上述方程通常稱為非平穩(wěn)隨機(jī)過(guò)程條件下的維納-霍甫(Wiener-Kolmogorov)積分方程。特別當(dāng)x(t),s(t)均為廣義(或?qū)?平穩(wěn)隨機(jī)信號(hào)，而濾波器是線性時(shí)不變系統(tǒng)的情況下，x(t)與s(t)必為聯(lián)合平穩(wěn)，式(7)可寫(xiě)為入0口1咪(，-印? (8)令”(■工,則有《(切■口(砒y°旬3?喇叫(動(dòng)此處，“*”號(hào)表示卷積，對(duì)上式兩邊取Fourier變換，可得TOC\o"1-5"\h\zs上WWS) (io)3霸 (11)對(duì)于因果線性系統(tǒng)，有?(/)-fh(r-rX(r)dr (⑵采用完全相同的分析方法，推得因果平穩(wěn)維納-霍甫積分方程如下(jy)-rii(4)R.(y-4)(U (13)則郡 (14)其中$<?)?￡(?)￡(?