CN114779358A - 一种航空大地电磁数据姿态校正方法 - Google Patents

Abstract

本发明提供了一种航空大地电磁数据姿态校正方法，包括以下步骤：A、获取测区地面基站的磁场水平分量数据以及空中接收线圈的磁场垂直分量数据，然后对时域信号进行傅里叶变换转换到频率域，计算得到若干个观测频率点的倾子实虚部数据；B、利用空中接收线圈的三个GPS天线测量数据，计算空中接收线圈的空中姿态角度误差；C、逐个观测频率点计算校正后的倾子实、虚部值；D、绘制各观测频率点校正后的倾子实、虚部等值线断面图，若校正后倾子等值线断面图中出现线性条带，则调整校正系数。本发明在一定程度上消除了航空大地电磁测量过程中由于空中线圈姿态变化所引起的电磁噪声，从而提高了测量数据的信噪比，避免了在反演过程中带来虚假异常。

Description

一种航空大地电磁数据姿态校正方法

技术领域

本发明涉及电磁勘探技术领域，具体涉及一种航空大地电磁数据姿态校正方法。

背景技术

航空大地电磁法通过在空中观测天然电磁场的磁场垂直分量，同时在地面基站观测磁场水平分量，通过观测的磁场垂直和磁场水平分量计算得到倾子矢量，进而获取地下介质的电性结构信息。该方法兼具了天然场勘探深度大以及航空物探测量效率高、成本低、绿色环保、对地形影响小等优点，能够克服复杂地形条件下传统地面勘探手段无法有效开展且成本高的问题，是一种新兴的航空地球物理勘查技术手段。目前，已经在铜、镍、金、铀等多金属矿勘探和工程勘察中取得了较好的应用效果。

航空大地电磁数据采集过程中，空中接收线圈一般保持水平姿态。但是由于接收线圈与飞行平台之间通过软连接固定，在测量飞行过程中受飞行姿态以及风速等的影响，空中接收线圈的姿态会产生一定变化，使得测量的磁场垂直分量Hz的大小会与真实磁场值大小有偏差，导致在实测数据中产生一些电磁噪声。因此，需要对该姿态变化引起的误差进行校正，提高测量数据质量以及反演解释的准确性。

发明内容

本发明的目的就是提供一种航空大地电磁数据姿态校正方法，以解决现有技术中实测数据因姿态变化产生的电磁噪音问题。

本发明是这样实现的：一种航空大地电磁数据姿态校正方法，包括以下步骤：

A、在工作区内开展航空大地电磁测量，获取测区地面基站的磁场水平分量数据Hx、Hy 以及空中接收线圈的磁场垂直分量数据Hz，通过GPS时间同步将水平和磁场垂直分量数据合并，然后对时域信号进行傅里叶变换转换到频率域，计算得到若干个观测频率点的倾子实虚部数据；

B、利用空中接收线圈的三个GPS天线测量数据，计算空中接收线圈的空中姿态角度误差ex、ey；

C、利用步骤a中计算的观测频率点倾子实虚部数据以及步骤b计算得到的空中接收线圈姿态角度误差ex、ey，给定一初始校正系数，采用姿态校正公式，逐个观测频率点计算校正后的倾子实、虚部值；

D、绘制各观测频率点校正后的倾子实、虚部等值线断面图，查看姿态校正效果；若校正后倾子等值线断面图中出现线性条带，则调整校正系数，按步骤c重新计算校正值，直至线性条带消失，保存校正系数；

E、根据各观测频率点确定的校正系数，编写Geosoft Correction.gs插件，将步骤a计算得到的倾子实、虚部数据导入Geosoft软件中生成数据库，调用GeosoftCorrection.gs插件，生成姿态校正后的航空大地电磁数据，用于后续倾子资料的反演解释。

进一步地，本发明可以按如下技术方案实现：

在步骤a中，通过GPS时间同步，对采集的航空大地电磁磁场垂直和磁场水平分量数据进行同步合并，然后对时域信号进行傅里叶变换转换到频率域。

选取25Hz、37Hz、75Hz、150Hz、300Hz、600Hz 5个最佳观测频率点，根据倾子的计算表达式，计算得到各个观测频率点的倾子实、虚部数据。

在步骤b中，根据空中接收线圈上三个GPS接收机的位置坐标计算以及坐标转换，得到空中接收线圈上这三个点的空间几何关系，计算出空中接收线圈的空间姿态。

选取空中接收线圈上三台GPS接收机GPS1、GPS2和GPS3的位置A、B、C三点，设定A、B、C三点在空中接收线圈本体坐标系中的坐标分别为：A(x_A,y_A,z_A)、B(x_B,y_B,z_B)和 C(x_C,y_C,z_C)，A、B、C三点测量参数的WGS-84坐标系的坐标分别为：

和

其中，x、y、z代表地心坐标系下的直角坐标，B、L、H分别为经纬度和高程，w表示在WGS-84 坐标系下；通过用如下方法求解三维姿态角偏航角，俯仰角和滚动角：

计算转换矩阵F_m＝(L_m)(B_m)^T，其中

计算坐标差分：

令

则有：

计算姿态变换矩阵T_m＝(γ_AB)_m(α_AB)_m，其中，

令

则得到滚动角ψ_BC为：

令

则姿态角度误差为：

在步骤c中，对计算的倾子实、虚部值以及空中接收线圈姿态角度误差ex、ey，采用姿态校正公式，逐个计算各个观测频率点姿态校正后的倾子实、虚部数据。

姿态校正公式为：

Xipc＝a_x*X_ip-e_x+c_x*Y_ip

Xqdc＝Xqd+b_x*X_ip+d_x*Y_ip

Yipc＝a_y*Y_ip-e_y+c_y*X_ip

Xqdc＝Yqd+b_y*Y_ip+d_y*X_ip

式中，a_x、a_y、b_x、b_y、c_x、c_y、d_x、d_y均为校正系数；其中a_x、a_y、b_x、b_y的取值范围在0.6～1.4之间，c_x、c_y、d_x、d_y的取值范围在-0.2～0.2之间；姿态校正过程中，首先给定一初始的校正系数(a_x＝b_x＝a_y＝b_y＝1，c_x＝d_x＝c_y＝d_y＝0.1)，然后逐个计算观测频率点校正后的倾子实、虚部值。

在步骤d中，要求校正后的倾子实、虚部数据的等值线断面图中无线性条带。

本发明在一定程度上消除了航空大地电磁测量过程中由于空中线圈姿态变化所引起的电磁噪声，从而提高了测量数据的信噪比，避免了在反演过程中带来虚假异常。

附图说明

图1是本发明的航空大地电磁系统数据采集示意图。

图2为本发明姿态校正的流程图。

图3为本发明中航空大地电磁姿态校正后150Hz倾子等值线断面适当矫正图。

图4为本发明中航空大地电磁姿态校正后150Hz倾子等值线断面过度矫正图。

具体实施方式

下面结合附图对本发明的具体实施方式作列举性说明。

如图1和图2所示，本发明的航空大地电磁数据姿态校正方法，其特征在于，包括以下步骤：

A、在工作区内开展航空大地电磁测量，获取测区地面基站的磁场水平分量数据Hx、Hy 以及空中接收线圈的磁场垂直分量数据Hz，磁场水平分量Hx和Hy相互垂直，并通过数据处理获得地下介质的电阻率分布信息，地下介质的电阻率分布信息是通过倾子反演计算得到的，在反演之前需对测量的倾子数据进行姿态校正，本发明没有涉及反演部分的内容。通过 GPS时间同步，对采集的航空大地电磁磁场垂直和磁场水平分量数据进行同步合并，然后对时域信号进行傅里叶变换转换到频率域，计算得到若干个观测频率点的倾子实虚部数据。

在具体操作过程中，通过倾子计算得到25、37、75、150、300、600Hz 5个最佳观测频率点的倾子实虚部数据，选取这5个频率对应的趋肤深度大致呈均匀等间距分布，有利于后期数据成图及反演计算。根据倾子的计算表达式，计算得到各个观测频率点的倾子实、虚部数据。

根据麦克斯韦方程组：

得到在三维响应情况下，航空大地电磁倾子矢量的计算公式为：

通过在地面基站设置两个在水平面上呈90度的接收线圈，测量磁场水平分量Hx、Hy以及空中接收线圈测量的磁场垂直分量Hz，计算得到沿测线和垂直测线方向的倾子响应。该过程在航空大地电磁测量结束后，进行倾子计算，然后对计算的倾子数据进行姿态校正。从公式中可以看出，磁场的垂直分量与磁场水平分量之间存在着复系数关系，该关系通过倾子矢量T表示。Tzx和Tzy是沿测线和垂直测线方向的倾子值，校正公式中Xip和Xqd分别为Tzx 的实部和虚部，同理，Yip和Yqd为Tzy的实部和虚部。

其中，H_z为实测过程中，空中接收线圈沿测线飞行时测量的磁场垂直分量数据，H_x、H_y为地面固定基站测量的磁场水平分量数据，T_zx和T_zy分别为倾子沿测线和垂直测线方向的值。通过傅里叶变换，将时域的观测数据转换到频率域，在频率域下对倾子计算公式的两端分别乘以H_x(ω)、H_y(ω)的共轭复数，联立求解方程组，得到倾子的计算表达式为：

其中，时域信号是指Hx、Hy、Hz，频域信号是将时域信号做傅里叶变换转换到频率域，描述了频率变化和幅度变化的关系，然后选取了几个频率计算得到对应频率的Tzx和Tzy。

B、利用空中接收线圈的三个GPS天线测量数据，计算空中接收线圈的空中姿态角度误差ex、ey；即根据空中接收线圈上三个GPS接收机的位置坐标计算以及坐标转换，得到空中接收线圈上这三个点的空间几何关系，计算出空中接收线圈的空间姿态。

在具体操作过程中，空中接收线圈的三个GPS天线也就是三个GPS接收机，GPS1位于空中接收线圈的前方，GPS2和GPS的位置无要求，选取空中接收线圈上三台GPS接收机GPS1、GPS2和GPS3的位置A、B、C三点，设定A、B、C三点在空中接收线圈本体坐标系中的坐标分别为：A(x_A,y_A,z_A)、B(x_B,y_B,z_B)和C(x_C,y_C,z_C)，A、B、C三点测量参数 WGS-84坐标系的坐标分别为：

和

其中，x、y、z代表地心坐标系下的直角坐标，B、L、H分别为经纬度和高程，w表示在WGS-84 坐标系下；通过用如下方法求解三维姿态角偏航角，俯仰角和滚动角。

计算转换矩阵F_m＝(L_m)(B_m)^T，其中

计算坐标差分：

令

则有：

计算姿态变换矩阵T_m＝(γ_AB)_m(α_AB)_m，其中，

令

则得到滚动角ψ_BC为：

令

则姿态角度误差为：

C、利用步骤a中计算的观测频率点倾子实虚部数据以及步骤b计算得到的空中接收线圈姿态角度误差ex、ey，给定一初始校正系数，采用姿态校正公式，逐个观测频率点计算校正后的倾子实、虚部值。

对计算的倾子实、虚部值以及空中接收线圈姿态角度误差ex、ey，采用姿态校正公式，逐个计算各个观测频率点姿态校正后的倾子实、虚部数据。

姿态校正公式为：

Xipc＝a_x*X_ip-e_x+c_x*Y_ip

Xqdc＝Xqd+b_x*X_ip+d_x*Y_ip

Yipc＝a_y*Y_ip-e_y+c_y*X_ip

Xqdc＝Yqd+b_y*Y_ip+d_y*X_ip

式中，a_x、a_y、b_x、b_y、c_x、c_y、d_x、d_y均为校正系数；其中a_x、a_y、b_x、b_y的取值范围在0.6～1.4之间，c_x、c_y、d_x、d_y的取值范围在-0.2～0.2之间；姿态校正过程中，首先给定一初始的校正系数(a_x＝b_x＝a_y＝b_y＝1，c_x＝d_x＝c_y＝d_y＝0.1)，(后续会根据姿态校正效果对校正系数进行调整)，然后逐个计算观测频率点校正后的倾子实、虚部值。

姿态校正公式中Xip和Xqd分别为Tzx的实部和虚部，同理，Yip和Yqd为Tzy的实部和虚部。

D、绘制并对比各观测频率点校正后的倾子实、虚部等值线断面图，查看姿态校正效果；若校正后倾子等值线断面图中出现线性条带(见图4)，则调整校正系数，按步骤c重新计算倾子校正值，直至线性条带消失(见图3)，保存校正系数。

E、根据倾子实、虚部等值线断面图的形态确定各个频点的校正系数，编写GeosoftCorrection.gs插件，将步骤a计算得到的倾子实、虚部数据导入Geosoft软件中生成数据库，调用Geosoft Correction.gs插件，生成姿态校正后的航空大地电磁数据(生成校正后的航空大地电磁倾子实、虚部数据，)，用于后续倾子资料的反演解释。

具体地，根据各个频点确定的姿态校正系数，调用Geosoft软件中的math.gx，将每个频点的倾子实、虚部校正公式编写成一个gx文件，将所有频点对应的gx文件合并在一起，编制航空大地电磁数据姿态校正Geosoft Correction.gs插件。

Claims (8)

1.一种航空大地电磁数据姿态校正方法，其特征在于，包括以下步骤：

A、在工作区内开展航空大地电磁测量，获取测区地面基站的磁场水平分量数据Hx、Hy以及空中接收线圈的磁场垂直分量数据Hz，通过GPS时间同步将水平和磁场垂直分量数据合并，然后对时域信号进行傅里叶变换转换到频率域，计算得到若干个观测频率点的倾子实虚部数据；

B、利用空中接收线圈的三个GPS天线测量数据，计算空中接收线圈的空中姿态角度误差ex、ey；

C、利用步骤a中计算的观测频率点倾子实虚部数据以及步骤b计算得到的空中接收线圈姿态角度误差ex、ey，给定一初始校正系数，采用姿态校正公式，逐个观测频率点计算校正后的倾子实、虚部值；

D、绘制各观测频率点校正后的倾子实、虚部等值线断面图，查看姿态校正效果；若校正后倾子等值线断面图中出现线性条带，则调整校正系数，按步骤c重新计算校正值，直至线性条带消失，保存校正系数；

E、根据各观测频率点确定的校正系数，编写Geosoft Correction.gs插件，将步骤a计算得到的倾子实、虚部数据导入Geosoft软件中生成数据库，调用Geosoft Correction.gs插件，生成姿态校正后的航空大地电磁数据，用于后续倾子资料的反演解释。

2.根据权利要求1所述的一种航空大地电磁数据姿态校正方法，其特征在于，在步骤a中，通过GPS时间同步，对采集的航空大地电磁磁场垂直和磁场水平分量数据进行同步合并，然后对时域信号进行傅里叶变换转换到频率域。

3.根据权利要求2所述的一种航空大地电磁数据姿态校正方法，其特征在于，选取25Hz、37Hz、75Hz、150Hz、300Hz、600Hz 5个最佳观测频率点，根据倾子的计算表达式，计算得到各个观测频率点的倾子实、虚部数据。

4.根据权利要求1所述的一种航空大地电磁数据姿态校正方法，其特征在于，在步骤b中，根据空中接收线圈上三个GPS接收机的位置坐标计算以及坐标转换，得到空中接收线圈上这三个点的空间几何关系，计算出空中接收线圈的空间姿态。

5.根据权利要求4所述的一种航空大地电磁数据姿态校正方法，其特征在于，选取空中接收线圈上三台GPS接收机GPS1、GPS2和GPS3的位置A、B、C三点，设定A、B、C三点在空中接收线圈本体坐标系中的坐标分别为：A(x_A，y_A，z_A)、B(x_B，y_B，z_B)和C(x_C，y_C，z_C)，A、B、C三点测量参数的WGS-84坐标系的坐标分别为：

和

其中，x、y、z代表地心坐标系下的直角坐标，B、L、H分别为经纬度和高程，w表示在WGS-84坐标系下；通过用如下方法求解三维姿态角偏航角，俯仰角和滚动角：

计算转换矩阵F_m＝(L_m)(B_m)^T，其中

计算坐标差分：

令

则有：

计算姿态变换矩阵T_m＝(γ_AB)_m(α_AB)_m，其中，

令

则得到滚动角ψ_BC为：

令

则姿态角度误差为：

6.根据权利要求1所述的一种航空大地电磁数据姿态校正方法，其特征在于，在步骤c中，对计算的倾子实、虚部值以及空中接收线圈姿态角度误差ex、ey，采用姿态校正公式，逐个计算各个观测频率点姿态校正后的倾子实、虚部数据。

7.根据权利要求6所述的一种航空大地电磁数据姿态校正方法，其特征在于，姿态校正公式为：

Xipc＝a_x*X_ip-e_x+c_x*Y_ip

Xqdc＝Xqd+b_x*X_ip+d_x*Y_ip

Yipc＝a_y*Y_ip-e_y+c_y*X_ip

Xqdc＝Yqd+b_y*Y_ip+d_y*X_ip

式中，a_x、a_y、b_x、b_y、c_x、c_y、d_x、d_y均为校正系数；其中a_x、a_y、b_x、b_y的取值范围在0.6～1.4之间，c_x、c_y、d_x、d_y的取值范围在-0.2～0.2之间；姿态校正过程中，首先给定一初始的校正系数(a_x＝b_x＝a_y＝b_y＝1，c_x＝d_x＝c_y＝d_y＝0.1)，然后逐个计算观测频率点校正后的倾子实、虚部值。

8.根据权利要求1所述的一种航空大地电磁数据姿态校正方法，其特征在于，在步骤d中，要求校正后的倾子实、虚部数据的等值线断面图中无线性条带。

Images