空时自适应处理（STAP）原理
文章来源：微信公众号 EW Frontier
一、核心背景与优势
常规雷达处理将波束形成（空间域）、匹配滤波（快时间域）、多普勒处理（慢时间域）分离执行，仅能利用M+N个自由度（N为线性相控阵元数，M为多普勒处理脉冲数），易因信息割裂牺牲性能。空时自适应处理（STAP）的核心创新的是联合空间与时间维度处理，通过在每个天线元件处执行多普勒处理后再形成波束，可获得MN个自由度，能在角度-多普勒联合空间中精准区分目标与干扰（杂波、干扰机等），仅抑制含干扰的单元，避免常规分离处理中对非干扰单元的误抑制。
二、空间处理原理（STAP空间维度）
空间处理的核心目标是通过设计天线阵列权重向量，最大化输出端信号与干扰加噪声功率比（SINR），其数学基础为柯西-施瓦茨不等式。
1. 基础模型与数学表达
考虑N元线性阵列，目标位于角度θₛ，定义目标导向矢量s(θₛ)，包含各阵元接收信号的相位信息：
s(θs)=[1ej(2πdsinθs)/λ⋮ej(2π(N−1)dsinθs)/λ]s(\theta_{s})=\left[\begin{array}{c}1 \\ e^{j\left(2 \pi d sin \theta_{s}\right) / \lambda} \\ \vdots \\ e^{j\left(2 \pi(N-1) d sin \theta_{s}\right) / \lambda}\end{array}\right]s(θs)=1ej(2πdsinθs)/λ⋮ej(2π(N−1)dsinθs)/λ
（d为阵元间距，λ为波长）
阵列输出电压V(θₛ)可表示为权重向量w与导向矢量的共轭转置乘积：V(θs)=PswHs(θs)V\left(\theta_{s}\right)=\sqrt{P_{s}} w^{H} s\left(\theta_{s}\right)V(θs)=PswHs(θs)（Pₛ为目标信号功率，H表示共轭转置）。
2. 权重设计逻辑
仅考虑信号时：约束权重向量范数有限，通过柯西-施瓦茨不等式最大化|V(θₛ)|²，得到权重向量与目标导向矢量成正比（w=κs(θₛ)，κ为常数），对应阵列辐射方向图在θ=θₛ处形成峰值。
含噪声时：噪声功率用协方差矩阵Rₙ表示（Rₙ=PₙI，Pₙ为噪声功率，I为单位矩阵），输出信噪比（SNR）为：
SNR=Ps∣wHs(θs)∣2wHRnwSNR=\frac{P_{s}\left|w^{H} s\left(\theta_{s}\right)\right|^{2}}{w^{H} R_{n} w}SNR=wHRnwPs∣wHs(θs)∣2
最大化SNR的权重解与仅信号场景一致。
含相关干扰（杂波、干扰机）时：干扰信号用导向矢量d(φₙ)描述（φₙ为干扰角度），总干扰协方差矩阵Rᵢ=ΣPᵢₙd(φₙ)dᴴ(φₙ)（Pᵢₙ为第n个干扰功率）。总协方差矩阵R=Rₙ+Rᵢ，此时需最大化信号与干扰加噪声比（SINR）：
SINR=Ps∣wHs(θs)∣2wHRwSINR=\frac{P_{s}\left|w^{H} s\left(\theta_{s}\right)\right|^{2}}{w^{H} R w}SINR=wHRwPs∣wHs(θs)∣2
利用矩阵变换与柯西-施瓦茨不等式，最终得到最优权重：w=R−1s(θs)w=R^{-1} s\left(\theta_{s}\right)w=R−1s(θs)，该权重可将主波束对准目标，同时在干扰角度位置形成零陷（抑制干扰）。
三、时间处理原理（STAP时间维度）
时间处理对应多普勒处理，核心是通过设计脉冲序列的滤波权重，利用目标与干扰的多普勒频移差异区分信号。
1. 信号模型
发射信号为M个脉冲的串：vT(t)=ej2πfct∑l=0M−1p(t−lT)v_{T}(t)=e^{j 2 \pi f_{c} t} \sum_{l=0}^{M-1} p(t-l T)vT(t)=ej2πfct∑l=0M−1p(t−lT)（f_c为载波频率，p(t)为脉冲形式，T为脉冲重复周期）。
目标回波因运动产生多普勒频移f_d（f_d=-2Ṙ/c，Ṙ为目标径向速度，c为光速），经外差、匹配滤波后，采样输出为：Vl(fd)=ej2πfdlTmRCV_{l}\left(f_{d}\right)=e^{j 2 \pi f_{d} l T} m_{R C}Vl(fd)=ej2πfdlTmRC（m_RC为匹配滤波响应采样值）。
2. 权重设计逻辑
定义多普勒导向矢量s(fₛ)（fₛ为目标多普勒频移），多普勒处理器（M阶FIR滤波器）输出可表示为：V(fs)=PSωHs(fs)V\left(f_{s}\right)=\sqrt{P_{S}} \omega^{H} s\left(f_{s}\right)V(fs)=PSωHs(fs)（ω为时间维度权重向量）。
干扰与噪声的建模形式与空间处理一致，最大化SINR的权重设计逻辑相同，即权重与多普勒导向矢量、总协方差矩阵相关，最终实现对目标多普勒频移的聚焦和干扰多普勒的抑制。
四、空时联合处理原理（STAP核心）
空时联合处理将空间阵列权重与时间多普勒权重融合，形成角度-多普勒联合空间的处理，是STAP的核心实现。
1. 联合模型与数学表达
联合输出电压为空间与时间维度输出的乘积，扩展为MN项求和，定义空时联合导向矢量s(θₛ,fₛ)，联合权重向量w（维度MN×1），则联合输出为：V(θs,fs)=wHs(θs,fs)V\left(\theta_{s}, f_{s}\right)=w^{H} s\left(\theta_{s}, f_{s}\right)V(θs,fs)=wHs(θs,fs)。
总干扰加噪声协方差矩阵仍为R=Rₙ+Rᵢ，其中Rₙ=PₙI（MN阶单位矩阵），Rᵢ=ΣPᵢₙd(φₙ,fₙ)dᴴ(φₙ,fₙ)（d(φₙ,fₙ)为干扰的空时联合导向矢量）。
2. 核心目标与特性
权重设计目标：通过w=R−1s(θs,fs)w=R^{-1} s\left(\theta_{s}, f_{s}\right)w=R−1s(θs,fs)，在角度-多普勒联合空间中，将主波束精准放置于目标对应的(θₛ,fₛ)位置，同时在杂波、干扰机对应的(φₙ,fₙ)位置形成零陷，实现对目标信号的最大化提取和干扰的精准抑制。
维度特性：联合处理需计算MN个权重（远多于常规处理的M+N个），且需为每个距离单元单独计算权重，存在显著计算负担，这也是STAP后续研究中需优化的核心方向。
MATLAB代码
%%clear all close all N=16;% Number of elements in ULAM=16;% Number of pulses?lambda=2;% Operating wavelengthd=lambda/2;% Distance of elements in ULAv=d;% Speed of the aircraft? Example in textbook didn't specifyT=v/2;% Duration of the pulse?Nclutter=200;% How many clutter patches% Target doppler and AoAfd=0.25;AoA=deg2rad(0);% Resolution of clutter patches% The response of the beamformer is plotted against these angles and% dopplersangVec=d/lambda*sin(deg2rad(linspace(-90,90,Nclutter)));% norm anglebeta=2*v*T/d;fdVec=beta*angVec;% Spatial steering vectors[iGrid,angGrid]=ndgrid(0:N-1,angVec);A=exp(1i*2*pi.*iGrid.*angGrid);% Doppler steering vectors[iGrid,fdGrid]=ndgrid(0:M-1,fdVec);B=exp(1i*2*pi.*iGrid.*fdGrid);% Space-time steering vectorV=zeros(N*M,Nclutter);foriClutter=1:NclutterV(:,iClutter)=kron(A(:,iClutter),B(:,iClutter));end% Response to a target[~,iMin]=min(abs(angVec-AoA));a=A(:,iMin);[~,iMin]=min(abs(fdVec-fd));b=B(:,iMin);% Covariance matrixR=eye(N*M);% Clutter covariance matrixR_clutter=V*V';% Jammer covariance matrixAoA_jammer=[32126];R_jammer=zeros(N*M);forang=AoA_jammer[~,iMin]=min(abs(angVec-d/lambda*sin(deg2rad(ang))));bj=B(:,iMin);aj=A;s=kron(bj,aj);R_jammer=R_jammer+db2mag(50)^2*(s*s');endRinv=inv(R);% An inefficient, but perhaps more informative way to do the code below% Rinv = eye(N*M);% Y = zeros(length(angVec));% for i = 1:length(angVec)% for j = 1:length(fdVec)%% x = kron(A(:, i), B(:, j));% s = kron(a, b);% w = s;% % w = Rinv*s;%% Y(i, j) = w'*x;% end% endW=kron(A,B);X=kron(a,b);Y=W'*Rinv*X;Y=reshape(Y,[Nclutter Nclutter]);% No jammer or cluttersubplot(2,2,1)Rinv=inv(R);W=kron(A,B);X=kron(a,b);Y=W'*Rinv*X;Y=reshape(Y,[Nclutter Nclutter]);s=surf(angVec,fdVec,mag2db(abs(Y')/(N*M)));s.EdgeColor='none';view(-90,90)colorbar% clim([-60 0])ylabel('Angle')xlabel('Doppler')% Cluttersubplot(2,2,2)Rinv=inv(R+R_clutter);W=kron(A,B);X=kron(a,b);Y=W'*Rinv*X;Y=reshape(Y,[Nclutter Nclutter]);s=surf(angVec,fdVec,mag2db(abs(Y')/(N*M)));s.EdgeColor='none';view(-90,90)colorbar% clim([-60 0])ylabel('Angle')xlabel('Doppler')title('Clutter Suppression')% Clutter and jammersubplot(2,2,3)Rinv=inv(R+R_clutter+R_jammer);W=kron(A,B);X=kron(a,b);Y=W'*Rinv*X;Y=reshape(Y,[Nclutter Nclutter]);s=surf(angVec,fdVec,mag2db(abs(Y')/(N*M)));s.EdgeColor='none';view(-90,90)colorbar% clim([-60 0])ylabel('Angle')xlabel('Doppler')title('Clutter and Jammer Suppression')% Jammersubplot(2,2,4)Rinv=inv(R+R_jammer);W=kron(A,B);X=kron(a,b);Y=W'*Rinv*X;Y=reshape(Y,[Nclutter Nclutter]);s=surf(angVec,fdVec,mag2db(abs(Y')/(N*M)));s.EdgeColor='none';view(-90,90)colorbar% clim([-60 0])ylabel('Angle')xlabel('Doppler')title('Jammer Suppression')saveas(gcf,'stap.png')%% Space-time Chebyshev windowingfigure hold onsubplot(1,2,1)Rinv=inv(R);W=kron(A,B);X=kron(a,b);Y=W'*Rinv*X;Y=reshape(Y,[Nclutter Nclutter]);s=surf(angVec,fdVec,mag2db(abs(Y')/(N*M)));s.EdgeColor='none';view(-90,90)colorbar% clim([-60 0])ylabel('Angle')xlabel('Doppler')% Chebyshev windowingsubplot(1,2,2)Rinv=inv(R);W=kron(A,B);X=kron(a,b);Y=W'*Rinv*(X.*kron(chebwin(N,30),chebwin(N,30)));Y=reshape(Y,[Nclutter Nclutter]);s=surf(angVec,fdVec,mag2db(abs(Y')/(N*M)));s.EdgeColor='none';view(-90,90)colorbar% clim([-60 0])ylabel('Angle')xlabel('Doppler')
