引用本文
任芳, 李振春, 张敏, 杨国权, 陈飞旭. 基于角度域共成像点道集的微分相似优化法偏移速度分析[J]. 物探与化探, 2017,41(5): 872-880.
REN Fang, LI Zhen-Chun, ZHANG Min, YANG Guo-Quan, CHEN Fei-Xu. Differential semblance optimization migration velocity analysis based on angle domain CIGs[J]. Geophysical and Geochemical Exploration, 2017,41(5): 872-880.
Doi:10.11720/wtyht.2017.5.13  
Permissions
Copyright©2017, 物探与化探编辑部
《物探与化探》编辑部
基于角度域共成像点道集的微分相似优化法偏移速度分析
任芳, 李振春, 张敏, 杨国权, 陈飞旭
中国石油大学(华东) 地球科学与技术学院,山东 青岛 266580

作者简介: 任芳(1990-),女,在读硕士研究生,主要从事地震波偏移成像以及速度分析方面的研究工作。

收稿日期: 2016-08-01
基金: 山东省自然科学基金(ZR2013DL012); 高等学校博士学科点专项科研基金(20120133120012); 国家自然科学基金(41504100)
摘要

微分相似优化(Differential Semblance Optimization,DSO)法是一种可以有效避免非线性最小二乘反演过程中不聚焦问题的速度反演方法。基于DSO原理的波动方程偏移速度分析将成像道集的聚焦或者拉平程度作为速度判定的准则,其目标泛函具有良好的凸性,可以有效避免局部极值问题,目标泛函的梯度较为光滑,并且对地震数据中低频成分的要求较低,对速度模型中背景速度的估计较为准确。一般情况下,DSO方法利用偏移距域共成像点道集建立目标泛函,文中利用角度域共成像点道集建立目标泛函对速度的准确性进行判定。模型试算结果表明,通过角道集可以更加直观地对速度的准确性进行判定,并且准确反映速度与深度的对应关系。

关键词: 偏移速度分析; 微分相似优化法; 共成像点道集; 目标泛函; 梯度
中图分类号:P631.4 文献标志码:A 文章编号:1000-8918(2017)05-0872-09
Differential semblance optimization migration velocity analysis based on angle domain CIGs
REN Fang, LI Zhen-Chun, ZHANG Min, YANG Guo-Quan, CHEN Fei-Xu
School of Geoscience and Technology,China University of Petroleum (East China),Qingdao 266580,China
Abstract

Differential semblance optimization (DSO) is an approach to the inversion of the velocity which avoids the severe convergence associated with nonlinear least-squares inversion.DSO method-based wave-equation migration velocity analysis measures the focusing or flatness of image gathers and uses the deviation as the criterion to update the velocity.The DSO objective function has a decent global convexity property,therefore it can avoid the local minima problem,the gradient of objective function is smooth,it can inverse the background velocity accurately when there lacks low frequency data.Generally,ODCIGs is used to construct the objective function.In this paper,the authors use ADCIGs to construct the objective function,and the result of the test shows that DSO method with ADCIGs can estimate the correctness of velocity model directly and reflect the coupling relation between velocity and depth accurately.

Keyword: migration velocity analysis; differential semblance optimization; common image gathers; objective function; gradient
0 引言

随着地震勘探目标的复杂化以及计算机软、硬件技术的发展, 地震成像技术逐渐从叠后时间偏移发展到叠前深度偏移, 对于任何偏移方法, 其偏移结果的精度都取决于速度模型和偏移算子的准确程度。叠前深度偏移方法具有较高的处理精度, 较强的适应性和良好的稳定性, 但是深度与速度的强耦合性使得偏移成像结果对速度误差具有很强的敏感性。利用深度偏移对速度误差的敏感性, 可以对偏移速度进行估计, 并通过迭代更新得到更为准确的速度场。

偏移速度分析建立在偏移成像的基础上, 是在偏移过程后根据偏移结果来修正偏移速度模型的方法, 属于一种反问题[1]。该过程利用偏移结果对偏移速度模型的敏感性, 建立适当的速度误差判断准则以及速度更新方程, 通过判断准则判断速度是否需要进行调整, 通过更新方程计算速度模型的修正量, 将成像结果差异反投影到速度空间, 对速度模型进行更新直到得到满意的偏移速度结果[1-3]。

波动方程偏移速度分析(wave equation migration velocity analysis, WEMVA)基于波场传播理论, 在本质上比射线类的偏移速度分析方法更稳健, 避免了射线传播过程中由于遇到复杂构造或突变界面而表现出的不稳定现象。WEMVA中带限波场的传播分量对速度场的光滑变化也更为敏感。因此, WEMVA可以得到稳定、光滑的速度更新量, 是准确估算背景速度场的有效手段。Toldi首先采用褶积正演模型进行叠加能量最大化自动速度分析[4], Chavent和Jacewitz利用炮域逆时偏移及对应的层析算子采用相同准则重建速度场[5]。Soubaras和Gratacos利用面炮偏移在叠加能量最大准则下进行自动偏移速度分析[6]。Sava和Biondi利用单程波建立成像扰动与速度扰动之间的线性关系, 然后通过线性化的Stolt剩余偏移求取成像扰动[1], 从而线性反演得到速度扰动。

DSO方法将数据拟合项作为约束条件, 使得扩展模型的一致性达到最优所对应的速度场即为该方法所追寻的解; 其核心思想是将模型进行扩展, 通过定义对应的扩展正演方程, 建立数据拟合和正则化合并之后的目标泛函。Symes在1991年提出DSO的反演概念[7], 该方法在成像域有自动反演偏移速度的能力, 且目标函数的凸性好, 梯度光滑, 随后有很多研究者被这一方法所吸引[8-11]。Shen先后利用单程波偏移和逆时偏移, 实现了波动方程偏下的DSO偏移速度分析[12-14]。其中Shen和 Symes基于炮域的单程波偏移展示了DSO偏移速度分析在复杂构造情况下的适用性, 并对目标函数进行优化, 得到较好的反演效果。Fei和Williamson分析了DSO方法所特有的梯度假象, 并提出了改进方法, 使得速度更新方向得到修正[15], 反演更准确。Yang T和P Sava于2009年, Yang于2013年提出了基于时空扩展的CIP道集下的DSO偏移速度分析方法[16, 17], 该方法对高陡构造有更强的反演能力。Tang和Biondi利用面向目标的思想对数据进行重新映射[18], 通过振幅规则化减小了振幅不均匀导致的梯度奇异分布的问题, 并且利用优化后的方法对盐丘下目标层进行了DSO偏移速度分析。Shan在RTM理论下实现了偏移速度分析[19], 得到了很好的反演效果。Shen和Symes提出水平收缩算子作用于地下偏移距域应用或者水平拉伸算子作用于角度域共成像点道集(angle domain common image gathers, ADCIGs)近似计算目标图像[20], 较好地避免了周期跳跃问题。

微分相似度最优化方法在一定条件下与波形反演等价, 其优势在于该目标泛函具有较好的凸函数性质, 可以在一定程度上解决波形反演存在的局部极值问题。经典的DSO 方法存在一些问题, 如梯度噪声、由照明导致的不聚焦等, 且自动拾取的计算效率较低。笔者基于逆时偏移的叠前深度偏移方法, 利用偏移距域共成像点道集和角度域共成像点道集共同建立速度判别准则, 对偏移速度进行迭代更新。相比于传统DSO方法, 笔者利用RTM提高偏移结果的准确性, 进而可以得到更加准确的共成像点道集; 利用偏移距域共成像点道集(offset domain common image gathers, ODCIGs)和ADCIGs建立速度判别准确, 可以有效提高速度判别准确性, 且判别条件更为直观。

1 DSO法偏移速度分析基本理论

DSO是波动方程偏移速度分析方法中的一类, 该方法利用道集之间的差异性来构造目标函数。其理论主要包括两方面:一是对成像结果优劣的判断, 这方面涉及剩余成像的提取, 即目标函数的建立; 二是利用波动方程对速度模型进行更新, 这方面涉及到目标函数对速度模型参数梯度求解。

1.1 目标函数的建立

DSO偏移速度分析基于偏移成像扰动和慢度之间的对应关系, 通过优化背景慢度使得成像域中的剩余场最小。目标函数表示为

J=12‖PI‖2, (1)

其中I为共成像点道集(common image gathers, CIGs), P为DSO算子, 速度正确时, DSO算子对CIGs的能量进行消减; 速度不正确时, DSO算子可以消除CIGs中的部分奇异值。

DSO属于一种对CIGs在扩展轴上的相似性进行量化的反演方法, 不同CIGs的表达式不同, 对应的DSO算子表达式也不同:

1)对于HOCIGs I(x, h), DSO算子P=|h|, 目标函数Jh表示为

Jh=12‖hI(x, h)‖2。(2)

2)对于时移共成像点道集(time shift common image gathers, TSCIGs)I(x, τ ), DSO算子P=|τ |, 目标函数Jτ 表示为

Jτ=12‖τI(x, τ)‖2, (3)

速度模型正确时, CIGs在时间和空间上只聚焦到一点, 即为真实成像点。

3)对于ADCIGs I(x, θ ), I(x, θ )=Aθ , hI(x, h), DSO算子P=口/口(tanθ ), 此时目标函数Jθ 为

Jθ=12‖∂I(x, θ)∂(tanθ)‖2=12‖∂Aθ, hI(x, h)∂(tanθ)‖2, 　(4)

当速度模型正确时, ADCIGs对入射角正切值的导数应为零, 角道集是水平的。

1.2 速度判定准则

基于DSO偏移速度分析的速度判别准则可以分为两类:

1)成像道集聚焦准则。速度正确时, ODCIGs聚焦于零偏移距或者零时刻处; 速度不正确时, 道集能量会向非零偏移距或非零时刻处发散。

2)道集拉平准则。速度正确时, ADCIGs同相轴拉平, 即成像深度不随入射角度变化; 速度不正确时, 偏移深度随地下反射角变化, 具体表现为:速度偏大时道集向下弯; 速度偏小时道集向上翘。

2 DSO梯度的求取

利用目标函数对剩余成像进行提取后, 需要计算剩余波场并将剩余波场进行反向传播进而对速度模型进行更新, 这部分涉及到目标函数对速度模型梯度的求取。

以ODCIGs为例, 目标函数表达式:

J[s]=12‖PhIh‖2, (5)

根据方程(2)可知, 地下偏移局域微分相似算子P=|h|, 所以在微分相似算子的作用下, 偏移距域共成像点道集变为地下偏移距与成像道集的乘机形式:

Ph∶Ih→hIh(x, h)。(6)

成像道集Ih是慢度s的非线性函数, 表示为

Ih=f[s], (7)

则

δIh=(DIh)[s]δs=(∂Ih/∂s)δs。(8)

其中DIh是对慢度(速度)的微分偏移, 通过对成像函数f相对于慢度s求一阶导数来获得。目标函数的扰动可以用内积的形式表示:

δJ=12< δ(PhIh), PhIh> +12< PhIh, δ(PhIh)> =12< PhDIhδs, PhIh> +12< PhIh, PhDIhδs> =12< δs, DIh* Ph* PhIh> +12< DIh* Ph* PhIh, δs> =12< δs, DIh* Ph* PhIh> +12< δs, DIh* Ph* PhIh¯> =< δs, Re(DIh* Ph* PhIh)> 。(9)

则目标函数J对于慢度s的梯度为

∇sJ=Re{(DIh)* (Ph* PhIh)}。(10)

ODCIGs用Green函数表示为

Ih(x, h)=∫G* (x-h, xs, ω)G* (x+h, xr, ω)·d(xr, xs, ω)dxrdxsdω, (11)

其中:d(xr, xs, ω )表示震源xs激发, 接收点xr接收到的数据在圆频率域的表示形式; G(x, y, ω )表示震源y在x点处的Green函数。成像扰动Δ Ih=DIhΔ s可以用Green函数G* (x-h, xs, ω )和G* (x+h, xr, ω )的Born扰动来表示为:

对求和项取复共轭可得梯度的Green函数表达式:

将式中的Δ I用 Ph* PhIh代替, Δ s则变为DSO目标函数的对于慢度的梯度∇sJ。即

其中u(x, h)= Ph* PhIh(x, h)。

上式等号右端可以分为两部分来进行分析。首先定义震源和接收点的背景波场为

震源和接收点的剩余波场为

这样, 可以分别得到震源和接收点的梯度:

梯度则可以表示为以上两部分的和:

g(y)=gS(y)+gR(y)。(19)

具体流程如图1所示。

图1 DSO偏移速度分析算法流程

3 模型试算

针对以上研究内容, 首先使用1个简单的平层模型来验证DSO方法的有效性, 然后使用Marmousi模型进一步的试算。笔者采用有限差分波动方程数值模拟方法进行正演, 利用RTM方法进行偏移成像。其实现过程主要需要解决:①波动方程的数值离散化求解; ②地下介质数值模型的建立; ③震源的数值模拟(震源加载); ④有限空间模型的边界反射压制和消除; ⑤稳定性条件和数值频散的压制。

3.1 水平层状模型

平层模型如图2a所示, 其中上、下层速度分别为2 000 m/s和2 500 m/s, 横向和纵向采样点数分别为737、200, 采样间隔分别为10 m和8 m。正演炮记录正演记录共200炮, 炮点初始位置x=1.8 km, 炮间距为20 m。接收方式为采用排列移动接收, 每个排列有301个检波器, 相邻检波器的间距为10 m, 炮点位于排列的中点, 每一炮的偏移距范围是[-1.5 km, 1.5 km]。采用主频为25 Hz的雷克子波进行激发。

图2b为正确速度下的偏移成像结果, 采用1 700 m/s的均匀速度场作为初始模型进行速度分析, 图2c、2d分别为迭代10次和迭代20次后的偏移结果, 图2e、2f分别为迭代10次和迭代20次后CDP=369处的偏移距道集和角道集。由图可见, 由于初始速度模型速度偏小, 偏移距道集中的能量向下发散, 角道集向上弯曲。随着迭代次数的增加, 偏移结果中同相轴深度逐渐向真实界面深度靠近, 能量逐渐收敛, 由于偏移距道集发散而产生的交叉假象得到消除; 偏移距道集能量逐渐收敛, 角道集逐渐被拉平, 偏移速度逐渐接近真实速度。

图2 水平层状模型试算

与图2b的真实偏移成像相比, 迭代20次的偏移结果(图2d)的构造形态及深度信息准确, 但由于照明等因素的影响, 构造范围并没有得到准确刻画, 在保真度与成像范围上存在一定的误差, 尤其在边界处存在照明不足的情况, 导致成像结果存在误差; 为了提高成像结果的质量, 可以在正演过程中对炮点进行加密, 并对目标函数进行优化, 进而对深部能量进行有效的补偿, 同时可以对反演算法进行优化更新, 提高反演的精度进而提高成像结果的准确程度。

3.2 Marmousi模型试算

真实速度模型如图3a所示, 该模型是对经典的Marmousi模型进行深度截取, 并在顶部加入200 m的水层得到的。模型横向采样点数为737, 采样间隔10 m, 纵向采样点数为200, 采样间隔8 m。地震记录正演采用时间2阶空间10阶网格, 正演时间间隔为0.5 ms, 采用20 Hz雷克子波, 记录时间总长度为2.4 s。

图3 Marmousi模型试算

真实速度模型得到的偏移结果如图3d所示, 采用图3b中的速度场做为初始模型, 对应的偏移成像结果如图3e, 利用DSO方法进行速度更新, 图3c和3f分别为迭代36次后的速度场以及对应的偏移结果。随着迭代次数的增加, 偏移结果中同相轴更加清晰, 基本构造的轮廓逐渐准确。

图4a~4c分别为初始模型、迭代36次模型和真实模型下的偏移距道集, 由图可知, 偏移距道集中发散的能量逐渐向零偏移距处聚焦, 深处的偏移速度逐渐准确; 图4d~4f分别为初始模型、迭代36次模型和真实模型下的角道集, 随着迭代次数的增加, 角道集逐渐被拉平, 尤其是深部的改变更为明显。

在模型计算过程中, 基于ADCIGs的偏移速度分析过程所需的时间相对于ODCIGs速度分析要更长, 但是所需要的迭代次数要少于传统的基于ODCIGs偏移速度分析, 总的来讲, 计算效率有了一定程度的提高, 一定程度上节约了计算费用。

为了更好地验证偏移速度的准确性, 从迭代反演结果中抽取x=2.4、4.6、5.8 km处的垂向速度, 并对初始速度、迭代15次速度、迭代36次速度以及真实速度进行对比, 对比结果如图5所示, 可以明显看出, 随着迭代次数的增加, 偏移速度逐渐靠近真实速度, 但是细节上刻画的仍然不够理想。

在模型试算的过程中, 对于断层、断点处的梯度噪声, 可以通过一致性补偿进行噪声的消除, 但去噪效果仍有待提高, 如何更好地处理梯度中存在的噪声也是下一步需要继续研究的内容。

图4 x=4.6 km处提取的共成像点道集

图5 抽取不同位置处不同迭代次数下CDP速度曲线对比

4 结论

文中详细介绍了基于角度域共成像点道集的DSO偏移速度分析方法的基本原理, 并对模型进行了试算, 得出如下结论:

1) DSO方法通过计算相邻共成像点道集之间的差异来获取成像扰动, 可以实现自动拾取, 有效避免人为因素的影响, 降低了速度分析的人工成本; DSO目标函数凸性较好, 可以有效避免局部极值的出现; 基于角度域共成像点道集的DSO方法在速度判别的过程中更加直观, 可以对偏移速度的正确性进行较为直观、准确的判断。

2) DSO方法作为一种波动方程反演方法, 其计算效率受迭代精度以及计算机性能的影响较大, 在计算时间上要高于射线类反演方法, 其计算效率比层析方法要低, 但对成像结果中细节的刻画更加准确, 在保证成像精度的基础上, 计算效率要高于波形反演类方法; 目标函数的构造也会对计算精度及效率产生影响。

3) DSO方法的反演过程中, 在反射层出现断层、断点等奇异值时, 会引入梯度噪声, 这种噪声在一定程度上会影响到反演的收敛速度; 在梯度的处理过程中, 需要对深部能量进行合理的补偿, 进而提高收敛速度与反演效果。

4) 对于文中利用的角度域共成像点道集, 可以更为直观的对偏移速度的准确性进行判断分析, 相比于传统的偏移距域共成像点道集的约束条件, 角道集中同相轴形态对速度误差更为敏感, 可以对更为精细的优化过程提高可靠依据, 进而实现更好的偏移效果。

致谢:感谢中国石油大学(华东)地震波传播与成像课题组的各位老师和同学的指导与帮助, 感谢审稿人的批评与指导。

The authors have declared that no competing interests exist.

参考文献
[1] Sava P, Biondi B. Wave-equation migration velocity analysis . I. Theory[J]. Geophysical Prospecting, 2004, 52(6): 593-606. [本文引用:2]
[2] Shen P, Symes W W. Automatic velocity analysis via shot profile migration[J]. Geophysics, 2008, 73(5): VE49-VE59. [本文引用:1]
[3] Albertin U, Etgen J, Maharramov M, et al. Adjoint wave-equation velocity analysis[C]//New Orleans: Expand ed Abstracts of 76th SEG Meeting, 2006: 3009-3013. [本文引用:1]
[4] Toldi J L. Velocity analysis without picking[J]. Geophysics, 1989, 54(2): 191-199. [本文引用:1]
[5] Chavent G, Jacewitz C A. Determination of background velocities by multiple migration fitting[J]. Geophysics, 1995, 60(2): 476-490. [本文引用:1]
[6] Soubaras R, Gratacos B. Velocity model building by semblance maximization of modulated-shot gathers[J]. Geophysics, 2007, 72(5): U67-U73. [本文引用:1]
[7] Symes W W, Carazzone J J. Velocity inversion by differential semblance optimization[J]. Geophysics, 1991, 56(5): 654-663. [本文引用:1]
[8] Symes W, Versteeg R. Velocity model determination using differential semblance optimization[C]//Washington: Expand ed Abstracts of 63rd SEG Meeting, 1993: 696-699. [本文引用:1]
[9] Symes W W, Kern M. Inversion of reflection seismograms by differential semblance analysis: algorithm structure and synthetic examples1[J]. Geophysical Prospecting, 1994, 42(6): 565-614. [本文引用:1]
[10] Mulder W A, Ten Kroode A P E. Automatic velocity analysis by differential semblance optimization[J]. Geophysics, 2002, 67(4): 1184-1191. [本文引用:1]
[11] Symes W W. Migration velocity analysis and waveform inversion[J]. Geophysical prospecting, 2008, 56(6): 765-790. [本文引用:1]
[12] Shen P, Symes W W, Stolk C C. Differential semblance velocity analysis by wave-equation migration[C]//Dallas: Expand ed Abstracts of 73rd SEG Meeting, 2003: 2132-2135. [本文引用:1]
[13] Shen P, Symes W W, Morton S, et al. Differential semblance velocity analysis via shot profile migration[C]//SEG Technical Program Expand ed Abstracts, SEG, 2005: 2249-2252. [本文引用:1]
[14] Shen P. An RTM based automatic migration velocity analysis in image domain[C]//Las Vegas: Expand ed Abstracts of 82nd SEG Meeting, 2012: 1-5. [本文引用:1]
[15] Fei W H, Williamson P. On the gradient artifacts in migration velocity analysis based on differential semblance optimization[C]//SEG Technical Program Expand ed Abstracts, SEG, 2010: 4071-4076. [本文引用:1]
[16] Yang T N, Sava P. Wave-equation migration velocity analysis using extended images[C]//ouston: SEG Technical Program Expand ed Abstracts, SEG, 2009: 3715-3719. [本文引用:1]
[17] Yang T N. Wavefield tomography using extended images[D]. Colorado: Colorado School of Mines, 2013. [本文引用:1]
[18] Tang Y X, Biondi B. Target-oriented wavefield tomography using synthesized Born data[J]. Geophysics, 2011, 76(5): WB191-WB207. [本文引用:1]
[19] Shan G J, Wang Y. RTM based wave equation migration velocity analysis[C]//SEG Technical Program Expand ed Abstracts, 2013: 4726-4731. [本文引用:1]
[20] Shen P, Symes W W. Subsurface domain image warping by horizontal contraction and its application to wave-equation migration velocity analysis[C]//Houston: 2013 SEG Annual Meeting, SEG, 2013: 4715-4719. [本文引用:1]