基于Hough变换和优化K均值聚类的反舰导弹编队识别目标选择方法

2023-03-21黄隽吴鹏飞李晓宝刘玥

兵工学报 2023年2期

黄隽，吴鹏飞，李晓宝，刘玥

(海军工程大学兵器工程学院，湖北武汉 430033)

0 引言

反舰导弹的目标选择通常有特征识别和编队识别两类技术。特征识别依据指定目标的单一特征或者加权综合特征，但随着反舰导弹自控时间增长，末端编队态势变化超出特征识别的能力。编队识别是将导弹末制导雷达得到的编队态势，与目标指示阶段实时探测或者情报数据库的编队态势相匹配，排序并识别编队预定目标的过程。实现编队识别[1]主要有两类技术思路，现有方法侧重队形识别，而非目标识别。

一类编队识别方法根据队形相关参数在量测空间结合编队模板[2-7]或者规则[8-9]直接计算判断。如以阵群成员位置观测集合和先验队形模板的最大势匹配子集中心为基准点描述观测和模板队形，通过队形描述匹配识别观测阵群目标队形，该方法工程上可操作性较弱；根据队形模板与待识别编队各舰方位角的方位距离建立空间场景方向相似性判别函数进行队形判别，对不同队形区分度不高、看齐角对判断影响过大；用队列线数量和队列角建立编队线型队形模板，给出基准目标、其他目标相对队列方位和相似性量度的数学模型，判断目标是否位于同一队列线，聚类获得队列线数量，从而识别队形，相似性量度要求方位角偏差必须小于探测误差，对目标位于同一队列线的判断不够严谨。

另一类编队识别将目标从量测空间转换到参数空间参数[10-13]，结合编队模板进行匹配判断，如针对队形形状和群体成员分布规律建立特征模型和模板，基于模板滑动匹配识别队型，可识别线型队形，但聚类过程复杂，实时性较弱。现有K均值聚类数优化方法[14-15]代价函数采用D+L，D为类内差异度，L为类间差异度，L取聚类中心到全体样本均值的距离之和，优化效果较好，但存在关键聚类数段区分度不明显、对典型数据样本聚类出现多局部最优值和聚类优化数略高于理论值问题；或采用D/L，但具有递减特性。共性问题是未考虑检测区间对聚类效果的影响，聚类数优化和聚类迭代过程相互独立。

本文基于Hough变换和优化K均值聚类算法，实现对目标指示阶段之后未制导阶段之前队列要素发生变化的典型编队目标识别与选择。

1 编队识别目标选择流程

基于Hough变换和优化K均值聚类的反舰导弹编队识别目标选择流程可以概括为目标指示阶段进行目标预选和末制导阶段目标确认两个阶段任务。

从舰艇雷达目标态势下达到导弹发射阶段，完成1.1节～1.3节三个步骤，预选打击目标。从导弹发射到导弹结束自控飞行段进入末制导雷达开机搜索跟踪阶段，编队目标发生旋转、缩放，或进行干扰。导弹末制导雷达开机后，重复1.1～1.3节，匹配目标指示态势，从而确认下达目标。

1.1 编队目标Hough变换映射

获取目标队形成员量测空间坐标，Hough变换映射到参数空间。Hough变换检测直线的实质是建立量测空间A={R:R=(X,Z)}点与参数空间B={r:r=(θ,ρ)}曲线的对应关系，将量测空间直线检测问题转化为参数空间曲线交点聚类问题。其中，X、Z为直角坐标系横坐标和纵坐标，θ为量测空间原点到过该点的任意一条直线的垂线与X轴正向所成的夹角，正常检测区间[12,16]取[0 rad,π rad)，|ρ|为垂线长度，ρ可正可负。

对量测空间任意点R=(X,Z)，可移动检测区间的Hough变换函数为

(1)

式中：ρ表达式为量测空间任意点转换到参数空间的基本公式，3个约束分别描述了正常检测区间、移动检测区间范围、正常检测区间到移动检测区间的极坐标转换公式；(θ0,ρ0)为需要移动检测区间的Hough曲线交点极坐标；Δθ为分辨率，当存在Hough曲线交点分布在[0 rad,5·Δθ)∪[π rad-5·Δθ，π rad)时，采用移动检测区间。

1.2 优化K均值聚类检测队列线

K均值聚类是将数据点到类中心的某种距离(常用欧氏距离)之平方和作为目标函数，通过迭代运算使得目标函数值最小，从而得到聚类中心的最优分类方法。假设将数据集合X划分为k类，利用式(2)定义距离误差平方和最小准则进行迭代聚类和类中心更新：

(2)

式中：E为目标函数；k为聚类数；ra为类Wj的任一样本；cj为类Wj的聚类中心；nj为类Wj中数据样本数目。

采用式(3)计算最优聚类数k*值，使得类内差异度D最小、类间差异度L最大：

(3)

式中：kmax为聚类数上限[15]；⎣·」表示向下取整；F(D,L)为差异度代价函数。

优化K均值聚类算法流程如下：

步骤1启动优化聚类数k迭代，k=1,2,…,kmax。

步骤2针对一定聚类数k，启动聚类迭代，rep=1,2,…,n，从样本集X按均匀分布规律随机选取k个样本作为初始聚类中心，由式(4)计算各样本到各初始聚类中心的目标函数矩阵：

D=(dij)n×k，dij=‖ri-cj‖2,ri∈X

(4)

式中：dij表示矩阵D的第i+(j-1)×n个元素。由式(5)确定各样本到各聚类中心最小距离平方和组成的列向量d及对应聚类号列向量id，进行初始分类，并计算各类样本数量。

d=[minD(1,:),…,minD(i,:),…,minD(n,:)]′id=[id1,…,idi,…,idn]′

(5)

式中：min(·)表示取向量中最小的元素。然后展开聚类迭代，贯序分为多样本更新聚类和单样本更新聚类两个阶段。

步骤3多样本更新聚类阶段。

①根据聚类号列向量id和样本发生增减的聚类号dc(见式(8))分批更新变化的聚类中心和变化样本数量；然后由变化的聚类中心分批计算变化的目标函数矩阵D。

②处理变为空集的聚类。

③由式(2)计算目标函数值，判断目标函数值是否减小。若未减小，则后撤一步，返回该迭代步的聚类中心并计算各类样本数量，转入单样本更新聚类阶段；若减少，则暂存聚类号列向量id和目标函数值。

④由式(5)同时确定各样本到各聚类中心最小距离平方和组成的列向量d及对应聚类号列向量id。

⑤由式(6)，比较新旧聚类号列向量id、id′是否变化，确定哪些样本改变聚类，得移动列向量dm(聚类号变化的样本号)：

(6)

式中：[]表示没有元素，降低列向量行数；dmi表示列向量dm的第i个元素。

如果有所改变则由式(7)比较样本到新旧聚类中心距离平方和，确定是否减小，改变移动列向量dm：

(7)

然后用新聚类号取代样本的当前聚类号，并获得样本发生增减的聚类号向量dc：

(8)

式中：dci表示行向量dc的第i个元素。式(8)表示元素i如果属于行向量[id(dm);id′(dm)]的元素，则该元素在行向量dc保留，否则该元素在行向量dc不保留。

步骤4单样本更新聚类阶段。

①针对聚类号列向量id发生变化的样本，由式(9)计算具有单样本聚类(至少存在一个样本数量为1个的聚类)样本集合的目标函数矩阵Del：

(9)

式中：mj为j号聚类的样本数量，j∈dc表示元素j属于行向量dc。

②暂存聚类号列向量id和目标函数值，由式(5)同时确定Del各样本到各聚类中心最小距离平方和组成的列向量del及对应的聚类号列向量id。

③由式(6)比较新旧聚类号列向量id、id′是否变化，确定哪些样本改变聚类，如果有所改变则由式(7)比较样本到新旧聚类中心距离平方和，确定是否减小；若没有改变则判断单样本更新迭代收敛。

④按照式(10)得到单元素移动列向量dm。如果dm小于上一次迭代的单元素移动d′m，则继续迭代。

dm=min (dm)

(10)

计算移动次数，并把新dm赋值给上一次迭代的单元素d′m，由式(11)暂存新旧移动聚类号单元素列向量mid、mid′。

mid′=id(dm),mid=mid(dm)

(11)

⑤由式(12)更新目标函数值E、聚类号列向量id。

E=E+Del(i,j)-Del(i,j′),id(dm)=midi∈dm,j∈mid,j′∈mid′

(12)

式中：i∈dm表示元素i属于列向量dm；j∈mid表示元素j属于列向量mid；j′∈mid′表示元素j′属于列向量mid′。

由式(13)更新聚类的样本数量mj，式(14)更新聚类中心cj：

(13)

(14)

由式(15)获得样本发生增减的2元素聚类号向量dc：

(15)

重复步骤4，直至单样本更新收敛。

步骤5每一次聚类迭代后判断各类样本是否为空集，计算各样本到各聚类中心的目标函数矩阵D=(dij)n×k，根据式(16)由已知聚类号列向量id，得到目标函数矩阵中各样本到各聚类中心最小距离平方和组成的列向量d：

d=[D(1+(id1-1)·n),…,D(i+(idi-1)·n), …,D(n+(idn-1)·n)]′

(16)

式中：idi表示样本所属的聚类号1,2,…,k；D(i+(idi-1)·n)表示目标函数矩阵D中第i+(idi-1)·n个元素。

(17)

步骤6按照式(3)计算聚类数k的差异度代价函数，重复步骤1～步骤5，遍历所有聚类数k=1,2,…,kmax。比较聚类数k迭代差异度代价函数，得到最优聚类数、聚类中心及对应分类。

1.3 反Hough变换编队识别目标排序

将聚类中心反映射回原量测空间对应直线。当在雷达态势图中判断各点是否在一条直线上时，就等价于在对应参数空间内找到Hough曲线交点的聚类。若设聚类中心为(θ,ρ)，则由式(18)可得参数空间的点到原量测空间的反映射表达式为

(18)

式中：θ∈(0 rad,π/2 rad)∪(π/2 rad，π rad)对应的一般反映射为直线z=ax+b，直线斜率a和截距b的表达式为该直线的约束；θ=0 rad、θ=π rad对应的反映射为平行于x轴的直线，θ=π/2 rad对应的反映射为平行于z轴的直线；x、z分别为量测空间的横轴、纵轴坐标；ρ0、ρπ和ρπ/2分别表示θ=0 rad、π rad和π/2 rad对应的ρ。

根据反映射直线识别出编队队形轮廓，结合各舰在编队中的运动信息，判断出编队组成和基准舰，并进行编队目标排序。

2 编队目标生成和选择模型

如图1所示，编队目标建模时以目标指示阶段的队列航向为z轴正方向，垂直z轴右手方向为x轴正方向，建立平面直角坐标系。编队队形要素包括：基准舰J是指为队列提供相对位置参照的舰艇，通常位于队列中较特殊的几何位置处。队列线是指由基准舰开始，队列中各舰所在点的连线。队列航向γ是指基准舰艏艉线与正北方向N间的夹角。编队间距Df可分解为沿垂直基准舰艏艉线方向的横向间距Dx和沿基准舰艏艉线方向的纵向间距Dz。看齐角α是指其他各舰艏艉线与队列线之间夹角。队列角β=180°-α是指基准舰艏艉线与队列线之间夹角。

图1 V形编队及平行编队示意图Fig.1 V-shaped and parallel warship formations

为构建编队模型，假设如下：1)编队中不同舰艇的距离方向可被传感器分辨；2)编队中舰艇位置定位误差和传感器探测误差符合均值为0的正态分布；3)末制导阶段并未改变目标指示阶段的编队类型，或者由一种典型线型编队队形变化成另一种典型线型编队队形，编队目标未发生位置交换，即以基准舰为参考的排序不发生变化；4)忽略海况、气象、网络电磁空间环境条件对导弹的影响，编队无任何硬抗击行动；5)目标编队队形出现变化类型、变化程度，均与导弹自控飞行时间无关；6)Δθ=0.1 rad；7)若相对参数|θ12|=|θ1-θ2|与无噪声、无干扰的理想条件下的相对参数偏差不大于30%[12]，则认为队形识别正确；若相对参数偏差大于等于10%，则启动移动检测区间检测，最终选择相对参数较小的检测区间；8)若末制导阶段确认目标与目标指示阶段预选目标一致，则对编队目标的识别选择正确。

常见线型编队队形[17-20]有V形、平行和环形编队等。

2.1 V形编队模型

2.1.1 编队目标模拟

如图1(a)所示，以V形编队两翼连接处舰艇为基准舰，设队列角为β，队列航向γ=0°，可得各舰的量测空间坐标为(xi,zi)(i=1，2，3，…，T)，V形编队左、右翼各舰量测空间坐标如式(19)、式(20)所示：

(19)

(20)

2.1.2 编队识别目标排序

(21)

目标指示阶段，识别基准舰编号为B，假设预设目标编号为P，则预设目标相对编号为P-B。末制导阶段，识别基准舰编号为B、确认目标相对编号为P-B，则确认目标编号为B+(P-B)=P。

2.2 平行编队模型

2.2.1 编队目标模拟

平行编队可分为双横队和双纵队。几何特征如图1(b)所示，双横编队前后队列线、双纵编队左右队列线相互平行。

平行编队通常以指定队列与方向的排头舰艇为基准舰。当设队列航向γ=0°，看齐角α始终为90°时，双横队内舰艇总数为T，各舰量测空间坐标为(xi,zi)(i=1，2，3，…，T)：

(22)

当设队列航向γ=0°、看齐角α始终为0°时，双纵队内舰艇总数为T，各舰量测空间坐标为(xi,zi)(i=1，2，3，…，T)：

(23)

2.2.2 编队识别目标排序

2.3 环形编队模型

2.3.1 编队目标模拟

环形编队通常以最中间的舰艇为基准舰。几何特征如图2(a)所示，水平队列线与垂直队列线相互垂直。图2中，C1～C4为假目标，R为假目标与舰艇间距，γm为导弹来袭方向，假目标与舰艇的连线与导弹来袭方向呈固定角度，记为ε、δ。

图2 环形编队及假目标菱形布阵模式示意图Fig.2 Ring warship formation and rhombus array of fake targets

设队列航向γ=0°，编队内舰艇总数为T，可得各舰在量测空间中的坐标为(xi,zi)(i=1，2，3，…，T)，横队如式(24)所示，纵队分两段如式(25)所示：

(24)

(25)

2.3.2 编队识别目标排序

(26)

3 编队目标变化模型

编队在战斗中，会根据战术需求与战场环境，改变队列要素，采取旋转队形、缩放队形、进行冲淡式干扰等措施。

3.1 旋转模型

整个编队绕某点旋转，调整队列航向。若旋转中心点坐标为(x0,z0)，逆时针旋转角为q，旋转前舰艇坐标为(xi,zi)，则旋转结束后坐标(x′i,z′i)为

x′i=(xi-x0)·cosq-(zi-z0)·sinq+x0z′i=(xi-x0)·sinq+(zi-z0)·cosq+z0

(27)

3.2 缩放模型

缩放是编队航渡和作战过程中，在不破坏整体队形情况下，增大/缩小编队间距以应对战场环境变化的过程。设k1、k2分别为横、纵坐标缩放系数，缩放变化前后的编队间隔(Dx,zz)和(D′x,z′z)关系为

D′x=k1·Dx,z′z=k2·Dz

(28)

当k1、k2∈(0,1)时，队形间距缩小，k1、k2∈(1,∞)时队形间距放大；当k1=k2时队形间距不变。

3.3 冲淡式干扰模型

当探测到来袭导弹时，编队舰艇会在导弹末制导雷达开机前进行冲淡式干扰，降低被导弹捕获的概率。常用假目标空间布阵有正方形、矩形和菱形等模式。设基准舰坐标为(Jx,Jz)，导弹来袭方向为γm([-180°,0°]∪[0°,180°])，则4个假目标坐标为

(29)

(30)

式中：正方形布阵，ε=δ=45°；矩形布阵，ε=δ；菱形布阵，ε≠δ，如图2(b)所示。

3.4 编队队型变化模型

3.4.1 V形编队与环形编队之间变化

目标指示阶段如式(19)、式(20)所示V形编队，末制导阶段如式(24)、式(25)所示环形编队，可叠加式(27)所示旋转、式(28)所示缩放变化和式(29)、式(30)所示冲淡式干扰；或者目标指示阶段如式(24)、式(25)所示环形编队，末制导阶段如式(19)、式(20)所示V形编队，可叠加式(27)所示旋转、式(28)所示缩放变化。

目标指示阶段，识别基准舰编号为B，假设预设目标编号为P，则预设目标相对编号为P-B。末制导阶段，识别基准舰编号为B1、确认目标相对编号为P-B+B-B1，则确认目标编号为B1+(P-B+B-B1)=P。

3.4.2 双横队和双纵队之间变化

目标指示阶段如式(22)所示双横队，末制导阶段如式(23)所示双纵队，可叠加式(27)所示旋转、式(28)所示缩放变化；或者目标指示阶段如式(23)所示双纵队，末制导阶段如式(22)所示双横队，可叠加式(27)所示旋转、式(28)所示缩放变化。

4 仿真分析

针对编队位置均方差为0.3的情况，本文方法对队形的识别率均为100%，略优于文献[12]的效果。对采用旋转变化的V形编队目标的识别选择如案例1所示；对采用旋转变化与冲淡式干扰的环形编队目标进行识别选择如案例2所示；对目标指示阶段为典型V形编队、末制导阶段变化为环形编队目标的识别选择如案例3所示；对叠加不同位置噪声的V形编队目标、平行横队目标、采用旋转变化与冲淡式干扰的环形编队目标、典型V形编队变化为环形编队目标的识别选择如案例4所示。

4.1 案例1

想定预选攻击编队3号舰，目标指示阶段态势如图3所示，设γ=0°，Dx=0.5 nmile，Dz=1 nmile，α=26.6°。x轴、z轴坐标分别叠加均值为0、均方差0.1和0.3的正态分布噪声。为检验本文方法对旋转队形的识别效果，假定末制导阶段γ变为330°。

图3 目标指示阶段V形编队识别目标预选图Fig.3 V-shaped formation recognition and target pre-selection in target indication stage

如图4所示，分别用D/L方法[12]、已有L+D方法[14]、本文L+D方法进行K均值聚类。其中，D/L方法差异度代价函数体现一定递减特性，已有L+D方法在关键聚类数段(2～8)缓慢增长、区分度不明显，而本文L+D方法区分度明显，显示最优聚类数为2。左右翼参数P1(2.692 0,5)、P2(0.505 9,5)，|θ12|=0.953 9 rad，识别为V形编队，在目标指示态势中画出队列线，判断基准舰为5号舰。方形标注预选目标3号舰(见图3)。末制导阶段态势如图5所示，γ变为330°。

图4 各种聚类数优化方法比较Fig.4 Comparison of optimal methods for clustering number

聚类分析得左右翼参数P1(0.074 0,5)、P2(1.029 5,5)，|θ12|=0.955 5 rad，识别为V形编队，判断基准舰为5号舰，方形标注确认打击目标3号舰。

案例1仿真结果表明，该方法能较好地识别选择采用旋转变化的V形编队目标。

图5 末制导阶段V形编队识别目标选择图Fig.5 V-shaped formation recognition and target confirmation in terminal guidance stage

4.2 案例2

想定预选攻击编队3号舰，目标指示阶段，设γ=0°，Dx=1 nmile，Dz=1 nmile。坐标分别叠加均值0、均方差0.1的正态分布噪声；假定末制导阶段γ变为330°；为检验本文方法应对冲淡式干扰效果，设编队在3号舰附近以矩形布阵释放假目标。采用移动检测区间[5Δθ，π rad+5Δθ)，取Δθ=0.1 rad，进行优化K均值聚类，得最优聚类数为2，聚类结果如图6所示。横、纵队参数P1(1.541 4,5)、P2(3.091 8,5)，相对参数|θ12|=1.550 4≈π/2 rad，相对参数偏差1.2%，识别为环形编队，判定基准舰为3号舰，方形标注预选目标3号舰，如图7所示。末制导阶段态势如图8所示，γm=0°，而γ变为330°，且放置了4个假目标。

图6 移动检测区间环形编队目标指示态势聚类参数空间图Fig.6 Ring formation line detection on the moving detection interval in target indication stage

图7 移动检测区间环形编队目标指示态势识别图Fig.7 Ring formation recognition and target pre-selection on the moving detection interval in target indication stage

图8 环形编队末制导雷达态势队形识别图Fig.8 Ring formation recognition and target confirmation in terminal guidance stage

聚类分析后得到横、纵队参数P1(2.065 0,5)、P2(0.473 8,5)，相对参数|θ12|=1.591 2≈π/2 rad，相对参数偏差1.4%，识别为环形编队，判断基准舰为3号舰，方形标注确认打击目标3号舰。

图9 原检测区间环形编队目标指示态势聚类参数空间图Fig.9 Ring formation line detection on normal detection interval in target indication stage

如图9所示，目标指示阶段采用原检测区间时，聚类分析后得到横、纵队参数P1(1.240 6,5)、P2(2.994 4,5)，相对参数|θ12|=1.753 8 rad，相对参数偏差11.7%，能够识别为环形编队；如图10所示，虽然目标指示阶段预选3号舰，但由于判断8号舰为基准舰，导致末制导阶段识别3号舰错误。

图10 原检测区间环形编队目标指示态势识别图Fig.10 Ring formation recognition and target pre-selection on normal detection interval in target indication stage

案例2仿真结果表明，本文方法能较好地识别选择采用旋转变化与冲淡式干扰的环形编队目标。

4.3 案例3

想定预选攻击编队3号舰，目标指示阶段态势如图3所示，设γ=0°，Dx=0.5 nmile，Dz=1 nmile，α=26.6°。x、z坐标分别叠加均值0、均方差0.1和0.3的正态分布噪声。设末制导阶段态势如图8所示，γ变为330°，Dx=1 nmile，Dz=1 nmile，坐标分别叠加均方差0.1的正态分布噪声，编队进行冲淡式干扰，在3号舰附近以矩形布阵释放4个假目标，γm=0°，而γ变为330°。

如图3所示，在目标指示态势中画出队列线，识别为V形编队，判断基准舰为5号舰，方形标注预选目标3号舰。如图8所示，末制导阶段，识别为环形编队，判断基准舰为3号舰，方形标注确认打击目标3号舰。

案例3仿真结果表明，本文方法能较好地识别选择目标指示阶段为典型V形编队、末制导阶段变化为环形编队的目标。

4.4 位置误差影响

对V形编队目标位置叠加均值为0、不同均方差0.1～0.4范围内的噪声，对每种情况分别仿真1 000次，对编队目标的识别正确率如表1所示，从中可见随着位置均方差的增大，识别正确率下降，但保持在77%以上。

表1 V形编队位置误差对识别正确率影响Table 1 Effects of location error of V-shaped formation on recognition rate

对平行编队目标位置叠加均值0、均方差0.1～0.4范围内的噪声，对每种情况分别仿真1 000次，识别正确率如表2所示。由表2可见：低均方差时目标指示和末制导阶段的队列线无交点比例高，两阶段无交点时识别正确，两阶段都有交点时识别基本正确，但存在相当比例前阶段无交点、后阶段有交点，由于对基准舰的选择方式不一样而识别错误；均方差提高时，两阶段队列线有交点的比例提高，两阶段无交点和前阶段无交点、后阶段有交点的比例均降低，反而出现识别率上升的现象。对编队目标的识别正确率保持在70%以上。

表2 平行横队位置误差对识别正确率影响Table 2 Effects of location error ofparallel formation on recognition rate

对环形编队目标位置叠加均值为0、均方差0.1～0.3范围内的噪声，对每种情况分别仿真1 000次，识别正确率如表3所示，从中可见随着位置均方差增大，对编队目标的识别正确率迅速下降，但保持在60%以上。

表3 环形编队位置误差对识别正确率影响Table 3 Effects of location error ofring formation on recognition rate

目标指示阶段V形编队目标位置叠加均值为0、不同均方差0.1～0.3范围内的噪声，末制导阶段队形变为γ=330°的环形编队，3号舰附近以矩形布阵释放冲淡式干扰假目标，且目标位置叠加和目标指示阶段相同均方差的噪声，对每种情况分别仿真1 000次，识别正确率如表4所示，从中可见对编队目标的识别正确率介于末制导阶段保持V形编队叠加噪声和末制导阶段保持环形编队叠加噪声情形之间。