music演算法程序

Ⅰ 急，誰有空間譜估計MUSIC演算法MATALB

clc
clear all
tic
% ---------------------------參數初始化------------------------------------
M=8; % 陣元數目
L=1024; % 信號長度/快拍數
N=3; % 信號源
lam=0.15; % 信號波長
d=1/2*lam; % 陣元間距
w=[pi/6 pi/9 pi/3]'; % 信號頻率
angle=[pi/4 pi/6 pi/9]; % 信號入射角度
% ---------------------------信號模型--------------------------------------
A1=exp(-j*2*pi/lam*([0:M-1]*d)*sin(angle(1)))'; % 導向向量
A2=exp(-j*2*pi/lam*([0:M-1]*d)*sin(angle(2)))';
A3=exp(-j*2*pi/lam*([0:M-1]*d)*sin(angle(3)))';
A=[A1 A2 A3];
S=[1.3*cos(w(1)*(1:L));1*sin(w(2)*(1:L));1*sin(w(3)*(1:L))]; % 構造輸入信號矢量
X=A*S+randn(8,L)+j*randn(8,L); % 加雜訊
% ---------------------------角度估計--------------------------------------
R=X*X'/L; % 協方差矩陣
Rcum=(kron(X,conj(X))*(kron(X,conj(X))'))/Num-kron(X,conj(X))/Num*kron(X,conj(X))'/Num-kron(X*X'/Num,conj(X*X'/Num));
[V D]=eig(Rcum); % 奇異值分解
[lambda,index] = sort((diag(D))); % 特徵值排序
UU=V(:,index(1:N)); % 雜訊子空間
theta=-pi/2:.001:pi/2; % ULA估計角度變化的范圍和頻率選擇 小步進0.001
for i = 1:length(theta) % 角度估計
AA=exp(-j*2*pi/lam*([0:M-1]*d)*sin(theta(i)))'; % 方向矢量
WW=AA'*UU*UU'*AA;
Pmusic(i)=abs(1/WW); % 角譜
end
Pmusic = 10*log10(Pmusic/max(Pmusic));
figure(1)
plot(theta*180/pi,Pmusic);title('Music-1')
figure(2);
polar(theta,Pmusic );grid;title('Music-2')
toc

這是線陣的MUSIC譜估計演算法。希望能幫助你。

Ⅱ music演算法的詳細內容

MUSIC演算法是空間譜估計測向理論的重要基石。演算法原理 如下：
（1） 不管測向天線陣列形狀如何，也不管入射來波入射角的維數如何，假定陣列由M個陣元組成，則陣列輸出模型的矩陣形式都可以表示為：Y(t)=AX(t)+N(t)
其中，Y是觀測到的陣列輸出數據復向量；X是未知的空間信號復向量；N是陣列輸出向量中的加性雜訊；A是陣列的方向矩陣；此處，A矩陣表達式由圖冊表示。
MUSIC演算法的處理任務就是設法估計出入射到陣列的空間信號的個數D以及空間信號源的強度及其來波方向。
（2） 在實際處理中，Y得到的數據是有限時間段內的有限次數的樣本（也稱快拍或快攝），在這段時間內，假定來波方向不發生變化，且雜訊為與信號不相關的白雜訊，則定義陣列輸出信號的二階矩：Ry。
（3） MUSIC演算法的核心就是對Ry進行特徵值分解，利用特徵向量構建兩個正交的子空間，即信號子空間和雜訊子空間。對Ry進行特徵分解，即是使得圖冊中的公式成立。
（4） U是非負定的厄米特矩陣，所以特徵分解得到的特徵值均為非負實數，有D個大的特徵值和M-D個小的特徵值，大特徵值對應的特徵向量組成的空間Us為信號子空間，小特徵值對應的特徵向量組成的空間Un為雜訊子空間。
（5） 將雜訊特徵向量作為列向量，組成雜訊特徵矩陣 ，並張成M-D維的雜訊子空間Un，雜訊子空間與信號子空間正交。而Us的列空間向量恰與信號子空間重合，所以Us的列向量與雜訊子空間也是正交的，由此，可以構造空間譜函數。

（6） 在空間譜域求取譜函數最大值，其譜峰對應的角度即是來波方向角的估計值。