技术领域
[0001] 本
发明属于雷达技术领域，涉及提供一种机载非正侧阵雷达用自适应子空间的杂波抑制方法，具体地说是通过估计杂波二维分布曲线的未知构型参数来计算待检测距离单元的杂波子空间，再将数据向杂波子空间对应的
正交补空间进行投影来抑制杂波。
背景技术
[0002] 机载雷达主要处于下视工作状态，不可避免的要观测到许多地杂波，同时由于载机运动，杂波谱大大展宽，传统的脉冲多普勒技术不再有效。空时自适应处理（space time adaptive processing，STAP）技术可以有效的抑制地杂波，提高运动目标的检测性能。STAP技术通常选取待检测距离单元的邻近距离单元作为训练样本来计算自适应权值，当训练样本满足独立同分布（independent identically distributed，IID）的条件时，随着样本数目的增多自适应处理的性能逐渐收敛于最优。实际中机载雷达为了实现全方位搜索，通常在
机身同时架设多个阵列天线，每个天线负责一定的方位扫描，此时各个天线阵列的轴向和载机速度方向就会存在一定的夹
角，出现阵列为非正侧放置的情况。非正侧阵结构会导致杂波角度多普勒特性随距离变化，使训练样本不再满足独立同分布条件，严重影响了协方差矩阵的估计
精度，使空时自适应处理性能降低。
[0003] 对于杂波角度多普勒特性随距离变化的非平稳问题（即杂波距离非平稳性），国内外学者已提出多种解决方法。Zatman提出导数更新法（derivative-based updating，DBU），该方法假定权矢量是距离的一次函数，但当杂波距离非平稳性十分剧烈时，以上假设不再成立。Melvin和Davis提出角度多普勒补偿的方法（angle doppler compensation，ADW），该方法对主瓣杂波中心进行补偿，对旁瓣杂波的补偿则不是十分理想。Ries和Lapierre等提出基于配准的补偿方法（registration-based compensation，RBC），该方法基于杂波分布曲线，将不同距离单元的杂波
功率谱进行配准处理来解决杂波的距离非平稳性问题，但配准时不同距离单元曲线的映射关系难以确定。Friedlander等提出基于杂波子空间的
算法，该方法直接利用杂波二维分布曲线的构型参数计算待检测距离单元的杂波子空间，然后将回
波数据向杂波子空间对应的正交补空间进行投影来抑制掉杂波分量。该方法无需利用训练样本，避免了样本随距离变化的非平稳问题。但杂波子空间的计算需要利用杂波二维分布曲线的构型参数，而杂波二维分布曲线的构型参数和雷达系统参数、载机运动参数等有关，实际中这些参数不一定全部已知。
发明内容
[0004] 针对机载非正侧阵雷达的杂波距离非平稳性和子空间方法需要利用杂波二维分布曲线的构型参数这一问题，本发明的目的在于提供一种机载非正侧阵雷达用自适应子空间的杂波抑制方法，无需利用训练样本而且杂波抑制性能优异。
[0005] 为了达到上述目的，本发明采用以下技术方案予以实现。
[0006] 一种机载非正侧阵雷达用自适应子空间的杂波抑制方法，其特征在于，包括以下步骤：
[0007] 步骤1，建立机载非正侧阵雷达
信号模型，其中，雷达工作
波长为λ，相干处理间隔内发射M个脉冲，脉冲重复
频率为fr，天线阵列是N个阵元组成的等距线阵，阵元间距为d，载机高度为H，载机速度为v，载机速度方向与天线阵列轴向的夹角为ψ，地面杂波散射体的
俯仰角与方位角分别为θ和 ，载机与地面杂波散射体的斜距为R；
[0008] 其回波信号可以表示为
[0009] l=1,2…Nc
[0010] 式中n为噪声分量，xc为杂波分量，Nc为斜距为R的地面杂波环中独立的杂波
块数目，αl为杂波环第l个杂波块的幅度，且不同杂波块之间相互独立，ul为杂波环第l个杂波块的空时导向矢量；
[0011] 其杂波环的杂波二维分布曲线为
[0012]
[0013] 上式的简化函数表达式
[0014] C=Γ(η(λ,d,fr,v,ψ,θ))
[0015] 式中C为杂波曲线，Γ为通过参数矢量计算杂波二维分布曲线的算子，η为杂波二维分布曲线的参数矢量，其中杂波二维分布曲线的参数矢量η中，λ、d、fr为雷达系统已知参数，v、ψ、θ为未知量；
[0016] 步骤2，对回波信号进行角度多普勒平面离散化，再做加权
空域和时域的二维
傅立叶变换得到傅立叶谱Pb，然后对傅立叶谱Pb取绝对值的平方，得到对应的杂波功率谱；
[0017] 设置第一
门限，滤除杂波二维分布曲线外的网格单元，利用得到的杂波二维分布曲线上的网格单元构造搜索子矩阵 根据回波信号在角度多普勒平面对应网格单元上的复幅度γ具有稀疏特性，通过稀疏算法 求解杂波二维分布曲线上网格单元的复幅度γ，式中||||p表示Lp范数，ε表示稀疏算法求解的误差限，为设定值；
[0018] 根据杂波二维分布曲线上网格单元的复幅度γ，设置第二门限，滤除杂波二维分布曲线上较小幅值|γi|的次要网格单元，得到反映杂波的功率谱分布在对应空时二维平面分布的典型网格单元；
[0019] 步骤3，对杂波环的杂波二维分布曲线
[0020]
[0021] 进行变换
整理，写成矢量形式为y=aTβ，式中β=[1/v2 cosψ/v sin2ψcos2θ]T，a=[(frλwd/2)2 -λ2frwswd/d -1]T， 其中β与η中未知参数v、ψ、θ有关，a,y与η中已知参数λ、d、fr及杂波的空时二维频率有关；T表示转置操作。
[0022] 然后将步骤2得到的超过第二门限的P2个典型网格单元的空时二维频率及已知构型参数代入y=aTβ的a,y中，可以得到以下方程组
[0023]
[0024] 式中P2个网格单元对应y的实测值 为P2×1维的响应矢量，P2个典型网格单元对应的a的实测值 为P2×3维的量测矩阵，ξ为P2×1维的误差矢
量，为设定值；
[0025] 步骤4，采用最小截取二乘方法求解步骤3得到的方程组 得到β的估计值进而求解得到杂波二维分布曲线的构型参数v、ψ、θ对应的估计值
[0026] 步骤5，利用估计出来的杂波二维分布曲线的构型参数 与已知参数λ、d、fr，计算出杂波二维分布曲线的对应的杂波空时导向矢量矩阵 子空间正交投影算子和滤波权值向量W；
[0027] 步骤6，根据滤波权值向量W构造相应的空时
滤波器，抑制回波信号X中的杂波，得到当前杂波块杂波抑制后的最终多普勒谱：Z=WHX。
[0028] 具体地，利用第k个多普勒通道滤波权值向量Wk构造相应的空时滤波器，抑制回波信号X中的杂波，得到当前杂波块、第k个多普勒通道的输出为： k=1,2…K，K为多普勒通道数，H表示共轭转置操作；得到当前杂波块杂波抑制后的最终多普勒谱：Z=[Z1,Z2,…,ZK]T，T表示转置操作。
[0029] 本发明的技术方案的特点和进一步改进在于：
[0030] （1）步骤2的具体子步骤为：
[0031] 2a)将角度多普勒平面离散化，多普勒频率对应的网格单元数目为Nd，空域频率对应的网格单元数目为Ns；
[0032] 基于杂波块中的网格单元建立空时导向矢量搜索矩阵Φ，具体形式其中u(wd,i,ws,j)中i=1,…Nd;j=1,…Ns
表示坐标为(i,j)的网格单元的空时二维频率，即 式中 为
Kronecker积，其中ut(wd，i)和us(ws,j)分别为杂波块中坐标为(i,j)的网格单元的时域导向矢量和空域导向矢量；
[0033] 则回波信号X表示为X=Φγ， 为回波信号在角度多普勒平面对应网格单元上的复幅度，其中T表示转置操作；
[0034] 2b)对角度多普勒平面离散化回波信号做加权空域和时域的二维傅立叶变换，其形式为 式中Pb为傅立叶谱，tw为空域和时域的加权系数，为Hadamard积，得到傅立叶谱Pb；
[0035] 2c)对加权空域和时域二维傅立叶变换后得到的傅立叶谱Pb取绝对值的平方，得到对应的杂波功率谱，并对其进行第一门限检测，检测出功率超过第一门限值TH1的P1个网格单元，其对应的空时二维频率为 上标1表示通过第一门限TH1；并构建一个搜索子矩阵 其为搜索矩阵Φ的一个子集；
[0036] 2d)回波信号在角度多普勒平面对应网格单元上的复幅度γ具有稀疏特性，通过稀疏算法 求解估计单元上的复幅度γ，式中||||p表示Lp范数，ε表示稀疏算法求解的误差限，为设定值；
[0037] 2e)估计出网格单元上的复幅度γ的幅值后，按照幅值|γi|≥TH2，i=1,…NdNs，式中TH2为设定的第二门限值，滤除杂波二维分布曲线上较小幅值|γi|的次要网格单元，得到反映杂波的功率谱分布在对应空时二维平面分布的典型网格单元；
[0038] 2f)检测出的P2个典型网格单元对应的空时二维频率为
[0039] 上标2表示通过第二门限TH2；典型网格单元在空时二维平面的分布，反映了杂波的功率谱分布。
[0040] （2）步骤4的具体子步骤为：
[0041] 4a)采用最小截取二乘方法求解优化问题
[0042]
[0043] 式中zi表示标志位，zi=1表示该数据为正常点，zi=0表示该数据为干扰点，e表示全1矢量，Q为z中非零元素的个数，且Q
[0044] 4b)对上面的优化问题求解得到β的估计值 进而分别计算得到v、ψ、θ对应的估计值 其计算公式如下
[0045]
[0046]
[0047]
[0048] （3）步骤5具体子步骤为：
[0049] 5a)将得到的估计构型参数 和已知参数λ、d、fr，代入归一化多普勒频率和归一化空间频率 中，得到第l个杂波块对应的 估计的时域导向矢量 和空域导向矢量
则此时的第l个杂波块对应的估计的空时导向矢量为
估计的杂波空时导向矢量矩阵 l=1,2,…Nc，
式中 为Kronecker积；
[0050] 5b)按下式计算正交于杂波子空间的投影算子
[0051]
[0052] 式中μ是正数，以保证矩阵 可逆；⊥表示正交投影算子，Η表示共轭转置操作，()-1表示矩阵求逆，I表示单位阵；
[0053] 5c)利用目标空时导向矢量 即s=[s1,s2,…,sk,…sK]和步骤5b估计出来的正交于杂波子空间的投影算子 计算得到第k个多普勒通道的空时滤波器的滤波权值向量 k=1,2…K，K为多普勒通道数。滤波权值向量
[0054] 其中 为目标的时域导向矢量，为目标的空域导向矢量，其中归一化多普勒频率 归一化空间频率
v表示的是目标的速度；θ0表示的是主波束的俯仰角， 表示的
是主波束的方位角，ψ表示的是载机速度方向与天线阵列轴向的夹角。
[0055] 本发明的机载非正侧阵雷达用自适应子空间的杂波抑制方法，首先通过拟合回波数据估计的杂波功率谱与理论杂波曲线来估计杂波二维分布曲线的未知构型参数，在拟合的过程中采用LTS方法提高参数估计的稳健性，然后利用估计出的构型参数计算杂波子空间，再将数据向杂波子空间对应的正交补空间进行投影来抑制杂波。其仿真实验结果表明本发明的参数估计方法具有稳健性和准确性，利用估计参数计算的权值向量可以取得较好的杂波抑制效果。
[0056] 实际中，由于惯导的精度问题及地形高程知识的缺乏，将导致杂波二维分布曲线的构型参数不准确，此时就需要利用接收回波数据来估计这些未知构型参数。由于杂波的功率谱在角度多普勒平面的轨迹由杂波二维分布曲线决定，本发明通过拟合杂波功率谱与杂波二维分布曲线，来估计相应的构型参数。
[0057] 为了减小杂波的功率谱的扩散，本发明利用杂波在角度多普勒平面的稀疏性，采用稀疏算法估计回波数据的功率谱。实际上雷达回波数据中不仅有杂波信号，还可能有目标信号（主瓣目标或者旁瓣目标）。当目标信号
信噪比较大时，通过稀疏算法求解后，对应于目标的角度多普勒平面的网格单元也将通过门限检测。这些目标点为干扰点，将严重影响最小二乘（LS）方法的参数估计精度；为了解决这个问题，本发明采用最小截取二乘（least trimmed squares，LTS）方法代替LS方法来求解。
[0058] 本发明中，通过拟合回波数据估计的功率谱与理论杂波曲线来估计杂波二维分布曲线的未知构型参数，在拟合的过程中采用LTS方法，提高了参数估计的准确性和稳健性。
[0059] 本发明中，直接利用构型参数计算待检测单元的杂波子空间，然后将数据向杂波子空间对应的正交补空间进行投影来抑制掉杂波分量。该方法无需利用训练样本，避免了样本的非平稳问题。
附图说明
[0060] 下面结合附图及具体实施方式对本发明做进一步详细说明。
[0061] 图1是本发明的机载非正侧阵雷达杂波抑制
流程图。
[0062] 图2是本发明的机载雷达几何结构图。
[0063] 图3a是本发明在无干扰点时真实杂波曲线与最小二乘法（LS）和最小截取二乘（LTS）拟合结果示意图；横坐标表示归一化多普勒频率,纵坐标表示归一化空间频率。
[0064] 图3b是本发明在有干扰点时真实杂波曲线与最小二乘法（LS）和最小截取二乘（LTS）拟合结果示意图；横坐标表示归一化多普勒频率,纵坐标表示归一化空间频率。
[0065] 图4是最优处理、
采样协方差、本发明方法（自适应子空间）的信杂噪比损失曲线示意图；横坐标表示归一化多普勒频率,纵坐标表示信杂噪比损失（dB）。
具体实施方式
[0066] 参照图1，说明本发明的机载非正侧阵雷达用自适应子空间的杂波抑制方法，其具体步骤如下：
[0067] 步骤1，建立机载非正侧阵雷达信号模型，其杂波环的杂波二维分布曲线[0068]
[0069] 上式的简化函数表达式
[0070] C=Γ(η(λ,d,fr,v,ψ,θ))
[0071] 式中C为杂波曲线，Γ为通过参数矢量计算杂波二维分布曲线的算子，η为杂波二维分布曲线的参数矢量，其中杂波二维分布曲线的参数矢量η中，λ、d、fr为雷达系统已知参数，v、ψ、θ为未知量。
[0072] 具体举例说明如下：其中，机载非正侧阵雷达结构如附图2所示，雷达工作波长为λ，相干处理间隔内发射M个脉冲，脉冲重复频率为fr，天线阵列是N个阵元组成的等距线阵，阵元间距为d，载机高度为H，载机速度为v，载机速度方向与天线阵列轴向的夹角为ψ，地面杂波散射体的俯仰角与方位角分别为θ和 载机与地面杂波散射体的斜距为R。按照以上信号模型，回波信号可以表示为
[0073] l=1,2…Nc
[0074] 式中n为噪声分量，xc为杂波分量，Nc为斜距为R的地面杂波环中独立的杂波块数目，αl为杂波环第l个杂波块的幅度，且不同杂波块之间相互独立，ul为杂波环第l个杂波块的空时导向矢量。
[0075] 不失一般性，杂波块的空时导向矢量可以表示为
[0076]
[0077] 式中 为Kronecker积，ut(wd)和us(ws)分别为杂波块的时域导向矢量和空域导向矢量；
[0078]
[0079]
[0080] 式中 为归一化多普勒频率， 为归一化空间频率，T表示转置操作。
[0081] 利用wd与ws之间的关系进行变换整理可以得到杂波环在角度多普勒域的二维分布曲线为
[0082]
[0083] 此时任意一个杂波环的杂波二维分布曲线可以表示为
[0084] C=Γ(η(λ,d,fr,v,ψ,θ))
[0085] 式中C为杂波曲线，Γ为通过参数矢量计算杂波二维分布曲线的算子，η为杂波二维分布曲线的参数矢量，其中杂波二维分布曲线的参数矢量η中，λ、d、fr为雷达系统已知参数，v、ψ、θ为未知量。
[0086] 步骤2，对回波信号进行角度多普勒平面离散化，再做加权空域和时域的二维傅立叶变换得到傅立叶谱Pb，然后对傅立叶谱Pb取绝对值的平方，得到对应的杂波功率谱；
[0087] 设置第一门限，滤除杂波二维分布曲线外的网格单元，利用得到的杂波二维分布曲线上的网格单元构造搜索子矩阵 根据回波信号在角度多普勒平面对应网格单元上的复幅度γ具有稀疏特性，通过稀疏算法 求解杂波二维分布曲线上网格单元的复幅度γ，式中||||p表示Lp范数，ε表示稀疏算法求解的误差限，为设定值；
[0088] 根据杂波二维分布曲线上网格单元的复幅度γ，设置第二门限，滤除杂波二维分布曲线上较小幅值|γi|的次要网格单元，得到反映杂波的功率谱分布在对应空时二维平面分布的典型网格单元。
[0089] 具体子步骤说明如下：
[0090] 2a)将角度多普勒平面离散化，多普勒频率对应的网格单元数目为Nd，空域频率对应的网格单元数目为Ns。
[0091] 基于杂波块中的网格单元建立空时导向矢量搜索矩阵Φ，具体形式其中u(wd,i,ws,j)中i=1,…Nd;j=1,…Ns
表示坐标为(i,j)的网格单元的空时二维频率，即 式中 为
Kronecker积，其中ut(wd，i)和us(ws,j)分别为杂波块中坐标为(i,j)的网格单元的时域导向矢量和空域导向矢量。
[0092] 则回波信号X可以表示为X=Φγ， 为回波信号在角度多普勒平面对应网格单元上的复幅度，其中T表示转置。
[0093] 2b)对角度多普勒平面离散化回波信号做加权空域和时域的二维傅立叶变换，其形式为 式中Pb为傅立叶谱，tw为空域和时域的加权系数，为Hadamard积，得到傅立叶谱Pb。
[0094] 2c)对空域和时域二维傅立叶变换后得到的傅立叶谱Pb取绝对值的平方，得到对应的杂波功率谱，并对其进行第一门限检测，假定检测出共有P1个网格单元的功率超过第一门限值TH1，TH1为设定的第一门限值，一般设定为5-10dB。仿真实验中取8dB。
[0095] 这些网格单元对应的空时二维频率为 上标1表示通过第一门限TH1。
[0096] 根据这些网格单元构建一个新的搜索子矩阵 其为搜索矩阵Φ的一个子集。
[0097] 2d)由于杂波功率谱沿杂波二维分布曲线分布，且杂波分量相对于噪声分量在回波信号中占据主要成分。这就意味着估计的功率谱只在靠近杂波二维分布曲线附近的少数网格单元上幅值较大，而在其它
位置的网格单元上幅值很小并接近于零。
[0098] 即回波信号在角度多普勒平面对应网格单元上的复幅度γ具有稀疏特性，故可以通过稀疏算法 求解网格单元上的复幅度γ，式中||||p表示Lp范数，ε表示稀疏算法求解的误差限，为设定值。
[0099] 为了降低运算量，采用2c)得到的搜索子矩阵 代替稀疏算法中的Φ，即改进的稀疏算法为
[0100] 2e)估计出网格单元上的复幅度γ的幅值后，按照过门限准则来提取大的幅值分量|γi|≥TH2，i=1,…NdNs，式中TH2为设定的第二门限值，一般为10-20dB，仿真实验取16dB。
[0101] 这样，我们就滤除杂波二维分布曲线上较小幅值|γi|的次要网格单元，得到反映杂波的功率谱分布在对应空时二维平面分布的典型网格单元。
[0102] 2f)假定检测出共有P2个典型网格单元的幅值|γi|超过第二门限值TH2，这些典型网格单元对应的空时二维频率为 上标2表示通过第二门限TH2。这些典型网格单元在空时二维平面的分布，反映了杂波的功率谱分布。
[0103] 步骤3，对杂波环的杂波二维分布曲线
[0104]
[0105] 进行变换整理，写成矢量形式为y=aTβ，式中β=[1/v2 cosψ/v sin2ψcos2θ]T，a=[(frλwd/2)2 -λ2frwswd/d -1]T， 其中β与η中未知参数v、ψ、θ有关，a,y与η中已知参数λ、d、fr及杂波的空时二维频率有关；
[0106] 然后将步骤2得到的超过第二门限的P2个典型网格单元的空时二维频率及已知构型参数代入y=aTβ的a,y中，可以得到以下方程组
[0107]
[0108] 式中P2个网格单元对应y的实测值 为P2×1维的响应矢量，P2个典型网格单元对应的a的实测值 为P2×3维的量测矩阵，ξ为P2×1维的误差矢
量，为设定值。
[0109] 步骤4，采用最小截取二乘方法求解步骤3得到的方程组 得到β的估计值进而求解得到杂波二维分布曲线的构型参数v、ψ、θ对应的估计值 值。
[0110] 具体的子步骤说明如下：
[0111] 4a)LTS是一种稳健的回归方法，可以表示成以下优化问题
[0112]
[0113] 式中zi表示标志位，zi=1表示该数据为正常点，zi=0表示该数据为干扰点，e表示全1矢量，Q表示z中非零元素的个数,且Q
[0114] 4b)对上面的优化问题求解得到β的估计值 进而可以分别计算得到v、ψ、θ对应的估计值 其计算公式如下
[0115]
[0116]
[0117]
[0118] 步骤5，利用估计出来的杂波二维分布曲线的构型参数 与已知参数λ、d、fr，计算出杂波二维分布曲线的对应的杂波空时导向矢量矩阵 子空间正交投影算子和滤波权值向量W。
[0119] 具体的子步骤说明如下：
[0120] 5a)将得到的估计构型参数 和已知参数λ、d、fr，代入归一化多普勒频率和归一化空间频率 中，得到第l个杂波块对应的 估计的时 域导向 矢量 和空域导向 矢量
则此时的第l个杂波块对应的估计的空时导向矢量为
估计的杂波空时导向矢量矩阵 l=1,2,…Nc，
式中 为Kronecker积；
[0121] 5b)按下式计算正交于杂波子空间的投影算子
[0122]
[0123] 式中μ是正数，以保证矩阵 可逆；⊥表示正交投影算子，Η表示共轭转置操作，()-1表示矩阵求逆，I表示单位阵；
[0124] 5c)利用目标空时导向矢量 即s=[s1,s2,…,sk,…sK]和步骤5b估计出来的正交于杂波子空间的投影算子 计算得到第k个多普勒通道的空时滤波器的滤波权值向量 k=1,2…K，K为多普勒通道数。滤波权值向量
[0125] 其中 为目标的时域导向矢量，为目标的空域导向矢量，其中归一化多普勒频率 归一化空间频率
v表示的是目标的速度；θ0表示的是主波束的俯仰角， 表示的是
主波束的方位角，ψ表示的是载机速度方向与天线阵列轴向的夹角。
[0126] 步骤6，根据滤波权值向量W构造相应的空时滤波器，抑制回波信号X中的杂波，得到当前杂波块杂波抑制后的最终多普勒谱：Z=WHX。
[0127] 具体地，利用第k个多普勒通道滤波权值向量Wk构造相应的空时滤波器，抑制回波信号X中的杂波，得到当前杂波块、第k个多普勒通道的输出为： k=1,2…K，K为多普勒通道数，H表示共轭转置操作；得到当前杂波块杂波抑制后的最终多普勒谱：Z=[Z1,Z2,…,ZK]T，T表示转置操作。
[0128] 下面通过仿真实验，对本发明的杂波抑制方法的性能进行详细说明。
[0129] （1）仿真实验1验证本发明参数估计的稳健性
[0130] （1.1）仿真参数
[0131] 仿真参数如表1所示，待检测单元斜距为7km，杂噪比为50dB，加入两个干扰点，分别代表主瓣目标与旁瓣目标，对应的空时二维频率为(-0.40.0)、(0.2-0.25)，信噪比为20dB。
[0132] 表1机载雷达系统参数
[0133]
[0135] 附图3a和图3b分别给出了无干扰点和有干扰点时最小二乘（least square，LS）和LTS这两种方法的杂波曲线拟合结果。图中‘○’对应杂波点迹，‘＊’对应目标点迹。由图3a可以看出无干扰点时两种方法估计的杂波曲线与真实杂波曲线基本重合，表明了两种方法在无干扰点时均可以取得较高的参数估计精度。由图3b可以看出在存在干扰点时LS估计的杂波曲线明显偏离真实杂波曲线，这是因为LS方法通过最小化所有数据的剩余误差平方来求解系数矢量，未考虑干扰点的影响。LTS估计的杂波曲线仍然和真实杂波曲线基本重合，这是因为LTS假定数据中存在干扰点，通过优化选取正常数据进行计算，所以可以取得较好的参数估计效果。
[0136] （2）仿真实验2验证本发明参数估计的稳健性
[0137] （2.1）仿真参数
[0138] 仿真参数同仿真实验1，重复进行100次蒙特卡罗实验，得到LTS方法参数估计的均方根误差（root mean square error，RMSE）如表2所示。
[0139] 表2参数估计RMSE
[0140]
[0141] （2.2）仿真数据处理结果及分析
[0142] 由表2可以看出LTS方法估计的均方根误差较小，估计参数接近于真值。这是由于稀疏谱估计的杂波谱位置较为准确，由图3a和图3b可以看出杂波点迹基本位于真实杂波谱线附近，而且LTS方法通过加权抑制了干扰点的影响。
[0143] （3）仿真实验3验证本发明的杂波抑制性能
[0144] （3.1）实验方法
[0145] 以信杂噪比损失为准则验证算法的性能，信杂噪比损失定义为
[0146]
[0147] 式中 分别表示单阵元单脉冲的目标功率与噪声功率， 为目标的相关矩阵，R为理论的杂波加噪声协方差矩阵。
[0148] （3.2）仿真数据处理结果及分析
[0149] 附图4给出了最优处理、采样协方差、本文方法（自适应子空间）的信杂噪比损失曲线。由图4可以看出采样协方差相对于最优处理性能明显下降，这是因为采样协方差利用距离样本平均得到的协方差矩阵估计Ravg为R的有偏估计。本发明的自适应子空间方法性能明显优于采样协方差方法，这是因为子空间方法直接利用构型参数计算杂波子空间投影算子，避免了样本非平稳问题。
