CN104198998A - 非均匀背景下基于聚类处理的恒虚警检测方法 - Google Patents

Abstract

该发明公开了一种非均匀背景下基于聚类处理的恒虚警检测方法，属于雷达弱目标检测技术的领域，特别涉及了聚类分析和低信噪比下恒虚警(CFAR)检测技术。该方法首先进行地形分类，利用改进的K-means算法自适应的确定最优聚类数(不同地形的数目)，基于最优聚类数对雷达监测区域按幅度大小进行聚类处理，实现地形分类，使分类之后的地形呈现出多种不同地形；然后对分类完成的地形进行编号；在恒虚警检测阶段，利用编号的地形，筛选与待检测单元具有相同地形的分辨单元作为参考单元，估计检测门限，实现恒虚警检测。从而具有在非均匀检测背景下具有通用性强、分类精度高的效果。

Description

非均匀背景下基于聚类处理的恒虚警检测方法

技术领域

本发明属于雷达弱目标检测技术的领域，特别涉及了聚类分析和低信噪比下恒虚警(CFAR)检测技术领域。

背景技术

随着科技的发展，雷达分辨率的提高，雷达探测环境往往具有典型的非均匀性，特别对于海滨、城市以及其他复杂地形的陆地，这些区域往往具有地形快速起伏、多种不同地形交替变化等特点，这些特点导致CFAR参考窗内的数据不再满足独立同分布的条件，不能够准确的估计出干扰背景的功率，最终导致传统检测器虚警概率升高和目标的漏检。如果能够实现对雷达监测场景地形的有效正确分类，我们就可以利用分类的结果选取均匀参考单元，准确估计干扰背景功率，改善非均匀环境下雷达目标的检测性能，提升雷达的探测威力。因此，提高复杂环境下雷达对目标的检测能力具有重要的理论价值和实际意义。

恒虚警检测是在雷达自动检测系统中给检测策略提供检测阈值并且使杂波干扰对系统的虚警概率影响最小的一种自适应信号处理算法。针对非均匀背景，国内外的学者已经做了许多工作，2003年，美国雪城大学的Biao chen针对瑞利杂波背景提出了利用期望最大方法估计杂波边缘的位置，从而筛选均匀参考单元；2006年，意大利的De Maio提出利用地理信息系统来挑选均匀参考单元，改善检测性能。2007年，土耳其的ASELSAN公司的Do_yuran提出利用Anderson-Darling拟合优度检验来检测参考单元中杂波的均匀性，并估计杂波距离变化的位置来改善非均匀背景下CFAR检测性能；2012年，针对非均匀Weibull背景，Pourmottaghi.A提出利用最大似然估计来寻找杂波边缘的位置，挑选均匀参考单元；然而，上述工作只是针对存在一个杂波边缘的简单非均匀场景展开的，当雷达监测区域存在多个杂波边缘，即多种复杂地形时，上述所有方法将不再有效。2013年，普渡大学的Jeong Hun Kin针对复杂地形提出基于图像融合的参考窗自适应选择算法，但是该算法对于初始的小数据窗的大小的选择比较敏感，会导致结果的不准确。2013年，电子科技大学的彭馨仪结合杂波边缘检测算法和二元积累实现对地形的有效划分，但是此方法需要事先知道监测区域的杂波统计特性，当实际地形的统计特性与假设的模型失配时，检测性能大大下降。且对于分类数目是随机确定的，没有给出合理准确的选择准则。

发明内容

本发明的目的是针对背景技术存在的缺陷，研究设计一种非均匀背景下基于聚类处理的恒虚警检测方法，从而达到在非均匀检测背景下具有通用性强、分类精度高的目的。

本发明提出了一种非均匀背景下基于聚类处理的恒虚警检测方法，该方法分为地形分类和恒虚警检测两个阶段。首先，在地形分类阶段，利用改进的K-means算法自适应的确定最优聚类数(不同地形的数目)，基于最优聚类数对雷达监测区域按幅度大小进行聚类处理，实现地形分类，使分类之后的地形呈现出多种不同地形；然后对分类完成的地形进行编号，方便下一步检测参考单元的选取；在恒虚警检测阶段，利用编号的地形，筛选与待检测单元具有相同地形的分辨单元作为参考单元，估计检测门限，实现恒虚警检测。因而本发明包括以下步骤：

步骤1、初始化参数包括：误差平方和J，二维CFAR参考窗大小L_min×L_min和L_max×L_max；

步骤2、从雷达接收机中读取第i帧回波数据：

Z(i)＝{z_i(m,n)}

其中，1≤m≤N_r,1≤n≤N_a，m和n分别为距离向和方位向的量化状态，i为帧数；N_r是距离维量化的总单元个数，N_a是方位向量化的总单元个数；z_i(m,n)表示第i帧回波数据的量测单元(m,n)中的量测值，为回波数据的幅度；

步骤3、设置聚类数范围[k_min,k_max]，其中，k_min为设定的聚类数的最小值，k_max为设定的聚类数的最大值；

步骤4、初始各聚类中心值 $C_j\left(I\right)=(j\×\max \left(Z\right))/(\hat{k}+1),j=\mathrm{1,2},......\hat{k},$ 初始化

其中 $C_j\left(I\right),j=\mathrm{1,2},......\hat{k}$ 为第I次聚类的个聚类中心；

步骤5、计算第I次聚类时每个量测值与聚类中心的距离D(Z_i(m,n),C_j(I))＝|Z_i(m,n)-C_j(I)|，

其中，D(Z_i(m,n),C_j(I))表示第I次聚类时量测单元(m,n)中的量测值与第j个聚类中心的距离，如果满足 $D(Z_i(m,n),C_l\left(I\right))=\min \left\{D(Z_i(m,n),C_j\left(I\right),j=\mathrm{1,2},......\hat{k})\right\},$ 那么量测值Z_i(m,n)属于类中的第l类；

步骤6、计算第I次聚类的误差平方和准则函数：

其中误差平方和准则函数J(I)是聚类的评价标准，表示由步骤5判断的属于第l类的雷达量测值，n_l表示属于第l类的雷达量测值的总个数；

步骤7、如果|J(I)-J(I-1)|>ε,则I＝I+1，然后计算个新的聚类中心C_j，新的个聚类中心的取值为：本次聚类结果中每类(共类)数据幅度的平均值，然后返回步骤5，否则进入步骤8；

步骤8、当步骤7完成后，计算出该聚类结果的评价指标Sil，然后当返回步骤4；若则进入步骤9；

步骤9、确定最优聚类数k_opt以及最优聚类结果：

利用公式k_min+1,......k_max求取最优聚类数，然后对应下的聚类结果即为最优的聚类结果；

步骤10、将属于第l个聚类中心的量测值的地形编号l；

步骤11、选取参考单元数据。若待检测单元所在地形面积小于参考单元要求的最小面积，则首先选取与待检测单元具有相同地形的所有分辨单元作为参考单元的一部分,然后选取与待检测单元所在位置距离最近的分辨单元来补充参考单元；若待检测单元所在地形面积大于参考单元要求的最大面积，则选取与待检测单元所在位置距离最近的分辨单元作为参考单元；若待检测单元所在地形面积小于参考单元要求的最大面积且大于最小面积，则选取与待检测单元具有相同地形的所有分辨单元作为参考单元；

步骤12、利用选取的参考单元数据及与待检单元杂波统计特性匹配的传统CFAR检测器进行目标检测。

所述步骤8中评价指标 ${\mathrm{Sil}}_{\hat{k}}={\mathrm{\Σ}}_{l=1}^{\hat{k}}{\mathrm{\Σ}}_{p=1}^{n_l}((b_{\mathrm{lp}}-a_{\mathrm{lp}})/\max \{a_{\mathrm{lp}},b_{\mathrm{lp}}\left\}\right),$ a_lp是第l类的第p个样本到第l类内其他样本的平均距离，b_lp是第l类内的第p个样本到其他每个类中样本平均距离的最小值。

所述步骤11中首先假设待检单元所在地形的编号为f_cut，雷达回波数据平面所有单元(除待检单元外)的地形编号为f_j,j＝1,2,......k_opt，当f_j＝f_cut时，记A(t)＝1，否则A(t)＝0；

如果参考单元选取距离f_cut位置最近的A(t)＝1的L_max×L_max个分辨单元的数据，L_max×L_max为参考单元的最大尺寸；

如果则参考单元选取所有A(t)＝1的分辨单元的数据，其中L_min×L_min为参考单元的最小尺寸；

如果则参考单元选取所有A(t)＝1的分辨单元的数据，以及距离f_cut位置最近的 $L_{\min }\×L_{\min }-{\mathrm{\Σ}}_{t=1}^{N_r\×N_a-1}A\left(t\right)$ 个分辨单元的数据。

本发明非均匀背景下基于聚类处理的恒虚警检测方法，利用聚类分析方法，自适应的对盲的、非均匀雷达监测区域进行准确的地形分类，准确的筛选出与待检测单元同分布的参考单元，从而具有在非均匀检测背景下具有通用性强、分类精度高的效果。

附图说明

图1为本发明非均匀背景下基于聚类处理的恒虚警检测方法的流程图；

图2为对计算机仿真产生的GIS地形图和MSTAR二维杂波公开数据进行地形分类处理的结果对比；

图3为没有执行地形分类与执行了地形分类的仿真结果对比；

图4为对IPIX实测数据进行处理的结果。

图2中图2-1为由计算机仿真产生的某块真实地形的GIS地形图，数据大小为100×100，不同的颜色代表不同的地形；图2-2为利用聚类方法对仿真产生的GIS地形图进行地形分类并且编号的结果；图2-3是MSTAR二维杂波公开数据中HB06198文件的灰度图；图2-4是对HB06198进行地形分类的结果；

图3中图3-1为对GIS地形图右下角50×50的数据进行检测的仿真结果对比，图3-2为对GIS地形图右上角50×50的数据进行检测的仿真结果对比，图3-3为对GIS地形图左下角50×50的数据进行检测的仿真结果对比，图3-4为对GIS地形图中间50×50的数据进行检测的仿真结果对比；

图4中图4-1为IPIX回波幅度分布图；图4-2为对图4-1进行地形分类的结果；图4-3为IPIX实测数据下传统CAFR检测结果；图4-4为IPIX实测数据下本发明所提CAFR方法检测结果。

具体实施方式

本发明主要通过计算机仿真产生的某块真实地形的GIS地形图和IPIX实测数据来验证，GIS数据大小为100×100，GIS地形图中不同颜色代表不同地形。根据已有结论，地面杂波可用四种统计分布模型模拟，分别是瑞利分布、Log-normal分布、Weibull分布、K分布。GIS地图右下角主要由农田构成，近似服从瑞利分布，GIS地图右上角主要由城市、森林以及农田构成，近似服从Log-normal分布，GIS地图中间部分主要由城市和森林构成，近似服从Weibull分布，GIS地图左下角主要有城市和铁路构成，近似服从K分布，所有步骤、结论都在MATLAB-R2010b上验证确认。具体实施步骤如下：

步骤1、初始化参数包括：误差平方和J(0)＝0，二维CFAR参考窗大小L_min×L_min＝4×4和L_max×L_max＝8×8；

步骤2、从雷达接收机中读取第i帧回波数据：

Z(i)＝{z_i(m,n)}

其中，1≤m≤N_r,1≤n≤N_a，m和n分别为距离向和方位向的量化状态，i为帧数；N_r是距离维量化的总单元个数，N_a是方位向量化的总单元个数；z_i(m,n)表示第i帧回波数据的量测单元(m,n)中的量测值，为回波数据的幅度；

步骤3、设置聚类数范围[k_min＝2,k_max＝20]。

步骤4、设定初始聚类中心值。

开始循环，为第I次聚类的个聚类中心，初始化I＝1，根据步骤2中得到的幅度量测值Z，一般设定个初始聚类中心值为：以所有数据幅度的最小值为等差序列的第一项，为公差的个等差数列值。

步骤5、计算第I次聚类时每个量测值与聚类中心的距离。

每个量测值与聚类中心的距离为D(Z_i(m,n),C_j(I))＝|Z_i(m,n)-C_j(I)|，其中，D(Z_i(m,n),C_j(I))表示量测单元(m,n)中的量测值与第j个聚类中心的距离，如果满足 $D(Z_i(m,n),C_l\left(I\right))=\min \left\{D(Z_i(m,n),C_j\left(I\right),j=\mathrm{1,2},......\hat{k})\right\},$ 那么量测值Z_i(m,n)属于类中的第l类。

步骤6、计算第I次聚类的误差平方和准则函数：

误差平方和准则函数J(I)是聚类的评价标准，表示由步骤5判断的属于第l类的雷达量测值，n_l表示属于第l类的雷达量测值的总个数；

步骤7、如果|J(I)-J(I-1)|>ε,则I＝I+1，然后计算个新的聚类中心C_j，新的个聚类中心的取值为：本次聚类结果中每类(共类)数据幅度的平均值，然后返回步骤5，否则进入步骤8；

步骤8、利用评价指标Sil评价聚类数下，聚类结果的好坏：

a_lp是第l类的第p个样本到第l类内其他样本的平均距离，b_lp是第l类内的第p个样本到其他每个类中样本平均距离的最小值。然后当 $\hat{k}=\hat{k}+1<k_{\max },$ 返回步骤4；若则进入步骤9；

步骤9、确定最优聚类数k_opt以及最优聚类结果：

利用公式k_min+1,......k_max求取最优聚类数，然后对应下的聚类结果即为最优的聚类结果；

步骤10、将属于第l个聚类中心的量测值的地形编号为l；

步骤11、选取参考单元数据：

假设待检单元所在地形的编号为f_cut，雷达回波数据平面所有单元(去除待检单元外)的地形编号为f_j,j＝1,2,......k_opt，当f_j＝f_cut时，记A(t)＝1，否则A(t)＝0；

如果参考单元选取距离f_cut位置最近的A(t)＝1的64个分辨单元的数据做参考单元数据；

如果则参考单元选取所有A(t)＝1的分辨单元的数据；

如果则参考单元选取所有A(t)＝1的分辨单元的数据，以及距离fcut位置最近的 $16-{\mathrm{\&Sigma;}}_{t=1}^{N_r\&times;N_a-1}A\left(t\right)$ 个分辨单元的数据；

步骤12、利用选取的参考单元数据及与待检单元杂波统计特性匹配的传统CFAR检测器进行目标检测。

通过上面的步骤，就可以在盲背景下，自适应选择地形分类数，实现地形的准确划分，提升非均匀背景下的恒虚警检测能力。

对GIS地形图分类的结果如图2(b)所示，对MSTAR二维杂波公开数据的地形分类的结果如图2(d)所示。分别与原图对比可知，本发明所提的地形分类方法能够对雷达监测地形进行有效的地形划分。

在上述仿真中，没有执行地形分类步骤与执行了该步骤的仿真对比结果如图3所示。由图3可知，在未知非均匀背景下，本发明提出的基于地形分类的的恒虚警检测性能，在四种常见杂波分布类型下都有一定程度的提升。

由图4的结果可知，对这一帧IPIX实测数据而言，本发明所提CFAR检测算法虚警点个数341个，实际虚警率P′_fa＝1.2×10^-3，满足设定虚警概率P_fa＝1×10^-3的要求。而传统CFAR在目标信杂比虚警点个数为3979个，实际虚警率P′_fa＝1.4×10^-2，虚警概率大大超出设定值。即本发明所提地形分类CFAR算法实际的虚警点数相比传统CFAR减少了91.43％，更好的抑制了杂波虚警点。在目标信杂比(SCR)为12dB时，地形分类CFAR检测算法开始检测到目标，而传统二维CA-CFAR在目标信杂比(SCR)为26dB时，才开始检测到目标，可见，对这一帧IPIX实测数据而言，本发明所提出的地形分类CFAR检测算法与传统CFAR检测算法相比输入信杂比有14dB的改善。

通过本发明的具体实施可以看出，本发明在盲背景下，能够自适应的确定分类数目，实现对地形的有效正确划分，准确筛选均匀参考单元，实现在地形快速变化、复杂未知非均匀场景下目标的有效检测。

本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。本领域的普通技术人员可以根据本发明公开的这些技术启示做出各种不脱离本发明实质的其它各种具体变形和组合，这些变形和组合仍然在本发明的保护范围内。

Claims (3)

1.非均匀背景下基于聚类处理的恒虚警检测方法，该方法包括：

步骤1、初始化参数包括：误差平方和J，二维CFAR参考窗大小L_min×L_min和L_max×L_max；

步骤2、从雷达接收机中读取第i帧回波数据：

Z(i)＝{z_i(m,n)}

其中，1≤m≤N_r,1≤n≤N_a，m和n分别为距离向和方位向的量化状态，i为帧数；N_r是距离维量化的总单元个数，N_a是方位向量化的总单元个数；z_i(m,n)表示第i帧回波数据的量测单元(m,n)中的量测值，为回波数据的幅度；

步骤3、设置聚类数范围[k_min,k_max]，其中，k_min为设定的聚类数的最小值，k_max为设定的聚类数的最大值；

步骤4、初始各聚类中心值 $C_j\left(I\right)=(j\&times;\max \left(Z\right))/(\hat{k}+1),j=\mathrm{1,2},......\hat{k},$ 初始化

其中 $C_j\left(I\right),j=\mathrm{1,2},......\hat{k}$ 为第I次聚类的个聚类中心；

步骤5、计算第I次聚类时每个量测值与聚类中心的距离D(Z_i(m,n),C_j(I))＝|Z_i(m,n)-C_j(I)|，

其中，D(Z_i(m,n),C_j(I))表示第I次聚类时量测单元(m,n)中的量测值与第j个聚类中心的距离，如果满足 $D(Z_i(m,n),C_l\left(I\right))=\min \left\{D(Z_i(m,n),C_j\left(I\right),j=\mathrm{1,2},......\hat{k})\right\},$ 那么量测值Z_i(m,n)属于类中的第l类；

步骤6、计算第I次聚类的误差平方和准则函数：

其中误差平方和准则函数J(I)是聚类的评价标准，表示由步骤5判断的属于第l类的雷达量测值，n_l表示属于第l类的雷达量测值的总个数；

步骤7、如果|J(I)-J(I-1)|>ε,则I＝I+1，然后计算个新的聚类中心C_j，新的个聚类中心的取值为：本次聚类结果中每类(共类)数据幅度的平均值，然后返回步骤5，否则进入步骤8；

步骤8、当步骤7完成后，计算出该聚类结果的评价指标Sil，然后当返回步骤4；若则进入步骤9；

步骤9、确定最优聚类数k_opt以及最优聚类结果：

利用公式k_min+1,......k_max求取最优聚类数，然后对应下的聚类结果即为最优的聚类结果；

步骤10、将属于第l个聚类中心的量测值的地形编号l；

步骤11、选取参考单元数据。若待检测单元所在地形面积小于参考单元要求的最小面积，则首先选取与待检测单元具有相同地形的所有分辨单元作为参考单元的一部分,然后选取与待检测单元所在位置距离最近的分辨单元来补充参考单元；若待检测单元所在地形面积大于参考单元要求的最大面积，则选取与待检测单元所在位置距离最近的分辨单元作为参考单元；若待检测单元所在地形面积小于参考单元要求的最大面积且大于最小面积，则选取与待检测单元具有相同地形的所有分辨单元作为参考单元；

步骤12、利用选取的参考单元数据及与待检单元杂波统计特性匹配的传统CFAR检测器进行目标检测。

2.如权利要求1所述的非均匀背景下基于聚类处理的恒虚警检测方法，其特征在于步骤8中评价指标 ${\mathrm{Sil}}_{\hat{k}}={\mathrm{\&Sigma;}}_{l=1}^{\hat{k}}{\mathrm{\&Sigma;}}_{p=1}^{n_l}((b_{\mathrm{lp}}-a_{\mathrm{lp}})/\max \{a_{\mathrm{lp}},b_{\mathrm{lp}}\left\}\right),$ a_lp是第l类的第p个样本到第l类内其他样本的平均距离，b_lp是第l类内的第p个样本到其他每个类中样本平均距离的最小值。

3.如权利要求1所述的非均匀背景下基于聚类处理的恒虚警检测方法，其特征在于步骤11中首先假设待检单元所在地形的编号为f_cut，雷达回波数据平面所有单元(除待检单元外)的地形编号为f_j,j＝1,2,......k_opt，当f_j＝f_cut时，记A(t)＝1，否则A(t)＝0；

如果参考单元选取距离f_cut位置最近的A(t)＝1的L_max×L_max个分辨单元的数据，L_max×L_max为参考单元的最大尺寸；

如果则参考单元选取所有A(t)＝1的分辨单元的数据，其中L_min×L_min为参考单元的最小尺寸；

如果则参考单元选取所有A(t)＝1的分辨单元的数据，以及距离f_cut位置最近的 $L_{\min }\&times;L_{\min }-{\mathrm{\&Sigma;}}_{t=1}^{N_r\&times;N_a-1}A\left(t\right)$ 个分辨单元的数据。