实值波束空间ESPRIT算法

%RB_ESPRITALGORITHM%DOAESTIMATIONBYR_ESPRITclearall;closeall;clc;source_number=2;%ÐÅÔªÊýsensor_number=8;%ÕóÔªÊýN_x=1024;%Ðźų¤¶Èsnapshot_number=N_x;%¿ìÅÄÊýw=[pi/4pi/6].';%ÐźÅƵÂÊl=((2*pi*3e8)/w(1)+(2*pi*3e8)/w(2))/2;%ÐźŲ¨³¤d=0.5*l;%ÕóÔª¼ä¾àsnr=0;%ÐÅÔë±Èsource_doa=[4560];%Á½¸öÐźŵÄÈëÉä½Ç¶ÈA=[exp(-j*(0:sensor_number-1)*d*2*pi*sin(source_doa(1)*pi/180)/l);exp(-j*(0:sensor_number-1)*d*2*pi*sin(source_doa(2)*pi/180)/l)].';%ÕóÁÐÁ÷ÐÍs=sqrt(10.^(snr/10))*exp(j*w*[0:N_x-1]);%·ÂÕæÐźÅ%x=awgn(s,snr);x=A*s+(1/sqrt(2))*(randn(sensor_number,N_x)+j*randn(sensor_number,N_x));%¼ÓÁ˸ß˹°×ÔëÉùºóµÄÕóÁнÓÊÕÐźÅ%¹¹ÔìȨ¾ØÕóW%ac=zeros(sensor_number,sensor_number);%fork=0:(sensor_number-1)%gama(k+1)=k*2*pi/sensor_number;%ac(:,k+1)=exp(-j*((sensor_number-1)/2)*(gama(k+1)/pi))*exp(j*(0:sensor_number-1)'*(gama(k+1)/pi));%end%W=ac;W=[exp(-j*(sensor_number-1)/2*2*0*pi/sensor_number)*exp(j*(0:sensor_number-1)'*2*0*pi/sensor_number)exp(-j*(sensor_number-1)/2*2*1*pi/sensor_number)*exp(j*(0:sensor_number-1)'*2*1*pi/sensor_number)exp(-j*(sensor_number-1)/2*2*2*pi/sensor_number)*exp(j*(0:sensor_number-1)'*2*2*pi/sensor_number)exp(-j*(sensor_number-1)/2*2*3*pi/sensor_number)*exp(j*(0:sensor_number-1)'*2*3*pi/sensor_number)exp(-j*(sensor_number-1)/2*2*4*pi/sensor_number)*exp(j*(0:sensor_number-1)'*2*4*pi/sensor_number)exp(-j*(sensor_number-1)/2*2*5*pi/sensor_number)*exp(j*(0:sensor_number-1)'*2*5*pi/sensor_number)exp(-j*(sensor_number-1)/2*2*6*pi/sensor_number)*exp(j*(0:sensor_number-1)'*2*6*pi/sensor_number)exp(-j*(sensor_number-1)/2*2*7*pi/sensor_number)*exp(j*(0:sensor_number-1)'*2*7*pi/sensor_number)];disp(W);%½«ÕóÁнÓÊÕµÄÊý¾Ý´Ó¸´Êýת»»³ÉʵÊýY=W'*x;Y1=[real(Y)imag(Y)];disp(Y1);R1=Y1*Y1'/(2*snapshot_number);[U,S,V]=svd(R1);disp(U);Es=U(:,1:source_number);disp(Es);KS1=[1cos(pi/sensor_number)cos(2*pi/sensor_number)cos(3*pi/sensor_number)cos(4*pi/sensor_number)cos(5*pi/sensor_number)cos(6*pi/sensor_number)cos(7*pi/sensor_number)];KS1=diag(KS1);fori=1:sensor_number-1KS1(i,i+1)=cos(i*pi/sensor_number);endKS1(sensor_number,1)=(-1).^sensor_number;KS2=[0sin(pi/sensor_number)sin(2*pi/sensor_number)sin(3*pi/sensor_number)sin(4*pi/sensor_number)sin(5*pi/sensor_number)sin(6*pi/sensor_number)sin(7*pi/sensor_number)];KS2=diag(KS2);forii=1:sensor_number-1KS2(ii,ii+1)=sin(ii*pi/sensor_number);endM=pinv(KS1*Es)*(KS2*Es);%¶ÔµÃµ½µÄÐýת²»±ä¾ØÕó½øÐÐÌØÕ÷·Ö½â[Vm,Dm]=eig(M);disp(Dm);beta(1)=-2*atan(Dm(1,1));beta(2)=-2*atan(Dm(2,2));disp(beta);%·´½âµÃµ½Ðźŵ½´ï½Çdoa(1)=asin(beta(1)/pi)*180/pi;doa(2)=asin(beta(2)/pi)*180/pi;disp(doa);