esprit2.m

来自「阵列信号处理_空间谱估计原理的一些基本算法——MusicEspritMp」· M 代码 · 共 119 行

%  the  data of signal
clear;
format long;
c=3*10.^8;
L=6;
m=L-1;
sam=256;
N=256;
w1=1;
w2=1.2;
phase=1;
%lamta=c/w1;
%len=lamta/2;
len=1/2;
for nn=1:11;
    snr=3*(nn-1);
    Amp=sqrt(2*10^(snr/10));


P1=30*pi/180;
P2=50*pi/180;

for t=1:sam,
    s1(t)=Amp*(exp(j*2*pi*w1*0.1*t));
    s2(t)=Amp*(exp(j*2*pi*w2*0.1*t));
end
 
s=[s1(1:N);s2(1:N)];

  i=1:L;
 
  x1=exp(j*2*pi*(len*(i-1)*sin(P1)));
  x2=exp(j*2*pi*(len*(i-1)*sin(P2)));  

  a11=x1(:,1:5);
  a12=x2(:,1:5);
  a21=x1(:,2:6);
  a22=x2(:,2:6);
  
  aa1=[a11.',a12.'];
  aa2=[a21.',a22.'];

  a=[aa1;aa2];
 for n=1:100;
 noise=rand(2*m,N)+j*rand(2*m,N);
 
 z=a*s+noise;
 
 Rz=(z*z')/(N-1);
 
 [e,v]=eig(Rz);
 
 es=e(:,1:2);
 en=e(:,(3:(2*m)));

 es1=es(1:m,:);
 es2=es((m+1):(2*m),:);
 eezz=(es1'*es2)/(es1'*es1);
 [ee,vv]=eig(eezz);
 
 ang=angle(vv);
 ang1=ang(1,1);
 ang2=ang(2,2);
 est1=asin((ang1)/(2*pi*len));
 est2=asin((ang2)/(2*pi*len));

 agl1=est1*180/pi;
 agl2=est2*180/pi;

 aa(n)=agl1;
 bb(n)=agl2;
 
end
for n=1:100;
     if aa(n)>bb(n)
     mp=aa(n);
     aa(n)=bb(n);
     bb(n)=mp;
     end 
end
  
aav=mean(aa);
bbv=mean(bb);

aas=std(aa);
bbs=std(bb);

aavv(nn)=aav;
bbvv(nn)=bbv;
aass(nn)=aas;
bbss(nn)=bbs;
end
i=1:11;
snr=(i-1)*3;
plot(snr,aass,'r-^',snr,bbss,'b-^');
%text(21,0.2,'70度');
%text(21,0.05,'30度');
xlabel('信噪比SNR');
ylabel('估计角度的标准差');