矩阵束的ESPRIT算法

 

clear all;
theta=2*pi*[0.10,0.20,0.30];
theta1=theta/(2*pi);
polar(theta,[1 1 1],'ko')
hold on;
fais=[];

K=100;
N=10;
M=3;
SNR=10;
Ps=1;
Pn=Ps/(10^(SNR/10));
A=zeros(N,M);

c=3*(10^8);
for m=1:M
    A(:,m)=exp(-j*[0:N-1]*theta(m));
end;
for i=1:40
Rxx=zeros(N-1,N-1);
Rxy=Rxx;

for k=1:K
    Sk=Ps*(randn(M,1)+j*randn(M,1));
    Nk=Pn*(randn(N,1)+j*randn(N,1));
    Xk=A*Sk+Nk;  
    X=Xk(1:N-1);
    Y=Xk(2:N);   
    Rxx=Rxx+X*X'/K; 
    Rxy=Rxy+X*Y'/K;
end;

d1=eig(Rxx);
NoisePow1=sum(d1(1:N-1-M))/(N-1-M);
Cxx=Rxx-NoisePow1*eye(N-1);                 
d2=sort(eig(Rxy));
NoisePow2=sum(d2(1:N-1-M))/(N-1-M);
Z=zeros(N-1);
Z(2:N-1,1:N-2)=eye(N-2);
Cxy=Rxy-NoisePow2*Z;                     %去噪声

temp1=sort(eig(Cxx,Cxy));
temp2=temp1(N-M:N-1);
fai=-angle(temp2);
fai1=fai/(2*pi);
fais=[fais fai1'];
polar(fai',[1 1 1],'rx');
hold on;
end;
figure;
temp4=-80:0.1:80;
hist(fais*180/pi,temp4);
xlabel('方位角/deg');
ylabel('出现次数统计');
title('独立信源不使用空间平滑算法');
grid on;
save ESRIT;

D147