检索
E-mail
RSS

作者投稿 主编审稿 专家审稿 远程编辑

热带海洋学报 >

2020 , Vol. 39 >Issue 4: 80 - 90

DOI: https://doi.org/10.11978/2019087

海洋地球物理学

Kirchhoff向下延拓法在海洋合成多道地震走时反演中的应用*

  • 张茂传 , 1, 2, 3 ,
  • 徐敏 , 1, 3 ,
  • 赵旭 1, 2, 3 ,
  • 张佳政 1, 3 ,
  • 查财财 1, 2, 3 ,
  • VithanaM.V.P. 1, 2 ,
  • 狄会哲 1, 2, 3 ,
  • 曾信 1, 3
展开
  • 1.中国科学院边缘海与大洋地质重点实验室, 南海海洋研究所, 南海生态环境工程创新研究院, 广东 广州 510301
  • 2.中国科学院大学, 北京 100049
  • 3.南方海洋科学与工程广东省实验室(广州), 广东 广州 511458
徐敏。E-mail:

张茂传(1994—), 男, 江西赣州人, 硕士研究生, 主要从事海洋地球物理研究。E-mail:

Copy editor: 姚衍桃

收稿日期: 2019-09-16

  要求修回日期: 2020-01-17

  网络出版日期: 2020-07-27

基金资助

国家自然科学基金项目(41676044)

国家自然科学基金项目(91858207)

国家自然科学基金项目(41890813)

中国科学院战略性先导科技专项(XDA13010105)

南方海洋科学与工程广东省实验室(广州)人才团队引进重大专项(GML2019ZD0205)

版权

版权所有，未经授权，不得转载、摘编本刊文章，不得使用本刊的版式设计。
收起

Application of Kirchhoff downward-continued method to travel-time inversion of synthetic multi-channel seismic data

  • ZHANG Maochuan , 1, 2, 3 ,
  • XU Min , 1, 3 ,
  • ZHAO Xu 1, 2, 3 ,
  • ZHANG Jiazheng 1, 3 ,
  • ZHA Caicai 1, 2, 3 ,
  • Vithana M. V. P. 1, 2 ,
  • DI Huizhe 1, 2, 3 ,
  • ZENG Xin 1, 3
Expand
  • 1. Key Laboratory of Ocean and Marginal Sea Geology, South China Sea Institute of Oceanology, Innovation Academy of South China Sea Ecology and Environmental Engineering, Chinese Academy of Sciences, Guangzhou 510301, China
  • 2. University of Chinese Academy of Sciences, Beijing 100049, China
  • 3. Southern Marine Science and Engineering Guangdong Laboratory (Guangzhou), Guangzhou 511548, China
XU Min. E-mail:

Received date: 2019-09-16

  Request revised date: 2020-01-17

  Online published: 2020-07-27

Supported by

Foundation item: National Natural Science Foundation of China(41676044)

Foundation item: National Natural Science Foundation of China(91858207)

Foundation item: National Natural Science Foundation of China(41890813)

The Strategic Priority Research Program of the Chinese Academy of Sciences(XDA13010105)

Key Special Project for Introduced Talents Team of Southern Marine Science and Engineering Guangdong Laboratory(Guangzhou)(GML2019ZD0205)

Copyright

Copyright reserved © 2020. Office of Acta Agronomica Sinica All articles published represent the opinions of the authors, and do not reflect the official policy of the Chinese Medical Association or the Editorial Board, unless this is clearly specified.
Fold

摘要

海洋多道地震数据建模和成像是获取洋壳速度和构造信息的重要手段。海水层的存在使得多道地震拖缆接收到的折射波走时信息仅仅存在于较远的炮检距, 近炮检距被强振幅海底反射波覆盖, 制约了走时数据拾取和反演效果。本文基于波动方程的Kirchhoff积分法, 成功实现了多道地震数据向下延拓, 获取到了更大炮检距区间的初至折射波走时拾取, 并将其应用于洋中脊新生洋壳2A/2B层的合成多道地震数据走时反演。比较向下延拓前后的走时拾取范围及走时反演结果表明, 向下延拓法能够保持地震波场的运动学和动力学特征不变, 在共炮集数据的更大炮检距范围内进行初至折射走时拾取, 从而增加反演的数据选择和浅层射线覆盖, 反演结果能更加准确地分辨出洋壳2A/2B层界面, 并得到更高的分辨率和更准确的速度结构剖面。

关键词： 海洋多道地震; Kirchhoff积分; 向下延拓; 走时反演; 洋壳2A/2B层

本文引用格式

张茂传 , 徐敏 , 赵旭 , 张佳政 , 查财财 , VithanaM.V.P. , 狄会哲 , 曾信 . Kirchhoff向下延拓法在海洋合成多道地震走时反演中的应用*[J]. 热带海洋学报, 2020 , 39(4) : 80 -90 . DOI: 10.11978/2019087

Abstract

Velocity model building and seismic imaging from marine multi-channel seismic (MCS) data are essential to obtain the structural information of oceanic crust. Existence of the seawater layer makes first-arrival refraction phases recorded by streamer only existed in far offset, while most of them in near offset are masked by the strong amplitude seafloor reflection phase, affecting travel-time picking and inversion. Based on the Kirchhoff integral approach to wave equation, we successfully implement the downward continuation (DC) method to the MCS data, and obtain first-arrival travel-time picks within a wider offset range. The method is applied to the synthetic MCS data simulated based on layer 2A/2B structure of newborn crust along a mid-ocean ridge, and first-arrival picking range and inversion result are compared before and after the DC. The DC method can expand offset range for first-arrival refraction picking and increase the ray coverage in the shallow depth, without changing the kinematics and dynamic characteristics of seismic wave field. The inversion results can more precisely distinguish the layer 2A/2B interface, and obtain higher resolution and more accurate velocity profile.

Key words： marine multi-channel seismic; Kirchhoff integration; downward continuation; travel-time inversion; oceanic crust layer 2A/2B

*感谢J. P. Canales教授为本论文提供的意见和帮助; 感谢主编和三位审稿人对论文给出的建设性意见。
海洋多道地震(Multi-Channel Seismic, MCS)实验采集的高精度数据是获得海底地层速度和构造信息的重要来源(Kent et al, 2000; Singh et al, 2006)。MCS数据的采集方式不同于陆地, 科考船在海面上拖着气枪震源和水听器电缆航行, 由气枪激发地震波, 拖缆检波器接收返回的地震信号。无论是地震波到达海底地层还是返回检波器, 都需要穿过海水层。除去海水层内部传播的信号, 在近偏移距地震数据中首先记录到强振幅的海底反射波往往掩盖了来自岩层的有效折射波, 同时还会记录到一些特殊的干扰波, 如鸣震和交混回响, 以及与海底有关的衍射和散射干扰(Larkin et al, 1996)。尤其是在深海区(水深>2000m)采集的多道地震记录, 用于反演速度结构的初至折射走时信息仅仅存在于远炮检距, 有时甚至完全被海底反射覆盖, 造成有效数据的缺乏和速度结构反演偏差。因此, 消除海水层对地震波场所造成的影响, 挖掘地震波的有效折射信息, 更好地进行多道地震初至折射波走时拾取和速度结构反演, 进而获得高分辨率的海底地层速度结构信息, 是进行海洋地质构造调查和海洋油气勘探的关键性技术。
基于波动方程的向下延拓法(Downward Continuation, DC)是根据地震波在传播过程中的运动学(如时深关系)和动力学(如振幅、相位等)特征, 将波场向下延拓到实际的地质界面或虚拟的水平界面(Berryhill, 1979; Arnulf et al, 2012)。该方法应用于海洋主动源多道地震数据, 从近海面向下延拓到海底面或近海底面, 相当于在海底激发和接收, 使得洋壳折射震相走时早于海底反射。从而消除海水层和崎岖海底地形对地震记录造成的畸变, 使得地壳的折射波成为初至波并出现在更大范围的炮检距中, 提供更多的射线覆盖和数据选择, 进而反演得到高分辨率的速度结构。
自从提出用Kirchhoff积分法进行波场延拓以来, 为了得到更好的地震成像效果, 国内外很多学者将此方法应用在了各种地震数据中进行静校正和偏移处理, 比如海底地震仪(Ocean Bottom Seismometer, OBS) (王祥春 等, 2012)、海底电缆(Ocean Bottom Cable, OBC) (金丹 等, 2011)、模拟地震数据(耿建华 等, 1996)和复杂地表(Yilmaz et al, 1986)等地震数据中。近几年, 随着走时层析成像和全波形反演等地震速度结构反演方法以及计算机性能的迅速发展, 向下延拓法能在保持地震波动力学和运动学特征不变的情况下, 重新定义基准面, 挖掘地震数据中更多有效信息的作用变得越来越重要。在海洋多道地震资料中, Harding等(2007)首次提出在多道地震数据中实现向下延拓, 并在一系列实际运用过程中不断完善(Arnulf et al., 2011, 2012, 2014a, b), 现在已成为多道地震初至走时反演乃至常规多道地震数据处理中极其重要的一步。
本文将回顾Kirchhoff积分向下延拓法的原理, 并应用到基于洋壳2A/2B层速度模型合成的多道地震数据中, 比较延拓前与延拓后地震数据有效初至折射波的炮检距拾取范围, 通过走时层析成像得到速度反演结果, 建立多道地震数据向下延拓方法的处理流程。

1 向下延拓法原理和海洋多道数据处理流程

向下延拓法是由Berryhill(1979)首次提出的波动方程基准面延拓(Wave Equation Datuming)思想发展而来的, 他利用波动方程Kirchhoff积分法建立了2D基准面向下延拓算法, 该方法消除了由于海底地形起伏引起的地震数据畸变。随后, 很多学者相继发展了各种波动方程的基准面延拓算法, 比如有限差分法(McMechan et al, 1990)、共焦点法(Common Focal Point, CFP)(Bolte, 2001)、叠前相移法(Reshef, 1991)以及逆时偏移法(Whitmore, 1983)等。Bevc(1995)用多种波场延拓方法建立了延拓算子, 比较了各种方法的效果, 认为Kirchhoff积分延拓能保证足够的精度, 计算效率很高, 成像效果最好。此外, 由于其本身的灵活性, 输入和输出的基准面可以是非线性平面(Wiggins, 1984), 能够适应横向和纵向的速度变化(Shtivelman et al, 1988), 因而很容易推广到3D数据处理的应用中, 这对于缺乏地质构造和速度结构等先验信息的地震数据处理是非常有帮助的(Bevc, 1997)。

1.1 基本原理

齐次声波方程U在一个封闭的区域$\Omega $内, 满足以下方程:
$\frac{{{\partial }^{2}}U}{\partial {{x}^{2}}}+\frac{{{\partial }^{2}}U}{\partial {{y}^{2}}}+\frac{{{\partial }^{2}}U}{\partial {{z}^{2}}}=\frac{1}{{{v}^{2}}}\frac{{{\partial }^{2}}U}{\partial {{t}^{2}}}$
为求解$\Omega $内任意一点的波场值, 可将上式代入Green公式:
$\iiint_{\Omega }{\left( U\Delta {R}'-{R}'\Delta U \right)}d\Omega =\mathop{{\int\!\!\!\!\!\int}\mkern-21mu \bigcirc}\nolimits_S\left( \frac{U\Delta {R}'}{\partial n}-\frac{{R}'\Delta U}{\partial n} \right) dS$
从而得到Kirchhoff积分解:
$\begin{align}& U\left( {{x}_{0}},{{y}_{0}},{{z}_{0}},t \right)=\frac{1}{4\pi }\iint_{S}{\left[ \frac{1}{R}\frac{\partial }{\partial n}-\frac{\partial \left( 1/R \right)}{\partial n}+\frac{1}{vR}\frac{\partial R}{\partial n}\frac{\partial }{\partial t} \right]}\times \\& \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ U\left( x,y,z,t-\frac{v}{R} \right)dS \\\end{align}$
式中, S为$\Omega $的封闭曲面; $R=$ $\sqrt{{{\left( x-{{x}_{0}} \right)}^{2}}+{{\left( y-{{y}_{0}} \right)}^{2}}+{{\left( z-{{z}_{0}} \right)}^{2}}}$, 并且满足Laplace方程。由式(3)可知, 对封闭曲面上波场的积分可用来表示封闭区域内任意一点的波场值, 这样就可应用式(3)进行地震波的波场成像, 也就是地下任一点的波场可用地面观测到的波场积分来表示。但由于既不可能在封闭曲面上进行实际数据采集, 也不可能获取到式(3)中所需要的一个封闭曲面上的所有$\partial U/\partial n$数据, 所以用有限的观测数据求解式(3)是不适定的。为此必须增加一定的约束条件, 以获取式(3)的适定解。假定在半无限空间上进行地震探测, 再根据半无限空间的对称性消除$\partial U/\partial n$项, 可以得到半无限空间上的Kirchhoff积分解:
$\begin{align}& U\left( {{x}_{0}},{{y}_{0}},{{z}_{0}},t \right)=\frac{1}{2\pi }\iint_{S}{\left[ \frac{\partial \left( 1/R \right)}{\partial n}+\frac{1}{vR}\frac{\partial R}{\partial n}\frac{\partial }{\partial t} \right]}\times \\& \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ U\left( x,y,z,t-\frac{v}{R} \right)dS \\\end{align}$
式(4)是三维Kirchhoff积分解, 在实际多道地震数据测线中可以看作是二维剖面。引入单位阶跃函数, 将在空间y方向上的积分替换为在时间t方向上的卷积, 经过一系列推导, 使其变换成二维Kirchhoff积分解。最后为了便于在多道地震数据中通过程序实现, 简化成如图1的二维基准面延拓示意图, 其对应的离散型延拓公式见式(5)。
图1 波动方程基准面延拓图示[修改自Berryhill (1984)]

观测面上的倒三角表示检波器, 延拓面上的圆点表示向下延拓后的检波器位置

Fig. 1 Wave equation datuming, modified from Berryhill (1984).

The inverted triangles in the observation plane stand for hydrophones, and the hydrophones’ location after DC is indicated by dots

${{U}_{j}}\left( t \right)=\frac{1}{\pi }\sum\limits_{i}{\Delta {{x}_{i}}\cos \left( {{\theta }_{ij}} \right)}\frac{{{t}_{i}}}{{{r}_{ij}}}\left[ {{U}_{i}}\left( t-{{t}_{i}} \right)*{{F}_{i}} \right]$
式中, ${{U}_{j}}\left( t \right)$为延拓目标面j位置的波场; ${{U}_{i}}\left( t \right)$为初始输入的波场记录, (t–ti)表示对输入地震道有一个ti走时的延迟, 因此${{U}_{i}}\left( t-{{t}_{i}} \right)$表示从观测面传播的上行波向下延拓; $\Delta {{x}_{i}}$为道间距; rij为i、j两点间的距离; ${{\theta }_{ij}}$为i、j两点连线与原始观测面的法向夹角; ti为i、j两点间的走时; *为卷积符号; Fi是预测反褶积滤波因子, 对于波形和振幅的保持很重要。
Berryhill(1979)基准面延拓思想的关键在于用恒定的海水速度收敛地震波在水中传播的过程。波场延拓是将波场从一个容易观察得到数据的基准面延拓到一个新的难以得到数据但有地质目标的基准面, 如从海平面向下延拓到海底面。Kirchhoff积分法的准确实现需要确切了解观测面与延拓面之间的速度信息, 以及没有空间噪音和子波变形的良好地震数据(Larkin et al, 1996)。对于海洋多道地震数据, 海水层的横向和纵向速度基本不变, 可以视为常数, 利用时深关系将波场反向从海面向下延拓到海底, 从而消除海水层带来的影响。

1.2 海洋多道地震数据向下延拓实现流程

在海洋多道地震实验中, 气枪源激发的地震波穿过海水层, 一部分穿透地层带回有效的初至折射波走时信息, 另一部分直接反射回来被检波器接收。由惠更斯原理可知, 经过海底返回的波场将以球面波的形式扩散, 并传播到检波器中被部分接收(图2a), 在海水面上观测得到的气枪震源多道地震数据可以表示为$U\left( x,z=0,t \right)$。向下延拓的作用是利用接收到的所有波场, 去除地震波在海水层传播带来的影响, 将经过水层的波场收敛至目标界面的某个点上, 使能够反映出海底地层特征的相对高速有效折射波从强振幅的海底慢速反射干扰波中优先被记录下来, 并成为初至波(图2b)。通过波场向下延拓, 可以从观测记录$U\left( x,z=0,t \right)$得到新的剖面$U\left( x,z={{z}_{1}},t \right)$, 相当于震源和检波器都是在海底或近海底, 使得反映地层信息的有效折射波在大部分炮检距中早于反射波并被记录在检波器中。
图2 海洋多道地震向下延拓前(a)和延拓后(b)的示意图

红色五角星为气枪震源; 观测面上的黑色倒三角表示检波器, 延拓面上的圆点表示延拓后检波器的位置; 黄色线条为水听器拖缆;图a中绿色、红色射线分别示意第一道和最后一道的收敛情况; 图b中蓝色射线表示最后的收敛情况

Fig. 2 Schematic diagrams of MCS data acquisition before (a) and after (b) DC.

Red stars represent the air-gun sources, inverted triangles in observation plane present hydrophones, and the hydrophones’ locations after DC are indicated by dots. Marine streamer is shown by the yellow line; green, red and blue rays indicate the convergence of the first and last channels and the final convergence, respectively

在海洋多道地震数据向下延拓的实现流程中(图3), 第一步是将共炮点道集的波场向下延拓到实际的地质界面或虚拟的界面。由于Kirchhoff积分法操作的灵活性, 可以延拓到近海底的平面或者与海底起伏一致的曲面(Wiggins, 1984)。在延拓前, 由导航数据和水深数据建立观测面和目标延拓面的网格模型, 使每道地震信号在延拓前后的炮点和接收点的位置都落在网格之中。为了得到更好的效果, 将原始地震数据进行道均衡、道编辑、预测反褶积和带通滤波等常规处理, 并去除因拖缆长度和模型大小有限带来的盲区效应。波场向下延拓最关键的参数是波的传播速度场, 准确的速度场可以保证射线追踪得到准确的走时, 以及波场向下延拓过程中同相轴不会偏离真实位置。若所选速度低于真实速度, 会导致绕射双曲线收敛较差, 而高于真实速度则会使地下的倾斜同相轴与真实位置有较大偏差。海洋多道地震实验在海水中进行, 海水层具有近乎均匀且没有横向变化的波传播速度, 因此延拓后能够准确反映波场的真实传播过程。
第二步根据互易定理, 将地震数据转换到共接收点道集。为了避免由于炮间距和道间距之间的差异而产生波形交叉混叠, 在延拓的过程中根据前后两道的信号对共接收点道集进行线性插值。为了保证延拓最终效果, 舍去少于100道的共接收点道集, 再滤波、道编辑和去除盲区效应处理, 进行共接收点道集延拓。完全延拓后的地震数据再次转换到共炮点道集, 至此即完成向下延拓的所有过程。
总之, 海洋气枪源多道地震向下延拓法的实现只需要多道地震数据、测线导航文件、水深数据及水体速度结构(通常可假设为1.5km·s-1), 就能够将在海面激发和接收得到的波场转换为相当于海底激发和接收得到的波场, 从而使清晰的初至折射震相能够在大部分的炮检距中被拾取, 并用于后续走时成像。

2 向下延拓法的应用

向下延拓法在走时层析成像中被广泛应用, 如Arnulf等(2011, 2012)利用向下延拓法处理大西洋中脊Lucky Strike火山活动区的3D多道地震数据, 结合走时层析成像, 得到了高分辨率的洋壳速度结构, 从而对火山结构和热液过程有了新的认识; 而且结合OBS-MCS联合数据的走时层析成像技术, 还原了胡安德富卡洋脊热点影响的脊轴海山的速度结构、地震活动性和增生过程(Arnulf et al, 2014a)。此外, Henig等(2012)利用向下延拓法处理大西洋中脊Atlantis Massif大洋核杂岩的多道地震数据, 结合走时层析成像, 得到了大洋核杂岩明显的横向和纵向变化的速度结构, 该结果有助于理解慢速扩张洋中脊非常不均一的岩石圈形成过程。
除了在走时层析成像中的应用, 向下延拓法在全波形反演和逆时偏移以及综合反演中也被广泛应用。Arnulf等(2012, 2013)利用向下延拓法, 结合全波形反演, 得到了该区域上层洋壳的精细速度模型; 之后结合逆时偏移地震成像, 深度解剖了活动海底火山之下从岩浆房到喷发点的岩浆运移路径(Arnulf et al, 2014b); Ghosal等(2014)利用向下延拓法, 结合走时层析成像和叠前深度偏移地震成像技术, 高分辨率地描绘了苏门答腊西北海岸俯冲带的沉积序列和沉积物结构性质的变化; Qin等(2018)利用向下延拓法和走时层析成像、全波形反演以及逆时偏移技术相结合, 对苏门答腊俯冲带的增生楔和板块界面性质进行了高精度的反演成像; Harding等(2016)利用向下延拓和全波形反演方法相结合, 在Atlantis Massif的IODP钻孔U1309D处进行了高精度的速度结构成像; Marjanović等(2017)和Vithana等(2019)也利用向下延拓和全波形反演方法相结合, 对东太平洋海隆9°16′—9°56′N洋脊段区域的海底热液系统的流体通道进行了精细的地震成像。以上向下延拓法的成功应用表明, 该方法正在不断走向成熟, 并在与走时层析成像、全波形反演、逆时偏移地震成像等技术相结合的海洋多道地震数据处理上发挥着重要作用。本文将利用向下延拓法与走时反演成像技术对洋壳2A/2B层结构进行反演测试, 并与传统未延拓的处理方法作对比。

2.1 2A和2B层分界

用地震速度结构确定海洋岩石圈的分层结构是一种普遍做法。根据地震波速度, 通常可以将快速扩张洋中脊增生的洋壳分成3层: 层1为沉积层, 层2为玄武岩和辉绿岩, 层3为辉长岩(图4a)。对于年轻的洋壳, 层1一般缺失, 层2通常分为2A和2B两个子层, 速度分别介于2.2~5.0kms-1和5.0~6.0kms-1。其中, 枕状喷出岩层2A为玄武岩, 侵入的席状岩墙层2B为辉绿岩, 并具有较大的速度梯度(>3s-1), 从而产生强烈的地震反射图像(Harding et al, 1993; Vera et al, 1994)。由于利用多道地震折射波成像的深度有限, 本文主要聚焦在图4a红线表示的2A和2B层深度部分, 模拟在新生洋壳区域进行多道地震实验, 并利用走时成像的方法得到向下延拓前后的地震速度模型。
图4 典型快速扩张洋脊新生洋壳Penrose分层模型(a)[修改自Vithana等(2019)]和原始共炮点道集射线(b)、第一步共炮点道集延拓射线(c)、第二步共接收点道集延拓射线(d)示意图

红色五角星表示炮点位置; 绿色线段表示初至折射波拾取的范围; 白色射线表示有效折射波, 红色射线表示海底反射波; 黑色射线表示延拓之后远偏移距中不可靠的初至波, 对应水面缆记录不到的远偏移距震相

Fig. 4 (a) Layered Penrose model of young oceanic formed along a fast-spreading ridge, modified from Vithana et al (2019); schematic ray paths of (b) original common shot gather data, (c) first-step DC of common shot gather data, and (d) second-step DC of common receiver gather data.

Red star indicates the shot location, and green lines indicate the available picking range of first-arrival. Effective refractive wave, reflective wave and the unreliable first arrivals after DC are marked by white, red and black ray paths, respectively. The black rays correspond to the far offset seismic phases that cannot be recorded at the sea level

图4b所示白色射线是来自2A/2B层的有效折射波, 红色射线为海底反射波, 绿色线段表示折射波比海底反射波早到的区间, 也就是初至折射波的拾取范围。在原始多道地震数据中, 拾取的初至折射波走时信息仅仅位于较远的炮检距, 主要来自深度和速度相对较大的2B层, 相对较浅且低速的2A层信息因被海底反射波覆盖而无法拾取, 合成的原始共炮点道集如图5a所示。经过共炮点道集和共接收点道集延拓, 原始多道地震转换为在海底激发和海底接收的数据。根据图4c、d所示, 延拓后绿色线段表示的初至折射波拾取范围得到增加, 覆盖了包括来自2A层和2B层的所有折射信息; 而图4c、d中的黑色射线表示延拓之后远偏移距中不可靠的初至波, 对应水面缆记录不到的远偏移距震相, 亦即图5b、d中远偏移距部分(400道之后)的折射波。
图5 实际延拓后数据与合成地震图的共炮点道集与抽取数据道对比

a. 原始共炮点道集; b. 共炮点道集延拓结果; c. 共炮点道集延拓后的合成地震图; d. 共接收点道集延拓结果; e. 共接收点道集延拓后的合成地震图; f. 抽取道的数据对比。图a—e中的红虚线和绿虚线分别表示2A层和2B层的震相, 白色虚线代表抽取的道位置; 图f中的绿线和红线分别代表实际延拓后数据和相应的合成地震图数据

Fig. 5 Comparison of seismic shot gather and traces from DC and synthetic seismograms.

(a) Original common shot gather; (b) DC result of common shot gather; (c) synthetic seismogram of DC common shot gather; (d) DC result of common receiver gather; (e) synthetic seismogram of DC common receive gather; and (f) comparison of selected traces from DC result and synthetic seismogram. Seismic phases of layer 2A and 2B in (a~e) are indicated by red and green dash lines, respectively. Vertical white dashed lines indicate the locations of selected traces, and the enlarged real and synthetic traces are marked by green and red dashed lines in (f), respectively

2.2 合成多道地震数据处理

海洋多道地震由于探测深度和目的不同而选用不同长度的水听器拖缆, 本文模拟的是6km缆长多道地震数据, 记录道数480道, 记录长度6s, 数据采样率4ms, 道间距为12.5m, 炮间距为37.5m, 最小炮检距200m。合成地震记录采用的一维速度模型如图4a所示, 速度模型长25km, 深8km, 海水深度2km, 速度设为1.5kms-1, 海底下2km为洋壳2A和2B层(图4a)。将海面激发接收的多道地震数据延拓到深度2km处的海底, 用主频为10Hz的雷克子波模拟地震波传播, 利用有限差分方法程序模拟地震波场传播, 并生成合成记录。文中所采用的地震波传播模拟方法适用于任何复杂介质模型, 包括海底及地层内部的固液边界地震波传播模拟。所用程序基于美国伍兹霍尔海洋研究所海洋声学实验室公开程序修改(Xu et al, 2017), 具体方法参见Xu(2012) (https://dspace.mit.edu/handle/1721.1/70780)。
基于真实一维的2A/2B层速度模型(图4a)合成的原始共炮点道集如图5a所示, 图中绿虚线表示的可拾取初至走时信息为来自于2B层的折射波, 它仅仅存在于较远端的炮检距中, 其余大部分折射波则被海底反射波所掩盖。
为了更好地展示向下延拓处理的效果, 本文将实际延拓结果与基于延拓后的炮检位置的合成地震图进行对比。图5b是共炮点道集延拓的结果, 与之对应的合成结果是图5c, 两图分别表示炮点在水面、检波点在海底的延拓与合成的共炮点道集; 图5d是共检波点延拓之后的结果, 与之对应的合成结果是图5e, 两图分别表示炮点和检波点都在海底的延拓与合成的共炮点道集。延拓之前每炮的初至拾取只有184道(图5a), 延拓之后去除受盲区效应影响的远端偏移距部分, 每炮的初至拾取可达到354道(图5d)。图5f对比了从图5d、e中抽取的同样位置的地震道数据, 可以发现实际延拓数据和合成地震数据之间的走时、振幅、波形和相位特征在近偏移距是十分吻合的, 但延拓之后远偏移距(400道之后)中的振幅信息不再可信。

3 初至折射波走时反演

3.1 海洋多道数据走时反演

海洋多道数据走时反演(走时层析成像)利用多道拖缆记录的初至折射波震相的走时信息反演得到地层速度, 是研究洋壳精细结构的重要途径。本文利用FAST(First-Arrival Seismic Tomography)程序(Zelt et al, 1998)对延拓前后的初至折射波走时数据进行反演, 程序中正反演网格的横向、纵向大小均为0.05km。FAST软件使用有限差分方法求解程函方程, 进而求得初至波时间。其主要特点是针对高分辨地震折射走时资料的高分辨和高精度, 正演模型采用尺度很小的格点插值的参数化方法, 有利于小尺度异常体的分辨。正则化反演算法不仅考虑了数据对模型的约束, 也兼顾了程函方程是波动方程高频近似的简化影响, 即模型是光滑的。但值得注意的是, FAST只能处理初至折射波, 并不具备处理反射震相和后至折射波的功能, 所以无法对不连续界面有任何深度约束。同时, FAST模型不能设置速度跳跃界面, 只能用大梯度薄层近似表示。走时反演初始速度模型设为具有固定速度梯度0.8s-1, 速度从海底的5.2kms-1递增到地下6km处的8kms-1。利用FAST迭代4次之后, 延拓前后反演的走时方差分别由初始的297ms和193ms降为22ms和14ms, 其走时拟合如图6所示。
图6 延拓前(a)和延拓后(b)的走时拟合

Fig. 6 Travel-time curves of observed and uncertainty (red dots with error bars), initial model (blue line) and final model (green line) before (a) and after (b) DC

3.2 向下延拓前后反演结果的对比分析

走时反演结果如图7所示, 其中图7a是生成原始多道反射地震的初始模型, 图7b和图7d是延拓前后FAST走时反演的二维结果对比。由于原始反演射线路径覆盖较密, 在图中显示成每1条射线代表50条射线。可以发现, 图7b中并没有2A层的射线路径, 而延拓之后的图7d中2A层也有射线覆盖, 使得反演得到的速度更加准确可靠。合成多道数据从模型右边开始激发, 为了避免盲区效应, 在模型的两边空出适当距离。将直接拾取的共炮点道集初至折射波进行走时反演, 得到向下延拓前速度反演结果如图7b所示, 对应的一维速度模型为图7c中的蓝线。基于向下延拓后拾取的初至折射波进行走时反演, 得到向下延拓后速度反演结果如图7d所示, 对应的一维速度模型为图7c中的红线。对比分析延拓前后的反演结果以及一维速度与图7c中黑线表示的真实速度可知, 延拓后的走时反演得到的速度结果更加接近于真实速度模型, 尤其是在速度变化较大的水层与地层起跳点的位置。相比于延拓前的反演结果, 向下延拓方法的应用成功地反演出了2A层的速度结构, 使得多道反射地震数据的走时反演具有浅层高分辨的效果。
图7 合成多道地震数据的真实二维速度模型(a)、向下延拓前走时反演的二维速度模型(b)、一维速度模型对比(c)和向下延拓后走时反演的二维速度模型(d)

图a、b和d中的白色虚线分别表示图c中黑线、蓝线和红线一维速度模型所在位置

Fig. 7 (a) Real velocity model for original synthetic shot gathers; (b) 2D tomography velocity model before DC; (c) 1D velocity comparison of real model, initial model and tomography model before and after DC; and (d) 2D tomography velocity model after DC.

The dotted white lines in (a), (b) and (d) correspond to the one-dimensional velocity models of black, blue and red lines, in (c), respectively

由于延拓后可拾取的初至折射走时数据增加, 即浅层射线覆盖变密, 使得向下延拓之后的反演速度结果比延拓之前更加接近真实模型, 并能够更准确地区分2A和2B层, 从而真实地反映了浅层洋壳的速度结构。因此, 海洋多道地震向下延拓结合初至折射波走时层析成像对了解洋壳的浅层地质构造、沉积层厚度以及2A/2B层速度变化具有明显的促进作用。

4 分析和讨论

合成海洋多道地震数据原始共炮点道集的初至折射波震相只存在于较远的炮检距中, 大部分近炮检距的初至折射波被强振幅的海底反射波所覆盖(图5a), 导致来自浅层的折射波无法被拾取并用于速度结构的走时反演。经过向下延拓处理之后, 共炮点道集的初至走时信息能够在大范围的炮检距内被拾取(图5d)。本文6km缆长的合成多道实验中, 在单个共炮点道集上拾取的初至折射波走时由延拓前的184道增加到延拓后的354道。与相应的合成道集(图5e)相比, 除去受盲区效应影响的远偏移距部分, 延拓后大部分道集的走时、振幅、波形和相位都互相吻合, 说明此延拓方法保持了地震波在传播过程中的运动学和动力学特征。
经过延拓后的共炮点道集数据如同在海底激发和接收的地震实验, 而实际地震实验往往很难做到海底激发和接收。为了获取洋壳浅层高分辨率的速度信息, 早期美国科学家尝试在近海底开展爆炸源地震实验(Near Ocean Bottom Explosives Launcher, NOBEL)(Christeson et al, 1994), 通过在海底布放检波器和从海面释放炸药并在近海底定时爆炸的方式来实现。这种地震作业方式非常不经济和危险, 同
时容易造成海洋环境污染, 所以后来没有得到进一步发展。事实上, 这种方式也难以满足水平方向高分辨率成像的要求, 因为海底接收器间距过大, 约为1km(Collins et al, 2009)。多道地震向下延拓近似于合成海底地震实验, 走时拾取的数据密度和连续性比传统海底地震仪实验拾取的射线大1个数量级(Harding et al, 2007), 增加了在速度结构反演中穿过每个网格的射线覆盖, 提高了横向和纵向的分辨率, 反演结果也能够更加准确地反映地层情况。
延拓的目标参考面可以是具有不规则形态的实际地质界面, 也可以是虚拟的水平界面(Wiggins, 1984), 延拓方法的有效性依赖于延拓目标参考面的深度。本文将海洋多道数据延拓至近海底的平面后, 地震波场不再经由海水层介质扩散传播(图4d), 收敛了来自海底的散射和衍射(Arnulf et al, 2011), 提高了数据的信噪比。从合成数据的延拓结果中可以看出, 海底反射波收敛到近零炮检距, 来自浅层的折射信息被挖掘出来, 使得洋壳折射成为初至波并可以在更大的炮检距范围内被拾取。
向下延拓法的实现基于波动方程的Kirchhoff积分法, 通过离散化的数学积分方法运用到地震资料上。从波动方程的积分解公式(4)中看出, 解析解需要对整个半无限空间S面进行积分求解; 转换到离散型的延拓公式(5)中, 也需要两点连线与原始观测面的法向夹角(θij)满足半无限空间的条件。但由于缆长限制, 共炮点道集接收到的波场范围有限, 其最终结果会产生盲区效应。延拓之后远炮检距的数据往往会出现“画弧”现象(图5d), 导致振幅大小不一致(图5f), 从而变得不可靠, 故在实际应用中应当舍去(Wang et al, 2014)。同时, 因为测线前后的共接收点道集并不是满覆盖的, 延拓过程中需要舍去一些位于模型边缘的道集, 因而会造成数据部分减少。延拓过程对数据质量的要求也很高, 因为延拓后每一个有效信号和噪音都会被放大, 这就要求在延拓过程中需要精细处理好每一步的数据。

5 结论

将Kirchhoff向下延拓法应用到海洋多道地震数据处理中, 能够将海面激发和接收得到的数据延拓到海底, 相当于海底激发和海底接收。该方法收敛了来自海底的衍射、散射和反射波, 使得来自浅层洋壳的折射波早于反射波到达检波器成为初至走时, 并可以在大部分的炮检距中被拾取。从而增加了走时层析成像所需的射线覆盖和数据选择, 提高了走时反演的横向和纵向分辨率, 对详细了解浅层洋壳结构具有重大作用。基于波动方程偏移的向下延拓结果保持了地震数据原有的运动学和动力学特征, 为后续的全波形反演也提供了可靠的波形, 进而能够得到更加精细的速度结构和地震成像。

参考文献

原文顺序 | 文献年度倒序 | 文中引用次数倒序
[1]
耿建华, 黄海贵, 马在田, 1996. Kirchhoff积分波场延拓基准面静校正方法研究[J]. 同济大学学报, 24(6):665-669.

GENG JIANHUA, HUANG HAIGUI, MA ZAITIAN, 1996. Kirchhoff wavefield extrapolation datuming[J]. Journal of Tongji University, 24(6):665-669 (in Chinese with English abstract).

[2]
金丹, 阎贫, 唐群署, 等, 2011. Kirchhoff波场延拓在OBC记录海水层基准面校正中的应用[J]. 热带海洋学报, 30(6):84-89.

JIN DAN, YAN PIN, TANG QUNSHU, et al, 2011. Application of Kirchhoff integral wave field extrapolation to water layer datuming for OBC record[J]. Journal of Tropical Oceanography, 30(6):84-89 (in Chinese with English abstract).

[3]
王祥春, 王延峰, 夏常亮, 等, 2012. Kirchhoff积分法OBS数据地震波场延拓[J]. 现代地质, 26(6):1231-1236.

WANG XIANGCHUN, WANG YANFENG, XIA CHANGLIANG, et al, 2012. Extrapolating the OBS data using the Kirchhoff integral method[J]. Geoscience, 26(6):1231-1236 (in Chinese with English abstract).

[4]
ARNULF A F, HARDING A J, KENT G M, et al, 2014a. Anatomy of an active submarine volcano[J]. Geology, 42(8):655-658.

[5]
ARNULF A F, HARDING A J, KENT G M, et al, 2014b. Constraints on the shallow velocity structure of the Lucky Strike Volcano, Mid-Atlantic Ridge, from downward continued multichannel streamer data[J]. Journal of Geophysical Research: Solid Earth, 119(2):1119-1144.

[6]
ARNULF A F, HARDING A J, SINGH S C, et al, 2012. Fine-scale velocity structure of upper oceanic crust from full waveform inversion of downward continued seismic reflection data at the Lucky Strike Volcano, Mid-Atlantic ridge[J]. Geophysical Research Letters, 39(8):L08303.

[7]
ARNULF A F, HARDING A J, SINGH S C, et al, 2013. Nature of upper crust beneath the Lucky Strike volcano using elastic full waveform inversion of streamer data[J]. Geophysical Journal International, 196(3):1471-1491.

[8]
ARNULF A F, SINGH S C, HARDING A J, et al, 2011. Strong seismic heterogeneity in layer 2A near hydrothermal vents at the Mid-Atlantic ridge[J]. Geophysical Research Letters, 38(13):L13320.

[9]
BERRYHILL J R, 1979. Wave-equation datuming[J]. Geophysics, 44(8):1329-1344.

[10]
BERRYHILL J R, 1984. Wave-equation datuming before stack[J]. Geophysics, 49(11):2064-2066.

[11]
BEVC D, 1995. Imaging under rugged topography and complex velocity structure[M]. PhD book, Stanford University.

[12]
BEVC D, 1997. Flooding the topography: Wave-equation datuming of land data with rugged acquisition topography[J]. Geophysics, 62(5):1558-1569.

[13]
BOLTE J F B, 2001. Velocity independent CFP redatuming: A strategy for subsalt imaging[C]//2001 SEG Annual Meeting. San Antonio, Texas: Society of Exploration Geophysicists, 901-904.

[14]
CHRISTESON G L, PURDY G M, FRYER G J, 1994. Seismic constraints on shallow crustal emplacement processes at the fast spreading East Pacific Rise[J]. Journal of Geophysical Research: Solid Earth, 99(B9):17957-17973.

[15]
COLLINS J A, BLACKMAN D K, HARRIS A, et al, 2009. Seismic and drilling constraints on velocity structure and reflectivity near IODP Hole U1309D on the Central Dome of Atlantis Massif, Mid-Atlantic Ridge 30°N[J]. Geochemistry, Geophysics, Geosystems, 10(1):Q01010.

[16]
GHOSAL D, SINGH S C, MARTIN J, 2014. Shallow subsurface morphotectonics of the NW Sumatra subduction system using an integrated seismic imaging technique[J]. Geophysical Journal International, 198(3):1818-1831.

[17]
HARDING A J, ARNULF A F, BLACKMAN D K, 2016. Velocity structure near IODP Hole U1309D, Atlantis massif, from waveform inversion of streamer data and borehole measurements[J]. Geochemistry, Geophysics, Geosystems, 17(6):1990-2014.

[18]
HARDING A J, KENT G M, BLACKMAN D K, et al, 2007. A new method for MCS refraction data analysis of the uppermost section at a Mid-Atlantic Ridge core complex[C]// American Geophysical Union, fall meeting 2007: S12A-03.

[19]
HARDING A J, KENT G M, ORCUTT J A, 1993. A multichannel seismic investigation of upper crustal structure at 9°N on the East Pacific Rise: Implications for crustal accretion[J]. Journal of Geophysical Research: Solid Earth, 98(B8):13925-13944.

[20]
HENIG A S, BLACKMAN D K, HARDING A J, et al, 2012. Downward continued multichannel seismic refraction analysis of Atlantis Massif oceanic core complex, 30°N, Mid-Atlantic ridge[J]. Geochemistry, Geophysics, 13(5):Q0AG07.

[21]
KENT G M, SINGH S C, HARDING A J, et al, 2000. Evidence from three-dimensional seismic reflectivity images for enhanced melt supply beneath mid-ocean-ridge discontinuities[J]. Nature, 406(6796):614-618.

DOI PMID

[22]
LARKIN S P, LEVANDER A, 1996. Wave-equation datuming for improving deep crustal seismic images[J]. Tectonophysics, 264(1-4):371-379.

[23]
MARJANOVIĆ M, FUJI N, SINGH S C, et al, 2017. Seismic signatures of hydrothermal pathways along the east pacific rise between 9°16′ and 9°56′N[J]. Journal of Geophysical Research: Solid Earth, 122(12):10241-10262.

[24]
MCMECHAN G A, CHEN H W, 1990. Implicit static corrections in prestack migration of common-source data[J]. Geophysics, 55(6):757-760.

[25]
QIN YANFANG, SINGH S C, 2018. Insight into frontal seismogenic zone in the Mentawai locked region from seismic full waveform inversion of ultra-long offset streamer data[J]. Geochemistry, Geophysics, Geosystems, 19(11):4342-4365.

[26]
RESHEF M, 1991. Depth migration from irregular surfaces with depth extrapolation methods[J]. Geophysics, 56(1):119-122.

[27]
SHTIVELMAN V, CANNING A, 1988. Datum correction by wave-equation extrapolation[J]. Geophysics, 53(10):1311-1322.

[28]
SINGH S C, HARDING A J, KENT G M, et al, 2006. Seismic reflection images of the Moho underlying melt sills at the East Pacific Rise[J]. Nature, 442(7100):287-290.

PMID

[29]
VERA E E, DIEBOLD J B, 1994. Seismic imaging of oceanic layer 2A between 9°30′N and 10°N on the East Pacific Rise from two-ship wide-aperture profiles[J]. Journal of Geophysical Research: Solid Earth, 99(B2):3031-3041.

[30]
VITHANA M V P, XU MIN, ZHAN XU, et al, 2019. Geological and geophysical signatures of the East Pacific Rise 8°-10° N[J]. Solid Earth Sciences, 4(2):66-83.

DOI

[31]
WANG HAIYANG, SINGH S C, CALANDRA H, 2014. Integrated inversion using combined wave-equation tomography and full waveform inversion[J]. Geophysical Journal International, 198(1):430-446.

DOI

[32]
WHITMORE N D, 1983. Iterative depth migration by backward time propagation[M]// SEG Technical Program Expanded Abstracts. Society of Exploration Geophysicists, 1983: 382-385.

[33]
WIGGINS J W, 1984. Kirchhoff integral extrapolation and migration of nonplanar data[J] Geophysics, 49(8):1239-1248.

DOI

[34]
XU MIN, 2012. Advanced geophysical studies of accretion of oceanic lithosphere in Mid-Ocean Ridges characterized by contrasting tectono-magmatic settings[D]. Cambridge, MA: Massachusetts Institute of Technology.

[35]
XU MIN, STEPHEN R A, CANALES J P, 2017. Waveform modeling of the seismic response of a mid-ocean ridge axial melt sill[J]. Marine Geophysical Researches, 38(4):373-391.

DOI

[36]
YILMAZ O, LUCAS D, 1986. Prestack layer replacement[J]. Geophysics, 51(7):1355-1369.

DOI

[37]
ZELT C A, BARTON P J, 1998. Three-dimensional seismic refraction tomography: A comparison of two methods applied to data from the Faeroe Basin[J]. Journal of Geophysical Research: Solid Earth, 103(B4):7187-7210.

Options
文章导航

/

〈 〉