CN112114308A - 一种扇扫雷达空时联合目标跟踪方法 - Google Patents

Abstract

本发明涉及一种扇扫雷达空时联合目标跟踪方法，该方法包括：从观测雷达获得当前周期的量测信息；判断当前周期是否存在航迹，若存在则顺序执行步骤，若不存在则跳转执行最后一步骤；对当前存在的所有航迹进行一步预测，得到当前周期的航迹存在性概率预测值，以及状态的一步预测值和量测预测值；对每一条航迹建立相关波门，获取与航迹相关的点迹；利用与该航迹相关的点迹进行滤波，获得多个状态估计值；将其对应的多个状态估计值及状态的一步预测值进行概率互联，获得目标最终状态估计值及航迹存在性概率；对未被用来更新航迹的点迹执行航迹初始化，获得新航迹目标初始状态。本发明能够有效的提高扇扫雷达多目标跟踪精度及效率。

Description

一种扇扫雷达空时联合目标跟踪方法

技术领域

本发明涉及空间目标跟踪技术领域，尤其涉及一种扇扫雷达空时 联合目标跟踪方法。

背景技术

扇扫雷达由于雷达天线在一定的扇扫范围内进行往返运动，一周 内对同一目标进行两次扫描，使得扇扫雷达必然存在以下的问题：

1.雷达不但在方位上有扫描边界,而且对同一方位还存在着扫描方 向的差别,使得航迹预相关不能只受点迹与航迹的方位控制，也与扫描 方向有关，特别是位于扫描边界的目标，不同扫描方向雷达对目标的 采样间隔往往差距较大。

2.不同方位上的目标预测时间(得到某个目标点迹后,预测同一目标 下次出现的时间)不同，同一目标相邻点迹之间采样时间间隔不同,使得 滤波与预测不能在固定重访间隔的条件下进行。

3.目标与雷达本身之间存在相对运动，同时存在噪声等因素影响， 因此目标下一次被照射的时间存在不确定性。

传统多目标跟踪方法中，天线对目标的重访间隔往往是确定且已 知的，这对扇扫雷达显然是不适用的。目前已有一些研究成果对于扇 扫雷达目标重访间隔进行粗略计算，但均是建立在雷达天线方位角速 度大小不变的前提下，由于扇扫雷达天线扫描方向会发生改变，显然 这种前提是不合理的，存在较大的误差。

发明内容

(一)要解决的技术问题

本发明要解决的技术问题是解决传统目标跟踪方法对于扇扫雷达 存在模型不匹配、检测偏差大的问题。

(二)技术方案

为了解决上述技术问题，本发明提供了一种扇扫雷达空时联合目 标跟踪方法，对每个扇扫周期依次执行如下步骤：

S1、从观测雷达处获得当前周期的量测信息；

S2、依据上一周期更新后的跟踪结果判断当前周期是否存在暂时 航迹和/或真实航迹；若存在，则顺序执行步骤S3，若不存在，则跳转 执行步骤S7；

S3、对当前存在的所有航迹进行一步预测，根据其上一周期的航 迹存在性概率得到当前周期的航迹存在性概率预测值，并根据扇扫雷 达空时联合系统模型得到航迹目标状态的一步预测值和量测预测值；

S4、对每一条航迹，依据所述量测预测值建立相关波门，从步骤 S1中得到的量测信息中获取与航迹相关的点迹；

S5、对每一条航迹，利用与该航迹相关的点迹采用空时联合不敏 卡尔曼滤波器进行滤波，获得多个状态估计值；

S6、对于每一条航迹，将其对应的多个状态估计值及状态的一步 预测值进行概率互联，获得目标最终状态估计值及航迹存在性概率， 更新为当前周期跟踪结果，依据航迹存在性概率判断航迹状态，输出 确认航迹并删除终结航迹；

S7、对所有未与任一航迹相关的点迹进行滤波初始化，获得航迹 目标状态值并设置初始航迹存在性概率，更新为当前周期跟踪结果。

优选地，所述步骤S7中初始化后，各个航迹初始均为暂时航迹， 所述步骤S6中依据航迹存在性概率判断航迹状态时，若一暂时航迹的 航迹存在性概率大于预先设置的确认门限t_c，则该航迹被确认为真 实航迹，并且保持确认状态直到被终结；若一暂时航迹或真实航 迹的航迹存在性概率跌至小于预先设置的航迹终结门限t_t，则该航 迹被终结。

优选地，所述步骤S3中的扇扫雷达空时联合系统模型为线性加减 速扇扫模式的雷达跟踪系统模型或三角函数加减速扇扫模式的雷达跟 踪系统模型。

优选地，所述线性加减速扇扫模式的雷达跟踪系统模型中，天线 的扇扫角度范围为[-β，β]，匀速阶段天线扫描角速度为α，加速阶段 天线加速度为

天线加速阶段扫描的角度为

天 线扫描一周时间为T_radar；设当前周期为第k周期，目标第k-1次扫描到 目标的时间为t_k-1，第k周期内目标被扫到时的目标状态向量为

方位角为A，目标方位角与目标状态向量之间关 系为

目标状态方程组的表达式为：

其中，F(k)是状态转移矩阵；v(k)是过程噪声向量；Γ(k)是噪声分 布矩阵；g(t_k-1,x(k))表示目标重访间隔与目标状态的非线性关系：

其中，T_k为目标重访间隔；f与目标扫描方向相关，若第k周期目 标顺时针扫描则f＝-1，若目标逆时针扫描则f＝1。

优选地，所述三角函数加减速扇扫模式的雷达跟踪系统模型中， 天线扫描角速度加减速以三角函数形式进行，天线加速阶段或减速阶 段扫描的角度范围为α/2，加速阶段天线扫描角速度从0增加到α，加 速阶段步长为t＝1，速度曲线采用α(t)＝α/2(1-cos(πt))，天线扫描一周 时间为T_radar；设当前周期为第k周期，目标第k-1次扫描到目标的时间为t_k-1，第k周期内目标被扫到时的目标状态向量为

方位角为A，目标方位角与目标状态向量之间关 系为

当A|<β-Δa时，目标状态方程组为：

当A|>β-Δa时，目标状态方程组为：

其中，F(k)是状态转移矩阵；v(k)是过程噪声向量；Γ(k)是噪声分 布矩阵；目标重访间隔T_k＝t_k-t_k-1，Δt为天线在第k周期内天线扫过目 标到距离目标较近的扫描边界之间的角度所用时间：

优选地，所述线性加减速扇扫模式的雷达跟踪系统模型或三角函 数加减速扇扫模式的雷达跟踪系统模型中，量测方程为：

其中，z(k)为目标的量测向量，包括目标相对所述观测雷达坐标系 原点的距离量测r_k、方位角量测θ_k；h表示目标状态与量测之间的数学 关系；W(k)为零均值的高斯白噪声序列，表示k时刻量测噪声，包括

两个分量，分别为距离、方位角量测噪声；具有协方差矩阵 E[W(k)W′(k)]＝R(k)δ_kj；

量测噪声协方差矩阵为：

其中，R_k,rr、R_k,θθ分别表示k时刻各量测噪声分量的自协方差，其 值分别为

R_k,rθ分别表示k时刻各量测噪声分量的互协方差， 假设各量测之间是互不相关的，因此各分量互协方差均为0。

优选地，所述步骤S3中根据其上一周期的航迹存在性概率得到当 前周期的航迹存在性概率预测值时，通过马尔科夫链1阶模型表述航 迹存在或不存在两种状态的转移概率：

其中，P₁₁表示上一周期航迹存在，下一周期航迹依然存在的概率； P₁₂表示上一周期航迹存在，下一周期航迹不存在的概率；P₂₁表示上一 周期航迹不存在，下一周期航迹存在的概率；P₂₂表示上一周期航迹不 存在，下一周期航迹依然不存在的概率；

已知k-1周期航迹τ存在的概率为

表示k-1周 期航迹τ存在；定义与航迹存在性相关的两种状态：

表示目标在 第k周期航迹τ存在，

表示目标在第k周期航迹τ不存在；

对航迹存在性概率进行预测，得到k周期该航迹存在性概率 的预测值为：

其中，Z^τ,k-1表示直到k-1时刻落入航迹τ波门的所有量测的集 合，

Z^τ(k)表示k时刻落入航迹τ相关波门内的

个 量测的集合，

优选地，所述步骤S3中根据扇扫雷达空时联合系统模型得到状态 的一步预测值时，选择均值周围的采样点δ作为非线性变换的输入，对 输出结果求取统计特性，得到目标状态的一步预测：

设第k-1周期中目标τ的状态估计值为

协方差为 P^τ(k-1|k-1)；依据该状态和协方差，产生一组长度为2L+1的采样点δ， 各采样点δ均匀分布在第k-1周期状态估计值

的附近，其中L 是状态向量

的维数：

其中，

是(L+λ)P^τ(k-1|k-1)的矩阵平方根的第i 列；

将上述2L+1个采样点分别作为扇扫雷达空时联合系统模型中目 标状态方程组的输入，采用迭代法求解得到2L+1组解，获得2L+1个 重访间隔及状态预测值；

目标的重访间隔是2L+1个采样点重访间隔的加权总和：

对2L+1个状态预测值进行加权，获得目标最终状态预测值及状态 预测协方差：

其中，

Γ_i,k是依据

得到的噪声分布矩阵，状态和协方差的加权权重分别 是：

α和κ控制采样点δ的传播；β与x的分布有关；

根据扇扫雷达空时联合系统模型中的量测方程，预测测量的采样 点是：

量测的预测值和相应的协方差分别为：

其中，

状态预测值及量测预测值的互协方差矩阵是：

优选地，所述步骤S3中，将2L+1个采样点分别作为扇扫雷达空 时联合系统模型中目标状态方程组的输入，采用迭代法求解得到2L+1 组解时，包括如下步骤：

以第i个采样点

为输入，若在第k-1周期里目标状态更新时间 为t_k-1，则天线扫过目标到达扫描边界时间为T₁＝T_radar·(k-1)-t_k-1，迭代法 求解状态方程组中目标重访间隔的初始值为：

T′＝2·T₁＝2·(T_radar·(k-1)-t_k-1)；

以T′为目标重访间隔，依据目标状态方程组中的第二个方程获得 迭代初始的状态预测值x′_k|k-1为：

其中，状态转移矩阵F′(k)依据T′得到：

计算第k周期目标被天线扫描到时的目标方位角：

依据目标方位角，计算天线在第k周期里从天线扫描边界到扫过 目标的时间间隔T₂；

存在时间差为ΔT＝T′-(T₁+T₂)，设置时间差的最大阈值为σ；若 |ΔT|<σ，则认为T′与目标真实的重访间隔差距可忽略不计，重访间隔

目标状态预测值为

否则对T′修正为T′＝T′-0.5·ΔT，并 返回至以T′为目标重访间隔，依据目标状态方程组中的第二个方程获 得迭代初始的状态预测值x′_k|k-1步骤，直到|ΔT|<σ；最终获得以第i个采 样点为输入获得的目标状态的一步预测值

及重访间隔

其中，计算天线在第k周期里从天线扫描边界到扫过目标的时间 间隔T₂时，T₂的计算方式与目标方位角和天线扫描方向有关，对于所述 线性加减速扇扫模式的雷达跟踪系统模型，T₂的表达式为：

对于三角函数加减速扇扫模式的雷达跟踪系统模型，若目标位于 匀速区即|A|<β-α/2，则T₂＝1+(β-Δa+f·A)/α；

若目标位于加速或减速区，则通过求解

求取T₂，若目标位于加速区则T₂＝Δt，目标位于减速区则T₂＝T_radar-Δt。

优选地，所述步骤S5中利用与该航迹相关的点迹采用空时联合不 敏卡尔曼滤波器进行滤波时，设一航迹τ有m_k个相关的点迹，分别利 用m_k个点迹应用空时联合不敏卡尔曼滤波器进行滤波，获得该航迹τ 相关的m_k个状态估计值；

其中，滤波包括获得航迹目标状态的一步预测值和量测预测值后， 计算卡尔曼增益：

依据一点迹i得到的状态更新值是预测状态与被卡尔曼增益加权 的新息之和，得到：

x^τ,i(k|k)＝x^τ(k|k-1)+K(k)[z_i(k)-z^τ(k|k-1)]；

协方差更新值为：

P^τ,i(k|k)＝P^τ(k|k-1)-K(k)P_zzK′(k)。

(三)有益效果

本发明的上述技术方案具有如下优点：本发明提供了一种扇扫雷 达空时联合目标跟踪方法，该方法采用递推循环的方式实现多目标跟 踪，采用增加了描述时间的方程对状态方程进行约束的扇扫雷达空时 联合系统模型对扇扫雷达跟踪系统进行描述，并提出了相应的滤波方 法，有效解决了传统目标跟踪方法中系统模型与扇扫雷达存在模型不 匹配的问题，提高了目标跟踪在准确性以及真实航迹确认速度。

附图说明

图1示出了本发明实施例中扇扫雷达空时联合目标跟踪方法 (STJ-LMIPDA方法)步骤示意图；

图2示出了本发明实施例中迭代法求解目标状态方程组解的步骤 示意图；

图3示出了仿真实验中构造的笛卡尔坐标系下三个目标沿直线匀 速运动的目标轨迹；

图4示出了线性加减速扇扫机制下，对目标2运用STJ-LMIPDA 方法与区分奇偶周期以2倍天线扫描周期为重访间隔运用传统 LMIPDA方法进行比较，得到的位置分量估计的均方根误差结果；

图5示出线性加减速扇扫机制下两种方法对目标2速度估计均方 根误差对比结果；

图6示出了三角函数加减速扇扫机制下两种方法对目标2位置分 量估计的均方根误差结果；

图7示出了三角函数加减速扇扫机制下两种方法对目标2速度分 量估计的均方根误差结果；

图8显示了线性加减速扇扫机制下两种方法真实航迹随扫描次数 增加的关系图；

图9显示了线性加减速扇扫机制下两种方法虚假航迹随扫描次数 增加的关系图；

图10显示了三角函数加减速扇扫机制下两种方法真实航迹随扫描 次数增加的关系图；

图11显示了三角函数加减速扇扫机制下两种方法虚假航迹随扫描 次数增加的关系图。

具体实施方式

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

如图1所示，本发明实施例提供的一种扇扫雷达空时联合目标跟 踪方法，是递推循环的过程，设当前周期为第k周期，对每个扇扫周 期依次执行如下步骤：

S1、从观测雷达处获得当前周期的量测信息；

S2、依据上一周期更新后的跟踪结果判断当前周期是否存在暂时 航迹和/或真实航迹；若存在，则顺序执行步骤S3，若不存在，则跳转 执行步骤S7；

S3、对当前存在的所有航迹，即所有暂时航迹和/或真实航迹，进 行一步预测，根据其上一周期的航迹存在性概率得到当前周期的航迹 存在性概率预测值，并根据扇扫雷达空时联合系统模型得到航迹目标 状态的一步预测值和量测预测值；

S4、对每一条航迹，依据所述量测预测值建立相关波门，从步骤 S1中得到的量测信息中获取与航迹相关的点迹；

S5、对每一条航迹，利用与该航迹相关的点迹采用空时联合不敏 卡尔曼滤波器(STJ-UKF)进行滤波，获得多个状态估计值；

S6、对于每一条航迹，将其对应的多个状态估计值及状态的一步 预测值进行概率互联，获得目标最终状态估计值及航迹存在性概率， 更新为当前周期跟踪结果，依据航迹存在性概率判断航迹状态，输出 确认航迹并从内存里删除终结航迹；

S7、对所有未与任一航迹相关的点迹进行滤波初始化，获得航迹 目标状态值并设置初始航迹存在性概率，更新为当前周期跟踪结果。

由于航迹包括暂时航迹和真实航迹，优选地，在(上一周期)所 述步骤S7中进行初始化后，各个航迹初始均为暂时航迹，所述步骤S6 中依据航迹存在性概率判断航迹状态时，预先设置门限，若一暂时航 迹的航迹存在性概率大于预先设置的确认门限t_c，则该航迹被确认 为真实航迹，并且保持确认状态直到被终结；若一暂时航迹或真 实航迹的航迹存在性概率跌至小于预先设置的航迹终结门限t_t，则 该航迹被终结，对判断终结的航迹从内存中删除，减少占用空间。

优选地，所述步骤S3中根据其上一周期的航迹存在性概率得到当 前周期的航迹存在性概率预测值时，通过马尔科夫链1阶模型表述航 迹存在或航迹不存在这两种状态的转移概率：

其中，P₁₁表示上一周期航迹存在，下一周期航迹依然存在的概率； P₁₂表示上一周期航迹存在，下一周期航迹不存在的概率；P₂₁表示上一 周期航迹不存在，下一周期航迹存在的概率；P₂₂表示上一周期航迹不 存在，下一周期航迹依然不存在的概率；P₁₁、P₁₂、P₂₁和P₂₂的具体数值 可由工程师根据实际情况进行设置。

已知k-1周期航迹τ存在的概率为

表示k-1周 期航迹τ存在；定义与航迹存在性相关的两种状态：

表示目标在 第k周期航迹τ存在，

表示目标在第k周期航迹τ不存在；

对航迹存在性概率进行预测，得到k周期该航迹存在性概率 的预测值为：

其中，Z^τ,k-1表示直到k-1时刻落入航迹τ波门的所有量测的集 合，

Z^τ(k)表示k时刻落入航迹τ相关波门内的

个 量测的集合，

考虑到扇扫雷达的通用性，本发明针对线性加减速扇扫雷达及三 角函数加减速扇扫雷达进行了建模。优选地，所述步骤S3中的扇扫雷 达空时联合系统模型为线性加减速扇扫模式的雷达跟踪系统模型或三 角函数加减速扇扫模式的雷达跟踪系统模型。

优选地，所述线性加减速扇扫模式的雷达系统模型中，天线的扇 扫角度范围为[-β，β]，匀速阶段天线扫描角速度为α，加速阶段天线 加速度为

天线加速阶段扫描的角度为

天线扫 描一周时间为T_radar；设当前周期为第k周期，目标第k-1次扫描到目标 的时间为t_k-1，第k周期内目标被扫到时的目标状态向量为

方位角为A，目标方位角与目标状态向量之间关 系为

目标状态方程组的表达式为：

其中，F(k)是状态转移矩阵；v(k)是过程噪声向量；Γ(k)是噪声分 布矩阵；g(t_k-1,x(k))表示目标重访间隔与目标状态的非线性关系：

其中，T_k为目标重访间隔；

f与目标扫描方向相关， 若第k周期目标顺时针扫描则f＝-1，若目标逆时针扫描则f＝1。

优选地，所述三角函数加减速扇扫模式的雷达跟踪系统模型中， 天线扫描角速度加减速以三角函数形式进行，天线加速阶段或减速阶 段扫描的角度范围为α/2，加速阶段天线扫描角速度从0增加到α，加 速阶段步长为t＝1，速度曲线采用α(t)＝α/2(1-cos(πt))，天线扫描一周 时间为T_radar；设当前周期为第k周期，目标第k-1次扫描到目标的时间为t_k-1，第k周期内目标被扫到时的目标状态向量为

方位角为A，目标方位角与目标状态向量之间关 系为

当A|<β-Δa时，目标状态方程组为：

当A|>β-Δa时，目标状态方程组为：

其中，F(k)是状态转移矩阵；v(k)是过程噪声向量；Γ(k)是噪声分 布矩阵；目标重访间隔T_k＝t_k-t_k-1，Δt为天线在第k周期内天线扫过目 标到距离目标较近的扫描边界之间的角度所用时间：

优选地，线性加减速扇扫模式的雷达跟踪系统模型或三角函数加 减速扇扫模式的雷达跟踪系统模型中，两种扇扫方式的量测方程是一 致的，量测方程的表达式为：

其中，z(k)为目标的量测向量，包括目标相对所述观测雷达坐标系 原点的距离量测r_k、方位角量测θ_k；h表示目标状态与量测之间的数学 关系；W(k)为零均值的高斯白噪声序列，表示k时刻量测噪声，包括

两个分量，分别为距离、方位角量测噪声；具有协方差矩阵 E[W(k)W′(k)]＝R(k)δ_kj；

量测噪声协方差矩阵为：

其中，R_k,rr、R_k,θθ分别表示k时刻各量测噪声分量的自协方差，其 值分别为

R_k,rθ分别表示k时刻各量测噪声分量的互协方差， 假设各量测之间是互不相关的，因此各分量互协方差均为0。

优选地，所述步骤S3中根据扇扫雷达空时联合系统模型得到状态 的一步预测值时，依据UT变换，选择均值周围的采样点δ作为非线 性变换的输入，对输出结果求取统计特性，得到目标状态的一步预测。

设第k-1周期中目标τ的状态估计值为

协方差为 P^τ(k-1|k-1)；依据该状态和协方差，产生一组长度为2L+1的采样点δ， 各采样点δ均匀分布在第k-1周期，即上一周期，更新后的状态估计值

的附近，其中L是状态向量

的维数；

其中，

是(L+λ)P^τ(k-1|k-1)的矩阵平方根的第i 列；

将上述2L+1个采样点分别作为扇扫雷达空时联合系统模型中目 标状态方程组的输入，采用迭代法求解得到2L+1组解，获得2L+1个 重访间隔及状态预测值；

目标的重访间隔是2L+1个采样点重访间隔的加权总和：

对2L+1个状态预测值进行加权，获得目标最终状态预测值及状态 预测协方差：

其中，

Γ_i,k是依据

得到的噪声分布矩阵，状态和协方差的加权权重分别 是：

α和κ控制采样点δ的传播；β与x的分布有关；

根据扇扫雷达空时联合系统模型中的量测方程，预测测量的采样 点是：

量测的预测值和相应的协方差分别为：

其中，

状态预测值及量测预测值的互协方差矩阵是：

优选地，如图2所示，所述步骤S3中，将2L+1个采样点分别作 为扇扫雷达空时联合系统模型中目标状态方程组的输入，采用迭代法 求解得到2L+1组解时，对于线性加减速扇扫模式的雷达跟踪系统模 型，包括如下步骤：

1)以第i个采样点

为输入，若在第k-1周期里目标状态更新 时间为t_k-1，则天线扫过目标到达扫描边界时间为T₁＝T_radar·(k-1)-t_k-1，迭 代法求解状态方程组中目标重访间隔的初始值为：

T′＝2·T₁＝2·(T_radar·(k-1)-t_k-1)；

2)以T′为目标重访间隔，依据目标状态方程组中的第二个方程获 得迭代初始的状态预测值x′_k|k-1为：

其中，状态转移矩阵F′(k)依据T′得到：

3)计算第k周期目标被天线扫描到时的目标方位角：

其中，x′_k|k-1(1)、x′_k|k-1(3)分别表示4x1维向量x′_k|k-1的第1个值和第3 个值；

4)依据目标方位角，计算天线在第k周期里从天线扫描边界到扫 过目标的时间间隔T₂，T₂的计算方式与目标方位角和天线扫描方向有 关，表达式为：

5)存在时间差为ΔT＝T′-(T₁+T₂)，设置时间差的最大阈值为σ，若 |ΔT|<σ，则认为T′与目标真实的重访间隔差距非常小可以忽略不计， 重访间隔

目标状态预测值为

否则对T′修正为T′＝T′-0.5·ΔT，并返回至以T′为目标重访间隔，依 据目标状态方程组中的第二个方程获得迭代初始的状态预测值x′_k|k-1步 骤，即在修正T′后重新进行步骤2)到5)，直到|ΔT|<σ；最终获得以第 i个采样点为输入获得的目标状态的一步预测值

及重访间隔

优选地，所述步骤S3中，将2L+1个采样点分别作为扇扫雷达空 时联合系统模型中目标状态方程组的输入，采用迭代法求解得到2L+1 组解时，对于三角函数加减速扇扫模式的雷达跟踪系统模型，包括如 下步骤：

1)以第i个采样点

为输入，若在第k-1周期里目标状态更 新时间为t_k-1，则天线扫过目标到达扫描边界时间为T₁＝T_radar·(k-1)-t_k-1， 目标重访间隔的初始值为：

T′＝2·T₁＝2·(T_radar·(k-1)-t_k-1)；

2)以T′为目标重访间隔，依据状态方程组中的第二个方程获得迭 代初始的状态预测值x′_k|k-1，即为：

其中，状态转移矩阵F′(k)依据T′得到，即：

3)计算第k周期目标被天线扫描到时目标方位角：

4)判断目标所在的区域，依据目标方位角，计算在天线k周期内， 从天线扫描边缘扫到目标的时间T₂：

若目标位于匀速区即|A|<β-α/2，则T₂＝1+(β-Δa+f·A)/α；

若目标位于加速或减速区，求取T₂需要求解

鉴于该方程是非线性的，同样可利用迭代 法求解该方程的解，求解过程如下：

1、天线在加速或减速区域扫描的时间为1s，令迭代初始值T₂′＝1， 由于T₂′的值存在误差，存在角度误差

2、比较角度误差的阈值是σ₁。它可能是一个非常小的值，可按需 自行设定。如果角度误差小于阈值即|ε|＜σ₁，那么认为角度误差可以忽 略，即Δt＝Δt′。否则需要对Δt′进行调整，令

重新执行第1 步骤，计算角度误差并对时间Δt′进行调整，直至角度误差满足|ε|＜σ₁， 若目标位于加速区那么T₂＝Δt，目标位于减速区T₂＝T_radar-Δt。

5)依据目标方位角A获得时间T₂，求取时间差为ΔT＝T′-(T₁+T₂)， 设置时间差的最大阈值为σ₂。若|ΔT|<σ₂则认为T′与目标真实的重访间 隔差距非常小可以忽略不计，重访间隔

目标状态预测值为

否则对T′修正为T′＝T′-0.5·ΔT，并返回至以T′为目标重访间 隔，依据状态方程组中的第二个方程获得迭代初始的状态预测值x′_k|k-1步骤，即在修正T′后重新进行步骤2)到5)，直到|ΔT|<σ₂。至此我 们得到第i个采样点为输入获得的目标状态的一步预测值

及重访间 隔

优选地，步骤S4中对每一条航迹，依据所述量测预测值建立相关 波门，从步骤S1中得到的量测信息中获取与航迹相关的观测点时，对 k时刻雷达接收所有的点迹进行波门相关，相关波门依据被跟踪目标在 当前时刻所在位置的预测值为中心，确立一个量测向量可能出现的范 围的决策门限，所有落入该门限的量测值均考虑为候选的量测，而在 门限值之外的量测则认为是杂波。相关波门通常用来滤除杂波，找到 所有可能来自于目标的量测值。常用的相关波门有矩形波门、环形波 门、椭圆(球)波门、极坐标系下的扇形波门等，这里我们将使用椭 圆波门。

对于航迹对于传感器在第k周期内收到的任何点迹

若存在：

则认为该点迹落在目标τ的波门，那么存在一个假设即，点迹i为 目标τ的量测值，并且

将用于更新目标τ的状态估计值。给定目标的 波门概率γ。

一条航迹在k时刻落入波门内的点迹可能有多个，也可能没有。

优选地，若k时刻有m_k个点迹落入航迹τ的相关波门，对于其中 的第i个点迹z_i(k)(这里i与前文中的i意义不同，表征落入相关波门 内的某一点迹，前文中指代UT变换过程中某一采样点)，利用点迹z_i(k) 对航迹τ进行状态更新，得到该点迹对航迹τ的状态更新值

同时得到该点迹来自航迹τ所跟踪目标的似然函数

进一步优选地，所述步骤S5中对每一条航迹，利用与该航迹相关 的点迹进行滤波过程，获得多个状态估计值时，设一航迹τ有

个相 关的点迹，分别利用

个点迹应用空时联合不敏卡尔曼滤波器 (STJ-UKF)进行滤波，获得该航迹τ相关的

个状态估计值。STJ-UKF 方法前半部分内容即状态一步预测过程，这里不再重复描述。

滤波还包括获得航迹目标状态的一步预测值和量测预测值后，计 算卡尔曼增益：

依据某一点迹i，即z_i(k)，得到的状态更新值是预测状态与被卡尔 曼增益加权的新息之和，得到：

x^τ,i(k|k)＝x^τ(k|k-1)+K(k)[z_i(k)-z^τ(k|k-1)]；

相应的协方差更新值为：

P^τ,i(k|k)＝P^τ(k|k-1)-K(k)P_zzK′(k)。

优选地，步骤S6中对于每一条航迹，将其对应的多个状态估计值 及状态的一步预测值进行概率互联时，如果点迹i落在目标τ的波门内， 则表示点迹i与目标τ相关，点迹i可能是航迹τ的观测值，定义事件

为波门内的第i个点迹来自目标τ，其他点迹对该航迹为杂波。

点迹i来自航迹τ的先验概率为：

其中

分别表示目标被检测的概率及波门概率，可由工程师根 据需要进行设定。

为航迹存在性概率的预测值，前文已经介 绍。

N表示目标数量，ρ_k,i是目标的杂波密度，

是在z_k(i)是目标τ的 量测的假设下的预测测量概率密度函数，其由下式定义：

在多目标跟踪系统中，目标τ波门中的点迹i可能来自另一个目标。 在线性化方法下，通过排除目标τ的所有可能来源来修改量测i所在位 置的杂波密度。跟踪τ的验证门中的测量i的先验散射体测量密度，由 Ω^τ(i)表示，可以表示为：

结合综合概率数据互联算法(IPDA)，目标存在的后验概率被修改 为：

其中，

表示落入目标门的测量总数，量测i是目标τ(i＝0表示该目标没 有量测)的量测的后验概率由下式给出：

其中，

至此获得落入目标τ的波门内的各个量测属于该目标的后验概率， 以及航迹τ存在的后验概率值。

假设有

个点迹落在目标τ的波门内，则存在

种假设：目标 没有量测被探测到

个点迹均为杂波；这

个点迹分别是目标τ的量 测值，其他

个点迹为杂波。后验跟踪状态估计是：

其中

是以点迹

作为目标的量测值对目标进行状态估计得 到的估计值，i＝0表示目标没有量测被探测，目标状态的预测值，即

状态估计协方差矩阵由下式给出：

目标τ的更新时间是

至此，目标τ的状态更新已经完 成。

对所有航迹均进行以上航迹更新过程，依据航迹存在性概率对航 迹进行管理。某些点迹可能未落入任何一条航迹的相关波门内，则需 要进行航迹初始化。

考虑到未被用来更新航迹的点迹可能代表着出现新航迹，对其执 行航迹初始化，可获得新航迹目标初始状态。优选地，步骤S7中对所 有未与任一航迹相关的点迹，依据观测雷达当前周期和上一周期所获 得的未与航迹相关的点迹进行滤波初始化时，利用未与航迹相关的点 进行航迹初始化获得新航迹，以第一周期和第二周期的点迹为例，第 一次扫描雷达得到n₁个点迹，第二次扫描雷达得到n₂个点迹，利用第一 次扫描时刻的第i(1≤i≤n₁)个点迹和第二次扫描时刻的第j(1≤j≤n₂)个 点迹获得第τ(1≤τ≤n₁×n₂)条航迹在笛卡尔坐标系下第二次扫描时刻的 状态向量。

对应的状态协方差为：

其中下标1、2分别表示第一周期和第二周期对应的值。T表示雷 达扫描周期。

和

是笛卡尔坐标下的量测沿x，y方向的位置量测信 息，是通过无偏量测转换方法将雷达所获得量测(径向距离、方位角) 转换到笛卡尔坐标系下得到的转换位置量测，转换公式为：

其中r_k，θ_k是从雷达获取的距离，方位角量测；

是转换后 得到的沿x，y方向的笛卡尔坐标量测，

是转换后的量测向量；μ_θ是 去偏系数，可通过方位角测量噪声方差

求得：

对应的协方差矩阵为：

其中各分量的自协方差及互协方差：

上标“c”代表与转换量测相关的向量、矩阵和函数。初始化的状 态误差协方差为：

T可以通过两个量测的时间差得到。

航迹存在性概率初始值可由实际需要进行设定。

综上，本申请提供了一种扇扫雷达空时联合目标跟踪方法，该方 法包括如下优势：

(1)针对扇扫雷达重访间隔不确定导致系统模型不匹配的问题， 通过增加时间约束方程，分别针对线性加减速扇扫雷达与三角函数扇 扫机制雷达建立了适用于扇扫雷达的状态方程组，准确描述两种扫描 机制下目标状态随时间变化的关系。

(2)针对新建立的系统模型，研究了对应的空时联合状态估计滤 波方法STJ-UKF，通过迭代法求取非线性状态方程组的解，同时得到 目标状态预测值及重访间隔，该过程可以看作一个广义的非线性函数， 因此可应用不敏变换计算其输出变量的统计特性，最后获得目标状态 估计结果。

(3)将滤波方法应用于线性多目标联合概率密度互联(LMIPDA) 框架，研究了适应于的扇扫雷达的多目标跟踪方法STJ-LMIPDA方法

为验证本发明的效果，利用仿真数据进行蒙特卡洛实验。仿真实 验监视区域内设置三个目标，仿真试验中的目标沿直线匀速运动，运 动轨迹如图3所示。三个目标初始状态如表1所示。

表1 3个目标的初始状态

进行100次蒙脱卡洛仿真，每次模拟运行有100次扫描。监视区 域中的杂波密度为6.0×10^-6/scan/m²且杂波满足均匀的泊松分布。目标初 始位置和速度总结在表1中给出。过程噪声设置为ds＝0.01m/s²，范围 和方位角测量标准偏差分别为σ_r＝50m，σ_θ＝1deg。天线匀速扫描阶段扫 描角速度为α＝π/6rad/s，天线扫描范围为[-π/3,π/3]，迭代过程时间差 的阈值为σ₂＝0.001s。

线性加减速扇扫方式加速阶段天线角度加速度为α＝π/6rad/s²；三 角函数加减速扇扫方式加速时间为1s，角度误差的阈值为σ₁＝10^-8rad。

将区分奇偶扫描周期，以两倍天线扫描周期作为重访间隔的传统 LMIPDA和本文所提方法进行比较用于评估STJ-LMIPDA的性能。

图4示出了线性加减速扇扫方式两种方法对目标2位置滤波结果 对应的估计均方根误差对比，图5示出了线性加减速扇扫方式两种方 法对目标2速度滤波结果对应的估计均方根误差对比。图6示出了三 角函数加减速扇扫方式两种方法对目标2位置滤波结果对应的估计均 方根误差对比，图7示出了三角函数加减速扇扫方式两种方法对目标2 速度滤波结果对应的估计均方根误差对比。从图4到图6中可以看到 不管哪种扫描方式STJ-LMIPDA方法对目标估计精度有显著提高。这 是由于通过对扇扫雷达系统模型的准确建立，通过迭代法获得目标的 状态预测值及重访间隔较为精确，综合考虑了目标运动及目标方位角 对重访间隔的影响，使得数据互联较为准确。在精准建模的基础上， 状态滤波结果也较为准确。但本发明中示出的方法在提供状态估计结 果的同时，参照图8及图10可以看出本发明所提方法还能够显著提高 真实航迹确认速度；参照图9及图11可以看出除开始几个周期本文所 提方法虚警稍高外，后续周期显著降低了目标跟踪过程的虚警率。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而 非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领 域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技 术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修 改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方 案的精神和范围。

Claims (10)

1.一种扇扫雷达空时联合目标跟踪方法，其特征在于，对每个扇扫周期依次执行如下步骤：

S1、从观测雷达处获得当前周期的量测信息；

S2、依据上一周期更新后的跟踪结果判断当前周期是否存在暂时航迹和/或真实航迹；若存在，则顺序执行步骤S3，若不存在，则跳转执行步骤S7；

S3、对当前存在的所有航迹进行一步预测，根据其上一周期的航迹存在性概率得到当前周期的航迹存在性概率预测值，并根据扇扫雷达空时联合系统模型得到航迹目标状态的一步预测值和量测预测值；

S4、对每一条航迹，依据所述量测预测值建立相关波门，从步骤S1中得到的量测信息中获取与航迹相关的点迹；

S5、对每一条航迹，利用与该航迹相关的点迹采用空时联合不敏卡尔曼滤波器进行滤波，获得多个状态估计值；

S6、对于每一条航迹，将其对应的多个状态估计值及状态的一步预测值进行概率互联，获得目标最终状态估计值及航迹存在性概率，更新为当前周期跟踪结果，依据航迹存在性概率判断航迹状态，输出确认航迹并删除终结航迹；

S7、对所有未与任一航迹相关的点迹进行滤波初始化，获得航迹目标状态值并设置初始航迹存在性概率，更新为当前周期跟踪结果。

2.根据权利要求1所述的扇扫雷达空时联合目标跟踪方法，其特征在于：所述步骤S7中初始化后，各个航迹初始均为暂时航迹，所述步骤S6中依据航迹存在性概率判断航迹状态时，若一暂时航迹的航迹存在性概率大于预先设置的确认门限t_c，则该航迹被确认为真实航迹，并且保持确认状态直到被终结；若一暂时航迹或真实航迹的航迹存在性概率跌至小于预先设置的航迹终结门限t_t，则该航迹被终结。

3.根据权利要求1所述的扇扫雷达空时联合目标跟踪方法，其特征在于：所述步骤S3中的扇扫雷达空时联合系统模型为线性加减速扇扫模式的雷达跟踪系统模型或三角函数加减速扇扫模式的雷达跟踪系统模型。

4.根据权利要求3所述的扇扫雷达空时联合目标跟踪方法，其特征在于：

所述线性加减速扇扫模式的雷达跟踪系统模型中，天线的扇扫角度范围为[-β，β]，匀速阶段天线扫描角速度为α，加速阶段天线加速度为

天线加速阶段扫描的角度为

天线扫描一周时间为T_radar；设当前周期为第k周期，目标第k-1次扫描到目标的时间为t_k-1，第k周期内目标被扫到时的目标状态向量为

方位角为A，目标方位角与目标状态向量之间关系为

目标状态方程组的表达式为：

其中，F(k)是状态转移矩阵；v(k)是过程噪声向量；Γ(k)是噪声分布矩阵；g(t_k-1,x(k))表示目标重访间隔与目标状态的非线性关系：

其中，T_k为目标重访间隔；f与目标扫描方向相关，若第k周期目标顺时针扫描则f＝-1，若目标逆时针扫描则f＝1。

5.根据权利要求3所述的扇扫雷达空时联合目标跟踪方法，其特征在于：

所述三角函数加减速扇扫模式的雷达跟踪系统模型中，天线扫描角速度加减速以三角函数形式进行，天线加速阶段或减速阶段扫描的角度范围为α/2，加速阶段天线扫描角速度从0增加到α，加速阶段步长为t＝1，速度曲线采用α(t)＝α/2(1-cos(πt))，天线扫描一周时间为T_radar；设当前周期为第k周期，目标第k-1次扫描到目标的时间为t_k-1，第k周期内目标被扫到时的目标状态向量为

方位角为A，目标方位角与目标状态向量之间关系为

当A|<β-Δa时，目标状态方程组为：

当A|>β-Δa时，目标状态方程组为：

其中，F(k)是状态转移矩阵；v(k)是过程噪声向量；Γ(k)是噪声分布矩阵；目标重访间隔T_k＝t_k-t_k-1，Δt为天线在第k周期内天线扫过目标到距离目标较近的扫描边界之间的角度所用时间：

6.根据权利要求4或5所述的扇扫雷达空时联合目标跟踪方法，其特征在于：

所述线性加减速扇扫模式的雷达跟踪系统模型或三角函数加减速扇扫模式的雷达跟踪系统模型中，量测方程为：

其中，z(k)为目标的量测向量，包括目标相对所述观测雷达坐标系原点的距离量测r_k、方位角量测θ_k；h表示目标状态与量测之间的数学关系；W(k)为零均值的高斯白噪声序列，表示k时刻量测噪声，包括

两个分量，分别为距离、方位角量测噪声；具有协方差矩阵E[W(k)W′(k)]＝R(k)δ_kj；

量测噪声协方差矩阵为：

其中，R_k,rr、R_k,θθ分别表示k时刻各量测噪声分量的自协方差，其值分别为

R_k,rθ分别表示k时刻各量测噪声分量的互协方差，假设各量测之间是互不相关的，因此各分量互协方差均为0。

7.根据权利要求6所述的扇扫雷达空时联合目标跟踪方法，其特征在于：所述步骤S3中根据其上一周期的航迹存在性概率得到当前周期的航迹存在性概率预测值时，通过马尔科夫链1阶模型表述航迹存在或不存在两种状态的转移概率：

其中，P₁₁表示上一周期航迹存在，下一周期航迹依然存在的概率；P₁₂表示上一周期航迹存在，下一周期航迹不存在的概率；P₂₁表示上一周期航迹不存在，下一周期航迹存在的概率；P₂₂表示上一周期航迹不存在，下一周期航迹依然不存在的概率；

已知k-1周期航迹τ存在的概率为

表示k-1周期航迹τ存在；定义与航迹存在性相关的两种状态：

表示目标在第k周期航迹τ存在，

表示目标在第k周期航迹τ不存在；

对航迹存在性概率进行预测，得到k周期该航迹存在性概率的预测值为：

其中，Z^τ,k-1表示直到k-1时刻落入航迹τ波门的所有量测的集合，

Z^τ(k)表示k时刻落入航迹τ相关波门内的

个量测的集合，

8.根据权利要求7所述的扇扫雷达空时联合目标跟踪方法，其特征在于，所述步骤S3中根据扇扫雷达空时联合系统模型得到状态的一步预测值时，选择均值周围的采样点δ作为非线性变换的输入，对输出结果求取统计特性，得到目标状态的一步预测：

设第k-1周期中目标τ的状态估计值为

协方差为P^τ(k-1|k-1)；依据该状态和协方差，产生一组长度为2L+1的采样点δ，各采样点δ均匀分布在第k-1周期状态估计值

的附近，其中L是状态向量

的维数：

其中，

是(L+λ)P^τ(k-1|k-1)的矩阵平方根的第i列；

将上述2L+1个采样点分别作为扇扫雷达空时联合系统模型中目标状态方程组的输入，采用迭代法求解得到2L+1组解，获得2L+1个重访间隔及状态预测值；

目标的重访间隔是2L+1个采样点重访间隔的加权总和：

对2L+1个状态预测值进行加权，获得目标最终状态预测值及状态预测协方差：

其中，

Γ_i,k是依据

得到的噪声分布矩阵，状态和协方差的加权权重分别是：

α和κ控制采样点δ的传播；β与x的分布有关；

根据扇扫雷达空时联合系统模型中的量测方程，预测测量的采样点是：

量测的预测值和相应的协方差分别为：

其中，

状态预测值及量测预测值的互协方差矩阵是：

9.根据权利要求8所述的扇扫雷达空时联合目标跟踪方法，其特征在于，所述步骤S3中，将2L+1个采样点分别作为扇扫雷达空时联合系统模型中目标状态方程组的输入，采用迭代法求解得到2L+1组解时，包括如下步骤：

以第i个采样点

为输入，若在第k-1周期里目标状态更新时间为t_k-1，则天线扫过目标到达扫描边界时间为T₁＝T_radar·(k-1)-t_k-1，迭代法求解状态方程组中目标重访间隔的初始值为：

T′＝2·T₁＝2·(T_radar·(k-1)-t_k-1)；

以T′为目标重访间隔，依据目标状态方程组中的第二个方程获得迭代初始的状态预测值x′_k|k-1为：

其中，状态转移矩阵F′(k)依据T′得到：

计算第k周期目标被天线扫描到时的目标方位角：

依据目标方位角，计算天线在第k周期里从天线扫描边界到扫过目标的时间间隔T₂；

存在时间差为ΔT＝T′-(T₁+T₂)，设置时间差的最大阈值为σ；若|ΔT|<σ，则认为T′与目标真实的重访间隔差距可忽略不计，重访间隔

目标状态预测值为

否则对T′修正为T′＝T′-0.5·ΔT，并返回至以T′为目标重访间隔，依据目标状态方程组中的第二个方程获得迭代初始的状态预测值x′_k|k-1步骤，直到|ΔT|<σ；最终获得以第i个采样点为输入获得的目标状态的一步预测值

及重访间隔

其中，计算天线在第k周期里从天线扫描边界到扫过目标的时间间隔T₂时，T₂的计算方式与目标方位角和天线扫描方向有关，对于所述线性加减速扇扫模式的雷达跟踪系统模型，T₂的表达式为：

对于三角函数加减速扇扫模式的雷达跟踪系统模型，若目标位于匀速区即|A|<β-α/2，则T₂＝1+(β-Δa+f·A)/α；

若目标位于加速或减速区，则通过求解

求取T₂，若目标位于加速区则T₂＝Δt，目标位于减速区则T₂＝T_radar-Δt。

10.根据权利要求1所述的扇扫雷达空时联合目标跟踪方法，其特征在于：所述步骤S5中利用与该航迹相关的点迹采用空时联合不敏卡尔曼滤波器进行滤波时，设一航迹τ有m_k个相关的点迹，分别利用m_k个点迹应用空时联合不敏卡尔曼滤波器进行滤波，获得该航迹τ相关的m_k个状态估计值；

其中，滤波包括获得航迹目标状态的一步预测值和量测预测值后，计算卡尔曼增益：

依据一点迹i得到的状态更新值是预测状态与被卡尔曼增益加权的新息之和，得到：

x^τ,i(k|k)＝x^τ(k|k-1)+K(k)[z_i(k)-z^τ(k|k-1)]；

协方差更新值为：

P^τ,i(k|k)＝P^τ(k|k-1)-K(k)P_zzK′(k)。

Images