MATLAB Examples

# Figure 38. SINR loss for Element Space pre-Doppler. Zero intrinsic cluller motion.

## Contents

```clc; clear; close all; ```

```fo = 450e6; % Operating Frequency in Hz Pt = 200e3; % Peak Transmit Power 200 kW Gt = 22; % Transmit Gain in dB Gr = 10; % Column Receive Gain in dB B = 4e6; % Receiver Instantaneous Bandwidth in Hz Ls = 4; % System Losses in dB fr = 300; % PRF in Hz Tr = 1/fr; % PRI in sec. M = 18; % Number of Pulses per CPI. K1 = 2; % Number of Pulses per sub-CPI. K2 = 3; % Number of Pulses per sub-CPI. M1 = M - K1 + 1; % Number of sub-CPI's for K=2. M2 = M - K2 + 1; % Number of sub-CPI's for K=3. Tp = 200e-6; % Pulse Width in sec. N = 18; % Number of Array Antenna Elements Gel = 4; % Element Gain in dB be = -30; % Element Backlobe Level in db Nc = 360; % Number of clutter patches uniformly distributed in azimuth. c = 299792458; % Speed of Light in m/sec. lambda = c/fo; % Operating wavelength in meters. d = lambda/2; % Interelement Spacing % Azimuth angle in degrees: phi = -180:179; Lphi = length(phi); f = zeros(1,Lphi); AF = zeros(1,Lphi); % Array Factor pre-allocation. % Platform Parameters: beta = 1; % beta parameter. ha = 9e3; % Platform altitude in meters. ```

## Thermal Noise Power Computations

```k = 1.3806488e-23; % Boltzmann Constant in J/K. To = 290; % Standard room Temperature in Kelvin. F = 3; % Receiver Noise Figure in dB; Te = To*(10^(F/10)-1); % Effective Receiver Temperature in Kelvin. Lr = 2.68; % System Losses on receive in dB. Ts = 10^(Lr/10)*Te; % Reception System Noise Temperature in Kelvin. Nn = k*Ts; % Receiver Noise PSD in Watts/Hz. Pn = Nn*B; % Receiver Noise Power in Watts sigma2 = 1; % Normalized Noise Power in Watts. ```

## Clutter Patch Geometry computations

```Rcik = 130000; % (clutter) range of interest in meters. dphi = 2*pi/Nc; % Azimuth angle increment in rad. dR = c/2/B; % Radar Range Resolution in meters. Re = 6370000; % Earth Radius in meters. ae = 4/3*Re; % Effective Earth Radius in meters. psi = asin(ha/Rcik); % Grazing angle at the clutter patch in rad (flat earth model). theta = psi; % Elevation (look-down angle) in rad. Flat earth assumption. gamma = 10^(-3/10); % Terrain-dependent reflectivity factor. phia = 0; % Velocity Misalignment angle in degrees. ```

## Clutter-to-Noise Ratio (CNR) Calculation

```% Calculate the Voltage Element Pattern: for i =1:Lphi if abs(phi(i))<=90 f(i) = cos(phi(i)*pi/180); else f(i) = 10^(be/10)*cos(phi(i)*pi/180); end end % Calculate the Array Factor (AF) (Voltage): steering_angle = 0; % Angle of beam steering in degrees. for k=1:Lphi AF(k) = sum(exp(-1i*2*pi/lambda*d*(0:N-1)*(sin(phi(k)*pi/180) ... - sin(steering_angle*pi/180)))); end % Calculate the Full Array Transmit Power Gain: Gtgain = 10^(Gt/10)*abs(AF).^2; % Calculate the Element Receive Power Gain: grgain = 10^(Gel/10)*abs(f).^2; % Clutter Patch RCS Calculation: PatchArea = Rcik*dphi*dR*sec(psi); sigma0 = gamma*sin(psi); sigma = sigma0*PatchArea; % Calculate the Clutter to Noise Ratio (CNR) for each clutter patch: ksi = Pt*Gtgain.*grgain*10^(Gr/10)*lambda^2*sigma/((4*pi)^3*Pn*10^(Ls/10)*Rcik^4); Ksic = sigma2*diag(ksi); ```

## Clutter Covariance Matrix Computations

```% Platform Velocity for beta parameter value: va = round(beta*d*fr/2); Ita = d/lambda*cos(theta); % Calculate Spatial and Doppler Frequencies for k-th clutter patch. % Spatial frequency of the k-th clutter patch: fsp = Ita*sin(phi*pi/180); % Normalized Doppler Frequency of the k-th clutter patch: omegac = beta*Ita*sin(phi*pi/180 + phia*pi/180); % Clutter Steering Vector Pre-allocation for sub-CPI of size K: a = zeros(N,Nc); b = zeros(M,Nc); b1 = zeros(K1,Nc); b2 = zeros(K2,Nc); Vc = zeros(M*N,Nc); Vc1 = zeros(K1*N,Nc); Vc2 = zeros(K2*N,Nc); for k=1:Nc a(:,k) = exp(1i*2*pi*fsp(k)*(0:N-1)); % Spatial Steering Vector. b(:,k) = exp(1i*2*pi*omegac(k)*(0:M-1)); % Temporal Steering Vector. b1(:,k) = exp(1i*2*pi*omegac(k)*(0:K1-1)); % Temporal Steering Vector for K = 2. b2(:,k) = exp(1i*2*pi*omegac(k)*(0:K2-1)); % Temporal Steering Vector for K = 3. Vc(:,k) = kron(b(:,k),a(:,k)); % Space-Time Steering Vector. Vc1(:,k) = kron(b1(:,k),a(:,k)); % Space-Time Steering Vector for K = 2. Vc2(:,k) = kron(b2(:,k),a(:,k)); % Space-Time Steering Vector for K = 3. end Rc = Vc*Ksic*Vc'; Rcsub1 = Vc1*Ksic*Vc1'; Rcsub2 = Vc2*Ksic*Vc2'; Rn = sigma2*eye(M*N); Rnsub1 = sigma2*eye(K1*N); Rnsub2 = sigma2*eye(K2*N); ```

## Jammer Covariance Matrix Calculation

```J = 2; % Number of Jammers. thetaj = 0; phij = [-40 25]; % Jammer elevation and azimuth angles in degrees. R_j = [370 370]*1e3; Sj = 1e-3; % Jammer ERPD in Watts/Hz. fspj = d/lambda*cos(thetaj*pi/180)*sin(phij*pi/180); % Spatial frequency of the j-th jammer. Lrj = 1.92; % System Losses on Receive in dB. Aj = zeros(N,J); for j=1:J Aj(:,j) = exp(1i*2*pi*fspj(j)*(0:N-1)); % Jammer Spatial Steering Vector. end indices= zeros(1,J); for j=1:J indices(j) = find(phi == phij(j)); end grgn = grgain(indices); ksi_j = (Sj*grgn*lambda^2)./((4*pi)^2.*Nn*10^(Lrj/10).*R_j.^2); Ksi_j = sigma2*diag(ksi_j); Phi_j = Aj*Ksi_j*Aj'; % Eq. (47) Rj = kron(eye(M),Phi_j); Rjsub1 = kron(eye(K1),Phi_j); Rjsub2 = kron(eye(K2),Phi_j); ```

## Analytic Interference Covariance Matrix for Fully Adaptive STAP:

```Ru = Rc + Rj + Rn; ```

## Analytic Interference Covariance Matrix for first sub-CPI (K=2):

```Rusub1 = Rcsub1 + Rjsub1 + Rnsub1; ```

## Analytic Interference Covariance Matrix for second sub-CPI (K=3):

```Rusub2 = Rcsub2 + Rjsub2 + Rnsub2; ```

## 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); omegact = fdt/fr; at = exp(1i*2*pi*fspt*(0:N-1)).'; % Target Spatial Steering Vector. ta = chebwin(N,30); % 30 dB Chebychev Spatial Tapper. bt1 = [1; -1;]; bt2 = [1; -1; 1;]; tb1 = [1; 1;]; % Binomial Doppler Taper for K = 2. tb2 = [1; 2; 1;]; % Binomial Doppler Taper for K = 3. gt1 = kron(tb1.*bt1,ta.*at); % sub-CPI #1 desired response. gt2 = kron(tb2.*bt2,ta.*at); % sub-CPI #2 desired response. ```

## Tapered, Element-Space STAP Solution

```wsub1 = Rusub1\gt1; % Weight Vector for K = 2; wsub2 = Rusub2\gt2; % Weight Vector for K = 3; ```

## W1 matrix construction. Equation (180).

```W1 = zeros(M*N,M1); % Assume that weight vectors for each sub-CPI are equal. for k=1:M1 W1((k-1)*N+1:(k+1)*N,k) = wsub1; end W2 = zeros(M*N,M2); % W2 matrix construction. Equation (180). % Assume that weight vectors for each sub-CPI are equal,. for k=1:M2 W2((k-1)*N+1:(k+2)*N,k) = wsub2; end td1 = chebwin(M1,40); % 40 dB Chebychev Doppler Taper. td2 = chebwin(M2,40); % 40 dB Chebychev Doppler Taper. ```

## SINR calculation for Optimum and pre-Doppler Element Space STAP:

```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); % Target Spatial Frequency at = exp(1i*2*pi*fspt*(0:N-1)).'; % Target Spatial Steering Vector. fd = 0:.5:300; Lfd = length(fd); omegad = fd/fr; dopplerfilterbank1 = linspace(0,300,M1+1); dopplerfilterbank2 = linspace(0,300,M2+1); omegadopplerbank1 = dopplerfilterbank1/fr; omegadopplerbank2 = dopplerfilterbank2/fr; LSINRopt = zeros(1,Lfd); SINRsub1_mat = zeros(length(dopplerfilterbank1),Lfd); SINRsub2_mat = zeros(length(dopplerfilterbank2),Lfd); InvRu = inv(Ru); SNRo = M*N; ```

## LSINR Computation for Optimum Fully Adaptive Case:

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

## SINR Computations for sub-CPI with K=2.

```for m=1:length(dopplerfilterbank1) um = exp(1i*2*pi*omegadopplerbank1(m)*(0:M1-1)).'; % Doppler Filter Steering Vector fm1 = td1.*um; wm1 = W1*fm1; % m-th Doppler bin composite weight vector. for n=1:Lfd bt = exp(1i*2*pi*omegad(n)*(0:M-1)).'; % Dummy Target Doppler Steering Vector vt = kron(bt,at); SINRsub1_mat(m,n) = abs(wm1'*vt)^2/real(wm1'*Ru*wm1); % Eq. (185) end end SINRsub1 = max(abs(SINRsub1_mat)); % Eq. (186) LSINRsub1 = SINRsub1/SNRo; ```

## SINR Computations for sub-CPI with K=3.

```for m=1:length(dopplerfilterbank2) um = exp(1i*2*pi*omegadopplerbank2(m)*(0:M2-1)).'; % Doppler Filter Steering Vector fm2 = td2.*um; wm2 = W2*fm2; % m-th Doppler bin composite weight vector. for n=1:Lfd bt = exp(1i*2*pi*omegad(n)*(0:M-1)).'; % Dummy Target Doppler Steering Vector vt = kron(bt,at); SINRsub2_mat(m,n) = abs(wm2'*vt)^2/real(wm2'*Ru*wm2); % Eq. (185) end end SINRsub2 = max(abs(SINRsub2_mat)); % Eq. (186) LSINRsub2 = SINRsub2/SNRo; ```

## Plot the SINR Losses

```figure('NumberTitle', 'off','Name', ... 'Figure 38. SINR loss for element space pre-Doppler. Zero intrinsic clutter motion',... 'Position', [1 1 700 600]); plot(fd,10*log10(LSINRopt),'LineWidth',1.5) hold on; plot(fd,10*log10(LSINRsub1),'r','LineWidth',1.5) hold on; plot(fd,10*log10(LSINRsub2),'g','LineWidth',1.5) ylabel('SINR Loss (dB)'); xlabel('Target Doppler Frequency (Hz)'); ylim([-30 1]); xlim([-5 305]); legend('Optimum Fully Adaptive', 'Pre-Doppler K=2','Pre-Doppler K=3', 'Location','South') grid on; ```