一种适用于捷联式海洋重力仪的重力异常数据处理方法与流程

技术特征：

1.一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：包括如下步骤：

(1)、将捷联式海洋重力仪安装于船舱内，捷联式海洋重力仪实时记录船舶的速度增量和角速度增量，并投影至载体坐标系下得到速度增量f^b和角速度增量

(2)、进行捷联式海洋重力仪的动基座初始对准，获得从载体坐标系到实际数学平台坐标系的姿态转移矩阵所述姿态转移矩阵为对准时间段内最后时刻k时刻的姿态转移矩阵；

(3)、根据对准时间段内最后时刻k时刻的姿态转移矩阵获得船舶导航过程中当前时刻从载体坐标系到实际数学平台坐标系的姿态转移矩阵东速v_e(t_i)、北速v_n(t_i)、纬度lat(t_i)和经度lon(t_i)，其中t_i为当前时刻；

(4)、将GPS测量获得的船舶运动信息进行差分处理，获得差分GPS数据；

(5)、根据捷联惯性系统误差方程，选取状态向量，构建卡尔曼滤波器的系统状态方程，根据步骤(3)中当前时刻的东速v_e(t_i)、北速v_n(t_i)、纬度lat(t_i)和经度lon(t_i)，以及步骤(4)中差分GPS数据中的当前时刻的东速V_E，北速VN、经度λ和纬度L，得到相应东速的差值v_e(t_i)-V_E、北速的差值v_n(t_i)-V_N、经度的差值lon(t_i)-λ和纬度的差值lat(t_i)-L，作为卡尔曼滤波器的观测量进行东速误差、北速误差、经度误差、纬度误差和姿态误差的估计；

(6)、根据步骤(5)中得到的当前时刻的姿态误差修正东北天向比力值，得到修正后的东北天向比力值f^n'；

(7)、计算重力异常粗值δg，公式如下：

$\mathrm{\δ}g=g_b+(f^u-f_0^u)-a^u+{\mathrm{\δa}}_E+{\mathrm{\δa}}_F-{\mathrm{\γ}}_0+{\mathrm{\δg}}_{drift}$

其中：g_b为码头处的重力基准值；

f^u为f^n'中的天向比力值；

为码头处的天向比力初值；

a^u为天向运动加速度；

δa_E为厄特弗斯改正；

δa_F为自由空间改正；

γ₀为正常重力改正；

δg_drift为零点漂移改正。

2.根据权利要求1所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：获得船舶导航过程中当前时刻的姿态转移矩阵东速v_e(t_i)、北速v_n(t_i)、纬度lat(t_i)和经度lon(t_i)的具体方法如下：

根据导航过程中起始时刻，即第k+1时刻的速度增量f^b(t_k+1)、角速度增量和所述姿态转移矩阵获得导航过程中船舶在第k+1时刻的姿态转移矩阵东速v_e(t_k+1)、北速v_n(t_k+1)、纬度lat(t_k+1)和经度lon(t_k+1)，根据船舶在第k+1时刻的姿态转移矩阵东速v_e(t_k+1)、北速v_n(t_k+1)、纬度lat(t_k+1)和经度lon(t_k+1)，以及第k+2时刻的速度增量f^b(t_k+2)、角速度增量获得的船舶在第k+2时刻的姿态转移矩阵东速v_e(t_k+2)、北速v_n(t_k+2)、纬度lat(t_k+2)和经度lon(t_k+2)，依次类推，获得船舶在导航过程中当前时刻的姿态转移矩阵东速v_e(t_i)、北速v_n(t_i)、纬度lat(t_i)和经度lon(t_i)。

3.根据权利要求1所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：所述姿态转移矩阵通过如下方法获得：

其中：

$C_e^n=\left[\begin{array}{ccc}0& 1& 0\\ -\sin L& 0& \cos L\\ \cos L& 0& \sin L\end{array}\right],C_i^e\left(t_k\right)=\left[\begin{array}{ccc}{\mathrm{cos\ω}}_{ie}(t_k-t_0)& {\mathrm{sin\ω}}_{ie}(t_k-t_0)& 0\\ -{\mathrm{sin\ω}}_{ie}(t_k-t_0)& {\mathrm{cos\ω}}_{ie}(t_k-t_0)& 0\\ 0& 0& 1\end{array}\right],$

其中：为载体坐标系到实际数学平台坐标系的姿态转移矩阵，L代表大地纬度，ω_ie为地球自转角速度，t₀为对准时间段内的初始对准时刻，t_k为对准时间段内的任意时刻；

为载体坐标系到载体惯性坐标系的姿态转移矩阵，具体表达式为：

$C_b^{i_{b0}}\left(t_k\right)=\left[\begin{array}{ccc}{q}_0^2+q_1^2-q_2^2-q_3^2& 2(q_1{q}_2-q_0{q}_3)& 2(q_1{q}_3+q_0{q}_2)\\ 2(q_1{q}_2+q_0{q}_3)& q_0^2-q_1^2+q_2^2-q_3^2& 2(q_2{q}_3-q_0{q}_1)\\ 2(q_1{q}_3-q_0{q}_2)& 2(q_2{q}_3+q_0{q}_1)& q_0^2-q_1^2-q_2^2+q_3^2\end{array}\right];$

式中：q₀ q₁ q₂ q₃为对准数据段最后时刻k的四元素；

为载体惯性坐标系到经线地心惯性坐标系的姿态转移矩阵，具体表达式为：

$C_{i_{b0}}^i={\left[\begin{array}{c}{\left(V^i\right)}^T\left(t_{k1}\right)\\ {\left(V^i\right)}^T\left(t_{k2}\right)\\ \[\left(V^i\right)\left(t_{k1}\right)\×\left(V^i\right)\left(t_{k2}\right)\]\end{array}\right]}^{-1}\left[\begin{array}{c}{\left({\tilde{V}}^{i_{b0}}\right)}^T\left(t_{k1}\right)\\ {\left({\tilde{V}}^{i_{b0}}\right)}^T\left(t_{k2}\right)\\ {\[\left({\tilde{V}}^{i_{b0}}\right)\left(t_{k1}\right)\×\left({\tilde{V}}^{i_{b0}}\right)\left(t_{k2}\right)\]}^T\end{array}\right]$

其中：g为地球重力值，分别计算Vⁱ(t_k1)和Vⁱ(t_k2)的值，t_k1和t_k2分别是对准时段内的两个时刻；

Δt_k＝t_k-t₀，

为当前时刻从载体坐标系到实际数学平台坐标系的姿态转移矩阵；

f^b(t_i)为当前时刻的速度增量。

4.根据权利要求3所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：所述载体坐标系到载体惯性坐标系的姿态转移矩通过如下方法获得：

(3.1)、初始对准时刻t₀，载体坐标系到载体惯性坐标系的姿态转移矩阵表示如下：

$C_b^{i_{b0}}\left(t_0\right)=I,$

其中：I为3阶单位矩阵，其对应的初始时刻四元素为Q(t₀)＝[1 0 0 0]；

(3.2)、根据t₀时刻的四元素Q(t₀)和t₁时刻的角速度增量获得t₁时刻的四元素其中，Φ＝|Φ|；

(3.3)、根据t₁时刻的四元素Q(t₁)和t₂时刻的角速度增量获得t₂时刻的四元素Q(t₂)，依次类推，获得对准数据段最后时刻k时刻的四元素Q(k)＝[q₀ q₁ q₂ q₃]，根据Q(k)计算如下：

$C_b^{i_{b0}}\left(t_k\right)=\left[\begin{array}{ccc}{q}_0^2+q_1^2-q_2^2-q_3^2& 2(q_1{q}_2-q_0{q}_3)& 2(q_1{q}_3+q_0{q}_2)\\ 2(q_1{q}_2+q_0{q}_3)& q_0^2-q_1^2+q_2^2-q_3^2& 2(q_2{q}_3-q_0{q}_1)\\ 2(q_1{q}_3-q_0{q}_2)& 2(q_2{q}_3+q_0{q}_1)& q_0^2-q_1^2-q_2^2+q_3^2\end{array}\right].$

5.根据权利要求1～4之一所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：所述步骤(4)中差分GPS数据包括GPS时间、经度λ、纬度L、海拔高、大地高、东北天速度(V_E，V_N，V_U)、东北天加速度、卫星数、PDOP、HDOP、VDOP、质量数Q以及GPS周。

6.根据权利要求1～4之一所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：所述步骤(5)中根据捷联惯性系统误差方程，选取的状态向量X_INS为13阶，具体表示如下：

$X_{INS}={\left[\begin{array}{ccccccccccccc}\mathrm{\δ}L& \mathrm{\δ}\mathrm{\λ}& {\mathrm{\δv}}_e& {\mathrm{\δv}}_n& {\mathrm{\φ}}_e& {\mathrm{\φ}}_n& {\mathrm{\φ}}_u& {\mathrm{\ϵ}}_x& {\mathrm{\ϵ}}_y& {\mathrm{\ϵ}}_z& {\▿}_x& {\▿}_y& {\▿}_z\end{array}\right]}^T;$

其中：δL为纬度误差；

δλ为经度误差；

δv_e、δv_n分别为东速误差和北速误差；

φ_e、φ_n和φ_u分别为三个姿态误差角；

ε_x、ε_y和ε_z为激光陀螺的零位；

和为加速度计零位。

7.根据权利要求6所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：所述步骤(6)中修正后的东北天向比力值f^n'表示如下：

其中：φ×为反对称阵，

为当前时刻从载体坐标系到实际数学平台坐标系的姿态转移矩阵；

f^b(t_i)为当前时刻的速度增量；

ΔT为系统采样间隔时间。

8.根据权利要求1～4之一所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：对所述步骤(7)中计算的重力异常粗值δg采用数字滤波器进行滤波处理，以提高数据精度。

9.根据权利要求8所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：所述数字滤波采用FIR和IIR低通滤波器，截止频率小于0.01 Hz；或者采用正反卡尔曼滤波器。

10.根据权利要求1～4之一所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：若船舶在行驶过程中存在形成网格的交叉点，则进行测线网平差处理方法，消除仪器的系统误差和半系统误差，提高测量精度。

11.根据权利要求1～4之一所述的一种适用于捷联式海洋重力仪的重力异常数据处理方法，其特征在于：所述步骤(4)中采用PPP技术将GPS测量获得的船舶运动信息进行差分处理。