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;

————————————————