多种传统DOA估计算法（MUSIC，Capon，ESPRIT等以及相互之间的对比）

%TLS_ESPRIT ALGORITHM %DOA ESTIMATION BY TLS_ESPRIT ALOGRITHM clear all; %close all; clc; source_number=2;%信元数 sensor_number=8;%原阵元数 sub_sensor_number=7;%子阵元数 N_x=1024; %信号长度 snapshot_number=N_x;%快拍数 w=[pi/4 pi/6].';%信号频率 l=sum(2*pi*3e8./w)/length(source_number);%信号波长 d=0.5*l;%阵元间距 snr=0;%信噪比/dB source_doa=[-10 60];%两个信号的入射角度 A=[exp(-1j*(0:sensor_number-1)*d*2*pi*sin(source_doa(1)*pi/180)/l);exp(-1j*(0:sensor_number-1)*d*2*pi*sin(source_doa(2)*pi/180)/l)].';%阵列流型 for n=1:50 s=sqrt(10.^(snr/10))*exp(1j*w*[0:N_x-1]);%仿真信号 x=A*s+(1/sqrt(2))*(randn(sensor_number,N_x)+1j*randn(sensor_number,N_x));%加了高斯白噪声后的阵列接收信号 x1=x(1:sub_sensor_number,:);%子阵1接受的数据矢量 x2=x(2:(sub_sensor_number+1),:);%子阵2接受的数据矢量 %对两个子阵的模型进行合并 X=[x1;x2]; R=X*X'/snapshot_number; %对R进行奇异值分解 [U,S,V]=svd(R); R=R-S(2*sub_sensor_number,2*sub_sensor_number)*eye(2*sub_sensor_number); [U,S,V]=svd(R); Us=U(:,1:source_number); disp(Us); Us1=Us(1:sub_sensor_number,:); Us2=Us((sub_sensor_number+1):2*sub_sensor_number,:); %形成矩阵Us12 Us12=[Us1,Us2]; %对“Us12'*Us12”进行特征分解，得到矩阵E [F,Sa,Va]=svd(Us12'*Us12); disp('F'); disp(F); disp(Sa); %将E分解为四个小矩阵 F11=F(1:2,1:2); F12=F(1:2,3:4); F21=F(3:4,1:2); F22=F(3:4,3:4); %按照公式得到旋转不变矩阵M E=-(F12*(inv(F22))); disp('E'); disp(E); %对得到的旋转不变矩阵进行特征分解 [V,D]=eig(E); disp(D); D=(diag(D)).'; doa=-asin(angle(D)/pi)*180/pi; doa=sort(doa); plot(n,doa(1),'*',n,doa(2),'square'); axis([0 51 -100 100]); grid on; xlabel('Number of Experiments'); ylabel('DOA Estimation/degree'); title('DOA Estimation by TLS ESPRIT Algorithm with Monte Carlo'); hold on; drawnow; end