CN109001841B - 一种基于地球自转角速度的重力梯度仪标定方法 - Google Patents

Abstract

本发明公开一种基于地球自转角速度的重力梯度仪标定方法，利用重力梯度仪在不同姿态下，地球自转产生的离心梯度激励及环境物体产生的自梯度激励不同，标定重力梯度仪的标度系数及自梯度模型参数。标定时，使重力梯度仪经历不同的姿态角，并记录重力梯度仪及离心梯度检测单元在不同姿态下的输出，即可标定重力梯度仪的标度系数及自梯度模型参数。该标定方法，可由计算机自动运行，实现自标定，该方法仅依靠重力梯度仪本身，不需要依赖外部检测质量即可完成标定，该标定方法，操作简单易实施。

Description

一种基于地球自转角速度的重力梯度仪标定方法

技术领域

本发明涉及一种基于地球自转角速度的重力梯度仪标定方法，属于精密测量技术领域。

背景技术

重力梯度仪是一种高精度万有引力梯度检测仪，能够测量仪器所在位置点的万有引力梯度值。航空重力梯度勘探能高精度、高效率的收集重力梯度信息，测量的重力梯度数据可直接用于分析勘探区域的地质密度、地质结构、重建重力场等。目前国外已经投入商业应用的重力梯度仪主要为旋转加速度计重力梯度仪，在研的具有应用前景的重力梯度仪有冷原子重力梯度仪、超导重力梯度仪、MEMS重力梯度仪等，我国的重力梯度仪样机正处于研制中。

从现有的文献来看，重力梯度仪的标度系数的标定方法主要有两类，基于万有引力梯度的标定方法和基于离心梯度的标定方法。基于万有梯度的标定方法，它通过施加外部检测质量利用万有引力场来激励重力梯度仪，标定重力梯度仪。这类标定方法，需要外部检测质量，且需要人为改变检测质量体的质量和位置，标定过程复杂，标定效率低。基于离心梯度的标定方法，它通过旋转重力梯度仪施加离心梯度，标定重力梯度仪的标度系数，它的优点是不需要施加外部检测质量。重力梯度仪对角运动、线运动敏感，旋转法标定重力梯度仪标度系数，如果转速过快，会给重力梯度仪引入测量误差，转速过慢会降低标定效率。现有的某些文献，在标定时忽略了自梯度对标度系数标定、零偏标定的影响，不符合实际情况。

发明内容

技术问题：本发明提供一种在标定时不需要旋转重力梯度仪，在静止状态进行标定，能同时标定重力梯度仪的标度系数和自梯度模型参数的基于地球自转角速度的重力梯度仪标定方法。

技术方案：本发明的基于地球自转角速度的重力梯度仪标定方法，包括以下步骤：

1)改变重力梯度仪的姿态，保持重力梯度仪静止，记录重力梯度仪的姿态及该姿态下重力梯度仪的输出及离心梯度检测单元的输出，重复M次上述操作，总共记录M组不同的姿态，及相应的重力梯度仪输出、离心梯度检测单元的输出，记录的第k个姿态表示为(θ_x(k),θ_y(k),θ_z(k))，该姿态下重力梯度仪的inline通道及cross通道的输出分别为F_in(k)，F_cs(k)，该姿态下离心梯度检测单元的inline通道输出为

cross通道输出为

2)依次将M个姿态依次代入下式，计算M个姿态对应的M组姿态特征参数，重力梯度仪的每个姿态对应1组36个姿态特征参数，这36个姿态特征参数为a_1,1,…,a_1,18,b_1,1,…,b_1,18：

式中c,c²分别表示cos(),cos²()，s,s²分别表示sin(),sin²()，θ_x,θ_y,θ_z为姿态角，第k个姿态(θ_x(k),θ_y(k),θ_z(k))，计算得到的第k组姿态特征参数为a_1,1(k),…,a_1,18(k),b_1,1(k),…,b_1,18(k)；

3)将M个姿态对应的M组姿态特征参数及该姿态下重力梯度仪的输出、离心梯度检测单元的输出代入下式，计算重力梯度仪标度系数及自梯度模型参数的标定值：

上式中p₁,…,p₂₀是标定的20个自梯度模型参数，k_in，k_cs分别是标定的重力梯度仪的 inline通道及cross通道的标度系数，[·]⁺表示最小二乘逆，a_1,1(1),…,a_1,18(1)，b_1,1(1),…,b_1,18(1)是第1个姿态(θ_x(1),θ_y(1),θ_z(1))计算得到的第1组姿态特征参数；a_1,1(k),…,a_1,18(k)， b_1,1(k),…,b_1,18(k)是第k个姿态(θ_x(k),θ_y(k),θ_z(k))计算得到的第k组姿态特征参数； a_1,1(M),…,a_1,18(M)，b_1,1(M),…,b_1,18(M)是第M个姿态(θ_x(M),θ_y(M),θ_z(M))计算得到的第M 组姿态特征参数；F_in(1)是第1个姿态(θ_x(1),θ_y(1),θ_z(1))对应的重力梯度仪的inline通道的输出，F_cs(1)是第1个姿态(θ_x(1),θ_y(1),θ_z(1))对应的重力梯度仪的cross通道的输出；F_in(k) 是第k个姿态(θ_x(k),θ_y(k),θ_z(k))对应的重力梯度仪的inline通道的输出，F_cs(k)是第k个姿态(θ_x(k),θ_y(k),θ_z(k))对应的重力梯度仪的cross通道的输出；F_in(M)是第M个姿态 (θ_x(M),θ_y(M),θ_z(M))对应的重力梯度仪的inline通道的输出，F_cs(M)是第M个姿态 (θ_x(M),θ_y(M),θ_z(M))对应的重力梯度仪的cross通道的输出；是第1个姿态 (θ_x(1),θ_y(1),θ_z(1))对应的离心梯度检测单元的inline通道的输出，

是第1个姿态 (θ_x(1),θ_y(1),θ_z(1))对应的离心梯度检测单元的cross通道的输出；

是第k个姿态(θ_x(k),θ_y(k),θ_z(k))对应的离心梯度检测单元的inline通道的输出，

是第k个姿态(θ_x(k),θ_y(k),θ_z(k))对应的离心梯度检测单元的cross通道的输出；

是第M个姿态(θ_x(M),θ_y(M),θ_z(M))对应的离心梯度检测单元的inline通道的输出，

是第M个姿态(θ_x(M),θ_y(M),θ_z(M))对应的离心梯度检测单元的cross通道的输出。

进一步的，本发明方法中，所述步骤1)中按照以下改变重力梯度仪的姿态：以当地地理坐标系为参考坐标系，当地地理坐标系绕z轴旋转θ_z，再绕x轴旋转θ_x，再y轴旋转θ_y和重力梯度仪测量坐标系重合，姿态(θ_x,θ_y,θ_z)符号的正负服从右手法则。

进一步的，本发明方法中，所述步骤1)中M的取值要使以下矩阵A的秩等于19：

进一步的，本发明方法中，所述步骤1)中的离心梯度检测单元是一种能够实时测量重力梯度仪的离心梯度的装置。

进一步的，本发明方法中，步骤1)中重力梯度仪的inline通道是指测量重力梯度张量的inline分量Γ_xx-Γ_yy的通道，重力梯度仪的cross通道是指测量重力梯度张量的cross分量Γ_xy的通道，离心梯度检测单元的inline通道是指测量离心梯度张量的inline分量

的通道，离心梯度检测单元的cross通道是指测量离心梯度张量的cross分量-ω_xω_y的通道。

进一步的，本发明方法中使用的重力梯度仪标定模型为：

式中是重力梯度仪的inline通道的测量值，

是重力梯度仪的cross通道的测量值； k_in是重力梯度仪的inline通道的标度系数，k_cs是重力梯度仪的cross通道的标度系数；F_in是重力梯度仪inline通道的原始输出，F_cs是重力梯度仪cross通道的原始输出，

是离心梯度检测单元的inline通道的输出，

是离心梯度检测单元的cross通道的输出；

是重力梯度仪自梯度模型的inline通道的输出，

是重力梯度仪自梯度模型的cross通道的输出，

是自梯度及电路零偏叠加。

由下式进行计算：

p₁,…,p₂₀是自梯度模型参数，它和重力梯度仪的标度系数k_in，k_cs构成22个待标定的参数；a_1,1,…,a_1,18,b_1,1,…,b_1,18是自梯度模型的36个姿态特征参数，它和姿态(θ_x,θ_y,θ_z)的关系如下：

式中c，c²分别表示cos()，cos²()；s，s²分别表示sin()，sin²()；θ_x,θ_y,θ_z为姿态角。

重力梯度仪是一种测量目标物体在空间某点产生的万有引力梯度的仪器。重力梯度仪自身的一些部件比如稳定平台组件、伺服电机、载体等能够产生万有引力梯度干扰重力梯度仪对目标物体产生的万有引力梯度的测量。这部分由重力梯度仪自身部件产生的万有引力梯度称为自梯度。自梯度和重力梯度仪电路因素造成的常值零偏耦合，构成梯度仪的等效零偏。自梯度随载体姿态变化，是姿态的函数。本发明建立了重力梯度仪的等效零偏模型，等效零偏的主要成分是自梯度，称它为自梯度模型。重力梯度仪的标定包含了对重力梯度仪标度系数的标定及自梯度模型参数的标定。地球自转产生的离心梯度可以达到5.28Eo，能够满足现有精度重力梯度仪的标定需求。

本发明中，重力梯度仪及离心梯度检测单元的inline通道、cross通道的定义如下：重力梯度张量有6个分量，它们是Γ_xx，Γ_yy，Γ_zz，Γ_xy，Γ_xz，Γ_yz，其中Γ_xx，Γ_yy，Γ_zz称为inline分量，Γ_xy，Γ_xz，Γ_yz称为cross分量，重力梯度仪的inline通道测量的是inline 分量Γ_xx-Γ_yy，重力梯度仪的cross通道测量的是cross分量Γ_xy。离心梯度张量同样有6个分量，它们是

-ω_xω_y，-ω_xω_z，-ω_yω_z，其中，

是inline分量，-ω_xω_y，-ω_xω_z，-ω_yω_z是cross分量，离心梯度检测单元的inline通道测量的是inline分量

离心梯度检测单元的cross通道测量的是cross分量-ω_xω_y。

有益效果：本发明与现有技术相比，具有以下优点：

本发明利用地球自转产生的离心梯度激励重力梯度仪，通过改变重力梯度仪的姿态角改变重力梯度仪敏感的离心梯度。标定时，改变重力梯度仪的姿态，保持重力梯度仪静止，记录不同姿态下，重力梯度仪的输出及离心梯度检测单元的输出。根据记录的重力梯度仪姿态，重力梯度仪的输出，离心梯度检测单元的输出可以同时标定重力梯度仪的标度系数，自梯度模型的参数，该标定方法使重力梯度仪的标定更简单、快捷。该标定方法在静止状态下标定重力梯度仪，不会引入动态测量误差，标定精度高。

附图说明

图1重力梯度仪系统。

图2θ_x＝23°,θ_y＝15°，θ_z＝0°：285°标定的自梯度模型的输出与真实自梯度的对比。

具体实施方式

下面结合实施例和说明书附图对本发明作进一步地说明。

如图1所示，重力梯度仪安装在三轴稳定平台并置于载体中构成动基座重力梯度仪系统，重力梯度仪系统产生自梯度的组件主要有转台的内框，安装托板，中框，外框，轴承，电机，转台基座、载体等。重力梯度仪的标定模型为：

式中

是重力梯度仪的inline通道的测量值，

是重力梯度仪的cross通道的测量值； k_in是重力梯度仪的inline通道的标度系数，k_cs是重力梯度仪的cross通道的标度系数；F_in是重力梯度仪inline通道的原始输出，F_cs是重力梯度仪cross通道的原始输出，是离心梯度检测单元的inline通道的输出，是离心梯度检测单元的cross通道的输出；

是重力梯度仪自梯度模型的inline通道的输出，

是重力梯度仪自梯度模型的cross通道的输出，

是自梯度及电路零偏叠加。

由下式进行计算：

p₁,…,p₂₀是自梯度模型参数，它和重力梯度仪的标度系数k_in，k_cs构成22个待标定的参数；a_1,1,…,a_1,18,b_1,1,…,b_1,18是自梯度模型的36个姿态特征参数，它和姿态(θ_x,θ_y,θ_z)的关系如下：

式中c，c²分别表示cos()，cos²()；s，s²分别表示sin()，sin²()；θ_x,θ_y,θ_z为姿态角。重力梯度仪的姿态是重力梯度仪测量坐标系与当地地理坐标系之间的欧拉角，具体为当地地理坐标系绕z轴旋转θ_z，再绕x轴旋转θ_x，再y轴旋转θ_y和重力梯度仪测量坐标系重合，姿态(θ_x,θ_y,θ_z)的符号的正负服从右手法则。

本发明方法包括以下步骤：

1)改变重力梯度仪的姿态，保持重力梯度仪静止，记录重力梯度仪的姿态及该姿态下重力梯度仪的输出及离心梯度检测单元的输出，重复M次上述操作，总共记录M组不同的姿态，及相应的重力梯度仪输出、离心梯度检测单元的输出，记录的第k个姿态表示为(θ_x(k),θ_y(k),θ_z(k))，该姿态下重力梯度仪的inline通道及cross通道的输出分别为F_in(k)，F_cs(k)，该姿态下离心梯度检测单元的inline通道及cross通道输出分别为

2)依次将M个姿态依次代入下式，计算M个姿态对应的M组姿态特征参数，重力梯度仪的每个姿态对应1组36个姿态特征参数，这36个姿态特征参数为a_1,1,…,a_1,18,b_1,1,…,b_1,18：

式中c,c²分别表示cos(),cos²()，s,s²分别表示sin(),sin²()，θ_x,θ_y,θ_z为姿态角。第k个姿态(θ_x(k),θ_y(k),θ_z(k))，计算得到的第k组姿态特征参数为a_1,1(k),…,a_1,18(k),b_1,1(k),…,b_1,18(k)；步骤1)中M的取值要使以下矩阵A的秩等于19：

3)将M个姿态对应的M组姿态特征参数及该姿态下重力梯度仪的输出、离心梯度检测单元的输出代入下式，计算重力梯度仪标度系数及自梯度模型参数的标定值：

上式中p₁,…,p₂₀是标定的20个自梯度模型参数，k_in，k_cs分别是标定的重力梯度仪的inline 通道及cross通道的标度系数，[·]⁺表示最小二乘逆，a_1,1(1),…,a_1,18(1)，b_1,1(1),…,b_1,18(1)是第1 个姿态(θ_x(1),θ_y(1),θ_z(1))计算得到的第1组姿态特征参数；a_1,1(k),…,a_1,18(k)，b_1,1(k),…,b_1,18(k) 是第k个姿态(θ_x(k),θ_y(k),θ_z(k))计算得到的第k组姿态特征参数；a_1,1(M),…,a_1,18(M)， b_1,1(M),…,b_1,18(M)是第M个姿态(θ_x(M),θ_y(M),θ_z(M))计算得到的第M组姿态特征参数； F_in(1)是第1个姿态(θ_x(1),θ_y(1),θ_z(1))对应的重力梯度仪的inline通道的输出，F_cs(1)是第1 个姿态(θ_x(1),θ_y(1),θ_z(1))对应的重力梯度仪的cross通道的输出；F_in(k)是第k个姿态 (θ_x(k),θ_y(k),θ_z(k))对应的重力梯度仪的inline通道的输出，F_cs(k)是第k个姿态 (θ_x(k),θ_y(k),θ_z(k))对应的重力梯度仪的cross通道的输出；F_in(M)是第M个姿态 (θ_x(M),θ_y(M),θ_z(M))对应的重力梯度仪的inline通道的输出，F_cs(M)是第M个姿态 (θ_x(M),θ_y(M),θ_z(M))对应的重力梯度仪的cross通道的输出；

是第1个姿态 (θ_x(1),θ_y(1),θ_z(1))对应的离心梯度检测单元的inline通道的输出，

是第1个姿态 (θ_x(1),θ_y(1),θ_z(1))对应的离心梯度检测单元的cross通道的输出；

是第k个姿态 (θ_x(k),θ_y(k),θ_z(k))对应的离心梯度检测单元的inline通道的输出，

是第k个姿态 (θ_x(k),θ_y(k),θ_z(k))对应的离心梯度检测单元的cross通道的输出；

是第M个姿态 (θ_x(M),θ_y(M),θ_z(M))对应的离心梯度检测单元的inline通道的输出，

是第M个姿态(θ_x(M),θ_y(M),θ_z(M))对应的离心梯度检测单元的cross通道的输出。为了对本发明实施例中的方案进行验证，可以进行如下仿真分析：

用四个不同参数的点质量模拟四类不同的环境物体产生自梯度，这四类环境的仿真参数如下表：

在仿真实验中，重力梯度仪的标度系数的理论值k_in＝1.7,k_cs＝0.85，重力梯度仪的所在位置点的纬度等于35°N。

1)标定时，使重力梯度仪经历27个不同的姿态，记录重力梯度仪的姿态，以及对应姿态下重力梯度仪的输出、离心梯度检测单元的输出。记录的实验数据列出，如下表：

表中(θ_x(k),θ_y(k),θ_z(k))表示记录的第k个姿态；F_in(k)，F_cs(k)分别表示该姿态下重力梯度仪的inline通道，cross通道的输出；

分别表示该姿态下离心梯度检测单元的inline通道，cross通道的输出。

2)计算27个姿态对应的27组姿态特征参数。重力梯度仪的每个姿态对应1组36个姿态特征参数，这36个姿态特征参数为a_1,1,…,a_1,18,b_1,1,…,b_1,18；依次将27个姿态代入下式，计算27组姿态特征参数：

将27个姿态对应的27组姿态特征参数及该姿态下重力梯度仪的输出、离心梯度检测单元的输出代入下式，计算重力梯度仪标度系数及自梯度模型参数的标定值：

标定的标度系数为k_in＝1.70，k_cs＝0.85，标定的环境物体特征参数p₁,…,p₂₀如下表：

从标定结果可以看出，重力梯度仪标度系数的标定值和设置的理想值完全一致。为了验证标定的自梯度模型参数是否正确，使重力梯度仪经历一系列新的姿态，将姿态数据输入估计的自梯度模型，可以得到估计的自梯度模型的输出；同时根据万有引力梯度定理计算出环境物体质量产生的自梯度真值，对比估计的自梯度模型输出的自梯度和真实的自梯度，可以判断估计的自梯度模型是否正确。使重力梯度仪的姿态角θ_x＝23°,θ_y＝15°，并使重力梯度仪的姿态角θ_z＝0°步增到θ_z＝285°，步距为15°，图2是估计的自梯度模型输出的自梯度和基于万有引力定理计算的真实自梯度的对比图；从图2中可以看出，估计的自梯度模型输出的自梯度和真实自梯度完全一致。实验结果表明，本发明标定方法可以标定重力梯度仪的标度系数和自梯度模型参数。

Claims (5)

1.一种基于地球自转角速度的重力梯度仪标定方法，其特征在于，该方法包括以下步骤：

1)改变重力梯度仪的姿态，保持重力梯度仪静止，记录重力梯度仪的姿态及该姿态下重力梯度仪的输出及离心梯度检测单元的输出，重复M次上述操作，总共记录M组不同的姿态，及相应的重力梯度仪输出、离心梯度检测单元的输出，记录的第k个姿态表示为(θ_x(k),θ_y(k),θ_z(k))，该姿态下重力梯度仪的inline通道及cross通道的输出分别为F_in(k)，F_cs(k)，该姿态下离心梯度检测单元的inline通道输出为

cross通道输出为

2)依次将M个姿态依次代入下式，计算M个姿态对应的M组姿态特征参数，重力梯度仪的每个姿态对应1组36个姿态特征参数，这36个姿态特征参数为a_1,1,…,a_1,18,b_1,1,…,b_1,18：

式中c,c²分别表示cos(),cos²()，s,s²分别表示sin(),sin²()，θ_x,θ_y,θ_z为姿态角，第k个姿态(θ_x(k),θ_y(k),θ_z(k))，计算得到的第k组姿态特征参数为a_1,1(k),…,a_1,18(k),b_1,1(k),…,b_1,18(k)；

3)将M个姿态对应的M组姿态特征参数及该姿态下重力梯度仪的输出、离心梯度检测单元的输出代入下式，计算重力梯度仪标度系数及自梯度模型参数的标定值：

上式中p₁,…,p₂₀是标定的20个自梯度模型参数，k_in，k_cs分别是标定的重力梯度仪的inline通道及cross通道的标度系数，[·]⁺表示最小二乘逆，a_1,1(1),…,a_1,18(1)，b_1,1(1),…,b_1,18(1)是第1个姿态(θ_x(1),θ_y(1),θ_z(1))计算得到的第1组姿态特征参数；a_1,1(k),…,a_1,18(k)，b_1,1(k),…,b_1,18(k)是第k个姿态(θ_x(k),θ_y(k),θ_z(k))计算得到的第k组姿态特征参数；a_1,1(M),…,a_1,18(M)，b_1,1(M),…,b_1,18(M)是第M个姿态(θ_x(M),θ_y(M),θ_z(M))计算得到的第M组姿态特征参数；F_in(1)是第1个姿态(θ_x(1),θ_y(1),θ_z(1))对应的重力梯度仪的inline通道的输出，F_cs(1)是第1个姿态(θ_x(1),θ_y(1),θ_z(1))对应的重力梯度仪的cross通道的输出；F_in(k)是第k个姿态(θ_x(k),θ_y(k),θ_z(k))对应的重力梯度仪的inline通道的输出，F_cs(k)是第k个姿态(θ_x(k),θ_y(k),θ_z(k))对应的重力梯度仪的cross通道的输出；F_in(M)是第M个姿态(θ_x(M),θ_y(M),θ_z(M))对应的重力梯度仪的inline通道的输出，F_cs(M)是第M个姿态(θ_x(M),θ_y(M),θ_z(M))对应的重力梯度仪的cross通道的输出；

是第1个姿态(θ_x(1),θ_y(1),θ_z(1))对应的离心梯度检测单元的inline通道的输出，

是第1个姿态(θ_x(1),θ_y(1),θ_z(1))对应的离心梯度检测单元的cross通道的输出；

是第k个姿态(θ_x(k),θ_y(k),θ_z(k))对应的离心梯度检测单元的inline通道的输出，

是第k个姿态(θ_x(k),θ_y(k),θ_z(k))对应的离心梯度检测单元的cross通道的输出；

是第M个姿态(θ_x(M),θ_y(M),θ_z(M))对应的离心梯度检测单元的inline通道的输出，是第M个姿态(θ_x(M),θ_y(M),θ_z(M))对应的离心梯度检测单元的cross通道的输出。

2.根据权利要求1所述的基于地球自转角速度的重力梯度仪标定方法，其特征在于，所述步骤1)中按照以下改变重力梯度仪的姿态：以当地地理坐标系为参考坐标系，当地地理坐标系绕z轴旋转θ_z，再绕x轴旋转θ_x，再y轴旋转θ_y和重力梯度仪测量坐标系重合，姿态(θ_x,θ_y,θ_z)符号的正负服从右手法则。

3.根据权利要求1所述的基于地球自转角速度的重力梯度仪标定方法，其特征在于，所述步骤1)中M的取值要使以下矩阵A的秩等于19：

4.根据权利要求1、2或3所述的基于地球自转角速度的重力梯度仪标定方法，其特征在于，所述步骤1)中的离心梯度检测单元是一种能够实时测量重力梯度仪的离心梯度的装置。

5.根据权利要求1、2或3所述的基于地球自转角速度的重力梯度仪标定方法，其特征在于，所述步骤1)中重力梯度仪的inline通道是指测量重力梯度张量的inline分量Γ_xx-Γ_yy的通道，重力梯度仪的cross通道是指测量重力梯度张量的cross分量Γ_xy的通道，离心梯度检测单元的inline通道是指测量离心梯度张量的inline分量

的通道，离心梯度检测单元的cross通道是指测量离心梯度张量的cross分量-ω_xω_y的通道。

Images