CN111856597B - 拖曳式海洋电磁地层电阻率与接收站位置联合反演方法 - Google Patents

Abstract

本发明提供一种拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，包括：读取并转换参与反演的海洋可控源电磁场数据；设置联合反演执行参数；设置联合反演初始模型参数；构建各向异性电阻率和接收站位置参数联合反演目标函数；求取电磁场关于各向异性电阻率和接收站位置参数雅各比矩阵和海森矩阵；基于反演参数特性自适应计算正则化因子；求取模型更新量；计算反演迭代模型的目标函数拟合差；判断是否满足反演要求；输出最终反演模型。该方法同时考虑海底介质电阻率各向异性和接收站的位置对反演结果的影响，其适用性更广，能够适用于电阻率各向异性资料的反演解释；容错性更强，能够处理电磁接收站位置参数存在一定误差时的反演问题。

Description

拖曳式海洋电磁地层电阻率与接收站位置联合反演方法

技术领域

本发明涉及多参数联合反演的海洋地球物理技术领域，具体涉及一种拖曳式海洋电磁地层电阻率与接收站位置联合反演方法。

背景技术

海洋可控源电磁法(CSEM)是探测海底油气资源和矿产资源的一种海洋地球物理勘探方法，在外国已被列为石油资源开采前必须采用的勘探方法之一。拖曳式施工方式是海洋可控源电磁法海上作业常用的方式之一，该方式将发射电偶极源与多个接收站串联，科考船拖曳着发射源和接收站同步前进，在发射源激发电磁波的同时接收站接收来自地下的电磁场信号。该观测方式具有施工速度快、采集效率高和实时采集监测的特点。然而该施工方式也存在弊端：由于发射源和接收站都同步拖曳施工，在拖曳的过程中发射源和接收站受海水运动和船速不稳定等因素影响，使得发射源和接收站不能在同一条测线上平行移动，致使发射源和接收站的相对位置出现偏差，虽然在实际的海洋探测中发射源和接收站的位置参数均有记录，但实测数据表明突发的接收站位置变化很难被记录到，导致采集的电磁信号存在误差，不准确的数据对海洋电磁资料的处理和反演解释带来了极大的困难。为提高采集电磁信号的质量，需对接收站的相对位置进行准确计算。

海洋电磁勘探实例表明，世界上大约30％的油气资源赋存于电性各向异性地层中，并且海底介质的电性各向异性是获得海底介质正确的电性分布的重要影响因素，若在解释海洋电磁资料时忽略电阻率各向异性，则很有可能无法获得合理的海底地电模型。然而，目前应用广泛的反演解释方法仍为传统的地球物理反演方法，该类方法在反演海洋资源勘探获取的海洋电磁资料时，常常假定海底介质为电性各向同性，这样的反演方法可能给后期资料解释提供巨大误差的反演结果，这也可能将导致获得错误的解释结果。为此，解释海洋电磁资料时必须考虑海底介质的电阻率各向异性的性质。

Ramananjaona(2011)分析电阻率垂直各向异性灵敏度特征，并利用Occam反演方法解释层状电阻率垂直各向异性海洋CSEM资料；罗鸣(2016)利用高斯-牛顿方法实现了一维垂直各向异性海洋CSEM反演算法，并对分析了多发射源、多频率数据能够更好的重构海底电性结构。然而，以上方法均仅反演了海底各向异性电阻率，而为考虑接收站位置对反演的影响。Swidinsky等(2011)提出海洋可控源电磁拖曳式采集系统的导航参数与海底各向同性电阻率的联合反演，并利用雅克比矩阵的特征值分解来分析反演的可行性，然而该算法未考虑海底电阻率存在各向异性的情况。此外，在现有的反演解释方法中，鲜有考虑接收站位置参数对反演结果的影响，大多方法均是直接利用存在误差的接收站参数进行反演解释，其反演结果必定也存在误差。

发明内容

本发明的目的在于提供一种拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，该方法可应用于海洋电磁探测和反演的精确性。

为了实现以上目的，本发明提供如下技术方案：

1.一种拖曳式海洋电磁地层电阻率与接收站位置联合反演方

法，其特征在于，主要包括：

S1、读取并转换参与反演的海洋可控源电磁场数据，所述数据包括电磁场实部与虚部数据、振幅与相位数据、极化椭圆长轴与短轴参数；

S2、设置联合反演执行参数，所述参数包括反演最大迭代次数、目标拟合差、最大迭代步长、步长比例系数、惩罚函数类型、正则化衰减系数；

S3、设置联合反演初始模型参数，所述参数包括背景层电阻率参数、厚度参数、观测系统参数；

S4、构建各向异性电阻率和接收站位置参数联合反演目标函数；

S5、求取电磁场关于各向异性电阻率和接收站位置参数的雅各比矩阵和海森矩阵；

S6、基于反演参数特性自适应计算正则化因子；

S7、求取模型更新量；

S8、计算反演迭代模型的目标函数拟合差；

S9、判断是否满足反演要求，若满足则转向S10，若不满足则转向S5；

S10、输出最终反演模型。

2.如权利要求1所述的拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，其特征在于，构建各向异性电阻率和接收站位置参数联合反演目标函数为：

式中，φ为反演算法的目标函数；m为模型反演参数向量，其包括海底各向异性电阻率参数m_ρ和接收站的位置参数m_P，即m＝m_ρ+m_P；

为模型参数向量的梯度；||·||为标准差算子；d为反演使用的观测数据向量；W_d为数据加权矩阵；W_m为模型加权矩阵；F(m)表示模型m的正演响应算子；μ_ρ和μ_p分别为反演模型中海底各向异性电阻率参数m_ρ和接收站位置参数m_P的正则化因子。

3.如权利要求1所述的拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，其特征在于，所述自适应正则化因子按照如下公式确定：

式中，i表示第i次反演迭代；μ_i为正则化因子；Max|·|为求取矩阵绝对值最大的元素；a_mj为矩阵乘积[(W_dJ)^T(W_dJ)]的元素；M为矩阵乘积[(W_dJ)^T(W_dJ)]的维度；χ为衰减系数；λ为反演模型海底地层的横向电阻率ρ_h、垂向电阻率ρ_v和接收站位置参数(x,y)的加权因子,采用以下公式确定：

其中，α为比例系数，m为反演模型参数。

4.如权利要求1所述的拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，其特征在于，基于反演参数计算模型参数更新量的方法为：

其中，Δm为下一次迭代的模型参数更新量，i为第i次反演迭代，H_i为海森矩阵，g_i为目标函数的梯度。

5.如权利要求4所述的海底各向异性电阻率与接收站位置联合反演方法，其特征在于，求取目标函数的极小值的方法为：

其中，i为第i次反演迭代，J_i为雅各比矩阵。

6.如权利要求1所述的海底各向异性电阻率与接收站位置联合反演方法，其特征在于，求取雅各比矩阵J_i的方法为：

其中，i为第i次反演迭代；J_i为正演响应F(m)的雅各比矩阵；ρ＝(ρ_h,ρ_v)为地层的各向异性电阻率分布，P＝(x,y)为接收站位置参数。

7.如权利要求1所述的海底各向异性电阻率与接收站位置联合反演方法，其特征在于，求取海森矩阵H_i的方法为：

其中，i表示第i次反演迭代次数。

较现有技术相比，本发明一些实施例中，提供的方法的有益效果在于：

本发明主要针对海洋可控源探测中广泛存在的海底介质各向异性问题，以及拖曳式海洋可控源电磁施工中电磁接收站位置存在误差导致反演结果误差大的问题，提出了一种拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，该方法在反演过程中考虑海底介质电阻率各向异性的同时，也考虑接收站的位置参数对反演结果的影响，能够同时对海底介质各向异性电阻率，以及电磁接收站的位置参数进行反演，最终可以获得与实际数据拟合最佳的海底各向异性电阻率分布情况，以及接收站的真实位置参数信息。该反演方法不仅为复杂的海底介质电阻率各向异性电磁资料的解释提供了一种有效的反演方法，同时也为电磁接收站位置参数存在误差情况下的数据处理问题提供了一种可行的技术手段。较于传统的反演算法，该方法同时考虑海底介质电阻率各向异性和接收站的位置对反演结果的影响，联合二者作为反演参数同时进行反演，该方法能够有效解决接收站位置不准确给海洋可控源电磁资料解释带来的影响；所提出的联合反演方法适用性更广，能够适用于电阻率各向异性资料的反演解释；容错性更强，能够处理电磁接收站位置参数存在一定误差时的反演问题。

附图说明

为了更清楚地说明本发明实施例中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的部分实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。

图1为本发明方法的程序流程框图；

图2为一维电阻率各向异性模型结构示意图；

图3为反演迭代模型目标拟合差与反演迭代次数的变化图；

图4为接收站位置参数(x,y)的反演结果示意图；

图5为海底各向异性电阻率的反演示意图。

具体实施方式

为了使本发明所要解决的技术问题、技术方案及有益效果更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

本发明提供一种拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，该方法同时考虑海底介质电阻率各向异性和接收站的位置参数对反演结果的影响，联合二者作为反演参数同时进行反演，该方法能够有效解决接收站位置参数不准确给海洋可控源电磁资料解释带来的影响。

实施例1

参见图1，一种拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，主要包括以下步骤：

S1.读取并转换参与反演的海洋可控源电磁场数据，所述数据包括电磁场实部与虚部数据、振幅与相位数据、极化椭圆长轴与短轴参数。

根据用户需求，将参与联合反演的海洋可控源电磁数据文件转换成特定格式的数据文件，并将输入的电磁场数据根据用户要求转换成特定参量(电磁场实部与虚部、振幅与相位、极化椭圆长轴与短轴参数等)，完成反演数据的读入并启动联合反演。

S2.设置联合反演执行参数，所述参数包括反演最大迭代次数、目标拟合差、最大迭代步长、步长比例系数、惩罚函数类型、正则化衰减系数。

联合反演过程是一个复杂且庞大的运算体系，在实现反演迭代直至获得最终的联合反演模型的过程中，需要运用到大量参数。不同的反演执行参数能够获得不同的反演效果，用户可根据参与联合反演数据的特点以及需要用于的特定需求设置反演执行参数。主要涉及的反演参数包括：反演最大迭代次数、目标拟合差、最大迭代步长、步长比例系数、惩罚函数类型、参与反演的数据分量、正则化衰减系数等。

具体的，联合反演执行参数的具体功能如下。

S21.反演最大迭代次数：反演迭代次数的最大值。设置该参数能够使联合反演过程达到该最大迭代次数后退出联合反演程序，避免反演程序持续循环运行，提高程序计算效率。

S22.目标拟合差：反演模型与真实模型拟合程度的目标值。联合反演的最终目标是让反演模型不断逼近于真实模型，数据拟合差即是判断反演模型与真实模型的拟合标准。用户根据数据质量设定目标拟合差，当反演模型的拟合差达到或者低于目标拟合差，联合反演退出；若反演模型的拟合差没有达到目标拟合差，联合反演程序将持续进行。

S23.最大迭代步长：用户设置联合反演过程中模型参数变化量的最大值。在联合反演中，由于海洋可控源电磁场(或其他参量)的数量级跨度较大，所以通常都是取对数后再参与联合反演运算。若迭代的模型变化量过大，其结果变化非常大，例如模型参量(取对数后)变化量为1.0，即其真实值即提高了一个数量级，这样的变化值通常情况是不合理的，因此设计该参数限制模型参数的变化量大小。

S24.步长比例系数：用户设置联合反演过程中模型参数变换量的比例系数。联合反演涉及的反演参数较多，在反演迭代过程中，不同反演参数的变化程度不同，为使得整体的模型变化量在一个较为合理的范围，可将模型参数变化量乘以一个比率系数(即步长比例系数),以保证反演模型稳健收敛。

S25.惩罚函数类型：模型参量之间的制约方式。在联合反演过程中，所有的模型参量同时参与反演，不同反演参量之间存在一定的联系，通过选择不同的惩罚函数可实现模型参量之间的制约方式，可使得模型更加光滑或曲度增大等效果。

S26.正则化衰减系数：自适应选择正则化因子的衰减系数。在反演迭代过程中，为保证不同参数之间的比率关系，正则化因子需随反演次数的改变而做相应改变，该参数即是设置正则化因子随迭代次数的衰减系数。

S3.设置联合反演初始模型参数，所述参数包括背景层电阻率参数、厚度参数、观测系统参数。

反演初始模型是联合反演模型参数的初始值，通常情况下需根据输入数据的观测系统和已知的水深数据等资料设定。

S4.构建各向异性电阻率和接收站位置参数联合反演目标函数。

海洋可控源电磁海底各向异性电阻率与电磁接收站位置参数联合反演涉及多属性的反演参数，参数间关系复杂。鉴于海底介质不同电阻率各向异性参数，以及接收站位置参数之间差异性，本发明利用正则化约束稳定海洋可控源电磁海底各向异性电阻率与电磁接收站位置参数联合反演过程。所采用的目标函数为：

式中，φ为反演算法的目标函数；m为模型反演参数向量，其包括海底各向异性电阻率参数m_ρ和接收站的位置参数m_P，即m＝m_ρ+m_P；

为模型参数向量的梯度；||·||为标准差算子；d为反演使用的观测数据向量；W_d为数据加权矩阵；W_m为模型加权矩阵；F(m)表示模型m的正演响应算子；μ_ρ和μ_p分别为反演模型中海底各向异性电阻率参数m_ρ和接收站位置参数m_P的正则化因子。

反演模型的海底电阻率反演参数m_ρ有以下形式：

式中，ρ_h和ρ_v分别为海底地层的横向电阻率和垂向电阻率，M为海底地层的层数。

接收站位置参数m_p有以下形式：

m_P＝[log₁₀x₁…log₁₀x_t,log₁₀y₁…log₁₀y_t]^T；

式中，(x,y)为接收站位置参数，t为电偶极源的个数。

S5.求取电磁场关于各向异性电阻率和接收站位置参数雅各比矩阵和海森矩阵。

S51.雅各比矩阵J_i是电磁场关于反演参数的偏导数矩阵，不同的反演参数有不同的形式：

其中，为第i次反演迭代；J_i为正演响应F(m)的雅各比矩阵；ρ＝(ρ_h,ρ_v)为地层的电阻率分布(横向电阻率，垂向电阻率)，P＝(x,y)为接收站位置参数。

若将电阻率张量ρ为反演模型横向电阻率ρ_h和垂向电阻率ρ_v的函数，即ρ＝f(ρ_h,ρ_v)，因此，电磁场关于电阻率张量ρ的雅克比矩阵有如下形式

同理的，可得到电磁场关于接收站的位置反演参数m_p的雅各比矩阵为

式中，(x,y)为接收站位置参数，t为电偶极源的个数

S52.海森矩阵H_i为目标函数对反演参数的二阶导数，忽略二阶导数项和非对称项，海森矩阵H_i可简化为

其中，i为第i次反演迭代，J_i为正演响应F(m)的雅各比矩阵；m_ρ和m_P分别为海底电阻率参数和接收站的位置参数；

为模型参数向量的梯度；||·||为标准差算子；d为反演使用的观测数据向量；W_d为数据加权矩阵；W_m为模型加权矩阵；μ_ρ和μ_p分别为反演模型中海底电阻率参数m_ρ和接收站的位置m_P的正则化因子。

S6.基于反演参数特性自适应计算正则化因子。

正则化因子μ的选择方式是求解海洋可控源电磁海底各向异性电阻率与接收站位置参数联合反演问题合理与否的关键。在地球物理反演方法中，出现了许多正则化因子选择方法。在进行海洋可控源电磁海底各向异性电阻率与接收站位置参数联合反演时，本发明针对不同的反演参数选择不同的正则化因子。首先，选择反演过程中的关于雅各比矩阵的特征参数作为基数，然后，再结合反演参数的特征及参数之间的关系调节不同反演参数的正则化因子大小，从而实现正则化因子的自适应选择。

正则化因子μ_i可写成如下形式：

式中，i表示第i次反演迭代；μ_i为正则化因子；Max|·|为求取矩阵绝对值最大的元素；a_mj为矩阵乘积[(W_dJ)^T(W_dJ)]的元素；M为矩阵乘积[(W_dJ)^T(W_dJ)]的维度；χ为衰减系数；λ为反演模型横向电阻率ρ_h、垂向电阻率ρ_v和接收站位置参数(x,y)的加权因子,采用以下公式确定

其中，α为比例系数，m为反演模型参数，权重向量的初始值λ₁为单位向量。当反演模型地层中存在电阻率各向异性时，正则化因子可根据电阻率各向异性率自适应地进行调节，从而调节目标函数中各元素之间的权重，从而达到自适应调节反演过程的效果。

S7.求取模型更新量。

通过求取目标函数的极小值可计算模型更新量，即令目标函数梯度g_i＝0，可由下式确定：

其中，i为第i次反演迭代，J_i为正演响应F(m)的雅各比矩阵；φ为反演算法的目标函数；m为模型反演参数向量，其包括海底电阻率参数m_ρ、接收站的位置m_P，即m＝m_ρ+m_P；

为模型参数向量的梯度；||·||为标准差算子；d为反演使用的观测数据向量；W_d为数据加权矩阵；W_m为模型加权矩阵；F(m)表示模型m的正演响应算子；μ_ρ和μ_p分别为反演模型中海底电阻率参数m_ρ和接收站的位置m_P的正则化因子。

计算得到下一次迭代的模型参量更新量为Δm，可表示为：

其中，Δm为下一次迭代的模型参数更新量，i为第i次反演迭代，H_i为海森矩阵，g_i为目标函数的梯度。

S8.计算反演迭代模型的目标函数拟合差。

在反演迭代过程中，反演迭代模型逐步收敛于真实模型，其目标函数拟合差ψ逐渐收敛于1.0。目标函数拟合差ψ有如下形式：

其中，i为反演迭代次数，d为反演使用的观测数据向量；F(m)表示模型m的正演响应算子；N为反演数据的个数；k为反演数据的编号；δ_k为第k个数据的标准差。

S9.判断是否满足反演要求，若满足则转向S10，若不满足则转向S5。

判断是否退出联合反演的标准有两个：1)是否满足联合反演目标拟合差；2)反演迭代次数是否达到最大迭代次数。若否，则转向S5，继续迭代求解模型更新量；若是，则转向S10。

S10.输出最终反演模型。

实施例2

参考图2，为一维典型地电模型原理图。为了验证海洋可控源电磁海底各向异性电阻率与接收站位置参数联合反演的有效性，以图2所示一维电阻率各向异性模型为例，利用合成数据进行海底各向异性电阻率与发射源参数的联合反演。假设发射源位于海底正上方50m处，布设26个接收站等间距拖曳于发射源后方0m-13km范围的海水中。设定发射电流为1安培。反演数据由电磁场各分量的实部和虚部组成，各数据加入了2％的随机高斯噪声。在反演算例中，反演初始模型为空气、海水和电阻率为1Ωm的均匀半空间，空气层和海水电阻率和深度固定。设定反演初始模型中接收站位置与真实测线位置交叉(存在一个约45°的夹角)。

参考图3，为反演迭代模型目标拟合差与反演迭代次数的关系示意图。由图可见，在联合反演过程中，目标函数拟合差持续减小并逐渐向真实模式靠近，并由最初较大的拟合差12.3，经过30次反演迭代后，最终收敛于1.0。

参考图4，为接收站位置参数反演结果。由图4可见，图中纵坐标和横坐标分别为(x,y)的位置，空心三角为电磁接收站的初始位置，实心三角为接收站真实位置，实心方块为联合反演获得的接收站位置，即使输入的接收站位置与真实位置差异较大(位置差异达100m)，联合反演得仍然能够得到相对准确的位置信息(反演得到的接收站位置与真实位置距离不超过20m)，由此说明，联合反演得到的接收站是准确的。

参考图5，为海底各向异性电阻率的反演结果图，浅灰色和深灰色线段为真实的海底介质横向和垂向电阻率分布情况，浅灰色和深灰色阶梯曲线为反演得到的海底介质横向电阻率和垂向电阻率曲线，由图可见，联合反演的围岩的各向异性电阻率与真实情况相符，高阻层的埋深、厚度和垂向电阻率值也得到了较准确的恢复。综上所述，本发明提出的海洋可控源电磁海底各向异性电阻率与接收站位置参数联合反演方法能够反演出海底介质的各向异性电阻率信息，以及电磁接收站的位置态参数。由此也说明，本发明提出的算法是有效的。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (4)

1.拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，其特征在于，包括：

S1、读取并转换参与反演的海洋可控源电磁场数据，所述数据包括电磁场实部与虚部数据、振幅与相位数据、极化椭圆长轴与短轴参数；

S2、设置联合反演执行参数，所述参数包括反演最大迭代次数、目标拟合差、最大迭代步长、步长比例系数、惩罚函数类型、正则化衰减系数；

S3、设置联合反演初始模型参数，所述参数包括背景层电阻率参数、厚度参数、观测系统参数；

S4、构建各向异性电阻率和接收站位置参数联合反演目标函数；

S5、求取电磁场关于各向异性电阻率和接收站位置参数的雅各比矩阵和海森矩阵；

S6、基于反演参数特性自适应计算正则化因子；

S7、求取模型更新量；

S8、计算反演迭代模型的目标函数拟合差；

S9、判断是否满足反演要求，若满足则转向S10，若不满足则转向S5；

S10、输出最终反演模型；

构建各向异性电阻率和接收站位置参数联合反演目标函数为：

式中，φ为反演算法的目标函数；m为模型反演参数向量，其包括海底各向异性电阻率参数m_ρ和接收站的位置参数m_P，即m＝m_ρ+m_P；

为模型参数向量的梯度；||·||为标准差算子；d为反演使用的观测数据向量；W_d为数据加权矩阵；W_m为模型加权矩阵；F(m)表示模型m的正演响应算子；μ_ρ和μ_p分别为反演模型中海底各向异性电阻率参数m_ρ和接收站位置参数m_P的正则化因子；

求取雅各比矩阵J_i的方法为：

其中，i为第i次反演迭代；J_i为正演响应F(m)的雅各比矩阵；ρ＝(ρ_h，ρ_v)为地层的各向异性电阻率分布，P＝(x，y)为接收站位置参数；

求取海森矩阵H_i的方法为：

其中，i表示第i次反演迭代次数。

2.如权利要求1所述的拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，其特征在于，所述自适应正则化因子按照如下公式确定：

式中，i表示第i次反演迭代；μ_i为正则化因子；Max|·|为求取矩阵绝对值最大的元素；a_mj为矩阵乘积[(W_dJ)^T(W_dJ)]的元素；M为矩阵乘积[(W_dJ)^T(W_dJ)]的维度；χ为衰减系数；λ为反演模型海底地层的横向电阻率ρ_h、垂向电阻率ρ_v和接收站位置参数(x，y)的加权因子，采用以下公式确定：

其中，α为比例系数，m为反演模型参数。

3.如权利要求1所述的拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，其特征在于，基于反演参数计算模型参数更新量的方法为：

其中，Δm为下一次迭代的模型参数更新量，i为第i次反演迭代，H_i为海森矩阵，g_i为目标函数的梯度。

4.如权利要求3所述的拖曳式海洋电磁地层电阻率与接收站位置联合反演方法，其特征在于，求取目标函数的极小值的方法为：

其中，i为第i次反演迭代，J_i为雅各比矩阵。

Images