CN110376580A - 一种异步多目标跟踪性能驱动的异构雷达网资源分配方法 - Google Patents

Abstract

本发明公开了一种异步多目标跟踪性能驱动的异构雷达网资源分配方法，首先建立异构雷达网，并建立目标运动模型和异步观测模型；然后基于融合时间间隔内的异步量测值集合将异构雷达网异质资源分配问题制定为凸问题，并根据其待优化资源变量可分为若干个互不重叠的块的特殊结构设计相应的块坐标下降方法搜索其最优解；基于资源分配的最优解和来自多个雷达对目标的异步测量值集合采用最大似然法估计目标的扩展状态，并将该状态作为卡尔曼滤波器的输入最终得到目标扩展状态的滤波值。该方法适用于异构雷达网有限异质发射资源条件下的资源分配，可以提高多目标跟踪精度。

Description

一种异步多目标跟踪性能驱动的异构雷达网资源分配方法

技术领域

本发明属于雷达技术领域，尤其涉及一种异步多目标跟踪性能驱动的异构雷达网资源分配方法。

背景技术

异构雷达网通常由不同工作模式和不同位置的雷达站组成，由于这些雷达站可以从不同视角观测目标，故相比单站雷达，异构雷达网可以获得更高的目标跟踪精度；通常，异构雷达网的目标跟踪精度与雷达数目，雷达带宽，信噪比有关；而信噪比又受雷达的发射资源、雷达目标散射截面积的影响；直观地来讲，通过提升异构雷达网的发射资源可以提高目标信噪比进而改善目标跟踪精度；然而，在实际应用中，异构雷达网中各节点的发射资源通常是有限的，比如机载、星载雷达装备的电池能量有限，另外融合中心的计算负载能力也是有限的；因此，在各雷达站资源有限的条件下为异构雷达网建立合理的资源分配机制对提升多目标跟踪精度具有重大意义。

现有的多目标跟踪性能驱动的雷达网资源分配方法(在各雷达站资源受限的条件下最大化目标跟踪精度)大多基于各雷达站的工作模式相同(同构)且同时地观测目标(同步)的假设，然而在实际应用中，由于各雷达站开机时间，采样频率的差异，其量测值通常是异步的，另外，同构的假设限制了可进行信息融合的雷达种类，造成了已有雷达的浪费，故很难满足实际目标跟踪精度需求；在异构雷达网中，由于各节点工作模式不同，所以待优化的发射资源变量也是不同的，例如，集中式MIMO雷达的发射功率可调但总发射功率有限，相控阵雷达的驻留时间可调但总驻留时间有限；这些不同类型的资源变量在跟踪性能度量函数(贝叶斯克拉美罗下界)中互相耦合，若采用传统方法按每个目标量测值的到达时刻顺序地计算其贝叶斯克拉美罗下界，所产生的资源分配问题是非凸的，因此只能获得次优解；综上，现有的资源分配方法不能直接应用于异构雷达网。

发明内容

为了解决现有技术中存在的上述问题，本发明提供了一种异步多目标跟踪性能驱动的异构雷达网资源分配方法。

本发明要解决的技术问题通过以下技术方案实现：

一种异步多目标跟踪性能驱动的异构雷达网资源分配方法，包括以下步骤：

步骤1：建立异构雷达网，所述异构雷达网包含个雷达站和一个融合中心；N、和均为大于0的正整数，其中，表示集中式MIMO雷达个数，表示相控阵雷达个数，表示机械扫描雷达个数；所述N个雷达站的探测区域中存在Q个目标，Q为大于零的正整数；

初始化：令t_k表示第k个融合时刻，k的初始值为1，k∈{1,2,…,K}，t_K为预先设定的最大融合时刻，K为大于0的正整数，(t_k,t_k+1)表示第k个融合时间间隔；

步骤2：设定所述异构雷达网中目标q为匀速直线运动，q∈{1,2,…,Q}，并设定t_k时刻目标q的扩展状态向量为

步骤3：异构雷达网中N个雷达异步的对目标q进行观测，并将第k个融合时间间隔内N个雷达站得到的所有异步量测值发送至融合中心；

步骤4：根据第k个融合时间间隔内N个雷达对目标q的异步观测集合计算目标q在第k+1个融合时刻的归一化贝叶斯克拉美罗下界并将其作为目标q的跟踪精度的指标函数，在每个雷达的发射总资源有限的约束条件下，最小化通过将发射资源变量分为个互不重叠的块,并采用块坐标下降方法对上述凸问题求解，进而获得第k+1个融合时刻的最优分配结果

步骤5：根据所述第k个融合时间间隔内异构雷达网中目标q的异步量测值集合Z_q,k和第k个融合时刻的最优资源分配结果采用最大似然法计算目标q的扩展状态在融合时刻t_k+1的估计值并计算其估计协方差矩阵

步骤6：将目标q的扩展状态的最大似然估计值及其对应的估计协方差矩阵作为卡尔曼滤波器的输入从而得到目标q的扩展状态在融合时刻t_k+1的最终状态估计值及其对应的协方差矩阵

步骤7：令k＝k+1，返回步骤2，直到得到第1个融合时刻的异质资源分配结果至第K个融合时刻的异质资源分配结果以及目标q在第1个融合时刻的状态估计值及其对应的协方差矩阵至第K个融合时刻的状态估计值及其对应的协方差矩阵q∈{1,2,…,Q}。

本发明的有益效果：

传统方法根据雷达量测值的到达时间顺序地计算贝叶斯克拉美罗下界(BCRLB)，在这种情况下产生的资源分配问题是非凸的，只能获得次优的资源分配结果；本发明通过引入融合时间间隔，并设定待优化资源变量在融合时刻可调，而在融合时间间隔内固定，然后基于异构雷达网在融合时间间隔内的异步量测值集合计算目标在融合时刻的贝叶斯克拉美罗下界；在此情况下，本发明将多目标跟踪性能驱动的异质资源分配问题制定为凸问题；通过将待优化的异质资源变量划分为若干个相互独立的块，并采用块坐标下降算法对上述凸问题进行求解可以快速获得资源分配的最优解。

以下将结合附图及实施例对本发明做进一步详细说明。

附图说明

图1是本发明实施例提供的一种异步多目标跟踪性能驱动的异构雷达网最优资源分配方法流程图；

图2是本发明实施例提供的异构雷达网探测区域内的目标运动场景示意图；

图3是本发明实施例提供的目标q在第k个融合时间间隔的异步采样示意图；

图4是本发明实施例提供的雷达2对目标1的雷达散射截面积；

图5是本发明实施例提供的异质资源分配方法的多目标跟踪性能；

图6是本发明实施例提供的目标q＝1的资源分配结果；

图7是本发明实施例提供的目标q＝2的资源分配结果；

图8是本发明实施例提供的每个目标的跟踪贝叶斯克拉美罗下界。

具体实施方式

下面结合具体实施例对本发明做进一步详细的描述，但本发明的实施方式不限于此。

实施例一

请参见图1，图1是本发明实施例提供的一种异步多目标跟踪性能驱动的异构雷达网最优资源分配方法流程图，包括以下步骤：

步骤1：建立异构雷达网，所述异构雷达网包含个雷达站和一个融合中心，N、和均为大于0的正整数，其中，表示集中式MIMO雷达个数，表示相控阵雷达个数，表示机械扫描雷达个数，N个雷达站的探测区域中存在Q个目标，Q为大于零的正整数；

初始化：令t_k表示第k个融合时刻，k的初始值为1，k∈{1,2,…,K}，t_K为预先设定的最大融合时刻，K为大于0的正整数，(t_k,t_k+1)表示第k个融合时间间隔；

具体地，设定雷达1的位置，并以雷达1正东方向70km，正南方向0.5km作为原点O，以正北方向为Y轴、正东方向为X轴建立平面直角坐标系。

请参见图2，图2是本发明实施例提供的异构雷达网探测区域内的目标运动场景示意图；异构雷达网包含个雷达站和一个融合中心，N、和均为大于0的正整数；其中，表示集中式MIMO雷达个数，表示相控阵雷达个数，表示机械扫描雷达个数；本实施例中N取值为8，取值为4；取值为3；取值为1；设定N个雷达站的探测区域中存在Q个目标，Q为大于零的正整数，本实施例中Q取值为2；令t_k表示第k个融合时刻，k表示融合时刻的索引，k的初始值为1，k∈{1,2,…,K}，t_K为预先设定的最大融合时刻，K为大于0的正整数，本实施例K取为20；(t_k,t_k+1)表示第k个融合时间间隔，T₀＝(t_k,t_k+1)，本实施例T₀＝6s，t_K＝121s。

步骤2：设定异构雷达网中目标q为匀速直线运动，q∈{1,2,…,Q}，并设定t_k时刻目标q的扩展状态向量为

具体地，目标q在融合时刻t_k的状态为：

上式中，表示目标q在融合时刻t_k的位置在X轴方向上的分量，表示目标q在融合时刻t_k的速度在X轴方向上的分量；表示目标q在融合时刻t_k的位置在Y轴方向上的分量，表示目标q在融合时刻t_k的速度在Y轴方向上的分量；本实施例目标1的初始状态为：目标2的初始状态为：

目标q的运动模型为：

式中，表示目标q在融合时刻t_k的状态转移矩阵，其中，I₂表示2维单位矩阵，表示直积操作，T₀表示融合时间间隔，本实施例T₀＝6s；式(2)中，表示目标q在融合时刻t_k+1的状态，表示目标q在融合时刻t_k的状态，表示目标q在融合时刻t_k的过程噪声向量，服从均值为零，协方差为的高斯分布。

式中表示目标q在融合时刻t_k的过程噪声的协方差，表示目标q的运动过程噪声强度，本实施例I₂表示2维单位矩阵，表示直积操作，T₀表示融合时间间隔；

目标q的雷达散射截面积转移模型为：

式中，表示通道状态向量，表示在融合时刻t_k雷达i对目标q的雷达散射截面积，上标T表示转置；本实施例中，目标1在融合时刻t₁的初始通道状态向量设置为目标2在融合时刻t₁的初始通道状态向量设置为表示目标q在融合时刻t_k的高斯白噪声，协方差为表示目标q在融合时刻t_k的过程噪声的协方差；(4)中，；将目标q在融合时刻t_k的状态和雷达散射截面积合并为一个向量，上标T表示转置，表示目标q在融合时刻t_k的扩展状态向量；的转移模型为：

上式，表示目标q的扩展状态向量在融合时刻t_k的转移矩阵，其中，表示目标q在融合时刻t_k的状态转移矩阵，I_N表示N维单位矩阵，blkdiag{·}表示矩阵对角化操作；(5)中，表示目标q在融合时刻t_k+1的扩展状态向量，表示目标q在融合时刻t_k的扩展过程噪声，服从均值为零，协方差为的高斯分布,其中，表示目标q在融合时刻t_k的运动过程噪声的协方差，表示目标q在融合时刻t_k的通道噪声的协方差。

步骤3：异构雷达网中N个雷达异步的对目标q进行观测，并将第k个融合时间间隔内N个雷达站得到的所有异步量测值发送至融合中心；

具体地，请参见图3，图3是本发明实施例提供的目标q在第k个融合时间间隔内的异步采样示意图，令M_i,q,k表示第k个融合时间间隔内雷达i对目标q的量测个数，在第k个融合时间间隔内雷达i对目标q的量测方程为：

上式中，表示第k个融合时间间隔内雷达i对目标q的第m个量测值，t_i,q,k(m)表示第k个融合时间间隔内雷达i对目标q的第m个量测值的采样时刻，本实施例中，雷达1至雷达8对目标q＝1的初始采样时刻分别为：2s，2.5s，3.1s，3.5s，4.1s，4.5s，5s，5.5s；雷达1至雷达8对目标q＝2的初始采样时刻分别为：2s，2.5s，3.2s，3.6s，4.2s，4.5s，5s，5.6s；雷达1至雷达8对目标q＝1和目标q＝2的采样时间间隔均为：2s，2s，2s，3s，2s，2s，2s，2s；表示在采样时刻t_i,q,k(m)目标q的扩展状态，式中，h(·)表示量测函数，

式(7)中，表示第k个融合时间间隔内雷达i对目标q的第m个距离量测值，表示第k个融合时间间隔内雷达i对目标q的第m个方位角量测值，表示在采样时刻t_i,q,k(m)雷达i对目标q的RCS量测值；表示目标q在采样时刻t_i,q,k(m)的位置在X轴方向上的分量，表示目标q在采样时刻t_i,q,k(m)的位置在Y轴方向上的分量；x(i)表示雷达i的位置在X轴方向上的分量，y(i)表示雷达i的位置在Y轴方向上的分量，本实施例雷达的位置见表1；arctan(·)表示反正切操作，表示n_x+N维的零向量的第n_x+i个元素为1，上标T表示转置。

表1

雷达	i＝1	i＝2	i＝3	i＝4	i＝5	i＝6	i＝7	i＝8
									位置/km	(-70,0.5)	(-40,-30)	(-40,30)	(40,30)	(70,0.5)	(40,-30)	(0,30)	(0,-30)

式(6)中，表示第k个融合时间间隔内雷达i对目标q的第m个量测值的量测高斯白噪声，噪声协方差为：

其中，表示第k个融合时间间隔内雷达i对目标q的第m个量测噪声的协方差矩阵，表示第k个融合时间间隔内雷达i对目标q的第m个距离量测值的量测方差，表示第k个融合时间间隔内雷达i对目标q的第m个方位角量测值的量测方差，表示第k个融合时间间隔内雷达i对目标q的第m个RCS量测值的量测方差；

其中，ζ_i表示雷达i的发射信号带宽，本实施例，ζ_i＝1MHz，B_i表示雷达i的3dB接收波束宽度，本实施例，B_i＝1°；P_i,q,k(m)表示第k个融合时间间隔内雷达i对目标q的第m次量测的发射功率资源，T_i,q,k(m)表示第k个融合时间间隔内雷达i对目标q的第m次量测的驻留时间资源，表示在采样时刻t_i,q,k(m)雷达i对目标q的雷达散射截面积量测值。提取式(9)中的公共项可得

其中，上标-1表示求逆运算，表示第k个融合时间间隔内雷达i对目标q的第m次量测的剩余参数组成的矩阵。

其中，Z_q,k表示第k个融合时间间隔内目标q的所有量测值构成的向量，表示第k个融合时间间隔内雷达i对目标q的第m个量测值，M_i,q,k表示第k个融合时间间隔内雷达i对目标q的量测值总个数，令M_q,k表示融合中心在第k个融合时间间隔内收到的所有量测值总个数，其中，N表示雷达的总数；在第k个融合时间间隔内，雷达将对目标的所有异步量测值发送到融合中心，用于更新目标在融合时刻t_k+1的状态及相应的估计协方差。

步骤4：根据第k个融合时间间隔内N个雷达对目标q的异步观测集合计算目标q在第k+1个融合时刻的归一化贝叶斯克拉美罗下界并将其作为目标q的跟踪精度的指标函数，在每个雷达的发射总资源有限的约束条件下，最小化通过将发射资源变量分为个互不重叠的块,并采用块坐标下降方法对上述凸问题求解，进而获得第k+1个融合时刻的最优分配结果

具体地，对于每个雷达发射总资源有限的资源分配问题，可以通过制定一个基于所有目标整体的归一化贝叶斯克拉美罗下界指标的优化问题来最小化多目标跟踪的均方误差，多目标跟踪性能驱动的异质资源分配方法的模型为：

其中，min表示最小化，s.t.表示约束条件，表示融合时刻t_k+1雷达i对目标q发射信号的能量，其中，表示在融合时刻t_k+1雷达i对目标q的发射功率，表示在融合时刻t_k+1雷达i对目标q的驻留时间；Q表示目标的总数，表示雷达i的总发射功率，表示雷达i的总驻留时间，本实施例i∈{1,2,...,N}；符号表示所有，符号∈表示属于，表示集中式MIMO雷达集合，表示相控阵雷达集合，表示在融合时刻t_k+1所有雷达的可调发射参数。

其中，表示集中式相控阵雷达的个数，表示相控阵雷达的个数，上标T表示转置，表示在融合时刻t_k+1的发射功率向量，表示在融合时刻t_k+1的驻留时间向量。

其中，表示在融合时刻t_k+1雷达i对目标q的发射功率，表示在融合时刻t_k+1雷达i对目标q的驻留时间。

其中，P_i,q,k+1(1)表示第k+1个融合时间间隔内雷达i对目标q的第1个量测的发射功率，P_i,q,k+1(M_i,q,k+1)表示第k+1个融合时间间隔内雷达i对目标q的第M_i,q,k+1个量测的发射功率，M_i,q,k+1表示第k+1个融合时间间隔内雷达i对目标q的总量测个数；T_i,q,k+1(1)表示第k+1个融合时间间隔内雷达i对目标q的第1个量测的驻留时间，T_i,q,k+1(M_i,q,k+1)表示第k+1个融合时间间隔内雷达i对目标q的第M_i,q,k+1个量测的驻留时间；式(12)中，表示融合时刻t_k+1目标q的归一化贝叶斯克拉美罗下界。

其中，Tr(·)表示求迹，上标-1表示求逆。

其中，blkdiag(·)表示矩阵块对角化操作，I₂表示2维单位矩阵，I_2N表示维数为2N的单位矩阵，表示矩阵直积，T₀表示融合时间间隔；式(16)中表示目标q的预测贝叶斯信息矩阵关于融合时刻t_k+1所有雷达的可调发射参数的函数的近似。

其中，表示目标q在融合时刻t_k+1的扩展状态，表示在融合时刻t_k+1雷达i对目标q的发射功率，表示在融合时刻t_k+1雷达i对目标q的驻留时间，M_i,q,k表示第k个融合时间间隔内雷达i对目标q的总量测个数，上标-1表示求逆运算，表示雷达i属于集中式MIMO雷达，表示雷达i属于相控阵雷达，表示第k个融合时间间隔内雷达i对目标q的第m个量测关于的雅克比矩阵的近似值，表示第k个融合时间间隔内雷达i对目标q的扩展状态的一步预测值，表示第k个融合时间间隔内雷达i对目标q的第m次量测的剩余参数组成的矩阵的近似，J_P(·)表示先验信息的费雪信息矩阵。

其中，表示目标q的扩展状态向量ξ在融合时刻t_k的过程噪声协方差，表示目标q的扩展状态向量ξ在融合时刻t_k的转移矩阵，表示目标q在融合时刻t_k的状态的贝叶斯信息矩阵，上标-1表示求逆，上标T表示转置，式(18)中，

其中，表示雷达i属于集中式MIMO雷达，表示在第k个融合时间间隔内集中式MIMO雷达i对目标q的第m个量测可调发射资源参数无关的参数组成的矩阵，表示雷达j属于相控阵雷达，表示在第k个融合时间间隔内相控阵雷达j对目标q的第m个量测可调发射资源参数无关的参数组成的矩阵，表示雷达l属于机械扫描雷达，表示在第k个融合时间间隔内机械扫描雷达l对目标q的第m个量测可调发射资源参数无关的参数组成的矩阵；将第k个融合时间间隔内的发射资源向量s_k分割成个不重叠的块，表示集中式MIMO雷达个数，表示相控阵雷达个数，

其中，表示在融合时刻t_k+1第i个雷达的发射资源向量，式(12)可以重新表述为：

其中，表示目标q的贝叶斯克拉美罗下界，min表示最小化，∑表示求和，表示在融合时刻t_k+1的发射资源向量，表示第1个雷达在融合时刻t_k+1的发射资源向量，表示第个雷达在融合时刻t_k+1的发射资源向量，符号∈表示属于，符号表示所有，表示在融合时刻t_k+1第i个雷达发射资源的凸集合。

其中，表示雷达i属于集中式MIMO雷达，表示雷达i属于相控阵雷达，1_Q表示维数为Q的列向量，表示在融合时刻t_k+1雷达i对目标q发射信号的能量，其中，表示在融合时刻t_k+1雷达i对目标q的发射功率，表示在融合时刻t_k+1雷达i对目标q的驻留时间，表示雷达i的总发射功率，表示雷达i的总驻留时间，符号表示所有，符号∈表示属于，M_i,k表示第k个融合时间间隔内雷达i对所有目标的量测个数组成的向量，M_i,k＝[M_i,1,k,...,M_i,Q,k]^T，其中，M_i,1,k表示第k个融合时间间隔内雷达i对目标1的量测总个数，M_i,Q,k表示第k个融合时间间隔内雷达i对目标Q的量测总个数，上标T表示转置。采用块坐标下降法结合CVX工具箱求解式(22)可以得到在融合时刻t_k+1资源分配的最优解

步骤5：根据第k个融合时间间隔内异构雷达网中目标q的异步量测值集合Z_q,k和融合时刻t_k的最优资源分配结果采用最大似然法计算目标q的扩展状态在融合时刻t_k+1的估计值并计算其估计协方差矩阵

具体地，根据第k个融合时间间隔内各个雷达对目标q的量测值相互独立，可得似然函数：

其中，Z_q,k表示第k个融合时间间隔内目标q的所有量测值构成的向量，表示目标q的扩展状态在融合时刻t_k+1的最大似然估计值，表示目标q的状态在融合时刻t_k+1的估计值，表示目标q的雷达散射截面积在融合时刻t_k+1的估计值；表示在已知的条件下Z_q,k的条件概率密度函数，上标T表示转置，Π表示求积操作，表示服从均值为方差为的高斯分布，表示第k个融合时间间隔内雷达i对目标q的第m个量测值；表示在采样时刻t_i,q,k(m)目标q的扩展状态，h(·)表示量测函数，表示第k个融合时间间隔内雷达i对目标q的第m个量测噪声的协方差矩阵，最大似然估计为

其中，表示当f(x)取得最大值时返回对应的参数x，ln(·)表示对数函数，根据迭代最小二乘法可以求具体地，在第r次迭代后的值为：

其中，表示在融合时刻t_k+1目标q的扩展状态在第r次迭代后的估计值，表示在融合时刻t_k+1目标q的扩展状态在第r次迭代后状态估计值对应的量测值；

表示雷达i在采样时刻t_i,q,k(m)对目标q的量测值，上标T表示转置；预测目标q在融合时刻t_k+1的扩展状态

其中，t_i,q,k(m)表示第k个融合时间间隔内雷达i对目标q的第m个量测值的采样时刻，t_k+1表示融合时刻，I₂表示2维单位矩阵，表示矩阵直积，表示在采样时刻t_i,q,k(m)的目标q的扩展状态；式(26)中Σ_q,k表示量测集合Z_q,k对应的噪声协方差矩阵。

其中，blkdiag(·)表示矩阵块对角化，表示第k个融合时间间隔内雷达i对目标q的第m个量测值的量测噪声协方差，式(26)中

其中，H_q,k,r表示第k个融合时间间隔内目标q的第r次迭代的雅克比矩阵，表示第k个融合时间间隔内雷达i对目标q的第m个量测值的雅克比矩阵。

其中，表示在采样时刻t_i,q,k(m)目标q的扩展状态，表示在采样时刻t_i,q,k(m)雷达i对目标q的量测函数，表示求偏导，表示对求偏导，表示将的值赋为为目标q在融合时刻t_k+1的混合量测的估计值，其对应的估计协方差可以用克拉美罗下界表示其中，Z_q,k表示第k个融合时间间隔内目标q的所有量测值构成的向量，表示目标q在融合时刻t_k+1的混合量测值的真实值，表示关于Z_q,k的费雪信息矩阵。

其中，上标N表示雷达的总个数，M_i,q,k表示第k个融合时间间隔内雷达i对目标q的总量测个数，P_i,q,k表示第k个融合时间间隔内雷达i对目标q的发射功率，T_i,q,k表示第k个融合时间间隔内雷达i对目标q的驻留时间，表示第k个融合时间间隔内雷达i对目标q的第m个量测值的雅克比矩阵，表示第k个融合时间间隔内雷达i对目标q的第m次量测的剩余参数组成的矩阵，上标-1表示求逆运算，上标T表示转置。

步骤6：将目标q的扩展状态的最大似然估计值及其对应的估计协方差矩阵作为卡尔曼滤波器的输入从而得到目标q的扩展状态在融合时刻t_k+1的最终状态估计值及其对应的协方差矩阵

具体地，根据目标q的扩展状态的最大似然估计值及其对应的估计协方差矩阵根据卡尔曼滤波可以得到目标q在融合时刻t_k+1的状态估计值及其对应的协方差矩阵

其中，表示目标q的扩展状态从融合时刻t_k到融合时刻t_k+1的预测值，表示目标q的扩展状态向量在融合时刻t_k的转移矩阵，表示目标q在融合时刻t_k的状态估计值；表示目标q的扩展状态从融合时刻t_k到融合时刻t_k+1的预测协方差，表示目标q在融合时刻t_k的状态估计估计协方差，上标T表示转置，式(33)中

其中，表示目标q的扩展状态在融合时刻t_k+1的新息协方差，表示目标q的扩展状态在融合时刻t_k+1的卡尔曼增益，Z_q,k表示第k个融合时间间隔内目标q的所有量测值构成的向量，表示目标q在融合时刻t_k+1的混合量测值的真实值，表示关于Z_q,k的费雪信息矩阵。

步骤7：令k＝k+1，返回步骤2，直到得到第1个融合时刻的异质资源分配结果至第K个融合时刻的异质资源分配结果以及目标q在第1个融合时刻的状态估计值及其对应的协方差矩阵至第K个融合时刻的状态估计值及其对应的协方差矩阵q∈{1,2,…,Q}。

实施例二

下面通过仿真实验对本发明效果作进一步验证说明。

(一)仿真条件：

本发明的仿真运行系统为Intel(R)Core(TM)i5-4590CPU@3.30GHz，64位Windows10操作系统，仿真软件采用MATLAB(R2016b)。

(二)仿真内容与结果分析：

请参见图2，图2是本发明实施例提供的异构雷达网探测区域内的目标运动场景示意图；在异构雷达网中雷达总个数为N＝8，其中，集中式MIMO雷达的个数为集中式MIMO雷达用正方形表示，相控阵雷达个数为相控阵雷达的位置用圆圈表示，机械扫描雷达的个数为机械扫描雷达用三角形表示；第i个雷达的发射信号有效带宽为ζ_i＝1MHz，第i个雷达的3dB波束宽度为B_i＝1°，第i个雷达的发射信号频率为f_i＝(1+0.1i)，i∈{1,2,...,N}；若第i个雷达为集中式MIMO雷达，则其总发射功率资源为第i个雷达为相控阵雷达，则其总驻留时间资源为设定当第k个融合时间间隔内第i个雷达与第q个目标的径向距离为并且在采样时刻t_i,q,k(m)雷达i对目标q的雷达散射截面积为时的信噪比为12dB；设定在异构雷达网的探测区域中存在的目标个数为Q＝2，目标q＝1的初始位置为(-40,0)km，初始速度为(50,0)m/s；目标q＝2的初始位置为(40,0)km，初始速度为(-50,0)m/s，假定目标均为匀速直线运动；假定每个集中式MIMO雷达采用同时多波束工作方式，因此其对应于不同目标的采样时刻相同；同时，每个机械扫描雷达以固定的转速连续照射多个目标，因此其对多个目标的重访间隔也是相同的；相反，相控阵雷达由于具有波束的灵活性，故其可以在不同的初始时刻以不同的重访时间间隔照亮多个目标；表2为每个雷达对每个目标的初始采样时刻和对应的重访时间间隔。

表2

(t<sub>i,q,1</sub>(1),T<sub>i,q</sub>)/s	i＝1	i＝2	i＝3	i＝4	i＝5	i＝6	i＝7	i＝8
									q＝1	(2,2)	(2.5,2)	(4.5,2)	(5,2)	(3.1,2)	(3.5,3)	(4.1,2)	(5.5,2)
q＝2	(2,2)	(2.5,2)	(4.5,2)	(5,2)	(3.2,2)	(3.6,3)	(4.2,2)	(5.6,2)

请参见图4，图4是本发明实施例提供的雷达2对目标1的雷达散射截面积示意图，而其它目标的雷达散射截面积参数都设置为1；在这种情况下，雷达2和雷达5不断照射目标1，但它们的雷达散射截面积不同，因此可以分析目标的雷达散射截面积对资源分配结果的影响。

请参见图5，图5是本发明实施例提供的异质资源分配方法的多目标跟踪性能；传统的异质资源分配方案不使用任何先验知识，为了方便起见，将每个雷达有限的资源统一分配给多个目标。将传统的异质资源分配方案作为基准，为了检验探测性能驱动的异质资源分配方法的最优性，将其归一化跟踪均方根误差和相应的贝叶斯克拉美罗下界作为度量标准，与基准进行比较。蒙特卡罗次数为100；虚线为基准的贝叶斯克拉美罗下界，正方形点表示基准的归一化跟踪均方根误差，实线为探测性能驱动的异质资源分配方法的贝叶斯克拉美罗下界，圆点表示探测性能驱动的异质资源分配方法的归一化跟踪均方根误差；结果表明，探测性能驱动的异质资源分配方法贝叶斯克拉美罗下界均小于基准的贝叶斯克拉美罗下界，所以本发明方法可以显著提高多目标跟踪精度。

请参见图6，图6是本发明实施例提供的目标q＝1的探测性能驱动的异质资源分配结果；在探测性能驱动的异质资源分配问题中，每个雷达需要将其总资源分配给多个目标，这些资源更有可能首选分配给对跟踪性能贡献大的目标。当雷达的资源丰富时，还可能将剩余资源分配给其他对跟踪性能贡献小的目标；在此，采用跟踪精度增长率(归一化克拉美罗下界对发射资源的负梯度)作为跟踪性能贡献的指标；例如，虽然雷达2对目标1的航迹条件相对较弱，但由于雷达2与目标2的距离非常大，所以雷达2对目标1的跟踪精度增长率大于目标2，故雷达2仍会被用来跟踪目标1；相似地，请参见图7，图7是本发明实施例提供的目标q＝2的探测性能驱动的异质资源分配结果，虽然雷达6对目标2的航迹条件相对较弱，但由于雷达6与目标1的距离非常大，所以雷达6对目标2的跟踪精度增长率大于目标1，故雷达6仍会被用来跟踪目标2；

请参见图8，图8是本发明实施例提供的每个目标的跟踪贝叶斯克拉美罗下界；目标1的跟踪贝叶斯克拉美罗下界(BCRLB)高于目标2的跟踪贝叶斯克拉美罗下界，目标1和目标2的位置关于雷达4对称，雷达4对目标1和目标2的量测雅克比矩阵是相同的，但是目标1的先验信息较少与目标2，雷达4对目标1的跟踪精度增长率会高于雷达4对目标2的跟踪精度增长率，故雷达4会被用于照射目标1。

综上所述，仿真实验验证了本发明的正确性，有效性和可靠性。

以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。

Claims (1)

1.一种异步多目标跟踪性能驱动的异构雷达网资源分配方法，其特征在于，包括以下步骤：

步骤1：建立异构雷达网，所述异构雷达网包含个雷达站和一个融合中心；N、和均为大于0的正整数，其中，表示集中式MIMO雷达个数，表示相控阵雷达个数，表示机械扫描雷达个数；所述N个雷达站的探测区域中存在Q个目标，Q为大于零的正整数；

初始化：令t_k表示第k个融合时刻，k的初始值为1，k∈{1,2,…,K}，t_K为预先设定的最大融合时刻，K为大于0的正整数，(t_k,t_k+1)表示第k个融合时间间隔；

步骤2：设定所述异构雷达网中目标q为匀速直线运动，q∈{1,2,…,Q}，并设定t_k时刻目标q的扩展状态向量为

步骤3：异构雷达网中N个雷达异步的对目标q进行观测，并将第k个融合时间间隔内N个雷达站得到的所有异步量测值发送至融合中心；

步骤4：根据第k个融合时间间隔内N个雷达对目标q的异步观测集合计算目标q在第k+1个融合时刻的归一化贝叶斯克拉美罗下界并将其作为目标q的跟踪精度的指标函数，在每个雷达的发射总资源有限的约束条件下，最小化通过将发射资源变量分为个互不重叠的块,并采用块坐标下降方法对上述凸问题求解，进而获得第k+1个融合时刻的最优分配结果

步骤5：根据所述第k个融合时间间隔内异构雷达网中目标q的异步量测值集合Z_q,k和第k个融合时刻的最优资源分配结果采用最大似然法计算目标q的扩展状态在融合时刻t_k+1的估计值并计算其估计协方差矩阵

步骤6：将目标q的扩展状态的最大似然估计值及其对应的估计协方差矩阵作为卡尔曼滤波器的输入从而得到目标q的扩展状态在融合时刻t_k+1的最终状态估计值及其对应的协方差矩阵

步骤7：令k＝k+1，返回步骤2，直到得到第1个融合时刻的异质资源分配结果至第K个融合时刻的异质资源分配结果以及目标q在第1个融合时刻的状态估计值及其对应的协方差矩阵至第K个融合时刻的状态估计值及其对应的协方差矩阵