基于马氏距离重加权的非均匀杂波抑制方法与流程

本发明属于雷达信号处理技术领域，特别涉及基于马氏距离重加权的非均匀杂波抑制方法，用于非均匀杂波环境下的杂波抑制及慢速目标检测。

背景技术：

机载预警机雷达由于常处在下视工作状态，其工作性能必然会受到地面杂波的强烈影响。同时，由于机载雷达平台与地面的相对运动，使得机载雷达接收的杂波呈现出很强的空时耦合性，给弱目标、慢速目标的检测带来很大困难。在上世纪七十年代，brennan等人提出了空时自适应处理(space-timeadaptiveprocessor，stap)方法，通过在空域和时域的二维处理来进行杂波抑制达到弱目标和慢速目标检测的目的。经过这些年的研究和发展，stap方法成为一种有效的杂波抑制手段，并且可以有效提高雷达对弱小、慢速目标检测的能力。

通常，stap方法在理想情况下，即环境均匀情况下，可以取得较好的性能。stap方法在实施时需要利用训练样本来估计待检测单元的杂波协方差矩阵并计算自适应权矢量。理想状况下训练样本满足独立同分布且不包含目标信号时，stap方法可以取得较好的杂波抑制性能。然而，在真实情况中机载雷达面临的实际杂波环境一般是非均匀的，尤其是在地面慢速运动目标检测时，由于机载雷达波束照射范围内运动目标比较多，训练样本中往往包含有目标信号，从而导致训练样本不符合独立同分布状态，此时样本的均匀假设条件将不再成立。如果利用这些包含目标信号的非均匀的样本来估计待检测单元的杂波协方差矩阵并计算自适应权矢量，进而利用自适应权矢量进行滤波，将导致待检测样本中目标信号自相消，使得目标检测概率下降。

为此，有学者提出了非均匀检测器(nonhomogeneitydetector，nhd)，nhd能有效的检测并剔除被污染的训练样本，使其不参与权矢量的运算，因而能改善stap的性能。现有的nhd有：(1)广义内积检测器(generalizedinnerproduct，gip)，gip对于训练样本中存在明显目标信号时检测效果较好，然而由于广义内积法只提取协方差矩阵中的大特征值，故在多目标信号的情况下，其检测性能将会下降；自适应功率剩余(adaptivepowerresidue，apr)检测器，同样该检测器在训练样本中的目标信号数量较少、信噪比较低时可以获取较好的检测效果，当训练样本中的目标信号数量较多、信噪比较高时，训练样本所构成的干扰加噪声协方差矩阵会受到目标信号严重扰动，apr检测器性能下降；对角加载的广义内积检测器(ioadedgeneralizedinnerproduct，lgip)，该检测器通过对由训练样本估计出的协方差矩阵对角加载来降低目标信号对检测器的扰动，存在如何选择合适的加载系量以提高检测性能的难题；基于长球波函数的广义内积检测器(prolatespheroidalwavefunctionsgip，pswf-gip)，该检测器通过长球波函数的基向量直接估计待检测样本的杂波协方差矩阵，避免了训练样本中目标信号的影响，在理论上可以取得较好的性能，但在计算长球波函数的基向量时要求雷达构型参数及阵列流型精确已知，而实际情况中的阵列误差以及其他因素会影响阵列流型及雷达构型参数，使得检测性能下降。

技术实现要素：

有鉴于此，本发明提供了基于马氏距离重加权的非均匀杂波抑制方法，能够在非均匀杂波情况下有效抑制杂波，从而提高检测性能。

为达到上述目的，本发明采用如下技术方案：

一种基于马氏距离重加权的非均匀杂波抑制方法，包括以下步骤：

步骤1，获取机载雷达天线阵列接收到的原始回波数据，天线阵列为包括n个阵元且阵元间距为半波长的均匀线阵，原始回波数据包括目标分量、杂波分量以及噪声分量；

步骤2，对原始回波数据进行下变频处理和脉冲压缩处理，得到空时快拍数据；

利用降维矩阵对空时快拍数据进行降维处理，得到降维后的空时快拍数据，降维后的空时快拍数据包含l个距离单元的数据；

获取目标空时导向矢量，利用降维矩阵对目标空时导向矢量进行降维处理，得到降维后的目标空时导向矢量；

步骤3，初始化：令待检测单元序号l＝1，即令l个距离单元中的第l个距离单元为当前待检测单元；

步骤4，假设第l个距离单元存在对应的kl个采样样本，对第l个距离单元的kl个采样样本进行求平均，得到第l个距离单元的基准样本，并计算第l个距离单元的基准样本对应的基准协方差矩阵；

利用基准协方差矩阵，计算得到kl个采样样本中每个采样样本与基准样本之间的马氏距离，进而利用各采样样本与基准样本之间的马氏距离，构造每个采样样本对应的剔除系数；

利用各采样样本对应的剔除系数以及预设的选取门限，确定kl个采样样本中的干扰样本，并剔除kl个采样样本中的干扰样本，得到kl′个数据样本；其中，干扰样本为剔除系数大于选取门限的采样样本；

步骤5，利用基准协方差矩阵，计算得到kl′个数据样本与基准样本间的马氏距离；

判断kl′是否大于或等于6n：若kl′大于或等于6n，则选取kl′个数据样本中与基准样本间马氏距离最小的6n个数据样本作为第l个距离单元的训练样本；若kl′小于6n，则选取kl′个数据样本作为第l个距离单元的训练样本；

利用每个训练样本与基准样本间马氏距离，构造该训练样本对应的重加权系数；

利用各训练样本及其对应的重加权系数进行重加权，构造得到第l个距离单元对应的采样协方差矩阵；

步骤6，利用降维后的目标空时导向矢量以及第l个距离单元对应的采样协方差矩阵，构造得到第l个距离单元对应的自适应权矢量；

利用第l个距离单元对应的自适应权矢量对第l个距离单元的数据进行杂波抑制，得到第l个距离单元对应的输出功率；

步骤7，令待检测单元序号l加1，重复执行步骤4至步骤6，直至l＝l，即完成对全部l个距离单元的杂波抑制。

基于本发明上述方案，利用非均匀杂波数据的特性，通过采样数据与基准数据间的马氏距离确定各采样数据对应的剔除系数，利用剔除系数及预设的选取门限剔除采样样本中的干扰样本得到数据样本，进而利用所得到的数据样本与基准样本之间的马氏距离，选取对应的训练样本，再进一步利用训练样本与基准样本之间的马氏距离构造对应的重加权系数进行重加权，构造得到对应的采样协方差矩阵，最后利用所得的协方差矩阵构造自适应权矢量，利用该自适应权矢量进行杂波抑制，得到对应的输出功率。因此，相比于传统广义内积法，本发明方法在采样样本数不同的情况下可以稳健的获得较好的杂波抑制性能。同时，本发明方法对于复杂杂波背景中的杂波非均匀抑制具有普适性和稳健性，能够有效保证stap在非均匀杂波环境中的杂波抑制性能。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明实施例提供的一种基于马氏距离重加权的非均匀杂波抑制方法的流程示意图；

图2为采样样本取2×3n时本发明方法和gip算法的杂波剩余曲线对比图，图中横坐标为多普勒通道号，纵坐标为功率大小，功率单位为db；

图3为采样样本取3×3n时本发明方法和gip算法的杂波剩余曲线对比图，图中横坐标为多普勒通道号，纵坐标为功率大小，功率单位为db；

图4为采样样本取4×3n时本发明方法和gip算法的杂波剩余曲线对比图，图中横坐标为多普勒通道号，纵坐标为功率大小，功率单位为db；

图5a为本发明方法在采样样本分别取2×3n、3×3n、4×3n情况下的杂波剩余曲线图；

图5b为gip算法在采样样本分别取2×3n、3×3n、4×3n情况下的杂波剩余曲线图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

图1所示为本发明实施例提供的一种基于马氏距离重加权的非均匀杂波抑制方法的流程示意图。

如图1所示，本发明实施例提供的基于马氏距离重加权的非均匀杂波抑制方法，包括以下步骤：

步骤1，获取机载雷达天线阵列接收到的原始回波数据。

其中，机载雷达天线阵列为包括n个阵元且阵元间距为半波长的均匀线阵，原始回波数据包括目标分量、杂波分量以及噪声分量。

为叙述方便，在本实施例中，假设机载雷达在一个相干处理间隔(coherentprocessinginterval，cpi)内发射m个脉冲，其脉冲重复频率为fr，雷达工作波长为λ，载机高度为h，载机速度为v。

步骤2，对原始回波数据进行下变频处理和脉冲压缩处理，得到空时快拍数据，进而利用降维矩阵对空时快拍数据进行降维处理，得到降维后的空时快拍数据，以及获取目标空时导向矢量，利用降维矩阵对目标空时导向矢量进行降维处理，得到降维后的目标空时导向矢量。

其中，降维后的空时快拍数据包含l个距离单元的数据。假设其中的第k个距离单元的空时快拍数据可表示为则该距离单元对应的二元假设检验可以表示为：

h0：xk＝ck+nk

h1：xk＝αss-t+xk/h0’

其中，h0假设表示该距离单元的回波数据不包含目标信号，h1假设表示该距离单元的回波数据包含目标信号；αss-t，ck，nk分别表示目标分量、杂波分量和噪声分量，其中α代表目标的复幅度，ss-t表示目标空时导向矢量，式中代表克罗内克积，ss表示时间导向矢量，st表示空间导向矢量。

其中，时间导向矢量ss和空间导向矢量st又可分别表示为：

式中，表示归一化多普勒频率，表示归一化空间频率，fd表示目标多普勒频率，θ表示方位角，表示俯仰角，d表示阵元间距。

具体的，步骤2中，利用降维矩阵对空时快拍数据进行降维处理，得到降维后的空时快拍数据，可以包括：

根据第三预设公式，利用降维矩阵对空时快拍数据进行降维处理，得到降维后的空时快拍数据及降维后的目标空时导向矢量。

其中，第三预设公式包括：

式中，x表示空时快拍数据，t表示降维矩阵，表示降维后的空时快拍数据，xl表示第l个距离单元的数据，l＝1，2，…，l，上标h表示共轭转置操作，的维数为3n×l，x的维数为mn×l，t的维数为mn×3n，m表示机载雷达在一个相干处理间隔内发射的脉冲个数。其中，降维矩阵t可表示为其中f表示m×3维的傅里叶矩阵加权，in表示n×n维的单位矩阵，表示克罗内克积。

具体的，步骤2中，利用降维矩阵对空时快拍数据进行降维处理，得到降维后的空时快拍数据，具体可以包括：

根据第四预设公式，利用降维矩阵对目标空时导向矢量进行降维处理，得到降维后的目标空时导向矢量。

其中，第四预设公式包括：

式中，ss-t表示目标空时导向矢量，表示降维后的目标空时导向矢量。

步骤3，初始化：令待检测单元序号l＝1，即令l个距离单元中的第l个距离单元为当前待检测单元。

步骤4，假设第l个距离单元存在对应的kl个采样样本，对第l个距离单元的kl个采样样本进行求平均，得到第l个距离单元的基准样本，并计算第l个距离单元的基准样本对应的基准协方差矩阵；利用基准协方差矩阵计算得到kl个采样样本中每个采样样本与基准样本之间的马氏距离，进而利用各采样样本与基准样本之间的马氏距离，构造每个采样样本对应的剔除系数；利用各采样样本对应的剔除系数以及预设的选取门限，确定kl个采样样本中的干扰样本，并剔除kl个采样样本中的干扰样本，得到kl′个数据样本。

其中，干扰样本为剔除系数大于选取门限的采样样本。需要说明的是，所述的选取门限具体可设置为0.3，即当某一采样样本对应的马氏距离占总体比重大于0.3时，即视其为干扰数据将其剔除。当然，本领域技术人员可以理解，以上所述仅为本发明实施例提供的一种可选的实现方式，在实际应用中可根据具体需求设置选取门限，本发明实施例对此不作具体限定。

具体的，步骤4中，对第l个距离单元的kl个采样样本进行求平均，得到第l个距离单元的基准样本，并计算第l个距离单元的基准样本对应的基准协方差矩阵，具体可以包括以下步骤：

对第l个距离单元的kl个采样样本进行求平均，得到第l个距离单元的基准样本

利用第l个距离单元的基准样本μ，计算第l个距离单元的基准样本对应的基准协方差矩阵ru＝μμ^t+σn²i。

其中，i＝1，2，…，kl，σn表示对角加载系数，σn取底噪水平，i表示单位矩阵，上标t表示转置操作。

进一步的，步骤4中，利用基准协方差矩阵，计算得到kl个采样样本中每个采样样本与基准样本之间的马氏距离，进而利用各采样样本与基准样本之间的马氏距离，构造每个采样样本对应的剔除系数，具体可以包括以下步骤：

利用基准协方差矩阵，计算得到kl个采样样本中每个采样样本与基准样本之间的马氏距离；

计算各采样样本与基准样本间马氏距离的倒数，进而对各采样样本与基准样本间马氏距离的倒数进行归一化处理；

将各采样样本对应的归一化处理后得到的数值作为该采样样本的剔除系数。

其中，第i个采样样本与基准样本之间的马氏距离为第i个采样样本的剔除系数为其中上标-1表示求逆操作。

步骤5，利用基准协方差矩阵，计算得到kl′个数据样本与基准样本间的马氏距离；判断kl′是否大于或等于6n：若kl′大于或等于6n，则选取kl′个数据样本中与基准样本间马氏距离最小的6n个数据样本作为第l个距离单元的训练样本；若kl′小于6n，则选取kl′个数据样本作为第l个距离单元的训练样本；利用每个训练样本与基准样本间马氏距离，构造该训练样本对应的重加权系数；利用各训练样本及其对应的重加权系数进行重加权，构造得到第l个距离单元对应的采样协方差矩阵。

本发明实施例的一种优选实现方式中，步骤5中，利用每个训练样本与基准样本间马氏距离，构造该训练样本对应的重加权系数，具体可以包括：

计算得到每个训练样本与基准样本间马氏距离的倒数，进而对各训练样本与基准样本间马氏距离的倒数进行归一化处理；

将各训练样本对应的归一化处理后的得到的数值作为该训练样本的重加权系数。

其中，第l个距离单元的第k个训练样本的剔除系数为dlk表示第k个训练样本与基准样本之间的马氏距离，k＝1，2，…，k，k表示第l个距离单元的训练样本个数。

可以理解，k的具体取值取决于kl′。具体来说，当kl′≥6n时，k＝6n；当kl′＜6n时，k＝kl′。

具体的，步骤5中，利用各训练样本及其对应的重加权系数进行重加权，构造得到第l个距离单元对应的采样协方差矩阵，具体包括：

利用各训练样本及其对应的重加权系数，根据第一预设公式进行重加权，构造得到第l个距离单元对应的采样协方差矩阵。

其中，第一预设公式包括：

式中，表示第l个距离单元对应的采样协方差矩阵，表示第l个距离单元的第k个训练样本，表示第l个距离单元的第k个训练样本对应的重加权系数，上标h表示共轭转置操作，∑表示求和符号，k＝1，2，…，k，k表示第l个距离单元的训练样本个数。

步骤6，利用降维后的目标空时导向矢量以及第l个距离单元对应的采样协方差矩阵，构造得到第l个距离单元对应的自适应权矢量；利用第l个距离单元对应的自适应权矢量对第l个距离单元的数据进行杂波抑制，得到第l个距离单元对应的输出功率。

具体的，步骤6中，利用降维后的目标空时导向矢量以及第l个距离单元对应的采样协方差矩阵，构造得到第l个距离单元对应的自适应权矢量，具体可以包括：

利用降维后的目标空时导向矢量以及第l个距离单元对应的采样协方差矩阵，根据第二预设公式，构造得到第l个距离单元对应的自适应权矢量。

其中，第二预设公式包括：

式中，w表示第l个距离单元对应的自适应权矢量，表示第l个距离单元对应的采样协方差矩阵，表示目标空时导向矢量，上标-1表示求逆操作，上标h表示共轭转置操作。

进一步的，步骤6中，利用第l个距离单元对应的自适应权矢量对第l个距离单元的数据进行杂波抑制，得到第l个距离单元对应的输出功率，具体包括：

根据第五预设公式，利用第l个距离单元对应的自适应权矢量对第l个距离单元的数据进行杂波抑制，得到第l个距离单元对应的输出功率。

其中，第五预设公式包括：

yl＝w^hxl，

式中，w表示第l个距离单元对应的自适应权矢量，yl表示第l个距离单元对应的输出功率。

步骤7，令待检测单元序号l加1，重复执行步骤4至步骤6，直至l＝l，即完成对全部l个距离单元的杂波抑制。

至此，即完成了对全部l个距离单元的杂波抑制，本发明实施例提供的基于马氏距离重加权的非均匀杂波抑制方法即结束。

基于本发明上述方案，利用非均匀杂波数据的特性，通过采样数据与基准数据间的马氏距离确定各采样数据对应的剔除系数，利用剔除系数及预设的选取门限剔除采样样本中的干扰样本得到数据样本，进而利用所得到的数据样本与基准样本之间的马氏距离，选取对应的训练样本，再进一步利用训练样本与基准样本之间的马氏距离构造对应的重加权系数进行重加权，构造得到对应的采样协方差矩阵，最后利用所得的协方差矩阵构造自适应权矢量，利用该自适应权矢量进行杂波抑制，得到对应的输出功率。因此，相比于传统广义内积法，本发明方法在采样样本数不同的情况下可以稳健的获得较好的杂波抑制性能。同时，本发明方法对于复杂杂波背景中的杂波非均匀抑制具有普适性和稳健性，能够有效保证stap在非均匀杂波环境中的杂波抑制性能。

以下通过仿真实验验证本发明方法的正确性及上述有益效果：

1、实测数据：

本发明通过采用机载雷达的实测数据验证算法的可行性，本发明采用的是mcarm数据，mcarm数据是在1996年由room实验室为研究stap及其相关技术而录取的一些列机载雷达数据。本仿真实验具体采用的是第5次飞行中录取序号为575的数据来进行实验分析，具体参数如表1所示。

表1

2、仿真内容及结果分析

(1)仿真内容

鉴于mcarm数据的前150个距离单元数据受污染，故本仿真实验从mcarm数据的第151个距离门开始采样，选取第151到第384个距离门的数据作为实验数据。利用本发明方法和传统的gip算法，分别在采样样本数为2×3n、3×3n、4×3n时进行杂波抑制，并对本发明方法和传统的gip算法的杂波抑制性能进行比较。其中，图3为采样样本取2×3n时本发明方法和gip算法的杂波剩余曲线对比图，图4为采样样本取3×3n时本发明方法和gip算法的杂波剩余曲线对比图，图5为采样样本取4×3n时本发明方法和gip算法的杂波剩余曲线对比图。其中，所述的杂波剩余曲线指的是距离向取平均后的功率曲线。

(2)仿真结果分析

从图2中可以看出，本发明方法和gip算法在采样数为2×3n时的杂波抑制性能相差不大，杂波功率曲线总体比较吻合。

从图3中可以看出，本发明方法和gip算法在采样数为3×3n时，本发明中提及的算法的杂波抑制性能优于gip的杂波抑制性能。

从图4中可以看出，本发明方法和gip算法在采样数为4×3n时，本发明中提及的算法的杂波抑制性能优于gip的杂波抑制性能。

从图5a、图5b中可以看出，本发明方法和gip算法在采样数分别为2×3n、3×3n、4×3n时，gip算法随着采样数的增多抑制杂波的性能逐渐下降，但是本发明提及算法在采样数增多时可以保持一个相对稳健的杂波抑制性能，总体上来说杂波抑制性能优于传统的gip算法。

本领域普通技术人员可以理解：实现上述方法实施例的全部或部分步骤可以通过程序指令相关的硬件来完成，前述的程序可以存储于一计算机可读取存储介质中，该程序在执行时，执行包括上述方法实施例的步骤；而前述的存储介质包括：rom、ram、磁碟或者光盘等各种可以存储程序代码的介质。

以上，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以权利要求的保护范围为准。