CN103455675A - 一种基于ckf的非线性异步多传感器信息融合方法 - Google Patents

Abstract

本发明公开了一种基于CKF的非线性异步多传感器信息融合方法。分别对各个传感器分别利用CKF估计出各自的状态变量，然后采用细分时间片方法将信息融合中心的时间间隔设定为各传感器间最高精度时间单位，在相应时刻对异步多传感器的估计结果进行判断和融合，得到更加精确的状态变量估计结果。本发明可以增强对异步多传感器信息的利用率，大幅提高多传感器***中状态变量的估计精度，增强***的生存能力。

Description

一种基于CKF的非线性异步多传感器信息融合方法

技术领域

本发明涉及的是一种非线性***中多传感器信息融合技术，特别是涉及一种在多传感器异步情况下的信息融合方法。

背景技术

由于多传感器信息融合技术在军事、国防及高科技领域的广泛应用，对它的研究一直未曾间断。在多传感器信息融合理论中，利用多个传感器量测中的冗余和互补信息，通过多源信息融合技术可以有效的提升***精度。但在解决实际问题过程中，常遇到所用传感器具有不同的采样频率，各传感器之间存在通信延迟和固有延迟不同等问题，统称为信息融合的异步问题。若不能有效地解决异步信息融合问题，提高***精度的目标很难实现。因此，对传感器异步信息融合问题的研究有较大的工程价值和实际意义。

尽管大多数传感器被假设成线性的，但是它们基本上或多或少都会有非线性特征，所以将***近似为线性必然会产生误差，所以多传感器信息融合***中状态估计应该选取具有非线性处理能力的滤波器。对于非线性***来说，通常使用的滤波方法为扩展卡尔曼和无迹卡尔曼滤波方法。虽然EKF实现简便、使用广泛，但是当***为强非线性时，很容易产生线性化误差，造成滤波器精度下降。虽然UKF滤波精度高于EKF，但是当处理高维数***时容易出现“维数灾难”。最近提出的容积卡尔曼滤波(Cubature Kalman filter，CKF)由于其估计精度高，不容易发散且计算量小的优点，迅速被应用到各种领域。因此本文研究的是基于CKF的非线性异步信息融合方法。

发明内容

本发明的目的是提供可以在非线性***中多传感器异步情况下大幅提高***状态变量估计精度的一种基于CKF的异步多传感器信息融合方法。

本发明的目的是这样是实现的：

(1)对多传感器***进行充分预热，并采集各个传感器的量测信息z_i(k_i)，其中i＝1，2，…，N(N为传感器的个数)，k_i为第i个传感器的采样时刻。

(2)建立非线性异步多传感器***的非线性状态方程及各传感器的量测方程。

(3)采用细分时间片方法，将融合中心的采样间隔设定为各传感器间最高精度时间单位，即取能同时被各传感器采样间隔整除的最大数为融合中心的采样间隔。

设任意两传感器的采样时间分别为T_i，T_j(i，j＝1，2，…，N，且i≠j)，N为传感器的总个数。且有：

则融合中心的采样间隔为：

$\mathrm{\ΔT}=\min \{{\mathrm{\ΔT}}_i=\frac{T_i}{N_i}|i=\mathrm{1,2},\·\·\·,N\}$

(4)根据融合中心的采样时刻，在融合中心对各传感器的估计状态进行判断融合：

a)k时刻，仅有一个传感器S有量测信息时，利用CKF对该传感器进行状态估计，得到该传感器的状态变量估计值

和协方差阵P_S(k)，则***的最优状态估计即为该传感器的状态变量估计值：

$\hat{x}\left(k\right)={\hat{x}}_S\left(k\right)---\left(1\right)$

b)k时刻，有m(m∈{2，3，…，N}，其中N为传感器的个数)个传感器S₁，S₂，…，S_m有量测信息时，利用CKF对这些传感器分别进行滤波估计，得到各传感器的状态变量估计值

和协方差阵P₁(k)，P₂(k)，…，P_m(k)。然后对各个传感器的估计值进行加权处理，融合后***的状态变量估计值为：

$\hat{x}\left(k\right)=D_1{\hat{x}}_1\left(k\right)+D_2{\hat{x}}_2\left(k\right)+\·\·\·+D_m{\hat{x}}_m\left(k\right)---\left(2\right)$

其中，D₁，D₂，…，D_m为相应传感器状态变量估计的权值。D₁，D₂，…，D_m可由式(3)～(5)获得：

$\left\{\begin{array}{c}{D}_1=S_1\mathrm{\Φ}(k+1,k)\\ D_2=S_2\mathrm{\Φ}(k+1,k)\\ \·\\ \·\\ \·\\ D_m=S_m\mathrm{\Φ}(k+1,k)\end{array}---\left(3\right)\right.$

$\left\{\begin{array}{c}{S}_1={(I+\underset{i=2}{\overset{m-1}{\mathrm{\Σ}}}\left(A_1{A}_i^{-1}\right)+A_1{A}_m^{-1})}^{-1}\\ S_2=S_1{A}_1{A}_2^{-1}\\ \·\\ \·\\ \·\\ S_m=S_1{A}_1{A}_m^{-1}\end{array}---\left(4\right)\right.$

$\left\{\begin{array}{c}{A}_1=P_1\left(k\right)+P_1^T\left(k\right)\\ A_2=P_2\left(k\right)+P_2^T\left(k\right)\\ \·\\ \·\\ \·\\ A_m=P_m\left(k\right)+P_m^T\left(k\right)\end{array}---\left(5\right)\right.$

其中，Ф(k+1，k)为异步多传感器***中k时刻的状态转移矩阵，P_m(k)为异步多传感器***中第m个传感器在k时刻的方差阵的估计值。

c)k时刻，各传感器均没有量测信息时，便利用前一时刻的估计值进行时间更新。

$\hat{x}\left(k\right)=\hat{x}\left(k\right|k-1)---\left(6\right)$

本发明具有如下优点：

应用多传感器信息融合技术，增大了量测信息的数据利用率，提高了传感器***状态变量的估计精度，增强了***的容错性和生存能力；应用新的非线性滤波方法CKF，提高了滤波精度和滤波收敛速度，降低了计算量，提高了异步信息融合的实时性，从而解决了非线性***中的滤波问题。所以本发明可以大幅提高非线性异步多传感器***中状态变量的估计精度。

附图说明

图1为本发明的流程图；

图2为异步多传感器信息融合采样原理图；

图3为利用本发明与利用单个传感器信息时位置变量估计结果的RMS对比曲线图；

图4为利用本发明与利用单个传感器信息时速度变量估计结果的RMS对比曲线图；

图5为利用本发明与利用传统算法时位置变量估计结果的RMS对比曲线图；

图6为利用本发明与利用传统算法时速度变量估计结果的RMS对比曲线图；

具体实施方式

以下结合具体实施例，对本发明进行详细说明。

(1)对多传感器***进行充分预热，并采集各个传感器的量测信息z_i(k_i)，其中i＝1，2，…，N(N为传感器的个数)，k_i为第i个传感器的采样时刻。

(2)建立非线性异步多传感器***的非线性状态方程及各传感器的量测方程。

(3)采用细分时间片方法，将融合中心的采样间隔设定为各传感器间最高精度时间单位，即取能同时被各传感器采样间隔整除的最大数为融合中心的采样间隔。

设任意两传感器的采样时间分别为T_i，T_j(i，j＝1，2，…，N，且i≠j)，N为传感器的总个数。且有：

则融合中心的采样间隔为：

$\mathrm{\ΔT}=\min \{{\mathrm{\ΔT}}_i=\frac{T_i}{N_i}|i=\mathrm{1,2},\·\·\·,N\}$

(4)根据融合中心的采样时刻，在融合中心对各传感器的估计状态进行判断融合：

a)k时刻，仅有一个传感器S有量测信息时，利用CKF对该传感器进行状态估计，得到该传感器的状态变量估计值

和协方差阵P_S(k)，则***的最优状态估计即为该传感器的状态变量估计值：

$\hat{x}\left(k\right)={\hat{x}}_S\left(k\right)---\left(1\right)$

b)k时刻，有m(m∈{2，3，…，N}，其中N为传感器的个数)个传感器S₁，S₂，…，S_m有量测信息时，利用CKF对这些传感器分别进行滤波估计，得到各传感器的状态变量估计值

和协方差阵P₁(k)，P₂(k)，…，P_m(k)。然后对各个传感器的估计值进行加权处理，融合后***的状态变量估计值为：

$\hat{x}\left(k\right)=D_1{\hat{x}}_1\left(k\right)+D_2{\hat{x}}_2\left(k\right)+\·\·\·+D_m{\hat{x}}_m\left(k\right)---\left(2\right)$

其中，D₁，D₂，…，D_m为相应传感器状态变量估计的权值。D₁，D₂，…，D_m可由式(3)～(5)获得：

$\left\{\begin{array}{c}{D}_1=S_1\mathrm{\Φ}(k+1,k)\\ D_2=S_2\mathrm{\Φ}(k+1,k)\\ \·\\ \·\\ \·\\ D_m=S_m\mathrm{\Φ}(k+1,k)\end{array}---\left(3\right)\right.$

$\left\{\begin{array}{c}{S}_1={(I+\underset{i=2}{\overset{m-1}{\mathrm{\Σ}}}\left(A_1{A}_i^{-1}\right)+A_1{A}_m^{-1})}^{-1}\\ S_2=S_1{A}_1{A}_2^{-1}\\ .\\ .\\ .\\ S_m=S_1{A}_1{A}_m^{-1}\end{array}\right.---\left(4\right)$

$\left\{\begin{array}{c}{A}_1=P_1\left(k\right)+P_1^T\left(k\right)\\ A_2=P_2\left(k\right)+P_2^T\left(k\right)\\ \·\\ \·\\ \·\\ A_m=P_m\left(k\right)+P_m^T\left(k\right)\end{array}\right.---\left(5\right)$

其中，Φ(k+1，k)为异步多传感器***中k时刻的状态转移矩阵，P_m(k)为异步多传感器***中第m个传感器在k时刻的方差阵的估计值。

c)k时刻，各传感器均没有量测信息时，便利用前一时刻的估计值进行时间更新。

$\hat{x}\left(k\right)=\hat{x}\left(k\right|k-1)---\left(6\right)$

Matlab仿真实验：

在以下的目标跟踪模型的仿真条件下，对本发明进行仿真实验：

4维目标跟踪的离散动态模型为：

$x\left(k\right)=\left[\begin{array}{cccc}1& \frac{\sin \mathrm{\Ω\Δ}}{\mathrm{\Ω}}& 0& \frac{\cos \mathrm{\Ω\Δ}-1}{\mathrm{\Ω}}\\ 0& \cos \mathrm{\Ω\Δ}& 0& -\sin \mathrm{\Ω\Δ}\\ 0& \frac{1-\cos \mathrm{\Ω\Δ}}{\mathrm{\Ω}}& 1& \frac{\sin \mathrm{\Ω\Δ}}{\mathrm{\Ω}}\\ 0& \sin \mathrm{\Ω\Δ}& 0& \cos \mathrm{\Ω\Δ}\end{array}\right]x(k-1)+w(k-1)---\left(7\right)$

其中，***的状态向量 $x={\left[\begin{array}{cccc}x& \stackrel{\·}{x}& y& \stackrel{\·}{y}\end{array}\right]}^T,$ x、y表示目标的位置，

表示目标的速度，Δ为两次量测的时间间隔，Ω为目标运动的角速度，且Ω=3°/s，***噪声w(k)～N(0，Q)。

为了简单起见，我们以两个传感器观测为例，传感器对运动目标进行测量，传感器与目标之间的斜距为r，方位角为θ，则观测方程为：

$\left[\begin{array}{c}{r}_i\left(k\right)\\ {\mathrm{\θ}}_i\left(k\right)\end{array}\right]=\left[\begin{array}{c}\sqrt{x_i{\left(k\right)}^2+y_i{\left(k\right)}^2}\\ {\tan }^{-1}\left(\frac{y_i\left(k\right)}{x_i\left(k\right)}\right)\end{array}\right]+v_i\left(k\right)---\left(8\right)$

其中，i＝1，２量测噪声v_i(k)～N(0，R_i)。

***的噪声方差阵：

$Q=0.1*\left[\begin{array}{cccc}{\mathrm{\Δ}}^3/3& {\mathrm{\Δ}}^2/2& 0& 0\\ {\mathrm{\Δ}}^2/2& \mathrm{\Δ}& 0& 0\\ 0& 0& {\mathrm{\Δ}}^3/3& {\mathrm{\Δ}}^2/2\\ 0& 0& {\mathrm{\Δ}}^2/2& \mathrm{\Δ}\end{array}\right]$

传感器1和传感器2的量测噪声方差阵：

$R_1=\left[\begin{array}{cc}{\mathrm{\σ}}_r/2& 0\\ 0& {\mathrm{\σ}}_{\mathrm{\θ}}/2\end{array}\right],R_2=\left[\begin{array}{cc}{\mathrm{\σ}}_r& 0\\ 0& {\mathrm{\σ}}_{\mathrm{\θ}}\end{array}\right]$

其中σ_r=10m， ${\mathrm{\σ}}_{\mathrm{\θ}}=\sqrt{10}\mathrm{mrad}.$

***状态向量的初始值：x₀=[100m 10m/s 100m 10m/s]^T；

状态向量的初始协方差：P₀=diag[100m² 10m²/s² 100m² 10m²/s²]；

传感器1的采样周期为：T₁=2s；

传感器2的采样周期为：T₂=3s；

融合中心的采样周期为：T₀=1s；

Monte Carlo仿真次数为：N=100；

每次仿直时间：T=100s。

采用均方根误差RMS评估滤波结果，k时刻位置、速度的RMS分别定义如下：

${\mathrm{RMS}}_{\mathrm{pos}}\left(k\right)=\sqrt{\frac{1}{N}\underset{n=1}{\overset{N}{\mathrm{\Σ}}}({(x_n\left(k\right)-x_n\left(k\right|k))}^2+{(y_n\left(k\right)-y_n\left(k\right|k))}^2)}$

${\mathrm{RMS}}_{\mathrm{vel}}\left(k\right)=\sqrt{\frac{1}{N}\underset{n-1}{\overset{N}{\mathrm{\Σ}}}({({\stackrel{\·}{x}}_n\left(k\right)-{\stackrel{\·}{x}}_n\left(k\right|k))}^2+{({\stackrel{\·}{y}}_n\left(k\right)-{\stackrel{\·}{y}}_n\left(k\right|k))}^2)}$

其中x_n(k)，y_n(k)，

为理论值，x_n(k|k)，y_n(k|k)，

为滤波估计值，N为Monte Carlo仿真次数。

利用本发明所述的方法对非线性异步多传感器***状态变量进行估计，并与利用单个传感器信息进行估计所得的结果进行了对比，比较结果如图2和图3所示。随后又将本发明与现有的基于UKF的多传感器信息融合算法进行了比较，比较结果如图4和图5所示。

通过图2到图5可以看出，在非线性***中多传感器异步的条件下，利用本发明可以获得很高的状态变量估计精度，而几估计曲线收敛速度很快，从而有效的满足非线性异步多传感器***中状态变量估计结果高精度与快速性的条件。

应当理解的足，对本领域普通技术人员来说，可以根据上述说明加以改进或变换，而所有这些改进和变换都应属于本发明所附权利要求的保护范围。

Claims (1)

1.一种基于CKF的非线性异步多传感器信息融合新方法，其特征是：

(1)对多传感器***进行预热，并采集各个传感器的量测信息；

(2)建立非线性异步多传感器***的非线性状态方程和量测方程；

(3)采用细分时间片方法，将融合中心的采样间隔设定为各传感器间最高精度时间单位，即取能同时被各传感器采样间隔整除的最大数为融合中心的采样间隔；

设任意两传感器的采样时间分别为T_i，T_j(i，j＝1，2，…，N，且i≠j)，N为传感器的总个数；且有：

则融合中心的采样间隔为：

$\mathrm{\ΔT}=\min \{{\mathrm{\ΔT}}_i=\frac{T_i}{N_i}|i=\mathrm{1,2},\·\·\·,N\}$

(4)根据异步多传感器***融合中心的采样时刻，在融合中心对各传感器的估计状态进行判断融合：

a)k时刻，仅有一个传感器S有量测信息时，利用CKF对该传感器进行状态估计，得到该传感器的状态变量估计值

和协方差阵P_S(k)，则***的最优状态估计即为该传感器的状态变量估计值：

$\hat{x}\left(k\right)={\hat{x}}_S\left(k\right)---\left(1\right)$

b)k时刻，有m(m∈{2，3，…，N}，其中N为传感器的个数)个传感器S₁，S₂，…，S_m有量测信息时，利用CKF对这些传感器分别进行滤波估计，得到各传感器的状态变量估计值

和协方差阵P₁(k)，P₂(k)，…，P_m(k)；然后对各个传感器的估计值进行加权处理，融合后***的状态变量估计值为：

$\hat{x}\left(k\right)=D_1{\hat{x}}_1\left(k\right)+D_2{\hat{x}}_2\left(k\right)+\·\·\·+D_m{\hat{x}}_m\left(k\right)---\left(2\right)$

其中，D₁，D₂，…，D_m为相应传感器状态变量估计的权值。D₁，D₂，…，D_m可由式(3)～(5)获得：

$\left\{\begin{array}{c}{D}_1=S_1\mathrm{\Φ}(k+1,k)\\ D_2=S_2\mathrm{\Φ}(k+1,k)\\ .\\ .\\ .\\ D_m=S_m\mathrm{\Φ}(k+1,k)\end{array}\right.---\left(3\right)$

$\left\{\begin{array}{c}{S}_1={(I+\underset{i=2}{\overset{m-1}{\mathrm{\Σ}}}\left(A_1{A}_i^{-1}\right)+A_1{A}_m^{-1})}^{-1}\\ S_2=S_1{A}_1{A}_2^{-1}\\ .\\ .\\ .\\ S_m=S_1{A}_1{A}_m^{-1}\end{array}\right.---\left(4\right)$

$\left\{\begin{array}{c}{A}_1=P_1\left(k\right)+P_1^T\left(k\right)\\ A_2=P_2\left(k\right)+P_2^T\left(k\right)\\ .\\ .\\ .\\ A_m=P_m\left(k\right)+P_m^T\left(k\right)\end{array}\right.---\left(5\right)$

其中，Φ(k+1，k)为异步多传感器***中k时刻的状态转移矩阵，P_m(k)为异步多传感器***中第m个传感器在k时刻的方差阵的估计值；

c)k时刻，各传感器均没有量测信息时，便利用前一时刻的估计值进行时间更新；

$\hat{x}\left(k\right)=\hat{x}\left(k\right|k-1)---\left(6\right).$

Images