2021, Vol. 64
随着油气勘探的不断深入，以构造圈闭为主的优质储层目标越来越少，勘探开发的对象逐步向岩性油气藏、构造-岩性复合油气藏拓展，对地震成像与储层预测提出了更高的精度要求.在这种情况下，以声压波场传播理论为基石的传统标量地震纵波勘探方法的局限性日益突出.多分量地震更全面地记录了纵波(P波)和横波(S波)的三维振动信息，联合P-S波的地震成像与反演除了可改善“气云区”和复杂地质体(如盐丘、火成岩)及其下部的构造成像，还有助于岩性区分、储层内部非均质性(包括裂缝系统)刻画以及含油气性检测，并在油气藏开发、二氧化碳填埋等过程的动态监控中体现出一定优势(Stewart et al., 2002; 张永刚等, 2004; 刘振武等, 2008; 何兵寿等, 2019).
近20年来，在现代地震观测与高性能计算技术的推动下，多分量地震方法在油气资源勘探领域愈发受到重视.在海上，多分量地震以海底电缆(OBC)、海底节点(OBN)和海底地震仪(OBS)等观测方式为主，由声压检波器(俗称水检)和三分量质点速度(或位移/加速度)检波器(俗称陆检)共同接收到达海底的矢量波场信号.相比于陆地三分量(3C)地震，海底四分量(4C)地震多记录一个声压场信号，就会受到海面多次波、水层鸣震以及涌浪噪声等干扰.国内外迄今采集了大量海底多分量地震资料，但资料的应用还基本处于联合利用水检和陆检信号压制接收端海水鸣震阶段，以此改善纵波成像质量.简便起见，有时也会直接将陆检垂直分量信号视为P波，将水平分量信号视为PS波.然而，忽视P波与S波在这些分量的模式耦合会不同程度地影响多分量数据处理的效果(何兵寿和张会星, 2006).理论与实践表明，多分量波场P/S波模式分离(或解耦)有利于后续地震数据处理，尤其是速度建模和偏移成像(陶春辉和何樵登, 1993; 许世勇等, 1999; 胡天跃等, 2004; Yan and Sava, 2008)，也能降低弹性波反演面临的非线性和多参数耦合效应(Wang et al., 2015, 2018; Wang and Cheng, 2017; Ren and Liu, 2016; 王腾飞, 2017; Xu et al., 2019).目前，联合水检和陆检信号实现上-下行P/S波分离，已被视为海底4C地震数据处理的一个最关键环节，有助于压制水层鸣震或海面多次波，提高构造成像和储层预测的可靠性(Edme and Singh, 2009; 李维新等, 2018; 芦俊等, 2018).
海底地震波场分离有两大类主流方法：第一类为声波场分离方法(只能实现上下行波分离)，它把分离界面置于紧挨海底的水层中，根据水听分量和陆检垂直分量上行波信号极性相同、下行波极性相反的原理，将二者匹配后叠加分离出上、下行波场，俗称PZ叠加或双检合并.该方法最先由Haggerty(1956)和White(1965)等提出，随后得到一系列扩展(Dragoset and Barr, 1994; Osen et al., 1999).这类声波理论基础上的波场分离方法仅适用于P波数据，无法区分P波与S波，因此难以支撑转换波(如PS波)速度建模与成像处理.
第二类为弹性波场分离方法(实现纵横波解耦和上下行波分离)，它把分离界面置于紧挨海底的固体介质中，具体实现要么依据海底波场极化信息(Cho et al., 1992)，要么采用弹性波特征分解方法(Schalkwijk et al., 1999).前者属于数据驱动方法，主要依据到达海底检波器的P波与S波极化方向的差异，但该方法过程繁琐且易受到噪声的影响，也不能区分上、下行波.后者基于弹性波动方程的平面波解，依据特征值分解(Frasier, 1970; Aki and Richards, 2002)或海底界面P/S波能量分配关系(Wang et al., 2002; Edme and Singh, 2009)构建依赖于海底介质参数的上-下行P/S波解耦方法.它一般在射线参数域实现，需要适当拓展才能处理海底介质的横向非均匀性(Amundsen et al., 2000).更重要的是，海底上-下行P/S波分离还受到水检与陆检在海底的耦合性和仪器响应差异、海底固体介质弹性参数的影响.于是，Schalkwijk等(1999, 2003)提出了一种有较强实用性的多步分离技术流程.该流程需在炮集数据上人为选取特定反射波时窗，人工交互导致工作量大，且在浅水多次波发育和信噪比较低情况下选择时窗存在困难.为此，Muijs等(2004, 2007)重新建立了一种面向浅水多分量地震数据的波场分离技术流程.总之，海底多分量地震数据上-下行P/S波分离是一个复杂的系统工程，步骤多，对工区条件、资料品质的依赖性很强，几乎每步都需针对特定的数据成分或震相求解相应的反问题.因此，现有方法与流程在有些实际地震资料处理中成效并不明显.
针对我国东海油气田绳系OBS四分量地震数据波场分离与偏移成像遇到的挑战，本文结合工区水深环境与资料信噪比特征，构建了一个海底多分量地震记录上-下行P/S波分离的新流程.通过二维模型数据与实际资料分离实验，分析几个关键步骤采用的方法以及整个流程的有效性.
1 方法原理与常规技术流程海底多分量地震数据采集时，通常认为水检处于紧邻海底的水层中，陆检位于紧邻海底的固体介质中(如图 1a，Δh接近于零).海底多分量地震数据上-下行P/S波分离的基本原理可由弹性波表示定理推导而来(Amundsen et al., 2000).对于复杂的海底介质情况，分离过程需要求解一系列复杂的积分方程.当海底介质横向均匀或变化缓慢时，积分方程退化为由频率域弹性波动方程特征分解(Frasier, 1970; Aki and Richards, 2002)建立的波场分离公式(Schalkwijk et al., 1999).下文对二维三分量数据的波场分离问题进行详细阐述和分析.
|
图 1 海底多分量地震记录波场分离示意图 (a)海底地震仪检波器位置；(b)声波场上-下行波分离；(c)弹性波场上-下行P/S波分离.图中水平虚线代表波场分离界面，蓝色箭头细线表示经海底反射立即在分离界面接收到的波场. Fig. 1 Schematic diagram of wavefield decomposition for ocean bottom multi-component seismic recordings (a) Detector positions of ocean bottom seismograph (OBS); (b) Up/down-decomposition of acoustic wavefield; (c) Up/down P- and S-wave decomposition of elastic wavefield. The horizontally dashed lines represent the interface of wavefield decomposition, and the thin arrows represent the wavefields that reflect at seabed and immediately arrive at the given interface. |
首先，海底多分量数据的声波场上-下行波分离发生在水检所在的界面(如图 1b中黑色水平虚线所示).由图可见，分离界面处的下行P波是经自由表面反射下行到海底的，而上行到分离界面的P波包含两种成分：一是在海底固体介质传播的上行P波(包括由S到P波上行转换波)，二是下行P波经海底反射并立即到达分离界面的上行干扰P波，通常情况下这种干扰波的能量会比较强.基于二维频率域声波方程的特征分解，可以推导出海底上-下行波分离公式(Amundsen, 1993)：
|
(1) |
其中P和vz分别表示水检分量和陆检垂直分量.公式左边项上标“±”分别表示下行波和上行波(如图 1所示)，公式右边项“±”单纯表示正负号，分别对应下行波和上行波计算所用符号.ρ0代表海水密度，α0为海水速度，p=kx/ω是水平慢度(或射线参数)，q0=(α0-2-p2)1/2表示垂向慢度，kx为水平波数，ω为圆频率.基于该公式的波场分离在ω-p域、ω-k域或τ-p域均可实现.因为海底界面满足质点振动速度连续性条件，公式(1)右端vz分量可用陆检记录垂直分量来代替，如果陆检分量记录的是位移或者加速度则需要先将其转换成速度.
其次，海底多分量数据的弹性波场上-下行P/S波分离发生在陆检所在的界面(图 1c黑色水平虚线表示).可见，在分离界面处可明确区分上行P波和S波，但是分离出的下行P波除了从海水透射的下行P波，还包含分别由上行P和S波在海底界面产生的下行反射和转换P波场.同理下行S波在分离界面处也由三种成分构成.当海底存在松软沉积层时，自由表面反射波到达海底透射到分界面处的下行P波(图 1c黑色带箭头粗实线)是下行P波的主导成分，而由上行S波经海底反射到分离界面的下行S波(图 1c浅蓝色带箭头粗线)在所有下行S波成分中占主导地位，其他海底相关的干扰体波成分会相对较弱.
基于二维频率域弹性波方程的特征分解，并考虑到海底切应力为零，可以导出海底上-下行P/S波分离公式(Schalkwijk et al., 2003)：
|
(2) |
以及，
|
(3) |
其中vx为陆检速度水平分量，τxz、τzz分别为切应力和正应力，α、β、ρ分别代表海底固体介质的纵、横波速度与密度，qP=(α-2-p2)1/2和qS=(β-2-p2)1/2分别为纵波和横波垂向慢度，F=4β2p2qPqS+β2(qS2-p2)2.φ±表示分离出的下、上行标量纵波，ϕ±表示分离出的下、上行标量横波.考虑到海底界面附近法向应力连续，海底多分量数据只能接收到水检分量，因此公式(2)右端τzz项用水检分量P代替，即-τzz=P.
假设海底多分量检波器真实地记录了质点振动矢量场信息，按上述方法原理可以分离出海底上-下行P/S波记录.然而，在实际海底地形与洋流条件下，某些地震仪与海底介质的耦合性得不到保障，导致声压分量与陆检分量在相位和振幅方面不符合理论上的数值关系.此外通常还会面临真实海底参数未知的问题.因此很难直接利用公式(1)-(3)对实际海底多分量数据进行数据分离.
随着转换波在海洋油气勘探中愈发受到重视，地球物理学者也在不断探索基于弹性波理论的上-下行P/S分离技术流程.Schalkwijk等(1999, 2003)经过多年的持续研究，提出了一种多步分离策略：首先，通过最小化一次反射波时窗中下行波的能量求取vz分量的校正因子.然后，通过最小化直达波时窗中上行波能量依次估计海底介质参数并求取vx分量的校正因子.最后，根据特征分解理论导出的波场分离公式计算上-下行P/S波记录.不过，该流程在耦合校正和参数估计时，需要在炮集数据上人工选定合适的时窗，除了工作量大，还容易受多次波及其他噪声干扰.为此，Muijs等(2007)提出了一种面向浅水多分量地震数据P/S模式分离的新流程，它基于海底接收的上行波与下行波总能量相等这一理论假设建立的目标函数来反演vz分量的校正因子，通过建立互相关目标函数来估计海底介质参数和反演vx分量校正因子.此外，有学者指出，依据海底接收的水听分量和陆检垂直分量中直达波的AVO关系，也可以估计海底介质参数(Amundsen and Reitan, 1995; Muijs et al., 2003).
2 浅水区海底多分量地震数据波场分离新流程海底地震仪在重力作用下降落至海底，与海底接触存在不可控性，且海底地形、环境与洋流状况复杂，记录中环境噪声和多次波发育.尤其是在浅水条件下，环境噪声的分布不均匀性远高于深水区(刘丽华等, 2012)，且多次波周期相对较短，对一次波的影响更为复杂.尽管全球范围内已经采集了许多海底4C地震资料，但在偏移成像与储层预测过程中真正利用好转换波(如PS波)数据的实例并不多，其中一个重要原因就是海底波场分离效果不理想(李维新等, 2018).
图 2显示了国内于2017年首次在东海油气区将某型国产海底地震仪以绳系方式投放后，在一条二维测线上某地震仪采集的三分量地震记录(已完成涌浪噪声压制、地震仪重定位与旋转校正等预处理).由于前文提到的浅水环境噪声与多次波、仪器和施工因素影响，声压分量的直达波受到严重污染，波形受到明显破坏, 浅层(0.7 s以内)反射波几乎被噪声掩盖(图 2a黑色箭头指示区域)，陆检垂直分量的近偏移距处存在很强的波型泄漏干扰(图 2b黑色箭头指示同相轴).由于很难挑选出信噪比有保证的直达波时窗以及多次波污染小的一次反射波时窗，前人针对垂直分量的校正方法无法使用.基于海底直达波的AVO关系估计海底介质参数也比较困难.在理想情况下，初至折射波数据是上行波信号(图 2红色箭头指示同相轴)，因此它为仪器校正因子的求取提供了可以凭借的依据.
|
图 2 东海QY探区海底多分量共接收点记录 (a)水检P分量; (b)陆检Z分量数据; (c)陆检X分量数据. Fig. 2 OBS common-receiver gathers from QY area in East China Sea (a) P- component of hydrophone; (b) Z- component of geophone; (c) X- component of geophone. |
因此，本文针对类似的海底多分量地震数据，依据弹性波场分离基本原理，提出由如下步骤组成的方法流程：
(1) P/vz校正：依据折射波时窗内下行波能量最小化准则反演vz分量相对于P分量的校正因子，进而完成陆检垂直分量的匹配滤波与耦合校正；并通过公式(1)分离获得纯模式下行波场；
(2) 估计海底介质参数：利用到达海底的上行波与下行波零延迟互相关最小化准则估计海底介质参数；
(3) P/vx校正：利用到达海底的上行波与下行波零延迟互相关最小化准则反演vx分量的校正因子，进而完成陆检水平分量的匹配滤波与耦合校正；
(4) 海底上-下行P/S数据分离：按前文方法原理，由三分量记录分离获得纯模式的上行P波、上行S波.
在流程处理之前需要对实际数据进行预处理，包括检波器重定位、方位角和倾斜校正、去噪、插值等工作.
2.1 P/vz校正通常水检比陆检的耦合性更好，故以水检记录的声压信号为标准首先对陆检vz信号进行校正.在认为海水声波速度和密度已知条件下，通过P/vz校正之后的vz分量能够比较真实地代表海底质点振动的垂直速度，那么在折射波时窗(中长偏移距初至折射波同相轴所在时窗)内下行波能量接近于零.设频率域的校正因子为Cpz(ω)，声波场信号的上、下行波分离遵循公式(1).当校正因子足够准确时, 时窗内下行波能量达到极小值，因此建立如下ω-k域的目标函数：
|
(4) |
其中上标*表示复共轭.利用最速下降法迭代最小化下行波能量可估计出校正因子Cpz(ω)，进而完成对vz分量的校正(原始数据在频率域与校正因子相乘).将校正后的vz分量代入公式(1)分离出声波情况下的纯模式下行波场(P+/vz+).
2.2 估计海底介质参数当海底附近没有明显分层结构时，上行波与下行波不会同时到达海底.当海底上行波和下行波记录分离得比较彻底时，它们的零延迟互相关应当是一个极小值.弹性波场分离的方法依赖于海底介质参数，故可通过最小化分离后上行波和下行波记录的互相关系数估计海底介质参数(Muijs et al., 2007).
在前面P/vz校正之后按声波波场分离原理分离海底P+波记录的过程与海底介质参数无关，故以它作为分离后的下行波.以公式(2)分离出的海底上行正应力τzz-记录代表分离后的上行波，式中vz分量是校正后的陆检垂直分量.据此建立τ-p域的目标函数：
|
(5) |
其中b(p)与海底介质参数α、β、ρ有关，即：
|
(6) |
首先使EP最小化，估计出随射线参数变化的b(p)曲线.然后根据(p)公式(6)采用全局寻优方法拟合b(p)曲线反演获得海底介质参数.该步时窗选择无具体要求，以信噪比较高的反射波时窗为宜.
2.3 P/vx校正vx分量仅参与海底弹性波场分离过程，因此其耦合校正与匹配滤波只能与这个过程结合.本文仍然采用Muijs等(2004)提出的上行波与下行波不同时到达海底的假设，在海底介质参数已知情况下，通过最小化上行波与下行波的零延迟互相关系数求取针对陆检vx分量的校正因子Cpx(ω).这里下行波同样采用声波场分离出来的海底P+记录，上行波则选用上行到海底的切应力τxz-信号.由公式(3)建立τxz-与vx分量之间的关系，进而建立如下τ-p域的目标函数：
|
(7) |
式中R代表取复信号的实部.对目标函数的迭代最小化可以估计出陆检水平分量的校正因子.
2.4 上-下行P/S波分离按照上述方法获得校正之后的陆检垂直与水平分量以及海底介质参数，代入公式(2)分离出上-下行应力波场记录τxz±、τzz±，随后代入公式(3)分离出海底弹性波场分离界面上的上-下行P/S波信号.其中上行P波与上行S波地震记录将用于后续的多次波压制、速度分析、偏移成像与波形反演等处理.
3 理论模型数值实验首先在一个含有100 m深水层的水平层状介质模型(如图 3)上开展方法与流程测试.利用一阶速度-应力方程交错网格有限差分算法合成海底二维三分量(2D3C)炮记录.P波震源设置在模型(1800 m, 50 m)处，采用主频为40 Hz的雷克子波，水检和陆检以4 m间隔均匀布放在海底，以0.5 ms间隔记录时长1.8 s.利用镜像法仿真水平自由表面效应.如图 4，2D3C共炮记录中一次反射波被大量的自由表面多次波所污染.
|
图 3 浅水各向同性层状介质模型 (a)纵波速度；(b)横波速度；(c)密度. Fig. 3 Layered isotropic medium model in shallow water (a) P-wave velocity; (b) S-wave velocity; (c) Density. |
|
图 4 二维三分量共炮点道集 (a) P分量；(b) vz分量；(c) vx分量. Fig. 4 Synthetic two-dimensional three-component (2D3C) common-shot gathers (a) P component; (b) vz component; (c) vx component. |
利用合成的海底多分量记录，验证海底介质参数的准确性对波场分离结果的影响.对任意两个参数组合，在真实参数(用黑色实线标识)附近一定范围扰动，将分离结果的误差平方和与真实结果平方和的比值作为衡量误差的标准，来评估P波或者S波分离对不同海底参数的敏感性，公式表示为
|
(8) |
式中da表示使用精确海底介质参数获得的分离结果，dw为使用试验参数获得的分离结果，t和r分别对应时间和检波器采样序号.图 5a-c误差图谱表明，P波分离受纵波速度影响大，受横波速度和密度影响相对较小；图 5d-f表明S波分离主要受横波速度影响，受密度的影响较小，且不受纵波速度的影响.总体而言，三个参数中纵、横波速度的估计比较关键.
|
图 5 波场分离对海底介质参数的敏感性 第1行:上行P波分离，第2行:上行S波分离.黑线表示真实参数值所在位置. Fig. 5 Sensitivity of wavefield decomposition to seabed medium parameters Top: decomposition of up-going P-waves; Bottom: decomposition of up-going S-waves. Black line denotes the actual value of the corresponding parameter. |
为了仿真海底仪器响应差异与耦合问题，对陆检vx和vz分量在频率域施加如图 6所示的频率畸变因子C(ω)，使其振幅特性发生改变(暂时不改变相位).首先采用前文所述第一步进行P/vz校正，该过程使用图 4箭头指示的折射波同相轴(中长偏移距初至折射波同相轴所在时窗).如图 7，从真实校正因子(黑色虚线)与反演的校正因子(红色实线)对比可知，在有效频带范围内得到了校正因子的合理估计，图中真实校正因子由畸变因子计算获得(真实校正因子的振幅是畸变因子振幅的倒数，相位与畸变因子符号相反).然后，采用流程第二步反演海底介质参数.将P分量和校正后的vz分量变换到τ-p域，在图 8矩形框所示的反射波时窗内，通过最小化公式(5)定义的目标函数反演获得b(p)曲线(如图 9红色实线所示).而后基于公式(6)在合理范围内遍历参数组合拟合b(p)曲线，估计出海底介质参数.本例中α搜寻范围1.5~2.5 km·s-1，β搜寻范围0.1~1.0 km·s-1，ρ搜寻范围1.5~2.3 g·cm-3，搜寻步长都为0.1.最优拟合b(p)曲线对应的三参数为α=1.6 km·s-1，β=0.4 km·s-1，ρ=1.85 g·cm-3.与真实海底介质参数α=1.60 km·s-1，β=0.40 km·s-1，ρ=1.80 g·cm-3对比发现，纵横波速度估计比较准确，密度参数估计误差较大.这与前面密度敏感性最低的认识一致.
|
图 6 频率域畸变因子 Fig. 6 Distortion factor in the frequency domain |
|
图 7 vz分量频率域校正因子 Fig. 7 Calibration operator in the frequency domain for vz component |
|
图 8 τ-p域vz分量数据 Fig. 8 vz component of seismogram in τ-p domain |
|
图 9 估计的b(p)曲线 黑线为真实值，红线为反演值. Fig. 9 Estimated b(p) curve Black: true value, Red: inverted value. |
接下来在同一反射波时窗按流程第三步反演vx分量的校正因子Cpx(ω).结果如图 10所示，与真实校正因子对比可知，在有效频带范围内同样合理地估计出了针对陆检水平分量的校正因子.在完成陆检信号的校正与海底参数的反演之后，就可按公式(2)和(3)完成海底波场分离.图 11a、b、c展示了不进行P/vz与P/vx校正，但采用精确的海底介质参数的情况下，按声波场分离原理分离出的下行P波记录，以及按弹性波场分离原理直接分离出的上行P波与S波记录.可见，忽略陆检信号校正环节会影响分离效果.下行P波记录仍然残留折射波信号(图 11a箭头指示同相轴)，上行P波与S波记录均含有大量波型泄漏(图 11b、c箭头指示同相轴).图 11d、e、f展示了依次进行P/vz校正、海底参数估计以及P/vx校正，然后实施波场分离取得的结果.可见，P/vz校正保证分离出来的下行P波记录基本无折射波残存，无明显上行波干扰；基于校正后的vx分量和vz分量与估计的海底参数，按弹性波场分离原理分离出的上行P波和上行S波记录明显消除了检波点端下行波(海水鸣震)干扰，避免了P波与S波信号相互泄漏.从偏移距100 m处抽道对比(图 12)看出，严格按流程步骤处理的分离结果与利用合成数据得到的真实结果(利用精确海底参数对未施加畸变因子的模拟数据做分离)基本一致.可见，水检与陆检信号校正及海底参数估计对P波与S波的高质量分离至关重要.
|
图 10 vx分量频率域校正因子 Fig. 10 Calibration operator in the frequency domain for vx component |
|
图 11 检波器校正前后上-下行P/S模式分离结果对比 (a-c)校正前多分量数据分离的下行P波、上行P波和上行S波；(d-f)校正后多分量数据分离的下行P波、上行P波和上行S波. Fig. 11 Comparison of up-down P/S decomposition results before and after geophone calibration (a-c) Down- and up-going P-waves, and up-going S-waves from the seismograms without geophone calibration; (d-f) Down- and up-going P-waves, and up-going S-waves from the calibrated seismograms. |
|
图 12 偏移距100 m处抽道对比图 三幅图依次对应下行P波、上行P波与上行S波；黑色虚线为真实分离结果，蓝色虚线为未做校正数据分离结果，红色实线为校正后数据分离结果. Fig. 12 Comparison the seismic traces at the offset of 100 m The black dotted line, blue dotted line and red solid line denote the true signals, the decomposition results without and with geophone calibration, respectively. |
对数据施加统一强度的白噪声，通常会导致局部数据信噪比(如初至折射波时窗信噪比)与整体信噪比有较大差别.在东海实际数据中初至折射波时窗信噪比相对较高(图 2所示)，且折射波时窗主要受涌浪等噪声干扰，去噪相对简单.本文方法最终分离结果的准确性依赖于第一步vz分量校正因子估计的准确性，所以我们向合成2D3C数据加入两种不同强度的随机噪声，使得初至折射波时窗信噪比分别为4 dB和-10 dB(注意此时数据整体信噪比分别为20 dB和6 dB，如图 13所示，本文只展示vz分量)，以检验本文流程的抗噪性.
|
图 13 含随机噪声vz分量记录 (a)整体信噪比为20 dB；(b)整体信噪比为6 dB. Fig. 13 vz component recordings with random noise (a) 20 dB of SNR; (b) 6 dB of SNR. |
对折射波时窗信噪比为4 dB的vz分量和vx分量合成记录施加如图 6所示的频率畸变因子C(ω)，然后按本文流程进行模式分离.图 14a与14b分别为vz分量和vx分量校正因子求取结果，在有效频带范围内与真实值非常接近.图 15为反演的b(p)曲线，对应估计的海底介质参数为：α=1.6 km·s-1，β=0.4 km·s-1，ρ=1.92 g·cm-3，只有密度误差较大.由图 16可知，最终分离的上行P/S波记录与精确结果非常接近.这说明折射波时窗信噪比为4 dB时，本文方法流程仍然可以获得合理结果.
|
图 14 校正因子结果对比 (a) vz分量校正因子; (b) vx分量校正因子.黑色虚线是精确解，红色实线是反演结果. Fig. 14 Comparison of calibration operators (a-b) Calibration operators of vz and vx components (SNR=20 dB). Black dotted line and red solid line denote exact solution and inversion results, respectively. |
|
图 15 20 dB信噪比数据b(p)曲线 黑线为真实值，红线为反演值. Fig. 15 The estimated b(p) curves of SNR=20 dB data The black dashed lines denote the true values, and the red solid lines denote the inverted values. |
|
图 16 偏移距100 m处分离的上行P/S波地震道对比 (a)信噪比20 dB；(b)信噪比6 dB.黑线为真实值，红线为校正后数据分离结果.注意P波与S波记录为不同时窗范围. Fig. 16 Comparison of decomposed seismic traces of P- and S-waves at the offset of 100 m (a) SNR=20 dB; (b) SNR=6 dB.The black and red curves denote the real and inverted value. Note that different time windows are shown for P- and S-wave recordings. |
对折射波时窗信噪比为-10 dB的vz分量和vx分量数据施加同样的畸变因子，然后利用本文流程进行模式分离.由于折射波时窗信噪比过低，导致vz分量校正因子结果存在很大误差，进而导致参数估计和vx分量校正因子误差很大，本文不再给出具体结果.因此，只对vx分量施加畸变因子(vz分量没有畸变)，随后按照本文流程进行海底参数估计和反演vx分量校正因子(这两步选取的窗口信噪比与整体信噪比接近，为6 dB)，以验证这两步对噪声的敏感性.最终分离出的上行P/S波结果如图 16b所示，与精确结果非常接近.所以综合来看，本文所述方法流程具有较好抗噪性.
4 实际OBS资料处理东海QY工区位于西湖凹陷三潭深凹中部，平均水深约85 m，处在较为宽缓的大陆架缓坡地带，海底较为平缓，构造总体呈东南倾斜，具有典型的凹中隆构造特征.经中海油、中石化与同济大学三方协商，依托2017年二次三维地震采集航次机会与震源船设备，以绳系OBS方式沿一条10 km二维测线均匀布放100台国产OBS仪器，接收由空气枪震源激发产生的地震波.本文把OBS二维测线正上方炮数据分选出来开展海底多分量地震数据的波场分离实验.该测线共834炮，炮间距37.5 m，共有100个OBS节点，间距为100 m.原始数据最大偏移距为20 km，但本文聚焦反射波成像问题，故只截取和处理偏移距10 km以内的数据.
用于波场分离实验的二维三分量(2D3C)地震记录依次经过了检波器重定位、倾斜校正和方向旋转校正、时间重采样、带通滤波和点源到线源振幅校正等预处理.如前文图 2所示，由于仪器与施工因素以及浅水噪声尤其是多次波的影响，声压分量在0.7 s以内信噪比非常低，垂直速度分量在近偏移距存在严重的P与S波模式串扰，按现有基于直达波、反射波时窗的陆检信号校正和海底参数估计方法很难开展.同时，基于低信噪比和波形振幅受到明显破坏的直达波信号利用AVO反演估计海底参数也比较困难.
下面介绍本文方法流程的试验结果.从100个OBS中选取了数据质量满足要求的85个，按检波器道集开展海底上-下行P/S波模式分离.这里以OBS编号为33的2D3C共接收点道集进行展示.首先选定中偏移距信噪比有保障的折射波时窗，按第一步反演vz分量校正因子(如图 17)，进而完成vz分量校正.同时按声波场分离理论得到下行P波记录，结果如图 18所示，在完成上下行波分离之后，下行P波记录中基本不存在折射波信号，这表明P/vz校正效果明显.接着按第二步在2.0s附近选定信噪比有保障的反射波时窗估计出海底纵、横波速度与密度.图 19为这个过程中首先反演得到的b(p)曲线，由它进一步估计出海底参数分别为α=1.6 km·s-1，β=0.3 km·s-1，ρ=1.7 g·cm-3(由它们拟合出的b(p)曲线由图 19中红色曲线表示).为了保证算法的稳定性，在反演b(p)曲线过程中对目标函数的梯度(局部优化的下降方向)进行了正则化(光滑)处理.这样的高纵-横波速度比与工区内海底存在未固结沉积物的情况比较吻合.然后在同样的反射波时窗按第三步，求取vx分量的校正因子(如图 20)，完成vx分量校正.
|
图 17 vz分量频率域校正因子 Fig. 17 Calibration operator in frequency domain of z-component |
|
图 18 声波场共接收点记录分离结果 (a)声压分量(左)和垂直速度分量(右)；(b)分离的下行波场. Fig. 18 Acoustic wavefield decomposition of common-receiver gathers at the ocean bottom (a) Pressure component (left) and vertical particle velocity component(right); (b) Decomposed down-going waves. |
|
图 19 实际数据b(p)曲线 红实线为反演曲线，黑虚线为拟合曲线. Fig. 19 The estimated b(p) curves of the field data The red solid and black dotted lines denote the inverted and fitted results, respectively |
|
图 20 vx分量频率域校正因子 Fig. 20 Calibration operator in frequency domain for x-component |
最后对完成P/vz与P/vx校正的2D3C共接收点道集，采用估计出来的海底介质参数，按第四步进行海底弹性波场分离.图 21分别展示分离出的上行P波记录以及上行S波共接收点记录.与陆检垂直分量对比可知，上行P波基本消除了S波的串扰，而且压制了一些与检波器一端水体鸣震相关的噪声；同样，与陆检水平分量的对比可看出，分离出的上行S波信噪比也有一定提升.vx分量中与P波视速度接近的同相轴几乎看不到，说明其中主要能量是视速度低的S波.海水声波速度为α0=1.5 km·s-1，密度为ρ0=1.0 g·cm-3，海底纵波速度为α=1.6 km·s-1，横波速度为β=0.3 km·s-1，密度为ρ=1.7 g·cm-3，根据Zoeppritz公式(Aki and Richards, 2002)，计算出上、下行P波、S波在海底界面的反射和透射系数曲线如图 22所示，在当前参数配置情况下，海底界面处纵-横波相互转换并不强烈，同样意味着vx分量基本为上行S波为主，这也一定程度上解释了vx分量和分离出的上行S波记录非常接近的原因.
|
图 21 海底共接收点数据弹性波场分离结果 (a)垂直速度(左)与分离的上行P波记录(右)；(b)水平速度(左)与分离的上行S波记录(右). Fig. 21 Elastic wavefield decomposition of common-receiver gathers at the ocean bottom (a) Vertical particle velocity (left) and the decomposed up-going P-wave recording (right); (b) Horizontal particle velocity (left) and the decomposed up-going S-wave recordings (right). |
|
图 22 海底界面反射/透射系数分析 (a)入射上行P波；(b)入射上行S波；(c)入射下行P波.图中左侧为平面波入射/透射/反射示意图，右侧红线、黑线和绿线分别为上行P波、下行P波以及上行S波的振幅随射线参数的变化系数. Fig. 22 Reflection and transmission coefficients at the seafloor interface (a) Up-going P-wave incidence; (b) up-going S-wave incidence and (c) down-going P-wave incidence. The left row denotes plane wave incidence/transmission/reflection, and the red, black and green lines on the right row denote the reflection and transmission coefficients of various wave modes. |
针对85个OBS分离出的上行P波与S波记录，进行标量声波高斯束叠前深度偏移处理，其中P波速度模型由水检记录初至波走时层析构建，S波速度模型按一定的纵-横波速度比换算得到.图 23为原始声压记录、分离出的上行P波记录对应的偏移剖面.可见后者反射界面层次更清晰，同相轴连续性得到改善，海水鸣震以及一些S波污染的压制对提升成像分辨率与信噪比起到了作用.图 24展示了对分离后上行S波记录进行高斯束偏移获得的PS反射成像剖面.由于上述流程仅能压制检波器端的水层鸣震，所以分离出的上行P波与S波均受到大量残余浅水多次波的干扰，今后需要研究针对性的多次波压制与偏移处理方法.
|
图 23 反射PP波叠前深度偏移结果 (a)基于原始声压分量的PP成像；(b)基于分离的上行P波记录的PP成像. Fig. 23 Prestack depth migration results of reflected PP waves (a) PP image of the raw pressure component; (b) PP image of the separated up-going P-wave recordings. |
|
图 24 基于分离后上行S波记录的反射PS波叠前深度偏移结果 Fig. 24 Prestack depth migration result of the separated up-going S-wave recording |
本文提出了由四个核心步骤组成的上-下行P/S波分离新流程，在至关重要的第一步依据中偏移距折射波时窗内下行波能量最小化准则，取代常规基于小偏移距直达波时窗的P/vz校正方法，扩充了浅水低信噪比海底多分量地震记录波场分离的技术选项.理论合成数据实验检验了方法流程的可行性，参数敏感性分析表明纵波速度对P波分离效果影响最大，横波速度对S波分离效果影响最大，而密度对波场分离影响较小.抗噪性分析实验表明本文所述方法具有较好抗噪性.
分析东海油气田国产地震仪采集的多分量地震数据发现，所有水检直达波信号都受到严重的噪声污染，因此基于中偏移距折射波时窗的P/vz校正方法正好克服常规海底波场分离方法流程的不足.从最终分离的上行P波与S波记录及其叠前深度偏移结果来看，本文技术流程基本达到了海底波场分离的预期效果.虽然它能压制检波端的水体鸣震，但仍然残余的浅水多次波能量直接影响一次反射PP波与PS波的偏移成像质量.今后针对该海底多分量地震资料的精细处理，需要补充其他有针对性的浅水多次波压制方法，发展针对转换波记录的多次波压制方法与技术流程，进而开展更细致的纵、横波速度建模与偏移成像处理.
致谢 感谢国家重点研发计划(2017YFB0202903 & 2018YFC0310100)和国家自然科学基金(41630964，41674117)以及国家科技重大专项(2017ZX05005-004)的资助.感谢于鹏飞博士以及中海油、中石化上海分公司有关人员在YQ探区多分量地震数据采集工作中付出的努力.
Aki K, Richards P G. 2002. Quantitative Seismology. 2nd ed. Sausalito, California: University Science Books.
|
Amundsen L. 1993. Wavenumber-based filtering of marine point-source data. Geophysics, 58(9): 1335-1348. DOI:10.1190/1.1443516 |
Amundsen L, Reitan A. 1995. Estimation of sea-floor wave velocities and density from pressure and particle velocity by AVO analysis. Geophysics, 60(5): 1575-1578. DOI:10.1190/1.1443890 |
Amundsen L, Ikelle L T, Martin J. 2000. Multiple attenuation and P/S splitting of multicomponent OBC data at a heterogeneous sea floor. Wave Motion, 32(1): 67-78. DOI:10.1016/S0165-2125(99)00047-5 |
Cho W H, Spencer T W. 1992. Estimation of polarization and slowness in mixed wavefields. Geophysics, 57(6): 805-814. DOI:10.1190/1.1443294 |
Dragoset W, Barr F J. 1994. Ocean-bottom cable dual-sensor scaling.//64th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 857-860.
|
Edme P, Singh S C. 2009. Receiver function decomposition of OBC data:theory. Geophysical Journal International, 177(3): 966-977. DOI:10.1111/j.1365-246X.2009.04147.x |
Frasier C W. 1970. Discrete time solution of plane P-SV waves in a plane layered medium. Geophysics, 35(2): 197-219. DOI:10.1190/1.1440085 |
Haggerty P E. 1956-07-31. Method and apparatus for canceling reverberations in water layers. US, 2757356.
|
He B S, Zhang H X. 2006. Vector prestack depth migration of multi-component wavefield. Oil Geophysical Prospecting (in Chinese), 41(4): 369-374. |
He B S, Tang P W, Li K R, et al. 2019. Difficulties of multi-component seismic wave reverse-time migration. Periodical of Ocean University of China (in Chinese), 49(1): 77-84. |
Hu T Y, Zhang G J, Zhao W, et al. 2004. Decomposition of multicomponent seismic wavefileds. Chinese Journal of Geophysics (in Chinese), 47(3): 504-508. |
Li W X, Cui J C, Wu W L, et al. 2018. State of the art and applications of the Multi-component seismic processing technology in CNOOC. Chinese Journal of Geophysics (in Chinese), 61(3): 1136-1149. DOI:10.6038/cjg2018L0355 |
Liu L H, Lv C C, Hao T Y, et al. 2012. Data processing methods of OBS and its development tendency in detection of offshore oil and gas resources. Progress in Geophysics (in Chinese), 27(6): 2673-2684. DOI:10.6038/j.issn.1004-2903.2012.06.047 |
Liu Z W, Sa L M, Zhang M, et al. 2008. Progress of application of multi-wave seismic technology in part of gas-fields of China. Oil Geophysical Prospecting (in Chinese), 43(6): 668-672. |
Lu J, Wang Y, Ji Y X, et al. 2018. Imaging techniques of multi-component seismic data. Chinese Journal of Geophysics (in Chinese), 61(8): 3499-3514. DOI:10.6038/cjg2018L0031 |
Muijs R, Robertsson J O A, Curtis A, et al. 2003. Near-surface seismic properties for elastic wavefield decomposition:Estimates based on multicomponent land and seabed recordings. Geophysics, 68(6): 2073-2081. DOI:10.1190/1.1635061 |
Muijs R, Robertsson J O A, Holliger K. 2004. Data-driven adaptive decomposition of multicomponent seabed recordings. Geophysics, 69(5): 1329-1337. DOI:10.1190/1.1801949 |
Muijs R, Robertsson J O A, Holliger K. 2007. Data-driven adaptive decomposition of multicomponent seabed seismic recordings:application to shallow-water data from the North Sea. Geophysics, 72(6): V133-V142. DOI:10.1190/1.2778766 |
Osen A, Amundsen L, Reitan A. 1999. Removal of water-layer multiples from multi-component sea-bottom data. Geophysics, 64(3): 838-851. DOI:10.1190/1.1444594 |
Ren Z M, Liu Y. 2016. A hierarchical elastic full-waveform inversion scheme based on wavefield separation and the multistep-length approach. Geophysics, 81(3): R99-R123. DOI:10.1190/geo2015-0431.1 |
Schalkwijk K M, Wapenaar C P A, Verschuur D J. 2003. Adaptive decomposition of multicomponent ocean-bottom seismic data into downgoing and upgoing P- and S-waves. Geophysics, 68(3): 1091-1102. DOI:10.1190/1.1581081 |
Schalkwijk K W, Wapenaar C P A, Verschuur D J. 1999. Application of two-step decomposition to multicomponent ocean-bottom data:Theory and case study. Journal of Seismic Exploration, 8(3): 261-278. |
Stewart R R, Gaiser J E, Brown R J, et al. 2002. Converted wave seismic exploration:Methods. Geophysics, 67(5): 1348-1363. DOI:10.1190/1.1512781 |
Tao C H, He Q D. 1993. Separation of P- and S-waves in multi-component surface seismic data. Geophysical Prospecting for Petroleum (in Chinese), 32(2): 47-55. |
Wang T F, Cheng J B. 2015. Elastic wave mode decoupling for full waveform inversion.//85th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 1461-1466.
|
Wang T F, Cheng J B. 2017. Elastic full waveform inversion based on mode decomposition:the approach and mechanism. Geophysical Journal International, 209(2): 606-622. DOI:10.1093/gji/ggx038 |
Wang T F. 2017. Elastic full waveform inversion based on wave mode decomposition[Ph. D. thesis] (in Chinese). Shanghai: Tongji University.
|
Wang T F, Cheng J B, Guo Q, et al. 2018. Elastic wave-equation-based reflection kernel analysis and traveltime inversion using wave mode decomposition. Geophysical Journal International, 215(1): 450-470. DOI:10.1093/gji/ggy291 |
Wang Y B, Singh S C, Barton P J. 2002. Separation of P- and SV-wavefields from multi-component seismic data in the τ-p domain. Geophysical Journal International, 151(2): 663-672. DOI:10.1046/j.1365-246X.2002.01797.x |
White J E. 1965. Seismic Waves: Radiation, Transmission, and Attenuation. New York: McGraw-Hill.
|
Xu S Y, Li Y P, Ma Z T. 1999. Seperation of P- and S-wave fields via the τ-q transfrom. China Offshore Oil and Gas (Geology) (in Chinese), 13(5): 334-337. |
Xu W C, Wang T F, Cheng J B. 2019. Elastic model low-to intermediate-wavenumber inversion using reflection traveltime and waveform of multicomponent seismic data. Geophysics, 84(1): R123-R137. |
Yan J, Sava P. 2008. Isotropic angle-domain elastic reverse-time migration. Geophysics, 73(6): S229-S239. DOI:10.1190/1.2981241 |
Zhang Y G, Wang Y, Wang M Y. 2004. Some key problems in the multi-component seismic exploration. Chinese Journal of Geophysics (in Chinese), 47(1): 151-155. |
何兵寿, 张会星. 2006. 多分量波场的矢量法叠前深度偏移技术. 石油地球物理勘探, 41(3): 369-374. |
何兵寿, 唐朋威, 李凯瑞, 等. 2019. 多分量地震资料逆时偏移的难点问题. 中国海洋大学学报, 49(1): 77-84. |
胡天跃, 张广娟, 赵伟, 等. 2004. 多分量地震波波场分解研究. 地球物理学报, 47(3): 504-508. DOI:10.3321/j.issn:0001-5733.2004.03.021 |
李维新, 崔炯成, 武文来, 等. 2018. 中海油多分量地震处理技术的发展与应用实例. 地球物理学报, 61(3): 1136-1149. DOI:10.6038/cjg2018L0355 |
刘丽华, 吕川川, 郝天珧, 等. 2012. 海底地震仪数据处理方法及其在海洋油气资源探测中的发展趋势. 地球物理学进展, 27(6): 2673-2684. DOI:10.6038/j.issn.1004-2903.2012.06.047 |
刘振武, 撒利明, 张明, 等. 2008. 多波地震技术在中国部分气田的应用和进展. 石油地球物理勘探, 43(6): 668-672. DOI:10.3321/j.issn:1000-7210.2008.06.010 |
芦俊, 王赟, 季玉新, 等. 2018. 多分量地震数据的成像技术. 地球物理学报, 61(8): 3499-3514. DOI:10.6038/cjg2018L0031 |
陶春辉, 何樵登. 1993. 地面多分量地震资料纵-横波分离方法. 石油物探, 32(2): 47-55. |
王腾飞. 2017.弹性波模式解耦全波形反演方法[博士论文].上海: 同济大学.
|
许世勇, 李彦鹏, 马在田. 1999. τ-q变换法波场分离. 中国海上油气(地质), 13(5): 334-337. |
张永刚, 王赟, 王妙月. 2004. 目前多分量地震勘探中的几个关键问题. 地球物理学报, 47(1): 151-155. DOI:10.3321/j.issn:0001-5733.2004.01.023 |