modifiedmusic.m

前后向空间平滑ＭＵＳＩＣ算法程序

%改进music算法10阵元等距阵波达方向估计（d=λ/2）,2信源，前向-后向 空间平滑 5个子阵（6阵元） 最小范数
clear
w1=pi/4;
w2=pi/3;
kp=2048;
a1=[exp(-j*pi*sin(w1)*[0:9])].';
a2=[exp(-j*pi*sin(w2)*[0:9])].';
A=[a1,a2];
k=1:kp;
s=[cos(k*0.015);cos(k*0.015)];
n=randn(10,kp)+j*randn(10,kp);
x=A*s+n;
r=x*(x)'/kp;
J=rot90(eye(6));
for k=1:5
    b=sprintf('r%d',k);
    c=sprintf('rb%d',k);
    rk=sprintf('rk%d',k);
    F=[zeros(6,k-1) eye(6) zeros(6,5-k)];
    Q=[zeros(6,k-1) J zeros(6,5-k)];
    f=F';
    q=Q';
    eval([b,'=','F','*r*','f;']);
    eval([c,'=','Q','*r*','q;']);
    eval([rk,'=(',b,'+',c,')/2;']);
end
rx=(rk1+rk2+rk3+rk4+rk5)/5;
[V,D]=eig(rx);
[bn,index]=sort(diag(D));
Vn=V(:,index);
Un=Vn(:,[1:4]);
II=Un*(Un)';
w=[-pi/2:pi/180:pi/2];
for l=1:length(w)
    ao=[exp(-j*pi*sin(w(l))*[0:5])].';
    pm(l)=(ao)'*ao*(inv((ao)'*II*ao));
     p(l)=10*log10(abs(pm(l)));
end
t=w*180/pi;
plot(t,p);grid;
title('改进music算法');
ylabel('P(\theta) (dB)');
xlabel('\theta (deg)');
[Wn Pn]=ginput;
hold on
plot(Wn,Pn,'or');
for k=1:length(Wn)
    text(Wn(k),Pn(k),['  \theta=',num2str(Wn(k)),'P=',num2str(Pn(k))]);
end
hold off