Figure 58. Beamspace post-Doppler in a clutter-only environment. K = 4. Ks = Kt = 2.

Contents

Coded by Ilias Konsoulas, 16 Sept. 2018. Code provided for educational purposes only. All rights reserved.

clc; clear; close all;

Radar System Operational Parameters

radar_oper_params;

Thermal Noise Power Computation

thermal_noise_power;

Thermal Noise Covariance Matrix

Rn = sigma2*eye(M*N);

Clutter Patch RCS Computation

clutter_patch_rcs;

Calculate the Array Transmit and Element Receive Power Gains

Tx_Rx_power_gains;

Calculate the Clutter to Noise Ratio (CNR) for each azimuth angle

ksi = Pt*Gtgain.*Grec*lambda^2*sigma/((4*pi)^3*Pn*10^(Ls/10)*Rcik^4);   % Eq. (58)

Clutter Covariance Matrix Computation

beta = 1;   % beta parameter.
phia = 0;   % Velocity Misalignment Angle.

Rc = clutt_cov(ksi,beta);

Analytic Interference Covariance Matrix for Fully Adaptive STAP

Ru = Rc + Rn;
InvRu = inv(Ru);

Target Space-Time Steering Vector:

phit = 0; thetat = 0;                    % Target azimuth and elevation angles in degrees.
fdt = 100;                               % Target Doppler Frequency.
fspt = d/lambda*cos(thetat*pi/180)*sin(phit*pi/180);
omegat = fdt/fr;
at = exp(-1i*2*pi*fspt*(0:N-1)).';       % Target Spatial Steering Vector.
bt = exp(-1i*2*pi*fdt*(0:M-1)).';        % Target Doppler Steering Vector.
ta = chebwin(N,30);
tb = chebwin(M,30);
fd = 0:.5:300;   Lfd = length(fd);
omegad = fd/fr;
SNRo = M*N;

LSINR Computation for Optimum Fully Adaptive Case:

LSINRopt = zeros(1,Lfd);
for n=1:Lfd
    bt = exp(-1i*2*pi*omegad(n)*(0:M-1)).'; % Target Doppler Steering Vector
    vt = kron(bt,at);
    w = InvRu*vt; %#ok<MINV>
    LSINRopt(n) = real(w'*vt)/SNRo;
end

Displaced Beam pre-Doppler Calculations:

Kt = 2;                  % sub-CPI length (fixed).
Ks = 2;                  % sub-apperture length or number of beams used.
M1 = M - Kt +1;          % Number of possible sub-CPI's.
N1 = N - Ks +1;          % Number of possible sub-appertures.

Create Doppler Filterbank Matrix F for Displaced Filter Beamspace post-Doppler STAP

dopplerfilterbank = linspace(0,300,M+1);
omegadopplerbank = dopplerfilterbank/fr;
F0 = zeros(M1,M);
for m=1:M
    F0(:,m) = 1/sqrt(M)*exp(-1i*2*pi*omegadopplerbank(m)*(0:M1-1));
end

F40 = diag(chebwin(M1,40))*F0;

Beamformer G Matrix Construction for Displaced Filter Beamspace post-Doppler STAP

g0  = exp(-1i*2*pi*fspt*(0:N1-1)).';
g30 = chebwin(N1,30).*g0;           % 30-dB spatial taper.
G0  = toeplitz([g0;  zeros(Ks-1,1);],  [g0(1) zeros(1,Ks-1)]);  % N1 x Ks Beamformer matrix G.
G30 = toeplitz([g30; zeros(Ks-1,1);], [g30(1) zeros(1,Ks-1)]);  % N1 x Ks Beamformer matrix G.

SINR0df_mat  = zeros(M,Lfd);
SINR30df_mat = zeros(M,Lfd);
SINR0af_mat  = zeros(M,Lfd);
SINR40af_mat = zeros(M,Lfd);

Create Doppler Filterbank Matrix F for Adjacent Filter Beamspace post-Doppler STAP

U2 = zeros(M,M);
Pt = floor(Kt/2);
if  mod(Kt,2)    % If Kt is odd
    for m=1:M
        U2(:,m) = 1/sqrt(M)*exp(-1i*2*pi*omegadopplerbank(m)*(0:M-1));
    end
else            % if Kt is even:
    outomegadoppler = zeros(1,M);
    for m=1:M
        outomegadoppler(m) = (omegadopplerbank(m) + omegadopplerbank(m+1))/2;
        U2(:,m) = 1/sqrt(M)*exp(-1i*2*pi*outomegadoppler(m)*(0:M-1));
    end
end

td0af   = ones(M,1);                      % Uniform Doppler Taper.
td40af = chebwin(M,40);                   % 40-dB Chebyshev Doppler Taper.

% Create Doppler Filter Bank in Fab matrix for Adjacent Filter Beamspace post-Doppler method:
Faf0   = diag(td0af)*U2;
Faf40 = diag(td40af)*U2;

Fmaf0   = zeros(M,Kt,M);
Fmaf40 = zeros(M,Kt,M);

for m=1:M
    if mod(Kt,2)                  % if K is odd and >1.
        if (m-Pt>0) && (m+Pt<=M)
            Fmaf0(:,:,m)  = Faf0(:,m-Pt:m+Pt);                       % Eq. 231.
            Fmaf40(:,:,m) = Faf40(:,m-Pt:m+Pt);
        elseif (m-Pt<=0) && (m+Pt<=M)
            Fmaf0(:,:,m)  = [Faf0(:,M+(m-Pt):M)   Faf0(:,1:m+Pt)];   % Eq. 231.
            Fmaf40(:,:,m) = [Faf40(:,M+(m-Pt):M) Faf40(:,1:m+Pt)];
        elseif m+Pt>M
            Fmaf0(:,:,m)  = [Faf0(:,m-Pt:M)   Faf0(:,1:m+Pt-M)];     % Eq. 231.
            Fmaf40(:,:,m) = [Faf40(:,m-Pt:M) Faf40(:,1:m+Pt-M)];
        end

    else      % if K is even.

        if (m-Pt>0) && (m+Pt<=M+1)
            Fmaf0(:,:,m)  = Faf0(:,m-Pt:m+Pt-1);                     % Eq. 231.
            Fmaf40(:,:,m) = Faf40(:,m-Pt:m+Pt-1);
        elseif (m-Pt<=0) && (m+Pt<=M)
            Fmaf0(:,:,m)  = [Faf0(:,M+(m-Pt):M)  Faf0(:,1:m+Pt-1)];  % Eq. 231.
            Fmaf40(:,:,m) = [Faf40(:,M+(m-Pt):M) Faf40(:,1:m+Pt-1)];
        elseif m+Pt>M+1
            Fmaf0(:,:,m)  = [Faf0(:,m-Pt:M)     Faf0(:,1:m-M+Pt-1)]; % Eq. 231.
            Fmaf40(:,:,m) = [Faf40(:,m-Pt:M)   Faf40(:,1:m-M+Pt-1)];
        end
    end

end

Create Beamformer Matrix G for Adjacent Filter Beamspace post-Doppler STAP

selbeamdist = 6.35;   % Selected Beam Angular Distance in degrees.
beamangles1 = [-9:-1  0:8]*selbeamdist; thetabeam = theta;

beamangleinc = beamangles1(2) - beamangles1(1);
beamangles2 = beamangles1 + beamangleinc/2;

beamfreqs1 = d/lambda*cos(thetabeam*pi/180)*sin(beamangles1*pi/180);
beamfreqs2 = d/lambda*cos(thetabeam*pi/180)*sin(beamangles2*pi/180);

Godd   = zeros(N,N);
Geven = zeros(N,N);

for n=1:N
    Godd(:,n)  = 1/sqrt(N)*exp(-1i*2*pi*beamfreqs1(n)*(0:N-1));
    Geven(:,n) = 1/sqrt(N)*exp(-1i*2*pi*beamfreqs2(n)*(0:N-1));
end

Create Beam Selector Matrix Jsm for Adjacent Filter Beamspace post-Doppler STAP

Ps = floor(Ks/2);
if mod(Ks,2)  % if Ks is odd
    Jsm = [zeros(N/2-Ps,Ks); eye(Ks); zeros(N/2-Ps-1,Ks)];    % Beam Selector Matrix Jsm.
    G0af = Godd;                                              % N x Ks Beamformer matrix.
    G30af = diag(chebwin(N,30))*G0af;
else
    Jsm = [zeros(N/2-Ps,Ks); eye(Ks); zeros(N/2-Ps,Ks)];      % Beam Selector Matrix Jsm.
    G0af = Geven;                                             % N x Ks Beamformer matrix.
    G30af = diag(chebwin(N,30))*G0af;
end

G0maf = G0af*Jsm;
G30maf = G30af*Jsm;

A. SINR Computations for Displaced Filter Beamspace post-Doppler STAP

Solve a separate adaptive problem in each Doppler bin m:

for m=1:M

    Fm0  = toeplitz([F0(:,m);  zeros(Kt-1,1)],[F0(1,m)  zeros(1,Kt-1)]);
    Fm40 = toeplitz([F40(:,m); zeros(Kt-1,1)],[F40(1,m) zeros(1,Kt-1)]);

    Tm0  = kron(Fm0,G0);
    Tm40 = kron(Fm40,G30);

    Rum0  = Tm0'*Ru*Tm0;
    Rum40 = Tm40'*Ru*Tm40;

    bdfb = exp(-1i*2*pi*omegadopplerbank(m)*(0:M-1)).';
    gt = kron(bdfb,at);

    utm0  = Tm0'*gt;
    utm40 = Tm40'*gt;

    wm0  = Rum0\utm0;
    wm40 = Rum40\utm40;

    w0  = Tm0*wm0;
    w40 = Tm40*wm40;

    for n=1:Lfd
        bt = exp(-1i*2*pi*omegad(n)*(0:M-1)).'; % Dummy Target Doppler Steering Vector
        vt = kron(bt,at);
        SINR0df_mat(m,n)  = abs(w0'*vt)^2/real(w0'*Ru*w0);
        SINR30df_mat(m,n) = abs(w40'*vt)^2/real(w40'*Ru*w40);
    end
end

SINR0df  = max(abs(SINR0df_mat));
SINR40df = max(abs(SINR30df_mat));

LSINR0df  = SINR0df/SNRo;
LSINR40df = SINR40df/SNRo;

B. SINR Computations for Adjacent Filter Beamspace post-Doppler STAP

for m=1:M

    F0maf  = Fmaf0(:,:,m);
    F40maf = Fmaf40(:,:,m);

    Tm0af = kron(F0maf,G0maf);
    Tm40af = kron(F40maf,G30maf);

    Rum0af  = Tm0af'*Ru*Tm0af;
    Rum40af = Tm40af'*Ru*Tm40af;

    bdfb = exp(-1i*2*pi*omegadopplerbank(m)*(0:M-1)).';
    gt = kron(bdfb,at);

    utm0af  = Tm0af'*gt;
    utm40af = Tm40af'*gt;

    wm0af  = Rum0af\utm0af;
    wm40af = Rum40af\utm40af;

    w0af  = Tm0af*wm0af;
    w40af = Tm40af*wm40af;

    for n=1:Lfd
        bt = exp(-1i*2*pi*omegad(n)*(0:M-1)).'; % Dummy Target Doppler Steering Vector
        vt = kron(bt,at);
        SINR0af_mat(m,n)  = abs(w0af'*vt)^2/real(w0af'*Ru*w0af);
        SINR40af_mat(m,n) = abs(w40af'*vt)^2/real(w40af'*Ru*w40af);
    end

end

SINR0af  = max(abs(SINR0af_mat));
SINR40af = max(abs(SINR40af_mat));

LSINR0af  = SINR0af/SNRo;
LSINR40af = SINR40af/SNRo;

Plot the SINR Losses

figure('NumberTitle', 'off','Name', 'Figure 58. Beamspace post-Doppler in a clutter only environment. K=4, Ksm = Ktm = 2.');

%  Displaced Filter Beamspace post-Doppler STAP Results:
plot(fd,10*log10(LSINRopt),'k','LineWidth',1.5)
hold on;
plot(fd,10*log10(LSINR0df),'LineWidth',1.5)
plot(fd,10*log10(LSINR40df),'g','LineWidth',1.5)
grid on;
ylabel('SINR Loss (dB)');
xlabel('Target Doppler Frequency (Hz)');
ylim([-30 1]);
xlim([-5 305]);

%  Adjacent Filter Beamspace post-Doppler STAP Results:
plot(fd,10*log10(LSINR0af),'r','LineWidth',1.5)
hold on;
plot(fd,10*log10(LSINR40af),'m','LineWidth',1.5)
grid on;
ylabel('SINR Loss (dB)');
xlabel('Target Doppler Frequency (Hz)');
ylim([-30 1]);
xlim([-5 305]);
legend('Optimum','Displaced-Uniform', 'Displaced-Tapered', 'Adjacent-Uniform', 'Adjacent-Tapered','Location','South');