[发明专利]一种改进SRCKF强跟踪滤波的机动目标跟踪方法

权利要求书

1.一种改进SRCKF强跟踪滤波的机动目标跟踪方法，其特征在于，该方法包括以下步骤：

S1根据目标运动情况和卫星的观测情况建立运动模型和观测模型；

S1.1根据目标运动情况建立运动模型，当目标为自由段运动的导弹目标时，在不机动的条件下，可以视作仅受重力影响，使用J₂摄动模型对目标运动进行建模；

S1.2根据卫星的观测情况建立观测模型，卫星为被动式观测，仅仅能够测得目标的方位角和俯仰角，在ECI坐标系下，假设k时刻卫星的状态向量为：

Sa(k)＝[x_s(k),y_s(k),z_s(k),v_sx(k),v_sy(k),v_sz(k)]^T

其中x_s(k),y_s(k),z_s(k)、v_sx(k),v_sy(k),v_sz(k)分别代表卫星在k时刻的位置和速度；

以卫星速度方向为X_b轴，以卫星到地球球心连线为Z_b轴，按右手螺旋定则给出Y_b轴，定义X_b轴，Y_b轴，Z_b轴的轴矢量分别为SX_b，SY_b，SZ_b，给出星体坐标系，α_s代表目标方位角，β_s代表目标俯仰角，SX代表卫星视线，此时可由如下方式计算α_s、β_s：

SX_b＝[v_x(k),v_y(k),v_z(k)]^T

SZ_b＝-[x(k),y(k),z(k)]^T

SX＝[x_s(k)-x(k),y_s(k)-y(k),z_s(k)-z(k)]^T

为了与滤波模型切合，一般将上式写为：

其中Z(k)代表卫星k时刻的量测，dβ_s和dα_s代表传感器测量噪声，r(k)为零均值高斯矢量；

S2：滤波初始化，滤波初始化过程为：

其中，代表对目标初始时刻的状态估计，E(·)代表对矩阵求期望，P(0)代表目标初始时刻状态估计值的协方差矩阵，S(0)代表目标初始时刻协方差矩阵的平方根，chol(·)代表对矩阵进行Cholesky分解得到下三角矩阵；

同时由于S5中需要统计渐消因子变化次数，初始时刻，定义变化次数flag＝0；

S3：进行滤波预测，也即在已知k时刻的系统状态时，对k+1时刻的系统状态进行外推的过程，仍按照SRCKF方法的框架，有：

对于k∈{0,…,∞}：

计算容积点：

其中代表k时刻容积点矩阵υ_k的第i列，N代表状态向量的维数，代表第i个容积点，代表k时刻的系统状态估计值，S(k)代表k时刻系统状态估计值协方差矩阵的平方根，[1]_i代表矩阵[1]中第i个元素，矩阵[1]如下所示：

预测过程如下：

υ_k+1|k＝f(υ_k)

S^e(k+1|k)＝qr([V_k+1|k S_Q]^T)

其中代表利用k时刻容积点值预测出k+1时刻值容积点矩阵υ_k+1|k的第i列，V_k+1|k代表容积点矩阵υ_k+1|k与k+1时刻的预测值的偏差，qr(·)即表示对矩阵进行QR分解得到的上三角矩阵求转置的函数，S_Q＝chol(Q(k))，代表预测过程中运动模型噪声的协方差矩阵的平方根，Q(k)＝q(k)q(k)^T代表k时刻预测过程中运动模型噪声的协方差矩阵；

即得到了k+1时刻未带有渐消因子的预测值及其协方差的平方根S^e(k+1|k)；

S4：进行第一次量测更新，具体步骤如下：

传播容积点：

其中代表利用k+1时刻的预测值和协方差平方根得到的容积点矩阵的第i列，表示对容积点矩阵按照S1.2的观测方程进行观测，代表对进行观测的容积点矩阵ψ_k+1|k的第i列；

未加渐消因子时量测更新过程如下：

S_zz^e(k+1)＝qr([Z_k+1S_R]^T)

P_xz^e(k+1)＝X_k+1Z_k+1^T

其中Z_k+1代表k+1时刻容积点矩阵ψ_k+1|k与预测量测值的偏差，X_k+1代表k+1时刻容积点矩阵与k+1时刻的预测值的偏差，S_R＝chol(R(k+1))，代表k+1时刻观测方程噪声所产生的协方差矩阵的平方根，R(k+1)＝r(k+1)r(k+1)^T即是k+1时刻观测方程噪声所产生的协方差矩阵，z(k+1)代表k+1时刻的真实量测值，Γ(k+1)代表对新息增量协方差的估计值，ρ为遗忘因子；

即得到了未加渐消因子时的预测量测值及其协方差的平方根S_zz^e(k+1)以及预测值和预测量测值的互协方差矩阵P_xz^e(k+1)，以及新息增量ε(k+1)；

S5：计算带有时间因子的渐消因子并进行门限判决，此时需要利用S3和S4得出的S^e(k+1|k)、S_zz^e(k+1)和P_xz^e(k+1)，此时有：

h_x(k+1)＝P_xz^e(k+1)^T/S^e(k+1|k)/S^e(k+1|k)^T

其中Norm(·)代表求矩阵各个元素的平方和，h_x(k+1)代表k+1时刻对观测过程中对预测过程中运动模型噪声的协方差矩阵平方根的增益矩阵，l(k+1)即代表k+1时刻的时间因子项，式中t_b代表观测时的采样间隔，T_s代表卫星观测时的收敛时间，λ₀代表渐消因子的计算初值，λ(k+1)代表k+1时刻的渐消因子终值；

此时若λ(k+1)＞1，则对变化次数flag加1，并将flag与预先设定的门限值n_λ比较，若flag大于等于n_λ则令k＝k-n_λ，即是向前回溯n_λ个时刻，将变化次数flag重置为0，将滤波过程中所有含有k的量都变为新的时刻所代表的值，并令若flag小于n_λ，则令S(k+1|k)＝S^e(k+1|k)，继续正常滤波过程；若λ(k+1)＝1，则直接将变化次数flag置为0，同时令S(k+1|k)＝S^e(k+1|k)，继续正常滤波过程；

S6：根据S5的门限判决结果再次进行量测更新；

根据S5得出的S(k+1|k)的值，对k+1时刻得到的真实量测进行更新：

传播容积点：

量测更新过程如下：

S_zz(k+1)＝qr([Z_k+1S_R]^T)

P_xz(k+1)＝X_k+1Z_k+1^T

W(k+1)＝P_xz(k+1)/S_zz(k+1)/S_zz(k+1)^T

S(k+1)＝qr([X_k+1-W(k+1)Z_k+1W(k+1)S_R]^T)

其中代表k+1时刻的预测量测值，S_zz(k+1)代表协方差的平方根，P_xz(k+1)代表k+1时刻的预测值和预测量测值的互协方差矩阵，W(k+1)表示k+1时刻的滤波增益矩阵；

即得到了k+1时刻的状态估计向量及其协方差矩阵的平方根S(k+1)，从而可以继续迭代，推导出所有时刻的状态估计值。