CN110231620A - 一种噪声相关系统跟踪滤波方法 - Google Patents

Abstract

本发明涉及一种噪声相关系统跟踪滤波方法，包括：从观测雷达处获取k时刻目标量测数据并建立噪声相关系统的状态空间模型，其中量测数据包括距离量测数据和方位角量测数据；将k时刻距离量测数据和方位角量测数据转换为直角坐标系下目标在x方向和y方向的位置量测数据，计算量测转换中的偏差项以及位置量测协方差矩阵，构建转换后的量测方程，并计算k‑1时刻和k时刻目标位置量测噪声的互协方差；在最小均方误差的准则下构建滤波器，利用k‑1时刻和k时刻量测噪声的互协方差补偿两个时刻噪声之间的相关性，对转换后的位置量测数据进行滤波处理，更新k时刻目标状态估计和状态估计误差协方差，完成目标跟踪。该方法可用于解决噪声相关问题，跟踪结果好，性能稳定。

Description

一种噪声相关系统跟踪滤波方法

技术领域

本发明涉及空间目标跟踪技术领域，尤其涉及一种噪声相关系统跟踪滤波方法。

背景技术

空间目标跟踪时，由于噪声干扰，在信噪比较低的环境下运动目标的检测和跟踪精度不高，因此在某些信噪比较低的场合中，需要通过时间累积来增强信号强度，减弱噪声的影响，以此来提高信噪比，完成后面的信号检测、跟踪、滤波等一系列步骤。为了在提高目标检测性能的同时，保证较高的探测数据率，一些雷达系统会采用滑窗积累的方式，即相邻时刻会用到同一部分数据进行积累处理，这导致相邻积累时刻的观测噪声存在相关性问题。也就是说，虽然量测噪声仍然是均值为零的高斯随机变量，但是由于不同时刻间的量测噪声存在相关问题，因此量测噪声不再是高斯白噪声。而以卡尔曼滤波为基础的现代跟踪滤波方法，均假设观测噪声为高斯白噪声，当观测噪声存在自相关性时，现有滤波方法无法实现有效的状态估计，因此不能直接应用于以上噪声相关系统。

目前，主要有两种方法解决时间相关量测误差问题，第一种方法是将量测矢量作为增强状态变量中的一组约束，将量测噪声相关转化为量测与状态向量相关。但是增强的状态空间方程导致量测协方差矩阵存在奇异，这可能导致系统变得发散。第二种方法通过使用量测差分步骤来构造一个新的量测方程，新构造的量测方程受高斯白噪声影响并且没有时间相关量测部分。但是在存在多个量测的情况下，该方法要求每个量测应该同时时间相关。

发明内容

本发明的目的是针对上述现有技术的至少一部分缺陷，提供一种目标跟踪滤波方法，实现噪声相关条件下的精确状态估计。

为了实现上述目的，本发明提供了一种噪声相关系统跟踪滤波方法，包括如下步骤：

S1、从观测雷达处获取k时刻目标量测数据并建立噪声相关系统的状态空间模型，其中量测数据包括距离量测数据和方位角量测数据；

S2、将k时刻距离量测数据和方位角量测数据转换为直角坐标系下目标在x方向和y方向的位置量测数据，计算量测转换中的偏差项以及位置量测协方差矩阵，构建转换后的量测方程，并计算k-1时刻和k时刻目标位置量测噪声的互协方差；

S3、在最小均方误差的准则下构建滤波器，利用k-1时刻和k时刻量测噪声的互协方差补偿两个时刻噪声之间的相关性，对转换后的位置量测数据进行滤波处理，更新k时刻目标状态估计和状态估计误差协方差，完成目标跟踪。

优选地，所述步骤S3包括如下步骤：

S3-1、利用由状态空间模型得到的状态方程计算状态一步预测及状态一步预测协方差，同时利用转换后的量测方程计算量测一步预测，结合位置量测协方差矩阵和相邻两时刻位置量测噪声的互协方差，计算量测预测协方差；

S3-2、利用状态一步预测、量测一步预测，并结合相邻两时刻位置量测噪声的互协方差，计算状态预测与量测预测间的互协方差；

S3-3、利用状态一步预测、状态预测协方差、量测一步预测、量测预测协方差，以及状态预测和量测预测之间的互协方差，更新目标状态估计和状态估计误差协方差，实现当前时刻的目标跟踪。

优选地，所述步骤S1中，建立噪声相关系统的状态空间模型为：

x_k＝Fx_k-1+w_k-1；

其中，x_k＝[p_x,k,p_y,k,v_x,k,v_y,k]^T是离散时刻k∈{1,2,3,...}目标的状态，(p_x,k,p_y,k)表示目标在x方向和y方向的位置，(v_x,k,v_y,k)表示目标在x方向和y方向的速度，F是状态x_k的转移矩阵，过程噪声w_k-1是均值为0、协方差为Q_k-1的白色高斯随机变量，n_k是量测噪声，z_k是极坐标下的目标量测数据，和分别表示距离和方位角量测数据；

状态与量测数据之间满足关系式：

其中，表示距离量测噪声和方位角量测噪声，是均值为0的高斯随机变量，且k-1时刻和k时刻之间的量测噪声存在相关性，对应的协方差分别表示为：

其中，ρ为相邻两时刻量测噪声相关系数，σ_r、σ_θ分别表示距离和方位角量测噪声的标准差，量测噪声n_k严格独立于过程噪声w_k；

匀速运动模型中，状态转移矩阵F和过程噪声协方差矩阵Q_k-1的表达式分别为：

其中，T是采样间隔，q_x和q_y分别代表x方向和y方向的过程噪声功率谱密度。

优选地，所述步骤S2中将k时刻距离量测数据和方位角量测数据转换为直角坐标系下目标在x方向和y方向的位置量测数据时，表达式为：

其中，表示经过转换且去除偏差项后的目标位置量测数据，分别表示利用距离和方位角量测数据计算出的目标在x方向、y方向的位置量测数据，分别表示量测转换中x方向、y方向上的偏差项。

优选地，所述步骤S2中计算量测转换中的偏差项以及位置量测协方差矩阵，构建转换后的量测方程，并计算k-1时刻和k时刻目标位置量测噪声的互协方差时，量测转换中的偏差项分别表示为：

位置量测协方差矩阵表示为：

构建转换后的量测方程，表达式为：

其中，H为量测矩阵，在匀速运动模型中表达式为：

表示转换后的量测噪声，是均值为零、协方差矩阵为的高斯随机变量，表达式为：

相邻两时刻目标位置量测噪声的互协方差的表达式为：

优选地，所述步骤S3-1中，利用由状态空间模型得到的状态方程计算状态一步预测及状态一步预测协方差时，状态一步预测表达式为：

状态一步预测协方差为：

P_k|k-1＝FP_k-1|k-1F^T+Q_k-1。

优选地，所述步骤S3-1中，利用转换后的量测方程计算量测一步预测，结合位置量测协方差矩阵和相邻两时刻位置量测噪声的互协方差，计算量测预测协方差时，量测的一步预测表达式为：

得到量测预测协方差的表达式为：

由于相邻两时刻的位置量测噪声相关，状态预测误差和量测噪声之间存在相关性，通过展开状态预测误差计算得到状态预测误差与量测噪声之间的协方差为：

进而得到量测预测协方差表达式为：

优选地，所述步骤S3-2中利用状态一步预测、量测一步预测，并结合相邻两时刻位置量测噪声的互协方差，计算状态预测与量测预测间的互协方差时，状态预测与量测预测间的互协方差表达式为：

优选地，所述步骤S3-3中利用状态一步预测、状态预测协方差、量测一步预测、量测预测协方差，以及状态预测和量测预测之间的互协方差，更新目标状态估计和状态估计误差协方差时，状态估计更新表达式为：

状态估计误差协方差表达式为：

本发明的上述技术方案具有如下优点：本发明针对由于滑窗处理、时间累积等步骤而导致的量测噪声相关问题，提供了一种噪声相关系统跟踪滤波方法，该方法将极坐标系中的相邻两时刻的距离和方位角量测噪声相关问题转化为直角坐标系中相邻两时刻的x方向和y方向的位置量测噪声存在相关性的问题。首先完成量测转换，其次计算转换过程中的偏差项和再次利用三角等式，计算相邻两时刻目标位置量测噪声的互协方差矩阵，最后基于最小均方误差准则构建滤波器，在滤波过程中对相邻时刻量测噪声的相关性进行补偿，实现目标的有效跟踪。该方法不受噪声相关系数小或大的制约，跟踪精度高、效果好，性能稳定。

附图说明

图1是本发明实施例中提供的噪声相关系统跟踪滤波方法步骤示意图；

图2示出了量测噪声相关系数为ρ＝0.6的噪声相关示意图；

图3示出了本发明实施例中噪声相关系统跟踪滤波方法(MC-MMSE)和转换位置量测卡尔曼滤波器(CPMKF)在ρ＝0.2时的位置RMSE性能；

图4示出了MC-MMSE和CPMKF在ρ＝0.2时的速度RMSE性能；

图5示出了MC-MMSE和CPMKF在ρ＝0.5时的位置RMSE性能；

图6示出了MC-MMSE和CPMKF在ρ＝0.5时的速度RMSE性能；

图7示出了MC-MMSE和CPMKF在ρ＝0.9时的位置RMSE性能；

图8示出了MC-MMSE和CPMKF在ρ＝0.9时的速度RMSE性能。

具体实施方式

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

如图1所示，本发明实施例提供的一种噪声相关系统跟踪滤波方法，包括如下步骤：

S1、从观测雷达处获取k时刻目标量测数据并建立噪声相关系统的状态空间模型，其中量测数据包括距离量测数据和方位角量测数据。

优选地，在笛卡尔直角坐标系中，建立噪声相关系统的状态空间模型为：

x_k＝Fx_k-1+w_k-1；

其中，x_k＝[p_x,k,p_y,k,v_x,k,v_y,k]^T是离散时刻k∈{1,2,3,...}目标的状态，(p_x,k,p_y,k)表示目标在x方向和y方向的位置，(v_x,k,v_y,k)表示目标在x方向和y方向的速度，F是状态x_k的转移矩阵。过程噪声w_k-1是均值为0、协方差为Q_k-1的白色高斯随机变量，n_k是量测噪声。z_k是极坐标下的k时刻目标量测数据，包括目标距离量测数据和目标方位角量测数据两种量测数据，和分别表示距离和方位角量测数据。

优选地，状态与量测数据之间满足关系式：

其中，表示距离量测噪声和方位角量测噪声，它们被假设为均值为0的高斯随机变量，且k-1时刻和k时刻之间的量测噪声存在相关性，对应的协方差分别表示为：

其中，ρ为相邻两时刻量测噪声相关系数，σ_r、σ_θ分别表示距离和方位角量测噪声的标准差。此外，量测噪声n_k严格独立于过程噪声w_k。

假设目标以几乎恒定的速度移动，即为匀速(constant velocity，CV)运动模型，因此状态转移矩阵F和过程噪声协方差矩阵Q_k-1的表达式分别为：

其中，T是采样间隔。而q_x和q_y分别代表了x方向和y方向的过程噪声功率谱密度。

S2、将k时刻距离量测数据和方位角量测数据转换为直角坐标系下目标在x方向和y方向的位置量测数据，计算量测转换中的偏差项以及位置量测协方差矩阵，构建转换后的量测方程，并计算k-1时刻和k时刻目标位置量测噪声的互协方差。

极坐标系下的距离量测数据和方位角量测数据可以被转换为直角坐标系下目标在x方向和y方向的位置量测数据，将k时刻距离量测数据和方位角量测数据转换为直角坐标系下目标在x方向和y方向的位置量测数据，并对转换后位置量测数据进行去偏处理，得到表达式为：

其中，表示经过转换且去除偏差项后的目标位置量测数据，分别表示利用距离和方位角量测数据计算出的目标在x方向、y方向的位置量测数据，分别表示量测转换中x方向、y方向上的偏差项。

所述步骤S2中计算量测转换中的偏差项以及位置量测协方差矩阵，构建转换后的量测方程，并计算k-1时刻和k时刻目标位置量测噪声的互协方差时，量测转换中的偏差项分别表示为：

转换后的位置量测协方差矩阵表示为：

其中，每一项为：

构建转换后的量测方程，具体表示为：

其中，H为量测矩阵，在匀速运动模型中表达式为：

表示转换后的量测噪声，是均值为零、协方差矩阵为的高斯随机变量，具体表达式为：

转换后的相邻两时刻目标位置量测噪声的互协方差的表达式为：

其中，每一项为：

S3、在最小均方误差(MMSE)的准则下构建滤波器，利用k-1时刻和k时刻量测噪声的互协方差补偿两个时刻噪声之间的相关性，对转换后的位置量测数据进行滤波处理，更新k时刻目标状态估计和状态估计误差协方差，完成目标跟踪。

对于每一时刻的滤波过程，需要重新计算量测转换中的偏差项、位置量测协方差矩阵和相邻两时刻目标位置量测噪声的互协方差矩阵对转换后的位置量测数据进行滤波处理，实现目标跟踪。

优选地，步骤S3具体包括如下步骤：

S3-1、利用由状态空间模型得到的状态方程计算状态一步预测及状态一步预测协方差，同时利用转换后的量测方程计算量测一步预测，结合位置量测协方差矩阵和相邻两时刻位置量测噪声的互协方差，计算量测预测协方差。

优选地，若k＝1或k＝2时，采用两点差分法进行状态和协方差矩阵的初始化。

k＝1时刻目标状态估计为：

k＝2时刻目标状态估计为：

其中，z(1,1)和z(2,1)表示k＝1时刻目标在x方向和y方向上的位置量测数据，z(1,2)和z(2,2)表示k＝2时刻目标在x方向和y方向上的位置量测数据；；

k＝1时刻状态协方差为：

k＝2时刻状态协方差为：

其中，为k＝2时刻的位置量测协方差矩阵；

优选地，初始化后，利用由状态空间模型得到的状态方程计算状态一步预测及状态一步预测协方差时，状态一步预测表达式为：

状态一步预测协方差为：

P_k|k-1＝FP_k-1|k-1F^T+Q_k-1；

利用转换后的量测方程计算量测一步预测，结合位置量测协方差矩阵和相邻两时刻位置量测噪声的互协方差，计算量测预测协方差时，量测的一步预测表达式为：

进一步得到，量测预测协方差的表达式为：

由于相邻两时刻的位置量测噪声相关，状态预测误差和量测噪声之间存在相关性，通过展开状态预测误差计算得到状态预测误差与量测噪声之间的协方差为：

进一步地，k时刻量测预测协方差可以重写为：

S3-2、利用状态一步预测、量测一步预测，并结合相邻两时刻位置量测噪声的互协方差，计算状态预测与量测预测间的互协方差。

优选地，状态预测与量测预测间的互协方差表达式为：

在步骤S3-1和步骤S3-2中，求解量测预测协方差以及状态预测和量测预测之间的互协方差时，状态预测误差与量测噪声存在相关性，对应的协方差和不再为零，需要根据相邻两时刻位置量测噪声的互协方差矩阵并结合最小均方误差准则进行计算，并在滤波过程中予以补偿，以此来提高目标状态估计的准确性。

S3-3、利用状态一步预测、状态预测协方差、量测一步预测、量测预测协方差，以及状态预测和量测预测之间的互协方差，更新目标状态估计和状态估计误差协方差，实现当前时刻的目标跟踪。

最后运用上述计算结果，根据最小均方误差准则完成本发明实施例中噪声相关系统跟踪滤波方法(MC-MMSE)的状态估计更新表达式为：

状态估计误差协方差表达式为：

如图2所示，滑窗处理、时间累积等步骤后，相邻时刻的量测噪声不再是高斯白噪声，而是存在一定相关性。本发明提供的噪声相关系统跟踪滤波方法通过转换量测的方式求解噪声相关互协方差矩阵，能够有效提高噪声相关数据的跟踪性能，并且适用于不同量测噪声相关系数。

为验证本发明所提供方法的效果，利用仿真数据进行蒙特卡洛实验。目标从位置[30,30]^Tkm开始以恒定速度[15,15]^Tm/s运动。观测雷达固定在原点位置，在采样间隔为T＝1s的规律时刻t_k＝kT,k∈{1,2,...,100}收集距离和方位角量测值。

k时刻的雷达的距离量测噪声σ_r＝0.1km，距离量测噪声之间的互协方差为：且k时刻雷达的方位角量测噪声σ_θ＝0.1°，方位角量测噪声之间的互协方差为：且

其中，为了确保协方差矩阵为半正定矩阵，相关系数的取值范围为：0≤ρ≤1，距离量测噪声与方位角量测噪声严格相互独立。过程噪声的功率谱密度设置为q_x＝q_y＝0.01m²/s³。仿真结果基于M＝100次的蒙特卡洛仿真实验。

图3-图8分别显示了量测噪声相关系数为ρ＝0.2,0.5,0.9时，目标位置或速度估计的最小均方根误差(RMSE)性能。从图3-图8中可以看出，在噪声相关系数比较小的情况下，传统方法转换位置量测卡尔曼滤波器(CPMKF)的RMSE性能与本发明所述新方法MC-MMSE的性能相近。随着相关系数的增大，传统的CPMKF方法的RMSE性能开始变差，在ρ＝0.5时，跟踪初期会有许多突出的奇点，虽然随着跟踪时间增大，跟踪效果有所改善，但是与新方法相比，其RMSE性能还是较差。这种情况在相关系数继续增大时变得更加明显，在ρ＝0.9时，CPMKF的RMSE性能极差。而本发明提出的新方法MC-MMSE，无论在量测噪声相关系数大或者小时都表现良好，滤波过程稳定，所得位置和速度估计的RMSE呈现出收敛的趋势，对量测噪声存在相关性的情况表现出良好的适应性。

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

Claims (9)

1.一种噪声相关系统跟踪滤波方法，其特征在于，包括如下步骤：

S1、从观测雷达处获取k时刻目标量测数据并建立噪声相关系统的状态空间模型，其中量测数据包括距离量测数据和方位角量测数据；

S2、将k时刻距离量测数据和方位角量测数据转换为直角坐标系下目标在x方向和y方向的位置量测数据，计算量测转换中的偏差项以及位置量测协方差矩阵，构建转换后的量测方程，并计算k-1时刻和k时刻目标位置量测噪声的互协方差；

S3、在最小均方误差的准则下构建滤波器，利用k-1时刻和k时刻量测噪声的互协方差补偿两个时刻噪声之间的相关性，对转换后的位置量测数据进行滤波处理，更新k时刻目标状态估计和状态估计误差协方差，完成目标跟踪。

2.根据权利要求1所述的噪声相关系统跟踪滤波方法，其特征在于，所述步骤S3包括如下步骤：

S3-1、利用由状态空间模型得到的状态方程计算状态一步预测及状态一步预测协方差，同时利用转换后的量测方程计算量测一步预测，结合位置量测协方差矩阵和相邻两时刻位置量测噪声的互协方差，计算量测预测协方差；

S3-2、利用状态一步预测、量测一步预测，并结合相邻两时刻位置量测噪声的互协方差，计算状态预测与量测预测间的互协方差；

S3-3、利用状态一步预测、状态预测协方差、量测一步预测、量测预测协方差，以及状态预测和量测预测之间的互协方差，更新目标状态估计和状态估计误差协方差，实现当前时刻的目标跟踪。

3.根据权利要求2所述的噪声相关系统跟踪滤波方法，其特征在于，所述步骤S1中，建立噪声相关系统的状态空间模型为：

x_k＝Fx_k-1+w_k-1；

其中，x_k＝[p_x,k,p_y,k,v_x,k,v_y,k]^T是离散时刻k∈{1,2,3,...}目标的状态，(p_x,k,p_y,k)表示目标在x方向和y方向的位置，(v_x,k,v_y,k)表示目标在x方向和y方向的速度，F是状态x_k的转移矩阵，过程噪声w_k-1是均值为0、协方差为Q_k-1的白色高斯随机变量，n_k是量测噪声，z_k是极坐标下的目标量测数据，和分别表示距离和方位角量测数据；

状态与量测数据之间满足关系式：

其中， 表示距离量测噪声和方位角量测噪声，是均值为0的高斯随机变量，且k-1时刻和k时刻之间的量测噪声存在相关性，对应的协方差分别表示为：

其中，ρ为相邻两时刻量测噪声相关系数，σ_r、σ_θ分别表示距离和方位角量测噪声的标准差，量测噪声n_k严格独立于过程噪声w_k；

匀速运动模型中，状态转移矩阵F和过程噪声协方差矩阵Q_k-1的表达式分别为：

其中，T是采样间隔，q_x和q_y分别代表x方向和y方向的过程噪声功率谱密度。

4.根据权利要求3所述的噪声相关系统跟踪滤波方法，其特征在于，所述步骤S2中将k时刻距离量测数据和方位角量测数据转换为直角坐标系下目标在x方向和y方向的位置量测数据时，表达式为：

其中，表示经过转换且去除偏差项后的目标位置量测数据，分别表示利用距离和方位角量测数据计算出的目标在x方向、y方向的位置量测数据，分别表示量测转换中x方向、y方向上的偏差项。

5.根据权利要求4所述的噪声相关系统跟踪滤波方法，其特征在于，所述步骤S2中计算量测转换中的偏差项以及位置量测协方差矩阵，构建转换后的量测方程，并计算k-1时刻和k时刻目标位置量测噪声的互协方差时，量测转换中的偏差项分别表示为：

位置量测协方差矩阵表示为：

构建转换后的量测方程，表达式为：

其中，H为量测矩阵，在匀速运动模型中表达式为：

表示转换后的量测噪声，是均值为零、协方差矩阵为的高斯随机变量，表达式为：

相邻两时刻目标位置量测噪声的互协方差的表达式为：

6.根据权利要求5所述的噪声相关系统跟踪滤波方法，其特征在于，所述步骤S3-1中，利用由状态空间模型得到的状态方程计算状态一步预测及状态一步预测协方差时，状态一步预测表达式为：

状态一步预测协方差为：

P_k|k-1＝FP_k-1|k-1F^T+Q_k-1。

7.根据权利要求6所述的噪声相关系统跟踪滤波方法，其特征在于，所述步骤S3-1中，利用转换后的量测方程计算量测一步预测，结合位置量测协方差矩阵和相邻两时刻位置量测噪声的互协方差，计算量测预测协方差时，量测的一步预测表达式为：

得到量测预测协方差的表达式为：

由于相邻两时刻的位置量测噪声相关，状态预测误差和量测噪声之间存在相关性，通过展开状态预测误差计算得到状态预测误差与量测噪声之间的协方差为：

进而得到量测预测协方差表达式为：

8.根据权利要求7所述的噪声相关系统跟踪滤波方法，其特征在于，所述步骤S3-2中利用状态一步预测、量测一步预测，并结合相邻两时刻位置量测噪声的互协方差，计算状态预测与量测预测间的互协方差时，状态预测与量测预测间的互协方差表达式为：

9.根据权利要求8所述的噪声相关系统跟踪滤波方法，其特征在于，所述步骤S3-3中利用状态一步预测、状态预测协方差、量测一步预测、量测预测协方差，以及状态预测和量测预测之间的互协方差，更新目标状态估计和状态估计误差协方差时，状态估计更新表达式为：

状态估计误差协方差表达式为：