music、esprit、mvdr算法的譜估計(jì)

課程(論文)題目：MUSIC、ESPRT、MVDR算法的譜估計(jì)內(nèi)容1算法原理MUSIC算法MUSIC算法利用信號(hào)子空間和噪聲子空間的正交性，構(gòu)造空間譜函數(shù)，通過(guò)譜峰

搜索，估計(jì)信號(hào)頻率。由APAhH=0，i=K+1,...,M 且矩陣AHA可逆得i(AhA)-1AhAPAh^=PAh#=0，i=K+1,...,M。又由于矩陣P為正定的對(duì)角矩陣，i i方程兩邊可再同時(shí)左乘P-1，推出aH佃)#=0，k=1,2,...,K，i=K+1,...,M。這就ki表明，信號(hào)頻率向量a他)與噪聲子空間的特征向量正交。信號(hào)角頻率的估計(jì)可以由掃k描函數(shù)P仙)的K個(gè)峰值位置確定。P (P (?)=MUSIC卜(?)G『=ah(?)GGa(?)5*匚兀，卻ESPRIT算法即基于旋轉(zhuǎn)不變技術(shù)的信號(hào)參數(shù)估計(jì)。連續(xù)MESPRIT算法即基于旋轉(zhuǎn)不變技術(shù)的信號(hào)參數(shù)估計(jì)。連續(xù)M個(gè)時(shí)刻的觀測(cè)值可表示為向量形式x(n)=As(n)+v(n)。定義隨機(jī)過(guò)程y(n)=x(n+1)，且向量y(n)和矩陣e分別為y(n)=[y(n)y(n一1) y(n一M+1)]T， diag^ej^1 ｝，則y(n)=Aes(n)+v(n+1)。向量x(n)的自相關(guān)矩陣為R=Exx相關(guān)矩陣為R=xy[x(n)xh(n)]APAh+◎21，向量x(n)和『(n)的互E「x(n)yh(n)]=AP?hAh+c2Z。對(duì)R進(jìn)行特征分解，找到R的v xx xx最小特征值九=九=c2(九'九'九)。定義矩陣：minM v1 2 MC=R-c21=R一尢I=APAh，C=R—c2Z=R—九Z=AP?hAhxx xx v xxmin xy xy v xyminxxxxxy可以通過(guò)求解方程式C-九C|=0來(lái)求得到矩陣(C,Ci的廣義特征值。當(dāng)

九=ej?k,k=1,…,K時(shí)，矩陣(C一九C)是奇異的；而九Hejq,k=1,…,K時(shí)，

(C-九C)是滿秩的。矩陣對(duì)｛C,C｝的廣義特征值恰為e?e陰…ej^K，這些根xxxxxyxx xy xxxy的相位即為信號(hào)的頻率估計(jì)。MVDR算法MVDR算法即最小方差無(wú)失真響應(yīng)算法，是有別于經(jīng)典功率譜估計(jì)和參數(shù)模型估計(jì)的另一類信號(hào)頻率估計(jì)方法。定義向量x(n)=[x(n)x(n-1)…x(n-M+1#，假定信號(hào)通過(guò)一個(gè)M抽頭的FIR濾波器W(z)=0w(k)Z-k，則輸出信號(hào)為：k=oHX*HX*XTw}=WHEL*XT}w—WHRwo處'x(n)無(wú)失真地通過(guò)'且piR-1a(?)MVDR aH(h)R-1； )ixx i1p— 丫 丫 aH(w)R-1a(w)ixx iy(n)=x(n)*a(n)=竣1w(k)x(n一k)=xtwy(n)的功率可以表示為p=E((n)2LE為求得濾波器的系數(shù)’需要滿足在對(duì)給定的某一頻率最小。此時(shí)’pjspjs)并不是真正意義上的功率譜'但則P —廠= 廠rMVDR aH(①)R-1a(①)xx它描述了信號(hào)真正譜的相對(duì)強(qiáng)度’可以由此估計(jì)信號(hào)頻率。2算法實(shí)現(xiàn)采用空間譜估計(jì)的典型代表MUSIC算法，ESPRIT算法和MVDR算法，對(duì)含有高斯白噪聲的復(fù)正弦信號(hào)進(jìn)行頻率估計(jì)。選取的信號(hào)數(shù)p=3，陣元數(shù)M=10，采樣數(shù)N=1000。待檢測(cè)信號(hào)的歸一化頻率為f=0.10,f=0.20,f=0.45，仿真的待檢測(cè)信號(hào)為：e-j2kf2e-j2kf3e-j2kf1x(M-1)e-j2kf2x(M-1)e-j2kf3x(M-1)采樣數(shù)為N,且彼此之間相互獨(dú)假設(shè)a,a,a采樣數(shù)為N,且彼此之間相互獨(dú)12 3立，則陣列響應(yīng)矩陣為e-j2ka1e-j2ka2e-j2ka3假定噪聲v,ie[1,M]為零均值，方差為1的高斯白噪聲，采樣數(shù)為N,則待檢測(cè)i信號(hào)引入的噪聲為v1v2于是仿真信號(hào)為：X=S*A+vo2.1MUSIC算法利用MUSIC算法進(jìn)行頻率估計(jì)時(shí)，首先求出仿真信號(hào)的自相關(guān)矩陣R=X*X'/N,然后對(duì)自相關(guān)矩陣進(jìn)行特征值分解得到特征值空間U和特征向量空間D，然后對(duì)特征值空間U進(jìn)行升序排序，取較小的m-P個(gè)特征值所對(duì)應(yīng)的特征向量組成噪聲子空間G,然后求得PMUSIW3E[°°5ImUSIC算法程序段如下所示：R=X*X'/N;[UD]=eig(R);[d,index]=sort(diag(D));U=U(:,index);fori=-N:Na=exp(-j*2*pi*[0:M-l]'*(i*0.5/N));Pmusic(i+l+N)=abs(l/(a'*U(:,l:end-p)*U(:,l:end-p)'*a));endplot(omg,10*log10(Pmusic/max(Pmusic)));%畫出MUSIC算法圖形2.2ESPRIT算法ESPRIT算法要先構(gòu)造相關(guān)矩陣R和R，然后對(duì)R進(jìn)行特征值分解得到最小特xx xy xx征值九.即為噪聲的方差b2,通過(guò)對(duì)矩陣對(duì)1C,Ci進(jìn)行廣義特征值分解，最接近min v xxxy單位圓的K個(gè)特征值相位即為信號(hào)的頻率估計(jì)。ESPRIT算法程序段如下所示：S1=U(1:end-1,8:10);S2=U(2:end,8:10);S12=[S1S2];[UuDd]=eig((S12'*S12));[dd,ind]=sort(diag(Dd));Uu=Uu(:,ind);dt=-Uu(1:3,1:3)*inv(Uu(4:6,1:3));dd=eig(dt);t_esprit=asin(-angle(dd)/(2*pi));plot(t_esprit,[000],'*','color','red')%畫出ESPRIT算法圖形，紅色表示MVDR算法MVDR算法要先求自相關(guān)矩陣R，然后將峰搜索矩陣a?)帶入最小方差譜估計(jì)xx i公式Po=MVDR1Po=MVDR1ah(①)R-1a(3)[0,0.5]通過(guò)觀察譜峰得到信號(hào)的頻率估計(jì)值。MVDR算法程序段如下所示:fori=-N:Na=exp(-j*2*pi*[0:M-1]'*(i*0.5/N));Pmvdr(i+N+1)=abs(1/(a'*inv(R)*a));endplot(omg,10*log10(Pmvdr/max(Pmvdr)),'--','color','black')2.4程序流程圖各算法程序流程圖如下所示:圖2.4.1MUSIC算法圖2.4.2ESPRIT算法圖2.4.3MVDR算法3仿真結(jié)果使用MATLAB軟件進(jìn)行仿真，得出MUSIC算法、ESPRIT算法、MVDR算法的仿真圖形分別如圖3.1、圖3.2、圖3.3所示。為了將這三種算法的譜估計(jì)效果更加直觀地比較，將三種算法的仿真結(jié)果圖形放在一張圖中進(jìn)行比較，如圖3.4所示。另外，MUSIC算法得出的頻率估計(jì)為f二0.10,f二0.20,f二0.45，運(yùn)算時(shí)間time=1.415；ESPRIT算TOC\o"1-5"\h\z1 2 3法的頻率估計(jì)為f二0.10,f二0.20,f二0.47，運(yùn)算時(shí)間time=0.0251；MVDR算法的頻1 2 3率估計(jì)為f二0.10,f二0.20,f二0.45，運(yùn)算時(shí)間time=0.0910。1 2 3QF：^ur.1FkEe￡|QF：^ur.1/r