CN111679270B - 一种反射点不确定场景下多路径融合目标检测算法 - Google Patents

Abstract

本发明公开了一种反射点不确定场景下多路径融合目标检测算法，属于雷达技术领域。本发明以迭代的方式同时估计目标和反射点的状态，用一批数据使用滑动窗口方式用于航迹初始化。每个迭代由两个步骤组成：预测和更新。在预测阶段，通过上个时刻最后一次迭代的估计状态来预测反射点的状态，在此基础上，通过多遍网格搜索全局最优解来估计目标的状态，在更新阶段，利用更新后的目标状态来重新估计反射点状态。重复上述步骤，迭代过程基于给定的最大迭代次数或判断是否过门限，如果迭代次数达到最大迭代次数或过门限则停止迭代，输出目标位置和反射点的位置。本发明提高了目标状态估计精度，对低空小目标进行有效的检测。

Description

一种反射点不确定场景下多路径融合目标检测算法

技术领域

本发明属于雷达技术领域，具体涉及一种反射点不确定场景下多路径融合目标检测算法。

背景技术

低小慢目标指的是飞行高度较低、体积相对较小、飞行速度慢的小型飞行器，典型的低小慢目标包括无人机、航模、滑翔伞等。最近几年随着我国民用无人机的发展，无人机性能达到国际一流水平，普通民众也能消费得起民用无人机。与此同时也带来了一些安全上的隐患，比如在民航机场经常会出现未被允许飞行的小型无人机，这对于机场的正常运行造成严重的影响，甚至会影响飞机的起飞与降落。另外，在一些活动场所中经常出现一些无人机，同样严重影响到了活动的安全保障。因此，如何有效管理和控制无人机，成为了当前雷达监测系统亟待解决的问题。

为了提高在低信噪比、高杂波场景下目标的检测与跟踪性能，在检测系统中应尽可能挖掘传感器接收量测中的目标信息。在一些应用场景中，由于电磁波传播存在多径效应，传感器在接收回波中同一帧中能产生多个目标量测。比如在城市环境中，高楼大厦随处可见，对于监测系统来说目标很容易被遮挡，存在视觉盲区，再者电磁波信号在传播的时候不仅可以通过目标直接反射回到接收机，还会通过建筑物反射回到接收机，这种情况下目标会有多条传播路径，如果不对其进行处理的话会影响目标的正常检测。换个角度来说，通过其他路径反射回来的目标信号也是目标信息的一部分，如果能将这一部分目标信息利用起来，势必会增加目标检测性能。由于目标在不停的运动，反射点的位置一般也会随着目标的运动而运动，通常情况下只知道建筑物的大概位置，但是并不知道反射点在建筑物的具体位置。

现有的基于多路径目标检测算法假定反射点已知，量测模式确定，如果在反射点未知场景下使用传统的多路径ML-PMHT(最大似然概率多假设跟踪)目标检测算法，会使得目标估计误差增大，而且容易产生虚假航迹，失跟率高。现有的多路径ML-PMHT目标检测算法，并没有将反射点不确定性建模到算法中，由于反射点初始位置不确定，导致量测误差增加，因此通过搜索优化得到的目标状态距离真实目标状态差距较大，并且容易将由杂波形成的对数似然比(LLR)误认为是目标形成的对数似然比(LLR)，从而造成虚假航迹。

发明内容

本发明的目的是克服ML-PMHT算法在反射点不确定多径环境中应用的不足，提供一种反射点不确定场景下高效多路径信息融合的最大似然概率多假设跟踪(ML-PMHT)的目标检测算法，具体来说是一种基于城市复杂环境的信号多路径特点，以及反射点不确定的传播特点，对低空小目标(比如小型飞行器)进行检测的方法。本发明所述方法通过对目标、多路径观测函数显式建模，多并对多路径观测信息进行融合处理，从而对低信噪比、高杂波、反射点不确定场景下的低空小目标的检测。

本发明所述方法的思路是，以迭代的方式同时估计目标和反射点的状态。用一批数据使用滑动窗口方式用于航迹初始化。在每次迭代中，使用相同的测量集。每个迭代由两个步骤组成：预测和更新。在预测阶段，通过上个时刻最后一次迭代的估计状态来预测反射点的状态，在此基础上，通过多遍网格搜索全局最优解来估计目标的状态，由于LLR的计算使用的是估计的反射点而不是使用真实反射点位置来计算的，因此LLR的全局最优并不能得到真实的目标状态。在更新阶段，利用更新后的目标状态来重新估计反射点状态。重复上述步骤，该迭代过程基于给定的最大迭代次数或判断是否过门限，如果迭代次数达到最大迭代次数或过门限则停止迭代，输出目标位置和反射点的位置。从而实现城市低信噪比、高杂波、反射点不确定场景下低空小目标的检测。

本发明所提出的技术问题是这样解决的：

一种反射点不确定场景下多路径融合目标检测算法，包括以下步骤：

步骤1.初始化多路径ML-PMHT算法环境参数：

初始化观测环境各项参数：监控空间V，反射点数L，反射点初始状态

表示均值为

协方差为

的高斯分布，1≤l≤L，每个反射点对应一个路径，量测通过传播路径l来源于目标n的先验概率为π_n，l，量测来源于杂波的先验概率π₀₀，当前迭代次数Iter＝1，最大迭代次数为Imax；

目标初始化场景下参数向量分别表示为：

其中，k为批处理的采样时刻数，1≤k≤N_w，N_w为批处理长度，X为N_w个时刻的目标状态参数，X_f为N_w个时刻反射点状态参数，Z为N_w个时刻的量测集合；

在k时刻目标的状态向量表示为：

X_k＝{x_1，k，...，x_n，k，...，x_N，k}n＝1，...，N

在k时刻反射点状态向量表示为：

在k时刻收到的量测集合表示为：

其中，x_n，k为k时刻第n个目标的状态，N为目标的个数，

表示k时刻第l个反射点的状态，L为反射点的个数，z_k，j表示传感器k时刻收到的第j个量测，m_k表示传感器k时刻收到的量测个数；

对目标的运动过程和对不同传播路径的观测方程分别表示为：

x_n，k＝F_n，kx_n，k-1+v_n，k

其中，F_n，k为目标运动矩阵，x_n，0为第n个目标运动的初始状态，h_l为x_n，k经路径l的观测函数；z_n，k，l为x_n，k通过路径l产生的目标量测；v_n，k和ω_n，l，k分别为目标过程噪声和观测噪声，都是均值为零的高斯白噪声，协方差矩阵分别为Q_n，k和R_n，l，k；

步骤2.预测反射点状态与方差：

其中，

表示k时刻第l个反射点预测状态，F_l，k表示反射点的运动矩阵，v_l，k表示反射点过程噪声，是零均值的高斯白噪声，协方差矩阵为

表示k时刻反射点l的预测状态协方差，

表示反射点l的初始预测状态协方差；

是测量模型对反射点状态的雅可比矩阵，其中

表示

关于

的偏导，上标T表示转置；

步骤3.根据预测的反射点状态，搜索所有可能的目标状态：

步骤3-1.令N′＝1，目标状态集合X′为空向量，

步骤3-2.搜索单目标LLR全局最大值；

其中，

p_l[z_k，j|x_n′，k]表示通过路径l以第n′个目标为中心的高斯概率；

其中，

表示高斯概率密度函数，高斯变量z_k，j的均值为

协方差为R_n′，l，k，且有：

步骤3-3.判断单目标LLR全局最大值

是否大于检测门限，若是，确认

为新目标，

即将

加入到向量X′中，转至步骤3-4；否则，结束检测过程，转至步骤4；

步骤3-4.从量测集合

中剔除与已确认为新目标关联的量测，剔除依据为目标后验关联概率最高的量测将被剔除；

ML-PMHT后验关联概率w_{j，n′，l，k}计算公式为：

步骤3-5.令N′＝N′+1，转至步骤3-2；此时搜索状态空间为4N′维，由于有(N′-1)个目标已经确认，在新的搜索过程无需搜索已经搜出的4(N′-1)维，所以新的搜索阶段只需要对未知的4维状态进行搜寻；

步骤4.更新反射点状态：

步骤4-1.计算k时刻的量测z_k，j通过第l条路径来源第n个目标的后验概率：

步骤4-2.计算每个目标在每条路径下的综合量测

及综合协方差

步骤4-3.不同路径下的雅克比矩阵、综合量测、综合量测协方差分别进行堆叠：

堆叠雅克比矩阵：

堆叠综合量测：

堆叠量测协方差：

其中，diag表示对角化矩阵；

步骤4-4.对反射点执行扩展卡尔曼平滑算法：

步骤4-4-1.前向滤波

对于k＝1∶N_w，

其中，

为由k-1时刻的反射点状态预测的k时刻的第l个反射点状态，F_f为反射点的运动矩阵，

为k-1时刻的第l个反射点状态；

其中，

为第n′个新目标通过第l个反射点得到的量测，

堆叠得到综合量测

其中，

为预测得到的第l个反射点协方差，

为k-1时刻的第l个反射点协方差；

其中，

为卡尔曼增益；

其中，

为卡尔曼滤波更新后的第l个反射点协方差，I为单位矩阵；

其中，

为卡尔曼滤波更新后的第l个反射点状态；

步骤4-4-2.后向平滑

对于k＝N_w-1∶1，

其中，

为平滑增益；

其中，

为k时刻平滑后的第l个反射点状态；

其中，

为k时刻平滑后的第l个反射点协方差；

步骤5.令Iter＝Iter+1，判断Iter是否大于Imax，若是，转至步骤6，否则转至步骤2；

步骤6.输出目标状态集合X′和反射点的状态和协方差。

本发明的有益效果是：

本发明对城市复杂环境的信号多路径以及反射点不确定的传播特点，同时利用ML-PMHT的算法特性有效降低量测与目标关联的计算复杂度，能有效利用多条路径的量测信息，提高了目标状态估计精度。对低空小目标(比如小型飞行器)进行有效的检测。

附图说明

图1为实施例中城市飞行器检测场景下目标与传感器的位置几何图；

图2为实施例中监控空间中2两个目标10个采样时刻的时延观测杂波图；

图3为实施例中监控空间中2两个目标10个采样时刻的多普勒观测杂波图；

图4为实施例中目标和反射点状态估计的迭代误差图；

图5为100次蒙特卡洛实验中目标状态估计航迹图。

具体实施方式

下面结合附图和实施例对本发明进行进一步的说明。

本实施例提供一种反射点不确定场景下多路径融合目标检测算法，包括以下步骤：

步骤1.初始化环境参数和算法参数，城市飞行器检测场景下，目标与传感器的位置几何图如图1所示，该场景下存在2个反射点、1个发射机、1个接收机和2个目标，即L＝2，N＝2， l＝1，2，n＝1，2，该场景下单个目标根据电磁波传播顺序会产生3条传播路径：

在城市小型无人机的检测场景中，接收机基站固定在[0m，0m]收集由建筑物和目标反射回来的信号，发射机基站固定在[-500m，1000m]，由于本发明属于外辐射源雷达，采用4G 基站做为发射基站。

多普勒观测空间V_γ＝500Hz，多普勒标准差σ_γ＝2Hz，时延观测空间V_r＝10μs，监控空间V＝V_γV_r，时延标准差σ_r＝0.03μs，航迹虚警率P_FA＝0.02，每条路径信噪比SNR＝7dB，假定杂波在单元内均匀分布，批处理长度为N_w＝10，且平均每个时刻的杂波数N_λ＝10。监控空间中2两个目标10个采样时刻的时延观测杂波图如图2所示，多普勒观测杂波图如图3 所示。

雷达系统采样时刻之间的间隔ΔT为2s，仿真总时长为10个时刻，每批处理中最大迭代次数Imax＝10次，当前迭代次数Iter＝1。

在采样过程中，2个目标的初始状态的运动向量分别为：

x_1，1＝[-104m 21m/s 1310m -9m/s]

x_2，1＝[-100m 20m/s 500m 9m/s]

检测空间中存在两个反射点，反射点1初始状态

其中

反射点2初始状态

其中

量测通过传播路径l来源于目标n的先验概率为π_n，l＝0.02，量测来源于杂波的先验概率π₀₀＝0.96。

步骤2.预测反射点状态与方差：

其中，

表示k时刻第l个反射点预测状态，F_l，k表示反射点的运动矩阵，v_l，k表示反射点过程噪声，是零均值的高斯白噪声，协方差矩阵为

表示k时刻反射点l的预测状态协方差，

表示反射点l的初始预测状态协方差；

是测量模型对反射点状态的雅可比矩阵，其中

表示

关于

的偏导，上标T表示转置；

步骤3.根据预测的反射点状态，搜索所有可能的目标状态，得到

为新目标的状态，n′＝1，2；

步骤4.更新反射点状态，得到k时刻平滑后的第l个反射点状态

k时刻平滑后的第l个反射点协方差

步骤5.令Iter＝Iter+1，判断Iter是否大于Imax，若是，转至步骤6，否则转至步骤2；

步骤6.输出目标状态集合X′和反射点的状态和协方差。

图4给出了目标和反射点状态估计的迭代误差图。可以看到随着迭代次数的增加，目标和反射点的估计误差越来越小，在迭代到4次之后，目标和反射点的估计误差基本收敛，后续迭代中误差在一个很小的范围内波动，因此迭代算法是有效的。这是因为在迭代的过程中初始反射点误差较大，因此估计出来的目标状态误差较大，随着迭代的进行，通过估计出来的目标状态来更新反射点的状态，再使用更新后的反射点状态带入似然公式，从而使得目标状态估计更加精确，通过几次迭代之后，估计出来的反射点位置不断接近真实反射点位置，从而达到提高估计精度的目的。

图5中蓝色实线为在低信噪比，高杂波、反射点不确定场景下对多目标检测100次目标初始化结果，且统计结果对两个目标的过门限次数均为100。其结果表明，初始化结果和目标真实航迹非常接近，解决了该场景下目标检测的问题。

最后说明的是，以上实施仅用以说明本发明的技术方案而非限制，凡依本发明申请专利范围所做的均等变化和修饰，皆应属本发明的涵盖范围。

Claims (1)

1.一种反射点不确定场景下多路径融合目标检测算法，其特征在于，包括以下步骤：

步骤1.初始化多路径ML-PMHT算法环境参数：

初始化观测环境各项参数：监控空间V，反射点数L，反射点初始状态

表示均值为

协方差为

的高斯分布，1≤l≤L，每个反射点对应一个路径，量测通过传播路径l来源于目标n的先验概率为π_n,l，量测来源于杂波的先验概率π₀₀，当前迭代次数Iter＝1，最大迭代次数为Imax；

目标初始化场景下参数向量分别表示为：

其中，k为批处理的采样时刻数，1≤k≤N_w，N_w为批处理长度，X为N_w个时刻的目标状态参数，X_f为N_w个时刻反射点状态参数，Z为N_w个时刻的量测集合；

在k时刻目标的状态向量表示为：

X_k＝{x_1,k,...,x_n,k,...,x_N,k}n＝1,...,N

在k时刻反射点状态向量表示为：

在k时刻收到的量测集合表示为：

其中，x_n,k为k时刻第n个目标的状态，N为目标的个数，

表示k时刻第l个反射点的状态，L为反射点的个数，z_k,j表示传感器k时刻收到的第j个量测，m_k表示传感器k时刻收到的量测个数；

对目标的运动过程和对不同传播路径的观测方程分别表示为：

x_n,k＝F_n,kx_n,k-1+ν_n,k

其中，F_n,k为目标运动矩阵，x_n,0为第n个目标运动的初始状态，h_l为x_n,k经路径l的观测函数；z_n,k,l为x_n,k通过路径l产生的目标量测；ν_n,k和ω_n,l,k分别为目标过程噪声和观测噪声，都是均值为零的高斯白噪声，协方差矩阵分别为Q_n,k和R_n,l,k；

步骤2.预测反射点状态与方差：

其中，

表示k时刻第l个反射点预测状态，F_l,k表示反射点的运动矩阵，v_l,k表示反射点过程噪声，是零均值的高斯白噪声，协方差矩阵为

表示k时刻反射点l的预测状态协方差，

表示反射点l的初始预测状态协方差；

是测量模型对反射点状态的雅可比矩阵，其中

表示

关于

的偏导，上标T表示转置；

步骤3.根据预测的反射点状态，搜索所有可能的目标状态：

步骤3-1.令N'＝1，目标状态集合X'为空向量，

步骤3-2.搜索单目标LLR全局最大值；

其中，

p_l[z_k,j|x_n',k]表示通过路径l以第n'个目标为中心的高斯概率；

其中，

表示高斯概率密度函数，高斯变量z_k,j的均值为

协方差为R_n',l,k，且有：

步骤3-3.判断单目标LLR全局最大值

是否大于检测门限，若是，确认

为新目标，

即将

加入到向量X'中，转至步骤3-4；否则，结束检测过程，转至步骤4；

步骤3-4.从量测集合

中剔除与已确认为新目标关联的量测，剔除依据为目标后验关联概率最高的量测将被剔除；

ML-PMHT后验关联概率w_j,n',l,k计算公式为：

步骤3-5.令N'＝N'+1，转至步骤3-2；

步骤4.更新反射点状态：

步骤4-1.计算k时刻的量测z_k,j通过第l条路径来源第n个目标的后验概率：

步骤4-2.计算每个目标在每条路径下的综合量测

及综合协方差

步骤4-3.不同路径下的雅克比矩阵、综合量测、综合量测协方差分别进行堆叠：

堆叠雅克比矩阵：

堆叠综合量测：

堆叠量测协方差：

其中，diag表示对角化矩阵；

步骤4-4.对反射点执行扩展卡尔曼平滑算法：

步骤4-4-1.前向滤波

对于k＝1:N_w，

其中，

为由k-1时刻的反射点状态预测的k时刻的第l个反射点状态，F_f为反射点的运动矩阵，

为k-1时刻的第l个反射点状态；

其中，

为第n'个新目标通过第l个反射点得到的量测，

堆叠得到综合量测

其中，

为预测得到的第l个反射点协方差，

为k-1时刻的第l个反射点协方差；

其中，

为卡尔曼增益；

其中，

为卡尔曼滤波更新后的第l个反射点协方差，I为单位矩阵；

其中，

为卡尔曼滤波更新后的第l个反射点状态；

步骤4-4-2.后向平滑

对于k＝N_w-1:1，

其中，

为平滑增益；

其中，

为k时刻平滑后的第l个反射点状态；

其中，

为k时刻平滑后的第l个反射点协方差；

步骤5.令Iter＝Iter+1，判断Iter是否大于Imax，若是，转至步骤6，否则转至步骤2；

步骤6.输出目标状态集合X'和反射点的状态和协方差。

Images