[发明专利]一种长时间导航系统全参数精度评估方法

摘要

本发明属于惯性技术领域，具体涉及一种长时间导航系统全参数精度评估方法，利用局部基准系统航行实验的原始数据及卫星信息，采用基于RTS平滑的离线处理方法，通过对存储的惯导数据及卫星数据的处理，得到基准系统高精度的位置、姿态及方位信息，作为基准系统的全参数精度评估的基准。在本发明技术方案中，k时刻的R-T-S平滑估计量是k时刻的正向Kalman滤波估计量和k+1时刻的平滑估计量的线性融合，而k+1时刻的平滑估计量有效利用了全局数据，所以对Kalman滤波进行R-T-S平滑比单纯采用Kalman滤波获得的导航精度高。

主权项

一种长时间导航系统全参数精度评估方法，其特征在于：利用局部基准系统航行实验的原始数据及卫星信息，采用基于RTS平滑的离线处理方法，通过对存储的惯导数据及卫星数据的处理，得到基准系统高精度的位置、姿态及方位信息，作为基准系统的全参数精度评估的基准；具体包括以下两个步骤以获取位置、速度及航姿数据：（1)利用惯性/卫星数据信息进行前向的闭环卡尔曼滤波过程；前向卡尔曼滤波算法采取GPS与惯性导航系统位置匹配模式，采用闭环校正方式进行；（1.1）确定算法模型和状态误差量算法采用18阶导航误差模型，选取18个状态误差量为：$X=\left[\begin{array}{cccccccccccccccccc}\mathrm{\ΔL}& \mathrm{\Δh}& \mathrm{\Δ\λ}& {\mathrm{\ΔV}}_n& {\mathrm{\ΔV}}_u& {\mathrm{\ΔV}}_e& {\mathrm{\Φ}}_n& {\mathrm{\Φ}}_u& {\mathrm{\Φ}}_e& {\▿}_x& {\▿}_y& {\▿}_z& {\mathrm{\ϵ}}_x& {\mathrm{\ϵ}}_y& {\mathrm{\ϵ}}_z& R_x^b& R_y^b& R_z^b\end{array}\right]$其中：ΔL、Δh、Δλ：惯导的纬度、经度、高度误差；ΔV_n、ΔV_u、ΔV_e：惯导北向、天向、东向速度误差；Φ_n、Φ_u、Φ_e：惯导北向、天向、东向失准角误差；惯导X、Y、Z轴加速度计零偏；ε_x、ε_y、ε_z：惯导X、Y、Z轴陀螺漂移；惯导与GPS间杆臂误差。（1.2）确定系统状态方程系统状态方程为：其中：$F=\left[\begin{array}{cccccc}{F}_{11}& F_{12}& 0_{3\×3}& 0_{3\×3}& 0_{3\×3}& 0_{3\×3}\\ F_{21}& F_{22}& F_{23}& C_b^n& 0_{3\×3}& 0_{3\×3}\\ F_{31}& F_{32}& F_{33}& 0_{3\×3}& {-C}_b^n& 0_{3\×3}\\ 0_{3\×3}& 0_{3\×3}& 0_{33}& 0_{3\×3}& 0_{3\×3}& 0_{3\×3}\\ 0_{3\×3}& 0_{3\×3}& 0_{3\×3}& 0_{3\×3}& 0_{3\×3}& 0_{3\×3}\end{array}\right]$$F_{11}=\left[\begin{array}{ccc}0& -\frac{V_n}{{(R_M+h)}^2}& 0\\ 0& 0& 0\\ \frac{V_e}{R_N+h}\tan L\sec L& -\frac{V_e\sec L}{{(R_N+h)}^2}& 0\end{array}\right];F_{12}=\left[\begin{array}{ccc}\frac{1}{R_M+h}& 0& 0\\ 0& 1& 0\\ 0& 0& \frac{\sec L}{R_N+h}\end{array}\right];$$F_{21}=\left[\begin{array}{ccc}-({2V}_e{\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{{V_e}^2{\sec }^{2}L}{R_N+h})& \frac{V_n{V}_u}{{(R_M+h)}^2}+\frac{{V_e}^2\tan L}{{(R_N+h)}^2}& 0\\ {-2V}_e{\mathrm{\ω}}_{\mathrm{ie}}\sin L& -(\frac{{V_n}^2}{{(R_M+h)}^2}+\frac{{V_e}^2}{{(R_N+h)}^2})& 0\\ {2\mathrm{\ω}}_{\mathrm{ie}}(V_u\sin L+V_n\cos L)+\frac{V_e{V}_n{\sec }^{2}L}{R_N+h}& \frac{V_e{V}_u-V_e{V}_n\tan L}{{(R_N+h)}^2}& 0\end{array}\right];$$F_{22}=\left[\begin{array}{ccc}-\frac{V_u}{R_M+h}& -\frac{V_n}{R_M+h}& -2({\mathrm{\ω}}_{\mathrm{ie}}\sin L+\frac{V_e\tan L}{R_N+h})\\ \frac{{2V}_n}{R_M+h}& 0& 2({\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_e}{R_N+h})\\ {2\mathrm{\ω}}_{\mathrm{ie}}\sin L+\frac{V_e\tan L}{R_N+h}& -({2\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_e}{R_N+h})& \frac{V_n\tan L-V_u}{R_N+h}\end{array}\right];$$F_{23}=\left[\begin{array}{ccc}0& {-f}_e& f_u\\ f_e& 0& {-f}_n\\ {-f}_u& f_n& 0\end{array}\right];$$F_{31}=\left[\begin{array}{ccc}{-\mathrm{\ω}}_{\mathrm{ie}}\sin L& -\frac{V_e}{{(R_N+h)}^2}& 0\\ {\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_e{\sec }^{2}L}{R_N+h}& -\frac{V_e\tan L}{{(R_N+h)}^2}& 0\\ 0& \frac{V_n}{{(R_M+h)}^2}& 0\end{array}\right];F_{32}=\left[\begin{array}{ccc}0& 0& \frac{1}{R_N+h}\\ 0& 0& \frac{\tan L}{R_N+h}\\ -\frac{1}{R_M+h}& 0& 0\end{array}\right];$$F_{33}=\left[\begin{array}{ccc}0& -\frac{V_n}{R_M+h}& -({\mathrm{\ω}}_{\mathrm{ie}}\sin L+\frac{V_e\tan L}{R_N+h})\\ \frac{V_n}{R_M+h}& 0& {\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_e}{R_N+h}\\ {\mathrm{\ω}}_{\mathrm{ie}}\sin L+\frac{V_e\tan L}{R_N+h}& -({\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_e}{R_N+h})& 0\end{array}\right]$其中：Vn、Vu、Ve表示惯导系统的北向、天向和东向速度；L表示纬度，h表示高度，R_M表示地球的子午圈半径，R_N表示地球的卯酉圈半径，ω_ie表示地球的自转角速度，表示惯导系统从载体坐标系到导航坐标系的姿态转换矩阵，f_n、f_u、f_e表示惯导系统加速度计敏感到的比力在导航坐标系内的北向、天向和东向表示值；（1.3）确定量测方程量测方程为：Z=HX+V其中：Z为量测量，即惯导系统的位置信息与GPS位置信息的差值，并考虑两者间杆臂误差：Z=[L‑L^GPS‑R_n h‑h^GPS‑R_u λ‑λ^GPS‑R_e]^T；h^GPS、λ^GPS、L^GPS分别表示GPS系统输出的高度信息、经度信息和纬度信息；R_n、R_u、R_e分别表示惯导系统的位置信息与GPS位置信息之间的北向杆臂误差、天向杆臂误差和东向杆臂误差；H为量测阵，$H=\left[\begin{array}{ccc}{I}_{3\×3}& 0_{3\×12}& {-C}_b^n\end{array}\right]$（1.4）确定卡尔曼滤波方程状态一步预测：其中：表示由k‑1时刻到k时刻的状态一步预测值，表示k‑1时刻的状态估计值，Φ_k,k‑1表示离散化的卡尔曼滤波状态转移矩阵；一步预测均方误差阵：其中：表示由k‑1时刻到k时刻的一步预测均方误差阵，表示k‑1时刻的估计均方误差阵，Q_k表示离散化后的系统噪声矩阵；滤波增益矩阵：其中：表示k时刻的滤波增益阵，H_k表示k时刻的系统量测阵，R_k表示离散化后的量测噪声矩阵；状态估计：其中：表示k时刻的状态估计值，Z_k表示系统观测量矩阵；估计均方误差阵：$P_k^F=(I-K_k^F{H}_k)P_{k,k-1}^F{(I-K_k^F{H}_k)}^T+K_k^F{R}_k{\left(K_k^F\right)}^T,$其中：表示k时刻的估计均方误差阵，F_i表示未离散化的系统状态转移矩阵对应元素值，I表示对应维数的单位矩阵，T_n表示惯导系统导航周期，T_f表示滤波周期；（2)在步骤（1）卡尔曼滤波的基础上，进行反向最优固定区间平滑；设整个导航时间区间为[0 N]，t表示此时间间隔中的任一时刻，0≤t≤N，固定区间平滑估值表示为（2.1）首先进行0→N的正向Kalman滤波，同时存储滤波估计出的状态转移矩阵Φ_k,k‑1、状态一步预测一步预测均方误差阵状态估计估计均方误差阵（2.2）滤波结束后，将正向滤波最后时刻N的估计值作为反向R‑T‑S平滑起始时刻的初始值，令k=N，进行N→0的R‑T‑S平滑过程，其中表示N时刻的平滑估计值，表示N时刻的卡尔曼滤波状态估计值，表示N时刻的平滑估计均方误差阵，表示k时刻的卡尔曼滤波估计均方误差阵，方程如下：平滑增益阵：其中，表示k时刻的平滑增益阵，Φ_k+1,k表示离散化的卡尔曼滤波状态转移矩阵，表示由k时刻到k+1时刻的一步预测均方误差阵；平滑的状态估计：${\hat{X}}_{k/N}^S={\hat{X}}_k^F+K_{k/N}^S({\hat{X}}_{k+1/N}^S-{\hat{X}}_{k+1,k}^F),$其中：表示第k时刻的平滑估计量，表示第k+1时刻的平滑估计量，表示由k时刻到k+1时刻的状态一步预测值；平滑的估计均方误差阵：$P_{k/N}^S=P_k^F+K_{k/N}^S(P_{k+1/N}^S-P_{k+1,k}^S){\left(K_{k/N}^S\right)}^T,$其中：表示k时刻的平滑估计均方误差阵，表示k+1时刻的平滑估计均方误差阵，表示由k时刻到k+1时刻的一步预测均方误差阵，表示k时刻的平滑增益阵；从公式可见，第k时刻的平滑估计量是在正向Kalman滤波估计量的基础上线性补偿了一个修正量得到的，该修正量是第k+1时刻的平滑估计量与正向一步预测估计量的差值；由上面的分析过程得到，k时刻的R‑T‑S平滑估计量是k时刻的正向Kalman滤波估计量和k+1时刻的平滑估计量的线性融合，k+1时刻的平滑估计量利用了全局数据。