ESPRIT-LIKE algorithm (incomplete row information utilization)+源码

 1.数学模型

2. ESPRIT-LIKE源码　　

%% readMe
% DOA estimation of ESPRIT-LIKE

% ULA + Toeplitze matrix + ESPRIT 参照以下文献
% [R2] Fang-Ming Han and Xian-Da Zhang, "An ESPRIT-like algorithm for coherent DOA estimation," 
%                                        in IEEE Antennas and Wireless Propagation Letters, vol. 4, pp. 443-446, 2005.
%%
clear all;
close all;
clc;
derad = pi/180;
radeg = 180/pi;
twpi = 2 * pi;
N = 11;  % ULA阵元数
M = floor(N / 2);

lamda = 1;
dd = lamda / 2;  %阵元间距
d = -M * dd : dd : M * dd;

theta = [-5 3 25];
K = length(theta);  % 信元数
A = exp(-1j * d.' * twpi * sin(theta * pi / 180) );    % 阵列流型
%% 信号产生
snr = 0;
Nsub = 500;   % 快拍数
signalFlag = 1;
if signalFlag == 0  % independent signal
    s0 = (randn(K, Nsub) + 1i * randn(K, Nsub)) / sqrt(2);
elseif signalFlag == 1  % coherent signal
    s1 = (randn(1, Nsub) + 1i * randn(1, Nsub)) / sqrt(2);
%     s2 = s1;
    s3 = (randn(1, Nsub) + 1i * randn(1, Nsub)) / sqrt(2);
    s0 = [repmat(s1, [2, 1]);s3];
else                    % corelated signal
    s1 = (randn(2, Nsub) + 1i * randn(2, Nsub)) / sqrt(2);
    s0 = [s1 sum(s1, 1)/2];
end
noiseFlag = 0;
if noiseFlag == 0
    y0 = A * s0;
    y = awgn(y0, snr); 
elseif noiseFlag == 1
    % 非均匀噪声+信号产生       
end
%% Construst the Toeplitz matrix
Rx = y * y'/ Nsub;
% Rx = y0 * y0'/ Nsub;
m = 0;
R = toeplitz(fliplr(Rx(m + M + 1,1 : M + 1)), Rx(m + M + 1,M + 1 : end));
doaHat = sort(tls_esprit(dd, R, K))

doaAve = [];
for m = -M : M
    R1 = toeplitz(fliplr(Rx(m + M + 1,1 : M + 1)), Rx(m + M + 1,M + 1 : end));
    doaAve = [doaAve;sort(tls_esprit(dd, R1, K))];
end
doaAve = mean(doaAve)

function estimate =  tls_esprit(dd,cr,Le)
%*******************************************************
% This function calculates TLS-ESPRIT estimator for
% uniform linear array.
%
% Inputs
%    dd       sensor separation in wavelength
%    cr(K,K)  array output covariance matrix
%    Le       estimated number of sources ==L=iwave
%   
% Output
%    estimate  estimated angles in degrees
%              estimated powers
%*******************************************************

twpi = 2.0*pi;
derad = pi / 180.0;
radeg = 180.0 / pi;

% eigen decomposition of cr
[K,KK] = size(cr);
[V,D]=eig(cr);
EVA = real(diag(D)');
[EVA,I] = sort(EVA);
%disp('Eigenvalues of correlation matrix in TLS-ESPRIT:')
EVA=fliplr(EVA);
EV=fliplr(V(:,I));

%  composition of E_{xy} and E_{xy}^H E_{xy} = E_xys
Exy = [EV(1:K-1,1:Le) EV(2:K,1:Le)];
E_xys = Exy'*Exy;

% eigen decomposition of E_xys
[V,D]=eig(E_xys);
EVA_xys = real(diag(D)');
[EVA_xys,I] = sort(EVA_xys);
EVA_xys=fliplr(EVA_xys);
EV_xys=fliplr(V(:,I));

% decomposition of eigenvectors
Gx = EV_xys(1:Le,Le+1:Le*2);
Gy = EV_xys(Le+1:Le*2,Le+1:Le*2);

% calculation of  Psi = - Gx [Gy]^{-1}
Psi = - Gx/Gy;

% eigen decomposition of Psi
[V,D]=eig(Psi);
EGS = diag(D).';
[EGS,I] = sort(EGS);
EGS=fliplr(EGS);
EVS=fliplr(V(:,I));
 
% DOA estimates
ephi = atan2(imag(EGS), real(EGS));
ange = - asin( ephi / twpi / dd ) * radeg;
estimate(1,:)=ange;
end
% End tls_esprit.m

单次运行的结果如下：

doaHat =

   -4.7073    3.1227   25.2112

doaAve =

  -12.2968    1.8439   32.4782

注意：如果只使用协方差矩阵某一行构造Tｏｅｐｌｉｔｚ矩阵求解的结果效果比利用协方差矩阵的所有行求平均的效果要好，不知道为什么？

有小伙伴如果知道原因，望告知哈🤭