CN117036400A - 基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法 - Google Patents
基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法 Download PDFInfo
- Publication number
- CN117036400A CN117036400A CN202310977257.1A CN202310977257A CN117036400A CN 117036400 A CN117036400 A CN 117036400A CN 202310977257 A CN202310977257 A CN 202310977257A CN 117036400 A CN117036400 A CN 117036400A
- Authority
- CN
- China
- Prior art keywords
- target group
- measurement
- target
- matrix
- group
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 125
- 239000000203 mixture Substances 0.000 title claims description 11
- 238000005259 measurement Methods 0.000 claims abstract description 212
- 230000008569 process Effects 0.000 claims abstract description 33
- 230000007704 transition Effects 0.000 claims abstract description 25
- 238000004364 calculation method Methods 0.000 claims abstract description 17
- 238000012883 sequential measurement Methods 0.000 claims abstract description 5
- 238000012544 monitoring process Methods 0.000 claims abstract description 4
- 239000011159 matrix material Substances 0.000 claims description 133
- 238000004422 calculation algorithm Methods 0.000 claims description 32
- 238000005070 sampling Methods 0.000 claims description 17
- 238000012545 processing Methods 0.000 claims description 12
- 239000000654 additive Substances 0.000 claims description 11
- 238000013178 mathematical model Methods 0.000 claims description 11
- 230000000996 additive effect Effects 0.000 claims description 10
- 238000001914 filtration Methods 0.000 claims description 4
- 230000001133 acceleration Effects 0.000 claims description 3
- 238000006243 chemical reaction Methods 0.000 claims description 3
- 238000001514 detection method Methods 0.000 claims description 3
- 238000009825 accumulation Methods 0.000 claims description 2
- 230000036461 convulsion Effects 0.000 claims description 2
- 230000007123 defense Effects 0.000 claims description 2
- 230000008030 elimination Effects 0.000 claims description 2
- 238000003379 elimination reaction Methods 0.000 claims description 2
- 238000010606 normalization Methods 0.000 claims description 2
- 239000007788 liquid Substances 0.000 claims 2
- 230000006835 compression Effects 0.000 claims 1
- 238000007906 compression Methods 0.000 claims 1
- 238000004088 simulation Methods 0.000 description 13
- 230000006870 function Effects 0.000 description 5
- 230000008859 change Effects 0.000 description 4
- 230000006399 behavior Effects 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000006467 substitution reaction Methods 0.000 description 2
- 241000284156 Clerodendrum quadriloculare Species 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000003416 augmentation Effects 0.000 description 1
- 238000012733 comparative method Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000004927 fusion Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000002093 peripheral effect Effects 0.000 description 1
- 230000035485 pulse pressure Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000017105 transposition Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/20—Analysis of motion
- G06T7/246—Analysis of motion using feature-based methods, e.g. the tracking of corners or segments
- G06T7/251—Analysis of motion using feature-based methods, e.g. the tracking of corners or segments involving models
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
- G01S13/66—Radar-tracking systems; Analogous systems
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/36—Means for anti-jamming, e.g. ECCM, i.e. electronic counter-counter measures
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
- G01S7/414—Discriminating targets with respect to background clutter
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N7/00—Computing arrangements based on specific mathematical models
- G06N7/01—Probabilistic graphical models, e.g. probabilistic networks
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/20—Analysis of motion
- G06T7/277—Analysis of motion involving stochastic approaches, e.g. using Kalman filters
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
- G06V10/00—Arrangements for image or video recognition or understanding
- G06V10/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
- G06V10/762—Arrangements for image or video recognition or understanding using pattern recognition or machine learning using clustering, e.g. of similar faces in social networks
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Remote Sensing (AREA)
- Radar, Positioning & Navigation (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Computer Networks & Wireless Communication (AREA)
- Multimedia (AREA)
- Artificial Intelligence (AREA)
- Computing Systems (AREA)
- Software Systems (AREA)
- Evolutionary Computation (AREA)
- Medical Informatics (AREA)
- Pure & Applied Mathematics (AREA)
- Algebra (AREA)
- Computational Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Probability & Statistics with Applications (AREA)
- General Engineering & Computer Science (AREA)
- Mathematical Physics (AREA)
- General Health & Medical Sciences (AREA)
- Databases & Information Systems (AREA)
- Health & Medical Sciences (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明公开了一种基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法,主要解决现有多目标群跟踪方法复杂度高、实用性差的问题,其实现方案包括:为各目标群建立状态参数模型、状态转移模型、量测空间分布模型;根据状态转移模型预测场景中多目标群的状态参数,并对场景杂波进行一次去除;根据量测空间分布模型计算量测与各种量测源的边际关联概率,以二次去除场景杂波;将杂波去除后的量测和各目标群进行数据关联,对多目标群状态参数进行贯序量测更新,并重复其状态参数的预测和更新过程,直至目标群离开雷达监视区域。本发明计算复杂度低,对群目标外形估计准确,在复杂多目标群跟踪场景中实用性强,可用于对杂波环境中多个目标群的跟踪处理。
Description
技术领域
本发明属于数据处理技术领域,特别涉及一种目标跟踪方法,可用于对杂波环境中多个目标群的跟踪处理。
背景技术
随着多目标跟踪技术的发展,适用于更大规模目标数量的群目标跟踪技术已成为当前的研究热点。一般而言,群目标按照传感器量测的分辨程度可分为可分辨的群目标和不可分辨群目标,其中,不可分辨群目标是指目标群产生量测的分辨率过低,以至于在传感器上无法将所有量测明显的区分开来,在这一情况下,对不可分辨群目标的跟踪通常是将其视为扩展目标。在跟踪过程中,不但要对其运动学参数实现跟踪,还要跟踪其扩展外形的变化情况,因此需要对目标群的外形进行建模。
根据目标群外形的复杂程度,目前对目标群外形的建模可大致可分为:①简单外形建模;②随机矩阵RMM建模及其改进方案,如乘性误差均匀采样MEM建模、期望最大非先验扩展演化EMNPM建模;③随机超曲面建模及其改进方案,如随机椭圆超曲面建模、随机高斯过程建模;④混合模型建模。
其中,最为常用的是随机矩阵RMM的建模方式,该方法将目标群的扩展外形描述为椭圆形,并用一个未知的半正定矩阵表示,基于RMM的建模方式实现简单且对计算量的需求较小。但由于该建模是通过一个遗忘因子标量来描述外形的不确定性,因而无法估计出椭圆长短半轴的长度以及朝向,对于需要精确估计外形大小的目标群跟踪而言,这是致命的。所述基于RMM模型改进的乘性误差均匀采样MEM模型对目标群进行建模,同RMM模型一样,也是将目标群的外形建模为椭圆形,但又同时将椭圆长短半轴的长度以及朝向作为需要估计的参数一同参与跟踪,该方法虽然能够获得比传统RMM方法更为精确的目标群外形估计,但它在多目标群跟踪领域中的应用仍处在初级阶段,无法适应复杂多变的多目标群跟踪场景。
杨世姗等人于2020年在第23届国际信息融合会议上,提出了一种将MEM模型与传统的联合强度概率数据关联JIPDA算法相结合的算法JIPDA-MEM,该算法引入强度函数表示目标的存在概率,并利用JIPDA进行数据关联,初步解决了目标群数目时变场景中的多目标群跟踪问题。
贾亚福杜尔等人于2023年在《Remote Sensing》期刊上,提出了一种将MEM与传统的联合概率数据关联JPDA算法相结合的PAJPDA-MEM算法,该算法使用JPDA算法进行数据关联,同时基于主轴卡尔曼滤波器PAKF完成对多目标群的状态估计,有效解决了利用海面雷达点云量测的多目标群跟踪问题。
但是,上述两种方法由于采用传统的数据关联过程,因而计算量过于庞大,无法对场景中的多目标群进行实时跟踪。
针对上述方法的不足,杨世姗等人于2018年在第10届国际传感器阵列和多通道信号处理会议上,提出了一种将MEM模型和近似的线性联合概率数据关联LJPDA算法相结合的LJPDA-MEM算法,该算法使用多对多匹配算法近似量测与目标群的边缘关联概率,在线性计算复杂度条件下解决了多目标群的数据关联问题;
李秦雷等人于2021年在《IET Signal Processing》期刊上,提出了一种将MEM模型与传统的截断式联合概率数据关联TJPDA算法相结合的TJPDA-MEM算法,该算法利用截断近似的方法削减传统JPDA中的关联事件数量,从而有效简化了使用JPDA算法进行群目标数据关联的复杂程度。
上述两种改进方法采用近似的数据关联方法,虽说能够在一定程度上减少算法计算量。但是,由于近似关联的误差,将导致其在实际工程应用时出现不小的弊端,如考虑跟踪时多个目标群之间时常发生合并与***行为,此时,多个目标群的外形之间将存在较长时间的重叠,导致对多个目标参数估计的混淆,最终发生严重的误跟现象。
发明内容
本发明的目的在于针对上述现有技术的不足,提出一种基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法,以降低多目标群跟踪的计算复杂度,避免在多目标群发生重叠现象跟踪场景中出现的误跟踪,提高跟踪精确度。
为实现上述目的,本发明的技术问题包括如下步骤:
(1)采用外形运动学解耦的椭圆扩展外形模型,建立单个不可分辨目标群的状态参数模型;
(2)根据目标群在场景中的运动模式,建立基于马尔科夫状态空间模型的目标群状态转移模型;
(3)在以雷达为传感器的跟踪场景中,将以二维极坐标形式表达的雷达量测转化为直角坐标表示;
(4)针对目标群产生量测的特点,采用乘性误差均匀采样MEM模型对群目标所产生量测的空间分布进行建模,并获得目标群量测的似然函数;
(5)利用步骤(2)建立的群目标状态转移模型,使用标准卡尔曼滤波预测公式对场景中多个目标群的状态参数进行一步预测;
(6)根据多个目标群预测状态,对场景中的量测使用椭圆外形的门限处理,实现对场景中的杂波的一次剔除;
(7)根据步骤(4)中建立的目标群量测空间分布模型,对背景杂波的数目与目标产生量测的数目进行建模,并基于量测与目标的边际关联概率在一次杂波剔除的基础上实现对场景中杂波的二次剔除;
(8)对两次杂波去除后剩余的有效量测,结合步骤(5)中所得多目标群的状态预测结果,采用具有线性计算复杂度的基于混合高斯模型迭代最大化期望模糊聚类的数据关联算法实现有效量测和目标群的数据关联,并获得关联概率矩阵;
(9)根据数据关联结果,基于扩展式卡尔曼滤波器的更新公式和概率数据关联法则,对每个目标群的状态参数使用有效量测进行基于关联概率的贯序量测更新;
(10)基于贝叶斯框架,将当前所得的多目标群状态参数估计结果作为下一次多目标群估计的先验信息;
(11)重复步骤(5)~步骤(10),直至目标群离开雷达监视区域完成跟踪过程。
本发明与现有技术相比具有如下有益效果:
本发明实施例提供的基于混合高斯模型迭代最大化期望模糊聚类数据关联的多目标群跟踪方法通过对传统基于乘性误差均匀采样模型的多目标群跟踪算法中的数据关联过程进行改进,有效解决了现有方法对群目标外形估计不准、算法复杂度高以及在目标群存在合并或***行为场景中实用性不足的问题。其与现有技术相比具有以下优点:
第一,本发明由于采用乘性误差均匀采样模型对目标群进行建模,在对目标运动状态跟踪的同时考虑对目标群的外形变化的跟踪,对比现有基于RMM的多目标群跟踪方法,本发明能够获得更加精确的目标群外形参数估计。
第二,本发明由于采用线性计算复杂度的模糊聚类算法进行数据关联,对比现有基于传统数据关联和乘性误差均匀采样模型的多目标群跟踪方法,本发明降低了计算复杂度,极大提升计算效率。
第三,本发明由于采用基于混合高斯模型迭代最大化期望模糊聚类进行数据关联,对比现有基于近似边缘概率数据关联和乘性误差均匀采样模型的多目标群跟踪方法,本发明提高了数据关联的精确度,避免了在发生合并或***行为场景中出现的多目标群误跟问题。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1是本发明的实现流程图;
图2是多群目标跟踪场景图;
图3是在图2场景下,分别用本发明与现有方法对场景中的多目标群进行跟踪仿真的目标群外形估计精确度比较图;
图4是在图2场景下,分别用本发明与现有方法对场景中的多目标群进行跟踪仿真的数据关联精确度比较图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
参照图1,本发明实施例提供一种基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法,包括如下步骤:
步骤1,在二维笛卡尔坐标系下,对单个不可分辨目标群的状态参数进行建模。
现有方法对单个不可分辨目标群的状态参数建模包括使用简单外形模型、星凸外形模型、复杂外形模型及外形运动学解耦的椭圆扩展外形模型。本实例使用但不限于使用外形运动学解耦的椭圆扩展外形模型对不可分辨目标群的目标群进行建模,即将目标群的空间分布建模为质心和椭圆扩展外形两部分组成,其中目标群质心用于描述目标群的行为模式,对应参数为运动学参数;椭圆扩展外形用于描述目标群的规模大小和***轮廓形状,对应参数为外形参数,其实现步骤包括如下:
1.1)对二维平面中的目标群运动学参数κ进行建模,其包含目标群质心的轴向位置信息,轴向速度信息或更高维度信息,如轴向加速度信息等,此时运动学参数κ可表示为:
其中,表示目标群质心的位置信息,px表示目标群质心在x轴向上的坐标,py表示目标群质心在y轴向上的坐标;/>表示群目标质心的速度信息,vx表示目标群质心在x轴向上的速度,vy表示目标群质心在y轴向上的速度;“...”表示省略了包含加速度、加加速度在内的更高维度的运动学信息;
1.2)对二维平面中以椭圆扩展外形近似的目标群外形参数进行建模,其包含椭圆外形的长短半轴长度以及椭圆倾斜角度三种信息,此时外形参数e可表示为:
其中,表示椭圆外形的倾斜角,即椭圆长轴与笛卡尔坐标系x轴正方向的夹角,la和lb分别是椭圆外形的长半轴长度和短半轴长度;
1.3)根据目标群运动学参数数学模型κ和目标群外形参数数学模型e,得到采用外形运动学解耦的椭圆扩展外形建模下目标群的状态参数模型ξ:
ξ=[κ,e]
由上式易知目标群状态参数被描述为增广状态,且由目标群运动学参数κ和目标群外形参数e两部分构成。
步骤2:在二维笛卡尔坐标系下,根据目标群的实际运动模式建立基于马尔科夫状态空间方程的目标群状态转移模型。
假设场景中的目标群均服从同一运动模型,现有方法包括采用非线性模型、线性模型对目标状态转移过程进行建模;本实例使用但不限于采用基于马尔科夫状态空间方程的线性模型对目标群状态转移过程进行建模,其实现步骤包括如下:
2.1)将k时刻下,目标群状态参数ξ中的运动学参数和外形参数的过程噪声均建模为高斯白噪声,表示为:
其中,为目标群运动学参数的过程噪声,/>为目标群外形参数的过程噪声;/>和/>分别表示目标群运动学参数和外形参数过程噪声的协方差矩阵;
2.2)建立在两个连续时刻k和k+1下的目标群状态转移的数学模型,由于目标群状态参数模型中两部分参数视为解耦,所以两部分参数的状态转移过程视为相互独立并分别进行数学建模,表示为:
其中,分别是根据群目标运动模式所获得的运动学参数和外形参数的状态转移矩阵,κk和ek分别表示k时刻的目标群运动学参数和外形参数,κk+1和ek+1分别表示k+1时刻的目标群运动学参数和外形参数。
步骤3,在二维笛卡尔坐标系下,将获得的以极坐标形式表达的雷达量测转化为直角坐标表示。
为了方便后续处理,需要将极坐标形式表达的雷达量测转化为直角坐标表示,其实现步骤包括如下:
3.1)设雷达位于二维笛卡尔坐标系的原点,对雷达接收机在一个积累时间里接收到的所有回波信号,依次进行脉压、MTD、恒虚警检测以及单脉冲测角处理,从而得到雷达在扫描时间k内接收到的所有回波信号中量测的总数Ik以及以极坐标形式表达的雷达量测:
其中,为k时刻获取的第i个雷达量测,/>为k时刻雷达量测构成的量测集合。
3.2)忽略雷达量测中包含的俯仰角信息,将单个雷达量测表示为极坐标下的组合形式:
其中,表示第i个雷达量测的测距信息,/>表示第i个雷达量测的测角信息,运算·T表示取向量转置;
3.3)根据坐标轴转换公式,将雷达量测转换为以直角坐标系表示量测/>
其中,和/>分别为直角坐标系量测/>在x轴向和y轴向上的坐标值,分别表示为:
3.4)用经过上述坐标轴转换的所有量测构成k时刻的量测集合Zk:
步骤4,在二维笛卡尔坐标系下,对群目标所产生量测的空间分布的数学模型进行建模,并获得目标群量测的似然函数。
现有方法对群目标所产生量测的空间分布的建模包括使用随机矩阵模型、随机超曲面模型以及高斯过程模型,考虑到目标群生成量测时具有的非单一性,本实例采用但不限于采用二维平面下的乘性误差均匀采样MEM模型对目标群所产生量测的空间分布进行建模,其实现步骤包括如下:
4.1)根据k时刻下群目标的运动学参数κk,求取群目标中心坐标sk:
sk=H·κk
其中,H=[I2 0]为观测矩阵或质心提取矩阵,其中I2为二阶单位矩阵,0为全0矩阵,κk表示k时刻的目标群运动学参数;
4.2)根据k时刻下群目标的外形参数,求取目标群外形旋转矩阵Γk:
其中,表示k时刻群目标外形的倾斜角;lk,a和lk,b分别是k时刻群目标外形的长半轴长度和短半轴长度;
4.3)根据目标群外形旋转矩阵Γk,建立k时刻的目标群量测空间分布的数学模型,该模型也称为目标群量测模型:
其中,表示k时刻下的第i个量测;/>为乘性误差向量,元素/> 的值满足区间[-1,1]上的随机抽样,其表示目标群中第i个点目标的位置;/>为第i个点目标产生量测的加性噪声;
4.4)根据目标群量测模型,获取目标群量测的似然函数,其服从高斯分布,表示为:
其中,ξk为目标群的状态参数,为乘性误差向量/>的协方差矩阵,∑δ为加性噪声/>的协方差矩阵。
步骤5,利用步骤2建立的群目标状态转移模型,对场景中多个目标群的状态参数进行一步预测,实现贝叶斯框架下的时间更新过程。
由于本实例的工作整体流程基于贝叶斯滤波框架,其考虑两个连续时刻间的目标状态的时间更新和量测更新,以下为表述方便,考虑两个连续时刻t和t+1,即当前时刻为t+1,前一时刻为t。则在本步骤中,需要利用前一t时刻的后验估计来对当前t+1时刻的状态进行预测,其实现步骤包括如下:
5.1)针对目标群数目固定的场景,设前一t时刻的场景中有Mt个目标,则多个由步骤1中状态参数模型建模的目标群将构成t时刻的状态参数集合Ξt:
其中,表示前一t时刻目标群m的状态参数,/>为目标群m在合并t-1时刻之前所有关联量测值和前一t时刻/>个关联量测值/>之后的状态估计结果:
5.2)根据步骤2所建立的目标群状态转移模型,利用前一t时刻的后验估计得到当前t+1时刻的预测状态,包括当前t+1时刻目标群m的运动学参数预测值及对应协方差矩阵/>当前t+1时刻目标群m的外形参数预测值/>及对应协方差矩阵/>其分别表示为如下线性形式:
其中,和/>分别表示前一t时刻对目标群m运动学参数的后验估计值以及对应协方差矩阵;/>和/>分别表示前一t时刻对目标群m外形参数的后验估计值以及对应协方差矩阵;/>和/>分别表示目标群运动学参数和外形参数过程噪声的协防差矩阵;/>和/>分别表示群目标运动学参数和外形参数的状态转移矩阵。
步骤6,根据步骤5中获得的目标群的预测状态,对步骤3中坐标转换后的当前时刻的量测使用门限处理,对场景中的杂波进行一次去除。
由于获取的量测中包含了场景中的噪声,这些噪声会影响并降低本实例后续对多目标群跟踪的准确度,故需要首先对场景中的杂波进行去除,保留由目标产生的量测,其实现步骤包括如下:
6.1)获取步骤5中所有群目标在当前t+1时刻的状态预测结果:
其中,和/>分别表示当前t+1时刻目标群m的运动学参数预测值和对应协方差矩阵,/>和/>分别表示当前t+1时刻目标群m的外形参数预测值和对应协方差矩阵;
6.2)以目标群m为例,获取其在当前t+1时刻预测的外形参数包括朝向角与长短半轴长度,表示为:
其中,表示当前t+1时刻目标群m的朝向角预测,/>和/>分别表示当前t+1时刻目标群m的长半轴长度和短半轴长度预测;
6.3)根据步骤6.2)中获取的群目标m的外形参数,计算目标群m在当前时刻的外形旋转矩阵的预测值
6.4)计算上述预测值矩阵的第一行和第二行对外形参数的雅可比矩阵,分别表示为/>和/>
6.5)根据步骤6.3)和步骤6.4)的计算结果,对目标群m求取其乘性误差与外形旋转矩阵的一阶协方差矩阵和二阶协方差矩阵/>
其中,运算tr(·)表示取矩阵的迹,不难看出和/>均为2×2的矩阵;
6.6)根据步骤6.5)的计算结果,求取群目标m的创新矩阵
其中,H为观测矩阵,表示加性噪声δ的协方差矩阵,/>表示当前t+1时刻运动学参数估计/>对应的协方差矩阵;
6.7)对当前t+1时刻场景中的所有量测依次计算其到各目标群的门限距离,该门限距离采用椭圆门限,即计算当前t+1时刻的量测到目标群m的椭圆门限距离为:
其中,表示新息向量,H为观测矩阵;
6.8)设门限阈值为常数gσ,将量测到各群目标的门限距离与门限阈值gσ进行比较:
若对某一目标群m,有则量测/>满足门限判断条件,将其认为是由目标群m产生,对在当前t+1时刻的量测集/>中进行保留,其中Ik+1表示当前t+1时刻下的量测总数;
若对任意目标群m,有则量测/>不满门限判断条件,将其视为背景杂波产生,直接从当前t+1时刻的量测集中进行一次剔除,
步骤7,根据步骤4中建立的目标群量测模型,对背景杂波的数目与目标产生量测的数目进行建模,计算量测与目标群或背景杂波关联的边缘概率,对一次剔除后的场景杂波进行二次去除。
采用门限处理只能对距离目标群较远的场景杂波进行去除,对距目标群较近的杂波往往去除效果不佳。为进一步去除场景杂波,提高本实例对多目标群的跟踪精确度,需对场景杂波进行二次去除,其实现步骤包括如下:
7.1)利用非同质泊松过程对背景杂波进行建模,即将杂波数目建模为服从参数λc的泊松分布,且在场景内均匀分布,杂波密度为ρc;
7.2)利用泊松点过程对当前时刻目标m生成的量测数量进行建模,即将目标产生量测的数目建模为服从参数λm的泊松分布,且在目标群的扩展外形内均匀分布;
7.3)根据步骤4.4)中获取的量测似然函数,计算每个量测到每个目标预测状态的似然
其中,目标群索引m∈{1,2,…,Mt+1},Mt+1为当前t+1时刻的目标群个数,量测索引i∈{1,2,…,It+1},It+1为当前t+1时刻的量测个数;表示目标群m在当前t+1时刻的外形旋转矩阵的预测值,/>为乘性误差向量协方差矩阵,∑δ为加性噪声协方差矩阵,/>表示前t+1时刻第i个量测,常矩阵H表示观测矩阵,/>表示对目标群在当前t+1时刻的运动学参数预测;
7.4)通过似然与步骤7.1)、步骤7.2)建立的量测数目模型,计算当前t+1时刻下量测/>分配到目标群m的边缘概率/>
其中,表示当前t+1时刻目标群m产生量测数目的泊松率,
为归一化因子,常数ρc表示背景杂波密度,当m=0时,/>意味着将量测分配给背景杂波的概率;
7.5)通过边缘分配概率对量测进行剔除,即判断量测/>到场景杂波的边缘分配概率/>是否满足条件/>
若满足,则视为量测来自于场景杂波,剔除该量测;
若不满足,则视为量测来自于目标群,进行保留;
7.6)将经过步骤7.5)剔除处理后剩余的量测均视为目标群所产生的量测,称为有效量测,表示为其中/>表示当前时刻的有效量测数量。
步骤8,对步骤6和步骤7中两次杂波去除处理后剩余的有效量测,结合步骤5中所得多目标群的预测状态,进行有效量测和目标群之间的数据关联,并获得关联概率矩阵。
现有方法对群目标的数据关联算法包括使用采用传统数据关联JPDA算法、近似数据关联LJPDA算法,本实例为获得更加精确的数据关联结果,采用但不限于采用基于混合高斯模型迭代最大化期望模糊聚类算法实现有效量测与目标群间的数据关联,其实现步骤包括如下:
8.1)初始化高斯混合模型中的多个高斯模型参数,每个高斯模型即代表一个目标群,高斯模型个数与目标群数目相等,即对于当前t+1时刻下的有效量测分布,可将其建模为混合高斯模型
其中,为高斯模型m的权重,/>为高斯模型m的均值、/>为高斯模型m的方差,其可分别初始化为:
其中,表示当前t+1时刻目标群m的运动学参数预测,
表示目标群m在当前t+1时刻的新息矩阵,
Mt+1为当前t+1时刻的目标群数目,
表示第m个高斯模型均值的初始值,
表示第m个高斯模型方差的初始值,
表示第m个高斯模型权重的初始值;
8.2)初始化目标群m与有效量测之间的关联概率/>
8.3)初始化迭代变量,包括初始迭代次数k=0,初始化增量阈值E和最大迭代次数K;
8.4)执行EM算法中的求期望步骤,即更新量测与目标的关联概率
其中,表示当前t+1时刻的第i个有效量测,/>表示第m个高斯模型均值在第k次迭代后的更新值,/>表示第m个高斯模型方差在第k次迭代后的更新值,/>表示第m个高斯模型权重在第k次迭代后的更新值;
8.5)执行EM算法中的最大化期望步骤:即更新各个高斯模型的参数,包括各高斯模型的权重均值/>和方差/>
其中,表示当前t+1时刻第i个有效量测,/>表示第m个高斯模型均值在第k次迭代后的更新值,/>表示第m个高斯模型方差在第k次迭代后的更新值,表示第m个高斯模型权重在第k次迭代后的更新值;
8.6)重复步骤8.4)和步骤8.5),直到达到最大迭代次数K,或两次迭代的间增量小于增量阈值E,即得到有效量测/>到目标群m的关联概率
8.7)根据步骤8.6)的关联概率得到当前时刻的关联概率矩阵:
步骤9,根据步骤8获得的关联概率矩阵,基于扩展式卡尔曼滤波器的更新公式和概率数据关联法则,对每个目标群的状态参数使用关联的有效量测进行基于关联概率的贯序量测更新。
9.1)对目标群m,根据利用当前t+1时刻第i-1个关联的有效量测参数估计后的运动学参数/>以及对应协方差矩阵/>求取该状态估计目标的质心/>作为第i个有效量测的期望值,并求取目标群运动学参数/>与量测/>的协方差矩阵
9.1.1)计算第i个有效量测的期望值
9.1.2)计算目标群运动学参数与量测/>的协方差矩阵/>
9.2)对目标群m,根据当前时刻第i-1个关联的有效量测下外形参数估计/>以及对应协方差矩阵/>计算新息矩阵/>
9.2.1)根据当前t+1时刻目标群m在第i-1个关联的有效量测估计的外形参数求取目标群m的旋转外形矩阵估计值/>
其中,表示有效量测/>估计下的目标群m外形的朝向,/>和/>分别表示有效量测/>估计下的目标群m外形的长半轴长度和短半轴长度;
9.2.2)求取上述旋转外形矩阵第一行和第二行对外形参数的雅可比矩阵,分别表示为和/>
9.2.3)根据步骤9.2.1)和步骤9.2.2)的计算结果,求取目标群m的乘性误差与外形旋转矩阵之间的一阶协方差矩阵和二阶协方差矩阵/>
其中,运算tr(·)表示取矩阵的迹,表示乘性误差向量/>的协方差矩阵,表示目标群m的外形参数估计/>对应的协方差矩阵;
9.2.4)求取目标群m的在有效量测下的新息矩阵/>
其中,∑δ表示加性噪声δ的协方差矩阵,表示采用当前t+1时刻目标群m在第i-1个关联的有效量测/>估计的运动学参数/>对应的协方差矩阵;
9.3)对目标群m,根据外形参数协方差矩阵和步骤9.2.4)计算所得的创新矩阵,以及步骤9.1.1)中获得的第i个有效量测的期望值/>构造当前时刻关联的第i个有效量测/>的对应伪量测/>并求取外形参数/>与伪量测之间/>的协方差矩阵/>及伪量测的方差矩阵/>/>
9.3.1)求取关联到目标群m的第i个有效量测对应的伪量测/>
其中,为当前t+1时刻关联到目标群m的第i个有效量测/>的期望值,
表示选择第二列信息的常数矩阵;
9.3.2)求取伪量测的期望值
其中,是目标群m在有效量测/>下的新息矩阵/>运算/>表示将矩阵列向量化的算子;
9.3.3)求取外形参数协方差和伪量测-外形参数协方差的之间的比例系数
其中,表示使用有效量测/>对目标m外形参数的估计下的旋转外形矩阵/>的第1行,/>表示该旋转外形矩阵的第2行,/>分别为该旋转外形矩阵第一行和第二行对目标外形参数的雅克比矩阵,/>表示乘性误差向量/>的协方差矩阵;
9.3.4)根据步骤9.3.3)的结果,求取伪量测与外形参数估计/>之间的协方差矩阵/>
其中,为外形参数估计/>对应的协方差矩阵/>
9.3.5)求取伪量测的方差
其中,为表示选择第三列信息的常数矩阵,运算/>表示克罗内克积。
9.4)对目标群m,使用当前t+1时刻的第i个关联的有效量测对目标群m的运动学参数/>及对应协方差矩阵/>进行基于概率数据关联的量测更新,得到目标群m更新后的运动学参数/>及对应协方差矩阵/>
其中, 表示当前t+1时刻第i个有效量测,/>表示当前t+1时刻有效量测/>与群目标m的关联概率;
9.5)对目标群m,使用当前t+1时刻的第i个关联的有效量测生成的伪量测对目标群m的外形参数/>以及对应协方差矩阵/>进行基于概率数据关联的量测更新,得到目标群m更新后的外形参数/>及对应协方差矩阵/>
其中, 表示当前t+1时刻第i个关联的有效量测生成的伪量测;
9.6)量测索引加1,即令i=i+1;
9.7)重复步骤9.1)到步骤9.6),直到当前时刻的每个目标都使用与其关联的全部量测对其状态参数进行估计,完成量测更新。
步骤10,基于贝叶斯框架,将当前时刻所得的多目标群状态更新结果作为下一时刻多目标群的先验状态信息,实现跟踪流程在时间步上的迭代。
10.1)将当前t+1时刻所得的目标群m,m∈{1,2,…,Mt+1}的外形参数更新结果作为下一t+2时刻该目标群m外形参数的先验信息,可表示为:
其中,是使用当前t+1时刻的全部有效量测对目标群m的外形参数的贯序估计结果,/>是对应的协方差矩阵;
10.2)将当t+1前所得的目标群m,m∈{1,2,…,Mt+1}运动学参数更新结果作为下一t+2时刻该目标群m运动学参数的先验信息,可表示为:
/>
其中,是使用当前t+1时刻的全部有效量测对目标群m的运动学参数的贯序估计结果,/>是对应的协方差矩阵;
10.3)令当前时刻为t+2时刻,前一时刻为t+1时刻,完成方法在时间步上的迭代。
步骤11,并重复步骤5~步骤10,直至目标群全部离开监视区域,跟踪过程结束。
以下通过仿真实验对本发明上述有益效果作进一步说明:
1.仿真条件:
1)仿真场景为二维场景,场景大小:[-10,980]×[-320,320],目标群中心位置信息表达采用直角坐标系。场景中设有两个大小不同的目标群,它们的运动模式相同且均保持线性运动,并随着时间步的推移发生合并与***的行为。两目标群在仿真中运行的真实轨迹如图2所示。
2)两目标群的起始状态如表1所示:
表1两目标群起始状态设定
3)两目标群的状态转移模型:
运动学参数的状态转移矩阵:
外形参数的状态转移矩阵:
运动学参数的过程噪声:
外形参数的过程噪声:
4)两目标群的量测模型:
量测矩阵:H=[I2 O2]
量测噪声(加性噪声):∑δ=diag[1 1]
其中:运算diag·表示取对角阵,In表示n阶单位矩阵,On表示n阶零矩阵;
5)仿真软件环境:Matlab R2016a;
6)仿真硬件环境:处理器:2.8GHz Intel Core i7 processor;内存:16G。
2.仿真内容
在上述仿真条件下,分别用本发明和现有两种多目标群跟踪方法GGIW-PMB、LJPDA-MEM对二维场景的两个目标群进行跟踪:在1s时仿真开始,两个大小不同的目标群保持相同的速率逐渐靠近,该过程从1s持续到39s;随着两个目标群于第40s发生交叠,此时,视为两个目标群发生合并;之后,两个目标群保持该种合并状态并行持续27s,最终于68s发生***并逐渐远离直至88s仿真结束。在整个仿真过程按照T=1s的间隔选取了88个采样时刻。
将本发明和现有两种多目标群跟踪方法GGIW-PMB、LJPDA-MEM进行上述场景多目标群跟踪后的高斯-瓦瑟斯坦距离GWD结果与作为基准的理论最优TDA-MEM方法进行对比,所得对比结果如图3所示;
将本发明和现有多目标群跟踪方法LJPDA-MEM进行上述场景多目标群跟踪后的数据关联正确率DACR结果与作为基准的理论最优TDA-MEM方法进行对比,所得对比结果如图4所示;
将本发明和现有两种多目标群跟踪方法GGIW-PMB、LJPDA-MEM进行上述场景多目标群跟踪后的广义最优子模式分配GOSPA结果和单步仿真时间CYC结果与作为基准的理论最优TDA-MEM方法进行五个指标的对比,所得对比结果如表2所示。
表2各方法的广义最优子模式分配GOSPA对比结果
GOS | NLE | FALE | MISE | CYC | |
TDA-MEM方法 | 13.8297 | 3.2235 | 0 | 0 | 0.29 |
LJPDA-MEM方法 | 19.6454 | 4.5864 | 0 | 0 | 0.73 |
GGIW-PMB方法 | 76.2798 | 6.6739 | 4.1599 | 4.7117 | 55.56 |
本发明方法 | 14.7183 | 3.4704 | 0 | 0 | 0.83 |
表2中用于对比方法的出处如下:
所述GGIW-PMB方法,其为Y.Xia,K.L.Svensson,M.Fatemi,/>F.García-Fernández and J.L.Williams等人于2022年在IEEE Transactions on Aerospace andElectronic Systems期刊上提出的一种基于多扩展目标泊松多伯努利近似的多目标群跟踪方法,即“Poisson Multi-Bernoulli Approximations for Multiple Extended ObjectFiltering”。
所述LJPDA-MEM方法,其为S.Yang,K.Thormann and M.Baum等人于2018年在IEEE10th Sensor Array and Multichannel Signal Processing Workshop国际会议上提出的一种基于多扩展目标线性联合概率数据关联实现的多目标群跟踪方法,即“Linear-TimeJoint Probabilistic Data Association for Multiple Extended Object Tracking”。
所述作为对比基准的理论最优数据关联TDA-MEM方法,其为理论上最优的方法,并基于乘性误差均匀采样MEM量测分布建模和人为先验提供的理论最佳数据关联实现多目标群跟踪;
所述的五个指标包括GOS、NLE、FALE、MISE、CYC,其中:
GOS指标表示整体误差,用于比较各方法的整体跟踪性能;
NLE指标表示定位误差,用于比较各方法的整体跟踪准确度;
FALE指标表示误跟踪误差,用于比较各方法对目标群误跟踪程度;
MISE指标表示丢失跟踪误差,用于比较各方法对目标群丢失跟踪程度;
由于FALE.指标和MISE.指标只适用于描述场景中目标群发生改变的情况,由于假定仿真场景中目标群数目不发生改变,因此表2中的FALE.指标和MISE.指标可忽略。
3.仿真结果分析:
从图3可以看出,在1s~39s两目标群逐渐靠近的阶段,本发明和现有的LJPDA-MEM方法与作为基准的TDA-MEM方法有着极为接近的表现,且远远优于现有的GGIW-PMB方法,表明本发明采用乘性误差均匀采样MEM模型对量测分布建模,较之于采用随机矩阵量测分布建模的GGIW-PMB方法,在目标群的外形估计上有更为精确的结果。随着时间增加,当两目标群于40s发生重叠时,现有LJPDA-MEM方法的高斯-瓦瑟斯坦距离GWD明显增大,并与此后逐渐发散,说明该方法出现了严重的误跟踪现象,而本发明方法的高斯-瓦瑟斯坦距离GWD一直和维持在和作为基准的TDA-MEM方法相近的水平线上,且在多目标群出现合并的40s和出现***的68s时不发生明显改变,说明使用本发明方法处理多目标群合并与行为时不会发生误跟现象。
从图4可以看出,在两个目标发生合并和并行的40s~67s时间段内,本发明与现有LJPDA-MEM方法的数据关联正确率DACR均出现了下降,但本发明下降的幅度较现有LJPDA-MEM方法更小,且在两目标群***后能够快速的恢复到正常水平,即恢复到90%以上,但现有LJPDA-MEM方法在两目标群***后无法照常恢复,这也体现了本发明在数据关联上的准确性。
从表2可以看出,在整体误差GOS和定位误差NLE这两个指标上,本发明与作为基准的TDA-MEM方法有着最为接近的表现,所以在整体跟踪精度上,本发明较其它现有方法也更为精确。同时可以看出,在衡量算法计算时间的CYC.指标上,本发明获得了与现有LJPDA-MEM方法非常近似的表现,说明本发明在提高方法整体新能的同时还保持了与现有LJPDA-MEM方法一样的计算复杂度。
综上所述,本发明基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法,在高斯-瓦瑟斯坦距离GWD指标、数据关联正确率DACR指标、广义最优子模式分配GOSPA指标较现有多目标群跟踪方法均存在优势,同时保持了和现有方法一样的计算复杂度,验证了本发明能够用于跟踪一定杂波背景下的多个目标群,且能够在线性计算复杂度内有效解决发生合并与***行为下的多目标群跟踪问题。
实现上述方法实施例的全部或部分步骤可以通过程序指令相关的硬件来完成,前述的程序可以存储于一计算机可读取存储介质中,该程序在执行时,执行包括上述方法实施例的步骤;而前述的存储介质包括:ROM、RAM、磁碟或者光盘等各种可以存储程序代码的介质。各步骤的标号仅仅是为了清楚的说明本发明的实现方案,其序号顺序不作限定。
以上所述,仅为本发明的一个具体实施方式,但本发明的保护范围并不局限于此,任何熟悉本技术领域的技术人员在本发明揭露的技术范围内,可轻易想到变化或替换,都应涵盖在本发明的保护范围之内。
Claims (13)
1.一种基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法,其特征在于,包括如下步骤:
(1)采用外形运动学解耦的椭圆扩展外形模型,建立单个不可分辨目标群的状态参数模型;
(2)根据目标群在场景中的运动模式,建立基于马尔科夫状态空间模型的目标群状态转移模型;
(3)在以雷达为传感器的跟踪场景中,将以二维极坐标形式表达的雷达量测转化为直角坐标表示;
(4)针对目标群产生量测的特点,采用乘性误差均匀采样MEM模型对群目标所产生量测的空间分布进行建模,并获得目标群量测的似然函数;
(5)利用步骤(2)建立的群目标状态转移模型,使用标准卡尔曼滤波预测公式对场景中多个目标群的状态参数进行一步预测;
(6)根据多个目标群预测状态,对场景中的量测使用椭圆外形的门限处理,实现对场景中的杂波的一次剔除;
(7)根据步骤(4)中建立的目标群量测空间分布模型,对背景杂波的数目与目标产生量测的数目进行建模,并基于量测与目标的边际关联概率在一次杂波剔除的基础上实现对场景中杂波的二次剔除;
(8)对两次杂波去除后剩余的有效量测,结合步骤(5)中所得多目标群的状态预测结果,采用具有线性计算复杂度的基于混合高斯模型迭代最大化期望模糊聚类的数据关联算法实现有效量测和目标群的数据关联,并获得关联概率矩阵;
(9)根据数据关联结果,基于扩展式卡尔曼滤波器的更新公式和概率数据关联法则,结合步骤(1)所建立的目标群状态参数模型,对每个目标群的状态参数使用有效量测进行基于关联概率的贯序量测更新;
(10)基于贝叶斯框架,将当前所得的多目标群状态参数估计结果作为下一次多目标群估计的先验信息;
(11)重复步骤(5)~步骤(10),直至目标群离开雷达监视区域完成跟踪过程。
2.根据权利要求1所述的方法,其特征在于,步骤(1)中采用外形运动学解耦的椭圆扩展外形模型,建立单个不可分辨目标群的状态参数模型,实现步骤包括如下:
(1a)建模包含目标群的轴向位置信息、轴向速度信息或更高维度信息在二维笛卡尔坐标系下的目标群的运动学参数数学模型κ:
其中,表示目标群质心的坐标,px、py表示目标群质心在直角坐标系中两个轴向上的位置信息,/>表示群目标质心的运动矢量,vx,vy表示目标群质心在直角坐标系中两个轴向上的速度信息,省略号省略了包含群目标质心的加速度、加加速度在内的更高维度的运动学信息。
(1b)建立包含椭圆的长短半轴长度及椭圆倾斜角度这三种信息在二维平面下的椭圆扩展外形近似的目标群外形参数数学模型e:
其中,表示椭圆外形的倾斜角,即椭圆长轴与X轴方向的夹角;la和lb分别是椭圆的长半轴长度和短半轴长度;
(1c)根据目标群运动学参数数学模型κ和目标群外形参数数学模型e,得到椭圆扩展外形近似的目标群的状态参数模型ξ:
ξ=[κ,e]。
3.根据权利要求1所述的方法,其特征在于,步骤(2)中根据目标群在场景中的运动模式,建立基于马尔科夫状态空间模型的目标群状态转移模型,实现步骤包括如下:
(2a)将状态参数模型ξ中的两部分参数的过程噪声建模为高斯白噪声,表示为:
其中,为目标群运动学参数的过程噪声,/>为目标外形参数的过程噪声;/>和/>分别表示目标群运动学参数和外形参数过程噪声的协方差矩阵;
(2b)根据(2a)的结果,建立在两个连续时刻k和k+1下的目标群状态转移的数学模型:
其中,分别是根据群目标运动模式所获得的运动学参数和外形参数的状态转移矩阵,κk和ek分别表示k时刻的目标群运动学参数和外形参数。
4.根据权利要求1所述的多目标群跟踪方法,其特征在于,步骤(3)中将以二维极坐标形式表达的雷达量测转化为直角坐标表示,实现步骤包括如下:
(3a)对雷达接收机在一个积累时间里接收到的所有回波信号,依次进行脉冲压缩、动目标检测MTD、恒虚警检测以及单脉冲测角处理,获取k时刻的雷达量测集合
其中Ik为获得的量测个数;为单个雷达目标量测,其使用极坐标表示:
表示第i个雷达量测的测距信息,/>表示第i个雷达量测的测角信息;
(3b)根据坐标轴转换公式,将雷达量测转换为直角坐标系量测/>
其中,和/>分别为直角坐标系量测/>在x轴向和y轴向上的坐标值,分别表示为:
5.根据权利要求1所述的方法,其特征在于,步骤(4)中采用乘性误差均匀采样MEM模型对群目标所产生量测的空间分布进行建模,并获得目标群量测的似然函数,实现步骤包括如下:
(4a)根据k时刻下群目标的运动学参数κk,求取群目标中心坐标sk:
sk=H·kk
其中,H=[I20]为观测矩阵或质心提取矩阵,其中I2为二阶单位矩阵,0为全0矩阵,κk表示k时刻的目标群运动学参数;
(4b)根据k时刻下群目标的外形参数,求取目标群外形旋转矩阵Γk:
其中,表示k时刻群目标外形的倾斜角;lk,a和lk,b分别是k时刻群目标外形的长半轴长度和短半轴长度;
(4c)根据(4a)和(4b)的结果,建立k时刻的目标群量测空间分布的数学模型:
其中,表示k时刻下的第i个量测;/>为乘性误差向量,元素/>的值满足区间[-1,1]上的随机抽样,其表示目标群中第i个点目标的位置;/>为第i个点目标产生量测的加性噪声;
(4d)根据(4c)中所建立的群目标量测模型,获取目标群量测的似然函数,其服从高斯分布,表示为:
其中,ξk为目标群的状态参数,为乘性误差向量/>的协方差矩阵,∑δ为加性噪声/>的协方差矩阵。
6.根据权利要求1所述的方法,其特征在于,步骤(5)中利用步骤(2)所建立的群目标状态转移模型,使用标准卡尔曼滤波预测公式对场景中多个目标群的状态参数进行一步预测。其实现步骤包括如下:
(5a)针对目标群数目固定的场景,设前一时刻t的场景中有Mt个目标,则多个由步骤(1)状态参数建模的目标群将构成t时刻的状态参数集合Ξt:
其中,表示t时刻目标群m的状态参数,/>为目标群m在合并t-1时刻之前所有关联量测值和t时刻/>个关联量测值/>之后的状态估计结果;
(5b)根据步骤(2)的目标群状态转移模型,利用前一t时刻的后验估计得到当前t+1时刻的预测状态,包括当前时刻目标群m的运动学参数预测值及对应协方差矩阵/>当前时刻目标群m的外形参数预测值/>及对应协方差矩阵/>表示为如下线性形式:
其中,和/>分别表示前一t时刻对目标群m运动学参数的后验估计值以及对应协方差矩阵;/>和/>分别表示前一t时刻对目标群m外形参数的后验估计值以及对应协方差矩阵;/>和/>分别表示目标群运动学参数和外形参数过程噪声的协防差矩阵;和/>分别表示群目标运动学参数和外形参数的状态转移矩阵。
7.根据权利要求1所述的方法,其特征在于,步骤(6)中根据多个目标群预测状态,对场景中的量测使用椭圆外形的门限处理,实现对场景中的杂波的一次剔除,其实现步骤包括如下:
(6a)根据步骤(5)中所求当前时刻的参数预测值以及对应参数预测协方差/>结合步骤(4)建立的目标群量测模型,求解每个目标群的量测预测协方差矩阵/>该矩阵也被称为目标群m的新息矩阵;
(6b)计算当前t+1时刻某一量测到目标m的椭圆门限距离为:
其中,表示新息向量,H为观测矩阵;
(6c)设门限阈值为常数gσ,将量测到各群目标的门限距离与门限阈值进行比较:
若且m∈{1,...,Mt+1},则量测/>满足门限判断条件,进行保留;
若且m∈{1,...,Mt+1},则量测/>不满门限判断条件,则视为背景杂波,直接从当前t+1时刻的量测集/>中进行剔除,
其中Ik+1表示当前t+1时刻下的量测总数,Mt+1为当前t+1时刻下目标群的总数。
8.根据权利要求7所述的方法,其特征在于,步骤(6a)中求解目标群m的量测预测协方差矩阵实现步骤包括如下:
(6a1)根据目标群m的状态参数预测值计算目标群m在当前t+1时刻的外形旋转矩阵的预测值/>
其中,表示对目标群m椭圆外形朝向角的预测值,/>分别表示对目标群m椭圆长半轴长度和短半轴长度的预测值;
(6a2)求取外形旋转矩阵的第一行和第二行对目标群外形参数的雅可比矩阵,分别表示为/>和/>
(6a3)根据步骤(6a1)和步骤(6a2)的计算结果,求取目标群m的乘性误差向量与外形旋转矩阵的一阶协方差矩阵和二阶协方差矩阵/>
其中,运算tr(·)表示取矩阵的迹,表示乘性误差向量/>的协方差矩阵;
(6a4)根据(6a3)的结果、观测矩阵H、目标群m在当前t+1时刻的运动学参数的预测协方差矩阵及加性噪声协方差矩阵∑δ,计算目标群m的量测预测协方差矩阵/>
9.根据权利要求1所述的方法,其特征在于,步骤(7)中对背景杂波的数目与目标产生量测的数目进行建模,通过计算量测与目标群或背景杂波的边际关联概率实现对场景中杂波的二次剔除,其实现步骤包括如下:
(7a)利用非同质泊松过程将杂波数目表示为服从参数λc的泊松分布,且在场景内均匀分布;
(7b)利用泊松点过程将当前时刻目标群m生成的量测数量表示为服从参数λm的泊松分布,且在目标群内均匀分布;
(7c)根据步骤(4d)中的量测似然函数,计算每个量测到每个目标预测状态的似然
其中,目标群索引m∈{1,2,…,Mt+1},Mt+1为当前t+1时刻的目标群个数,量测索引i∈{1,2,…,It+1},It+1为当前t+1时刻的量测个数;表示目标群m在当前t+1时刻的外形旋转矩阵的预测值,/>为乘性误差向量协方差矩阵,Σδ为加性噪声协方差矩阵,/>表示前t+1时刻第i个量测,常矩阵H表示观测矩阵,/>表示对目标群在当前t+1时刻的运动学参数预测;
(7d)通过似然与步骤(7a)、步骤(7b)建立的量测数目模型,计算量测分配到目标的边缘概率/>
其中,表示当前t+1时刻目标群m产生量测数目的泊松率,/> 为归一化因子,常数ρc表示背景杂波密度,当m=0意味着将量测分配给背景杂波;
(7e)通过边缘分配概率对量测进行剔除,即判断量测到场景杂波的边缘分配概率是否满足条件/>
若满足,则视为量测来自于场景杂波,剔除该量测;
若不满足,则视为量测来自于目标群,进行保留;
(7f)将步骤(7e)剔除处理后剩余的量测均视为目标群所产生的量测,称为有效量测,表示为其中/>表示当前时刻的有效量测数量。
10.根据权利要求1所述的方法,其特征在于,步骤(8)中采用基于混合高斯模型迭代最大化期望模糊聚类的数据关联算法对有效量测和目标群的进行数据关联,并获得关联概率矩阵,其实现步骤包括如下:
(8a)初始化高斯混合模型中的多个高斯模型参数,每个高斯模型代表一个目标群,高斯模型个数与目标群数目相等,高斯模型m的参数包括权重均值/>和方差/>
其中,表示当前t+1时刻目标群m的运动学参数预测,/>表示目标群m在当前t+1时刻的新息矩阵,Mt+1为当前t+1时刻的目标群数目,符号[·]表示最大化期望算法中的迭代次数,[0]表示各参数的初始值,即/>表示第m个高斯模型均值的初始值,/>表示第m个高斯模型方差的初始值,/>表示第m个高斯模型权重的初始值;
(8b)初始化关联概率
其中,目标群索引m∈{1,,2,…,Mt+1},有效量测索引 为当前t+1时刻的有效量测数目;
(8c)初始化迭代变量,包括初始迭代次数k=0,初始化增量阈值E和最大迭代次数K;
(8d)执行EM算法中的求期望步骤,即更新量测与目标的关联概率
其中,表示当前t+1时刻的第i个有效量测,符号[k]表示最大化期望算法中第k次迭代后参数的更新值,即/>表示第m个高斯模型均值在第k次迭代后的更新值,/>表示第m个高斯模型方差在第k次迭代后的更新值,/>表示第m个高斯模型权重在第k次迭代后的更新值;
(8e)执行EM算法中的最大化期望步骤:即更新各个高斯模型的参数,包括各高斯模型的权重均值/>和方差/>
其中,符号[k+1]表示最大化期望算法中第k+1次迭代后参数的更新值,表示第m个高斯模型均值在第k次迭代后的更新值,/>表示第m个高斯模型方差在第k次迭代后的更新值,/>表示第m个高斯模型权重在第k次迭代后的更新值;
(8f)重复步骤(8d)和步骤(8e)直到达到最大迭代次数K,或两次迭代的间增量小于增量阈值E,即得到有效量测/>到目标群m的关联概率/>
(8g)根据(8f)的结果得到当前时刻的关联概率矩阵:
11.根据权利要求1所述的方法,其特征在于,步骤(9)中对每个目标群的状态参数使用关联的有效量测进行基于关联概率的贯序量测更新,其实现步骤包括如下:
(9a)对目标群m,根据当前t+1时刻第i-1个关联的有效量测下的运动学参数估计/>以及对应协方差矩阵/>求取该状态估计目标的质心/>作为第i个有效量测的期望值,并求取目标群运动学参数/>与量测/>的协方差矩阵/>
(9b)对目标群m,根据外形参数估计以及对应协方差矩阵/>计算新息矩阵/>
(9c)对目标群m,根据其外形参数估计和对应协方差矩阵/>及步骤(9b)中的创新矩阵,以及步骤(9a)的第i个有效量测的期望值/>构造当前时刻目标群m关联的第i个有效量测/>的对应伪量测/>求取外形参数/>与伪量测之间/>的协方差矩阵/>及伪量测的方差矩阵/>
(9d)对目标群m,使用当前t+1时刻的第i个关联的有效量测对目标群m的运动学参数/>及对应协方差矩阵/>进行基于概率数据关联的量测更新,得到目标群m更新后的运动学参数/>及对应协方差矩阵/>
其中, 表示当前t+1时刻第i个有效量测,/>表示当前t+1时刻有效量测/>与群目标m的关联概率;
(9e)对目标群m,使用当前t+1时刻的第i个关联的有效量测生成的伪量测对目标群m的外形参数/>以及对应协方差矩阵/>进行基于概率数据关联的量测更新,得到目标群m更新后的外形参数/>及对应协方差矩阵/>
其中, 表示当前t+1时刻第i个关联的有效量测生成的伪量测;
(9f)对量测索引加1,即令i=i+1;
(9g)重复步骤(9a)到步骤(9f),直到当前时刻的每个目标都使用与其关联的全部量测对其状态参数进行估计,完成量测更新。
12.根据权利要求11所述的方法,其特征在于,步骤(9b)中根据目标群m的外形参数估计及对应协方差矩阵/>计算新息矩阵/>实现步骤包括如下:
(9b1)根据使用当前t+1时刻第i-1个有效量测对目标m外形参数的估计值计算目标群m的旋转外形矩阵/>
其中,为使用有效量测/>对目标群m椭圆外形朝向角的估计值;分别为使用有效量测/>对目标群m椭圆外形长半轴长度和短半轴长度的估计值;
(9b2)分别计算旋转外形矩阵第一行和第二行对外形参数的雅可比矩阵和/>
(9b3)根据步骤(9b1)和步骤(9b2)的计算结果,求取目标群m的乘性误差向量与外形旋转矩阵的一阶协方差矩阵和二阶协方差矩阵/>
其中,运算tr(·)表示取矩阵的迹,表示乘性误差向量/>的协方差矩阵,/>表示目标群m的外形参数估计/>对应的协方差矩阵;
(9b4)根据步骤(9b3)的计算结果,求取目标群m在有效量测下的新息矩阵/>
其中,∑δ表示加性噪声δ的协方差矩阵,表示运动学参数估计/>以及对应的协方差矩阵。
13.根据权利要求11所述的方法,其特征在于,步骤(9c)构造当前时刻目标群m关联的第i个有效量测的对应伪量测/>求取外形参数/>与伪量测之间/>的协方差矩阵/>及伪量测的方差矩阵/>实现步骤包括如下:
(9c1)求取关联到目标群m的第i个有效量测对应的伪量测/>
其中,为当前t+1时刻关联到目标群m的第i个有效量测/>的期望值,
表示选择第二列信息的常数矩阵;
(9c2)求取伪量测的期望值
其中,是目标群m在有效量测/>下的新息矩阵/>运算/>表示将矩阵列向量化的算子;
(9c3)求取外形参数协方差和伪量测-外形参数协方差的之间的比例系数
其中,表示使用有效量测/>对目标m外形参数的估计下的旋转外形矩阵的第1行,/>表示该旋转外形矩阵的第2行,/>分别为该旋转外形矩阵第一行和第二行对目标外形参数的雅克比矩阵,/>表示乘性误差向量/>的协方差矩阵;
(9c4)根据步骤(9c3)的结果,求取伪量测与外形参数估计/>之间的协方差矩阵/>
其中,为外形参数估计/>对应的协方差矩阵/>
(9c5)求取伪量测的方差矩阵
其中,为表示选择第三列信息的常数矩阵,运算/>表示克罗内克积。/>
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202310977257.1A CN117036400A (zh) | 2023-08-04 | 2023-08-04 | 基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202310977257.1A CN117036400A (zh) | 2023-08-04 | 2023-08-04 | 基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN117036400A true CN117036400A (zh) | 2023-11-10 |
Family
ID=88638428
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202310977257.1A Pending CN117036400A (zh) | 2023-08-04 | 2023-08-04 | 基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN117036400A (zh) |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117784115A (zh) * | 2023-12-26 | 2024-03-29 | 兰州理工大学 | 高斯过程回归模型多扩展目标pmbm跟踪方法 |
CN118194609A (zh) * | 2024-05-17 | 2024-06-14 | 北京理工大学前沿技术研究院 | 基于泛化多变量模型的群目标跟踪方法、装置及存储介质 |
-
2023
- 2023-08-04 CN CN202310977257.1A patent/CN117036400A/zh active Pending
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117784115A (zh) * | 2023-12-26 | 2024-03-29 | 兰州理工大学 | 高斯过程回归模型多扩展目标pmbm跟踪方法 |
CN118194609A (zh) * | 2024-05-17 | 2024-06-14 | 北京理工大学前沿技术研究院 | 基于泛化多变量模型的群目标跟踪方法、装置及存储介质 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN111985093B (zh) | 一种带噪声估计器的自适应无迹卡尔曼滤波状态估计方法 | |
CN117036400A (zh) | 基于混合高斯模型模糊聚类数据关联的多目标群跟踪方法 | |
CN111178385B (zh) | 一种鲁棒在线多传感器融合的目标跟踪方法 | |
CN104199022B (zh) | 一种基于目标模态估计的临近空间高超声速目标跟踪方法 | |
CN109460539B (zh) | 一种基于简化容积粒子滤波的目标定位方法 | |
CN108871365B (zh) | 一种航向约束下的状态估计方法及*** | |
CN111259332B (zh) | 一种杂波环境下的模糊数据关联方法及多目标跟踪方法 | |
CN115204212A (zh) | 一种基于stm-pmbm滤波算法的多目标跟踪方法 | |
CN113191427B (zh) | 一种多目标车辆跟踪方法及相关装置 | |
CN116047495B (zh) | 一种用于三坐标雷达的状态变换融合滤波跟踪方法 | |
CN113343436A (zh) | 一种高斯混合协方差演化的碰撞概率计算方法及*** | |
CN113030940A (zh) | 一种转弯机动下的多星凸型扩展目标跟踪方法 | |
CN116224320B (zh) | 一种极坐标系下处理多普勒量测的雷达目标跟踪方法 | |
CN116734860A (zh) | 一种基于因子图的多auv自适应协同定位方法及*** | |
CN111563918A (zh) | 一种多卡尔曼滤波器数据融合的目标跟踪方法 | |
CN115544425A (zh) | 一种基于目标信噪比特征估计的鲁棒多目标跟踪方法 | |
CN115169136A (zh) | 三维空间下快速uk-gmphd多目标跟踪方法 | |
CN114440875A (zh) | 一种基于概率神经网络的重力匹配方法 | |
CN107590509B (zh) | 基于最大期望近似的切尔诺夫融合方法 | |
CN117784114B (zh) | 异常噪声下基于混合熵的不规则扩展目标跟踪方法 | |
CN115906413B (zh) | 一种Dirichlet过程混合模型节点自定位方法 | |
CN114415157B (zh) | 一种基于水声传感器网络的水下目标多模型跟踪方法 | |
CN114494340B (zh) | 一种kl交互式多模型水下目标跟踪方法 | |
CN111504326B (zh) | 一种基于t分布的鲁棒glmb多目标跟踪方法 | |
CN112581616B (zh) | 基于序贯分块滤波的最近邻ukf-slam方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination |