CN106407989B - 基于快速自动密度聚类的可变尺度红外小目标的检测方法 - Google Patents

Abstract

基于快速自动密度聚类的可变尺度红外小目标的检测方法，涉及图像分析与图像理解。利用SURF算子解决尺度可变小目标的特征提取。解决特征对于环境的敏感问题；快速自动化密度聚类算法，针对问题域，直接从空域完成轨迹的分割、提取，克服了原有序列检测问题中数据融合带来的指数级复杂度问题。解决了聚类过程中产生的过分割以及欠分割问题，保证了轨迹提取的完整性和独立性以及聚类中心的自动选择。在后期的轨迹提取中，利用回溯算法寻找最优解，将轨迹具有的平滑不变约束融入剪枝函数的设计，快速剪掉无关的杂波分支，加快了搜索解的速度。结合鲁棒的特征检测算子和回溯策略，解决可变尺度红外小目标的检测问题，提高算法的实时性和鲁棒性。

Description

基于快速自动密度聚类的可变尺度红外小目标的检测方法

技术领域

本发明涉及图像分析与图像理解，尤其是涉及一种基于快速自动密度聚类的可变尺度红外小目标的检测方法。

背景技术

红外小目标的检测、跟踪、识别等系列问题，主要来自于红外搜索与跟踪系统(Infrared Search and Track,IRST)。其主要任务是对天空、海面等复杂场景下的飞机、船只等感兴趣目标进行有效地检测、识别以及跟踪。由于其作用距离较远，目标在图像平面成像一般只有几个像素大小，缺乏形状、纹理等可资利用的信息，同时，受到天空云层，海面杂波、反光带的影响，处于低信杂比(Signal-to-Clutter Ratio,SCR)的环境下，属于典型的弱小目标，长期以来，一直是目标自动识别领域(Automatic Target Recognition，ATR)的难点与焦点问题。

基于序列图像的边检测边跟踪方法(Track Before Detection,TBD)在早期成为该问题的主要技术框架。其将检测软判决(soft decision)推迟到稳定地跟踪目标一段时间之后，适用于目标运动轨迹具有时空一致性、相关性的低信杂比目标。围绕该技术路线，一些经典的算法相继被提出，如：动态规划(dynamic programming,DP)、多级假设检验(multiple–stage hypothesis testing,MHT)、高阶相关(high-order correlation)等，但是该类方法均需要构建初始的轨迹搜索空间，之后借助于数据融合(data association)利用目标的运动信息完成轨迹检测，伴随图像帧数的增加，计算复杂度呈指数级增长。与TBD方法不同，后期提出的基于单帧图像的先检测后跟踪(Detection before Track,DBT)方法，检测思路是来一帧检测一帧，具有更好的实时性。但是，该研究方法的有效性强度依赖于提取的目标特征的稳定性。近些年来，研究者在不同的角度(时域、空域)，不同变换域(频率域、小波变换域)等尝试挖掘更为鲁棒、有效的目标特征，诸如局部小区域的灰度差，对比度等等。然而此类特征对环境变化和噪声非常敏感，而且不满足尺度不变(scale-invariant)特性，即：目标的大小会伴随成像环境改变等。从而无法得到稳定的检测结果，算法具有局限性。

在该类目标的自动识别任务中，如何挖掘鲁棒的目标特征，使其满足尺度不变特性，适应成像距离以及光照等外部环境的影响是核心工作。同时，借助于序列图像处理低信杂比目标的优势，并且克服其较高的计算代价，设计出低复杂度的数据融合(dataassociation)方法也是亟待解决的问题。

作为计算机视觉领域中的一种局部特征检测算子以及描述算子(speeded-uprobust feature，SURF)是尺度不变特征转换算子(scale-invariant feature transform，SIFT)的加速版本。由于其具有尺度不变特性，经常被用于人脸等有形目标的识别、图像配准以及分类等任务，利用SURF完成感兴趣点的提取(interest point extraction)是进行匹配的前期工作，其与本发明的红外小目标的点状特征匹配。

聚类分析(clustering analysis)是无监督分类的主要方法，用于对未标注数据集进行分类。Alex等人在文献“Clustering by fast search and find of densitypeaks.In SCIENCE，2014”中提出一种快速的密度聚类方法，区别于已有的DBSCAN，K-means等方法，该方法只以数据集中各点的距离矩阵为输入参数完成检测，适用于非球形聚类，同时由于简洁、快速，成为非常实用的聚类方法。

分析小目标轨迹的特点，由于其作用距离远，运动缓慢，与随机性更强的噪声相比，其轨迹片段在图像平面会呈现一种局部密度聚集特性。这启发我们可以利用聚类分析对这些数据点进行自动聚类，直接从空域完成轨迹的提取。由于在二维图像平面进行操作，时间复杂度可以控制在多项式级别，避免使用传统数据融合方法而带来的计算代价过高问题。

密度聚类算法的主要缺陷是无法进行聚类中心的自动确定(尚需手工选取)以及无法克服聚类分割产生的两类错误问题：过分割(over-segmentation)、欠分割(under-segmentation)。应用到我们的数据集，会产生两种不理想的聚类结果，即：一条完整轨迹被聚集到不同的聚类(cluster)；多条轨迹被分配到一个聚类中。这两种错误聚类，将导致低检测率(detection probability)和高误检率(False detection)。

发明内容

本发明的目的在于针对上述方法的不足，提供一种结合鲁棒的特征检测算子和回溯策略，解决可变尺度红外小目标的检测问题，提高算法的实时性和鲁棒性的基于快速自动密度聚类的可变尺度红外小目标的检测方法。

本发明包括如下步骤：

1)采集真实云天或者海天背景，构建仿真小目标，利用matlab中的imnoise函数叠加高斯噪声，合成仿真图像序列F(t)(t＝1,...N)作为图像数据,其中，N表示图像序列的长度；

2)利用SURF加速鲁棒特征检测算子提取图像序列F(t)中每一帧图像的特征点；

3)筛选所有特征点中最强的K个点；

4)保存标记得到的所有候选目标点迹的图像坐标，作为后续聚类的数据集；

5)计算数据集合中每个点i与其它各点的距离，生成距离矩阵dist，距离函数选取欧式距离；

6)确定自适应参数S(S代表每个数据点平均邻居的个数)；

给出S的选取原则为：假定S等于截断距离dc，当S为序列长度N×帧间运动位移frame_off×1/2时，可以得到完整的轨迹片段。

7)由步骤6)中给出的S，确定截断距离dc；

8)通过高斯核计算每个数据点i的局部密度ρ_i以及δ_i，并计算得到数组γ_i＝ρ_i×δ_i；

9)针对每个γ_i，若γ_i＞TH，(此处参数TH选取190～200之间，可以满足轨迹的独立性)则确定为聚类中心；

10)将所有的数据点聚集到与自己最近的聚类中心中；

11)剔除离散点；

12)保存得到的聚类，标记其数量；

13)设置轨迹长度门限；

14)对于得到的所有聚类，逐个依据下述步骤进行处理；

15)取出第i个聚类中的一个点A，判断该点是否是聚类中倒数第1个点，若是，则处理下一个聚类，若不是，则继续下列步骤；

16)取出聚类中的下一个点B；

17)生成直线C：A->B，将两点存入直线集合E；

18)取出聚类中下一个点K，若线段AB和BK的夹角θ满足平滑性约束(smooth-keptconstraint),即：0≤θ≤π/2，则将其归入直线集合E，否则转步骤19)；

19)判断K是否是该聚类中的最后一个点，若是，则判断直线集合E中点的个数是否大于轨迹最小长度，若小于，则转步骤16)；若大于，则转步骤20)；若不是该聚类的最后一个点，则返回到步骤18)；

20)输出直线集合E，转步骤15)；

21)处理完所有聚类后，将得到的所有轨迹片段进行标记保存。

本发明具有如下优点：

1.本发明利用SURF算子解决尺度可变小目标的特征提取。由于算子具有优良的尺度不变特性，用于提取候选感兴趣的目标点，较好地解决了目前众多方法中特征对于环境的敏感问题，为方法的鲁棒性奠定了基础。

2.快速自动化密度聚类算法，针对问题域，直接从空域完成轨迹的分割、提取，克服了原有序列检测问题中数据融合带来的指数级复杂度问题。同时，通过设计自适应参数，解决了聚类过程中产生的过分割以及欠分割问题，保证了轨迹提取的完整性和独立性以及聚类中心的自动选择。

3.在后期的轨迹提取中，利用回溯算法寻找最优解，将轨迹具有的平滑不变约束融入剪枝函数的设计，快速剪掉无关的杂波分支，加快了搜索解的速度。

附图说明

图1是使用本发明与其它4种方法的比对结果。

在图1中，(a)原始图像序列中的一帧；(b)max-median；(c)top-hat；(d)背景剔除；(e)局部对比度映射；(f)本发明。

图2是本发明与其它3种方法的ROC曲线(真实序列1)。

图3是本发明与其它3种方法的ROC曲线(真实序列2)。

图4是本发明与其它3种方法的ROC曲线(真实序列3)。

具体实施方式

本发明的实施步骤主要包括：基于SURF的感兴趣点提取，快速自动密度聚类以及回溯算法。

第一步：合成仿真图像序列F(t),(t＝1,…N)，参照文献“Infrared patch-imagemodel for small target detection in a single image.IEEE Trans.On Imageprocessing,2013”提供的合成方法，其中，背景为真实采集的云天或者海天背景，嵌入的仿真小目标尺寸控制在3×3像素至10×10像素之间，噪声方差控制在(0.01，0.03)范围内。

第二步：对于仿真图像序列F(t),(t＝1,...N)，逐帧进行基于SURF的感兴趣点提取。

(1)将待检测的灰度图序列F(t)，保存在文件夹F中；

(2)用matlab的fopen函数在步骤(1)中的文件夹F建立一个名为‘图像坐标’的csv文件；

(3)用matlab的imread函数读入步骤(1)中的文件夹F中的图像F(T)(t的初始值为1)；

(4)若t＝1，则用matlab的imshow函数显示该图像，记为img。若t不等于1，则不做任何处理；

(5)用matlab的detectSURFFeatures函数寻找图像的特征点，存入结构变量points中；

(6)用matlab的selectStrongest函数选出特征点points中最强的30个点(结合时间效率以及杂波点的干扰，一般选取20～30为宜)，存入结构变量strongpoints，利用该结构变量的位置属性得到strongpoints中所有点在图像F(t)的位置坐标(x,y)；

(7)用matlab的plot函数将strongpoints中的点标记在步骤(4)中显示的图像img上；

(8)用matlab的fprintf函数将点的位置坐标写入图像坐标.csv文件中；

(9)若i等于N，则进入下一步。若i不等于N，i自增1，则返回到步骤(3)；

(10)用matlab的fclose函数关闭图像坐标.csv，得到图像序列F(t)所有特征点的坐标信息；

(11)用matlab的saves函数保存被标记后的图像img。

第三步：执行快速自动密度聚类。

(1)将存有特征点坐标信息的‘图像坐标.csv’文件保存在文件夹F中；

(2)用matlab的csvread函数将步骤(1)中的‘图像坐标.csv’文件载入到二维数组data中；

(3)计算步骤(2)中数组data所有特征点到其他点的距离，保存在数组dist中；

(4)确定截断距离dc。用matlab的sort函数对步骤(3)得到的距离进行升序排列。选取升序序列中下标为position的数组元素作为截断距离dc，S＝5％。

position＝1/2×N×(N-1)×S (1)

(5)通过高斯核计算每个数据点i的局部密度ρ_i，保存在数组rho，计算公式如式(2)所示：

其中，U代表所有的数据点构成的集合，d_ij表示数据点i,j按照欧氏距离计算得到的距离。

(6)生成数组delta(δ)和数组nneigh。其中数组delta(i)的含义是：若rho(i)(数组rho在步骤(5)中得到)具有最大局部密度时，则delta(i)表示与第i个点距离最大的数据点与第i个点的距离；否则，delta(i)表示在所有局部密度大于rho(i)的数据点中，与第i个点距离最小的那个数据点与第i个点的距离。数组nneigh(i)记录比rho(i)更大的数据点中与第i个点距离最小的数据点的编号；

(7)根据公式(4)，计算得到数组γ_i，并用matlab提供的sort函数进行降序排序，记作γ1_i；

γ_i＝ρ_i×δ_i (4)

(8)统计聚类个数，若γ1_i＞TH(TH＝200)，则聚类个数NCLUST自增1；

(9)对非聚类中心数据点进行归类。对非聚类中心数据点进行归类。遍历所有的数据点，若第i个数据点未被归类，则根据步骤(6)得到的数组nneigh，将第i个数据点归入下标为neigh(i)的聚类中。

(10)去除离散点，去除聚类中局部密度小于1的点以及尺寸过小的聚类(结合允许的漏检率，此处门限取10)，标记得到的最后聚类的数量是NCLUST_final。

(11)遍历所有数据点，将点的坐标按聚类分别存储在文件中，并使用不同颜色将其标注在图中，用matlab的saveas函数保存。

第四步：采用回溯算法进行杂波的剪枝以及轨迹的提取。

(1)在文件夹F中有文件夹data，文件夹data中的每个data_i.csv文件保存着(i＝1,2,..,NCLUST_final；NCLUST_final是聚类的个数)各个聚类中点的坐标。

(2)设置直线长度阈值Trace_Length＝8。

(3)用matlab的dir函数获得文件夹data下的所有’.csv’文件信息Files。

(4)对步骤(3)得到的Files(i)(i＝1,2,3,...,LengthFiles)指向的文件做以下操作。

(5)将文件中的点的坐标读入到数组points，计算点的个数num。建立1×num的标记数组surPoint，若surPoint(1,i)(i＝1,2,3,...,num)的值为1，则表示第i个点在所求直线上；若surPoint(1,i)的值为0，则表示第i个点不在所求直线上。

(6)取出聚类中第j(j＝1,2,3,...,num；j初始值为1)个点A。判断该点是否是聚类中的最后一个点，若是，处理下一个聚类，若不是，继续下列步骤。

(7)取出聚类中的下一个点B。

(8)由A和B计算直线的斜率与截距。若两点的纵坐标之差的绝对值小于pp(pp设为0.05)，则斜率设为0，截距设为A点的纵坐标；否则，若两点横坐标之差的绝对值小于zz(zz设为0.05)，则斜率设为9999999，截距设为A点的横坐标；否则，通过两点式求出斜率kk和截距vv。

(9)根据两点的纵坐标构建平滑不变约束标志yG，若B点的纵坐标大于A点的纵坐标，yG＝1；否则yG＝2。用步骤(5)得到的surPoint数组标记该两点在直线集E中。

(10)取出聚类中的下一个点K。

(11)判断K点是否满足步骤(9)得到的平滑不变约束：若K点的纵坐标大于B点的纵坐标，则平滑约束标记等于1，否则等于2。若该平滑约束标记与步骤(9)中的yG相等，则说明K点满足平滑不变约束，进行下一步。若不满足直接跳到步骤(13)。

(12)判断K点是否在步骤(8)得到的直线上：设k2为斜率的平方再加1的开方。设li是根据斜率kk，截距vv，K点的横坐标得到的纵坐标。若li/k2<10，则K点在由斜率kk，截距vv得到直线上，否则不在。若K点在直线上，在surPoint数组中标记K点属于直线集E，否则跳到步骤(13)。

(13)判断K点是否是聚类中的最后一个点：若是，判断标记在直线集内点的个数是否大于Trace_Length若大于，将该直线集内标记的点的坐标保存在单独的一个文件中。若小于则返回步骤(7)；若K点不是聚类中的最后一个点，返回到步骤(10)。

(14)j自增1，返回到步骤(6)。

(15)处理完所有聚类后，按照保存的直线集用不同颜色标记每条直线上的点，保存标记图。

本发明的优点可以通过以下仿真实验进一步证明。

1、仿真条件：CPU InterCeleron N2930；主频1.83GHZ；内存4G；操作系统win10；仿真平台Matlab2013b。

2、比对方法：

(1)、经典方法

Top-Hat(参考文献“Morphology-based algorithm for point targetdetection in infrared backgrounds,Proc.SPIE,Oct.1993”)。

Max-mean/Max-median(参考文献“Max-mean and max-median filters fordetection of small targets,Proceedings of the SPIE 3809，1999”)。

(2)、最新方法：IPI、Local Contrast Map以及Background elimination.

IPI(参考文献“Infrared patch-image model for small target detection ina single image.IEEE Trans.On Image processing,2013”)。

局部对比度映射(参考文献“A local contrast method for small infraredtarget detection.IEEE Trans.on Geo-science and remote sensing,2014”)。

背景剔除(参考文献“small target detection and tracking based on thebackground elimination and kalman filter.2015International symposium onartificial intelligence and signal processing”)。

3、实验内容及结果分析

实验1：与top-hat、max-median以及局部对比度映射、背景抑制比较实验(图1所示)。

图1给出了4个图像序列，分别是2个海天背景，2个云天背景。从图1可以看出，图像背景存在大量的云层边缘、结构化的目标、海天分界线等，属于典型的低信杂比环境。对于这种复杂背景下的目标检测，依赖于局部灰度、对比度等特征的局部对比度映射以及背景剔除方法由于无法提取稳定、突出的目标特征，导致算法性能很快衰减。比较而言，基于形态学的Top-hat以及max-median优于前述的两种方法，但是检测结果出现了大量的虚警。与前述这几个方法相比，本发明对于背景的适应性更好，均能够得到令人满意的检测结果，并且具有更少的虚警。

实验2：与性能较为稳定的IPI算法的比较结果(表1所示)。为了进一步展示本发明的性能，给出下述的度量指标：检测率Pd,虚警False alarm以及漏检missing detection。实验序列是8个海天，5个云天背景下的仿真序列。从表1的比较结果可以看出，对于大部分的实验序列，两种方法都给出了满意的检测率。但是从整体上看，算法IPI对于多个序列都产生了不同程度的虚警和漏检。本发明虽然也针对海天序列-2、4、7产生了虚警和漏检，但是明显低于IPI。因此，从性能上看，本发明更能适应不同背景下的目标检测，性能更为稳定。

本发明与IPI方法的比较结果参见表1。

表1

实验3：接收者操作特征曲线(Receiver Operating Characteristic,ROC)是公认的信号检测理论中测度算法性能的分析工具。为此针对三个真实的序列图像，图2～4给出了max-mean、top-hat、IPI与本发明的比较结果。对于序列1，只有IPI方法比本发明稍快一点达到检测率Pd为1，但是本发明具有更低的虚警。对于序列2和3，本发明明显优于其他算法，可以在同样的虚警下得到理想的检测率。

Claims (1)

1.基于快速自动密度聚类的可变尺度红外小目标的检测方法，其特征在于包括如下步骤：

1)采集真实云天或者海天背景，构建仿真小目标，利用matlab中的imnoise函数叠加高斯噪声，合成仿真图像序列F(t)作为图像数据,其中，t＝1,...N，N表示图像序列的长度；

2)利用SURF加速鲁棒特征检测算子提取图像序列F(t)中每一帧图像的特征点；

3)筛选所有特征点中最强的K个点；

4)保存标记得到的所有候选目标点迹的图像坐标，作为后续聚类的数据集；

5)计算数据集合中每个点i与其它各点的距离，生成距离矩阵dist，距离函数选取欧式距离；

6)确定自适应参数S，S代表每个数据点平均邻居的个数；

给出S的选取原则为：假定S等于截断距离dc，当S为序列长度N×帧间运动位移frame_off×1/2时，可以得到完整的轨迹片段；

7)由步骤6)中给出的S，确定截断距离dc；

8)通过高斯核计算每个数据点i的局部密度ρ_i以及δ_i，并计算得到数组γ_i＝ρ_i×δ_i；

9)针对每个γ_i，若γ_i＞TH，则确定为聚类中心，其中，TH选取190～200之间，可以满足轨迹的独立性；

10)将所有的数据点聚集到与自己最近的聚类中心中；

11)剔除离散点；

12)保存得到的聚类，标记其数量；

13)设置轨迹长度门限；

14)对于得到的所有聚类，逐个依据下述步骤进行处理；

15)取出第i个聚类中的一个点A，判断该点是否是聚类中倒数第1个点，若是，则处理下一个聚类，若不是，则继续下列步骤；

16)取出聚类中的下一个点B；

17)生成直线C：A->B，将两点存入直线集合E；

18)取出聚类中下一个点K，若线段AB和BK的夹角θ满足平滑性约束,即：0≤θ≤π/2，则将其归入直线集合E，否则转步骤19)；

19)判断K是否是该聚类中的最后一个点，若是，则判断直线集合E中点的个数是否大于轨迹最小长度，若小于，则转步骤16)；若大于，则转步骤20)；若不是该聚类的最后一个点，则返回到步骤18)；

20)输出直线集合E，转步骤15)；

21)处理完所有聚类后，将得到的所有轨迹片段进行标记保存。