【雷达学习】双基地MIMO-ESPRIT方法解互耦回忆

主要学习来源于这篇文章：《Joint DOD and DOA estimation of bistatic MIMO radar in the presence of unknown mutual coupling》

大致思路是：

1. 先回忆了下 ULA 的 ESPRIT 方法；
2. 然后是 MIMO 下的 ESPRIT 方法
3. 然后是看文章看不懂，然后 AI 帮我解析，先读出大致的框架，然后丰富细节，最后对自己的理解进行确认，然后画出代码仿真的关键点

下面在这个网盘下面提供 2 个 PDF 和一个代码方便之后回顾学习
https://a.siyouyun.ren:30597/OutSaveFile2/RadarXG/ESPRIT_Double_mimo

matlab 代码预览

clc,clear all,close all;
%% Parameters of the radar system
M = 7;                                                 % number of Tx
N = 7;                                                 % number of Rx
d = 0.5;                                                 % inter-element space
% theta = [15 -10 5];                                % DOA
% phi = [20 35 60];    
theta = [-50 20 60];   % DOD                             % DOA
phi = [20 40 60];   
K = length(theta);                                  % number of target                      % dopplor frequency shift
L = 100;                                                % number of sanpshot
SNR = 30;                                            % signal-to-noise ratio before matched filtering
Geo_Rx = [0:N-1];                                  % geometry of Rx
Geo_Tx = [0:M-1];                                  % geometry of Tx
At = exp(-j*2*pi*d*Geo_Tx.'*sind(phi));   % transmitting direction matrix
Ar = exp(-j*2*pi*d*Geo_Rx.'*sind(theta));  % receiving direction matrix
item = 10;     
%% 设置互耦系数;
mc_K = 2;
mc_nb = mc_K+1;
C_t = generateMCMmat(M,mc_nb);
C_r = generateMCMmat(N,mc_nb);
C = kron(C_t,C_r);
%% 设置选择Pt和Pt矩阵
Pt = [zeros(M-2*mc_K,mc_K),eye(M-2*mc_K),zeros(M-2*mc_K,mc_K)];
Pr = [zeros(N-2*mc_K,mc_K),eye(N-2*mc_K),zeros(N-2*mc_K,mc_K)];
M1 = M-2*mc_K;
N1 = N-2*mc_K;
P = kron(Pr,Pt);
% number of trials
%% Selective matrices
JM1 = [eye(M1-1),zeros(M1-1,1)];
JM2 = [zeros(M1-1,1),eye(M1-1)];
JN1 = [eye(N1-1),zeros(N1-1,1)];
JN2 = [zeros(N1-1,1),eye(N1-1)];
Jt1 = kron(JM1,eye(N1));
Jt2 = kron(JM2,eye(N1));
Jr1 = kron(eye(M1),JN1);
Jr2 = kron(eye(M1),JN2);

a_kronrao = zeros(M*N,K);
for i=1:K
    a_kronrao(:,i) = kron(At(:,i),Ar(:,i));
end

for item_num = 1:item
    disp(['SNR = ',num2str(SNR),' dB, ',num2str(item_num), ' # try : ']);
    %% Matched Filtering
    S = randn(K,L)+1i*randn(K,L);
    X_0 = C*a_kronrao * S;
    X = awgn(X_0,SNR,"measured","dB");
    Y = P * X;
    %% Eigen decomposition
    %Rx = (X*X')/L;                                          % Estimated covariance matrix
    Rx = (Y*Y')/L;
    [Es,D] = eigs(Rx,K,'LM');                          % Signal subspace
    %% Rotational invariant property
    Vt =  pinv(Jt2*Es)*Jt1*Es;
    %Vt =  pinv(Jt1*Es)*Jt2*Es;
    [T,Phit] = eig(Vt);
    
    Phir = inv(T)*pinv(Jr2*Es)*Jr1*Es*T;
    %Phir = inv(T)*pinv(Jr1*Es)*Jr2*Es*T;
    E_theta  = asind(angle(diag(Phir))/pi).';
    E_phi = asind(angle(diag(Phit))/pi).';
    plot(E_theta,E_phi,'k*');hold on;
end
H(1)=plot(E_theta,E_phi,'k*');hold on;
H(2)=plot(theta,phi,'rx','markersize',28);grid on;
xlabel('DOA'),ylabel('DOD');
legend([H(1),H(2)],'Esimated','Ture')

闲聊

没想到 matlab 还可以画出这样的图片，“赛博生命”，hhh
https://mp.weixin.qq.com/s/GS-jSHnasVxrrIn5C2E3WQ

figure('Position',[300,50,900,900], 'Color','k');
axes(gcf, 'NextPlot','add', 'Position',[0,0,1,1], 'Color','k');
axis([0, 400, 0, 400])
SHdl = scatter([], [], 2, 'filled','o','w', 'MarkerEdgeColor','none', 'MarkerFaceAlpha',.4);
t = 0;
i = 0:2e4;
x = mod(i, 100);
y = floor(i./100);
k = x./4 - 12.5;
e = y./9 + 5;
o = vecnorm([k; e])./9;
while true
    t = t + pi/90;
    q = x + 99 + tan(1./k) + o.*k.*(cos(e.*9)./4 + cos(y./2)).*sin(o.*4 - t);
    c = o.*e./30 - t./8;
    SHdl.XData = (q.*0.7.*sin(c)) + 9.*cos(y./19 + t) + 200; 
    SHdl.YData = 200 + (q./2.*cos(c));
    drawnow
end