基于速度模型融合与断控约束的三维层析反演方法
(1. 中国石化华东油气分公司勘探开发研究院, 南京 210000
2. 地球与空间科学学院, 北京大学, 北京 100871)
摘要: 常规的网格层析反演是基于反射波的走时与共成像点的剩余曲率来构建矢量方程从而更新速度, 然而, 在复杂地表区域的单炮信噪比低, 反射信息不足等因素导致近地表区域速度建模精度不足, 且缺少层位与断层的约束, 无法准确还原深层的实际构造。目前工业中广泛使用的方法是将初至波层析反演得到的浅层速度以向下平移一定深度的地表高程面作为融合面与反射波层析得到的中深层速度相融合, 但这种方法会在融合区出现明显的融合痕迹, 且融合处的速度值会发生突变。为此, 本文从层析反演的原理出发, 根据初至波射线密度与两模型差值来拾取融合面, 并构建符合工区地质背景的融合函数进行速度融合, 最大限度消除融合痕迹与融合处速度值突变的问题, 并在之后的层析反演过程中加入断层约束。本方法运用速度融合, 断控约束以及改善剩余曲率拾取方法来提高复杂山地区的速度建模与成像精度。复杂模型结果证明了正确性, 实际资料测试结果证明了适用性。
关键词: 速度融合, 层位约束, 三维角道集提取, 全局层析, 偏移成像
DOI: 10.48014/cpngr.20230425001
引用格式: 熊鹰杰, 王非翊, 孟庆利, 等. 基于速度模型融合与断控约束的三维层析反演方法[J]. 中国石油天然气研究, 2023, 2(1): 1-14.
文章类型: 研究性论文
收稿日期: 2023-03-25
接收日期: 2023-03-27
出版日期: 2023-03-28
0 引言
南川地区页岩气资源丰富，开发潜力巨大，该区基地为前震旦系板溪群浅变质岩，上覆盖层自下而上发育上震旦统、寒武系、奥陶系、中下志留统、二叠系、三叠系、侏罗系地层，厚度2000~6000m，区内页岩主要发育上奥陶统五峰组-下志留统龙马溪组，该区的数据处理中主要面临两个问题，首先是近地表结构复杂且浅层资料信噪比低，基于共成像点剩余曲率的网格层析无法精确的建立近地表速度模型。其次南川地区逆掩断层发育，常规网格层析[1-3]速度建模方法没有充分考虑断层对速度分布的控制因素，从而造成局部成像精度低的问题。针对这两个问题许多学者都进行过深入研究，传统的网格层析是利用偏移距域共成像点距道集[3,4]进行更新速度模型，该道集在浅层的信息少，剩余曲率精度不足，因此在更新浅层速度模型时，速度往往不准确，之后改进为角度域共成像点道集[5,6](ADCIGs)进行层析反演，精度虽有提高但无法满足工业需求。波动方程层析反演[7]基于Born近似与Rytov近似，分析对比理论波场与观测波场的最佳匹配，在复杂的地下构造中理论上可以反演出精度更高的速度场[8]，但计算效率大大降低，而且最优解问题还未解决。目前实际生产中应用最广泛的速度建模流程是使用初至波层析反演[9]得到近地表速度，并与反射波层析得到的中深层速度场相融合[10]。这种方法虽然可以改善近地表速度模型的精度，但融合过渡区域速度值的正确性[11]与逆掩断层发育区的速度值扰动[12]都没有妥善解决。因此，亟须一种更加合理且准确的复杂地区速度反演方法。本文提出了一种将浅中深层速度融合后加入断层约束的全局层析反演方法，在初始速度模型构建时，将初至波层析出的浅表层速度场与常规速度分析得到的中深层速度场相融合，首先确定两速度场融合的目标层位，之后对该区域附近进行加权计算选取的权值;对于速度空缺区域，通过重构补缺建立统一全速度函数，然后在速度模型上加入断层与层位构造进行约束，最后根据角度域共成像点道集剩余曲率计算出的走时残差构建并改进层析公式进行全局层析以得到一个高精度的速度场。
1 速度模型融合
初始速度模型的精度对层析反演来说至关重要，尤其是射线类层析反演中，一个精度较高的初始速度模型不仅可以减少迭代的次数，而且影响着后续反演的准确性，因此在建立初始度模型时，为了提高准确度，需要分别建立浅层速度模型与中深层速度模型，近地表速度模型由初至波层析反演获得[13]，中深层速度模型由常规速度分析获得[14]，将近地表速度模型与中深层速度模型按照一定的方法进行有机的融合，从而构建一个高精度的初始全速度模型。
1.1 融合界面选择
速度融合之前需要对速度模型进行处理，首先判断浅层速度模型与中深层速度模型的大地坐标是否在同一区域内，如果不同则将两者截取相同部分。之后检查二者数据道与横向网格点数是否保持一致，网格点不一致可重采样处理。
浅层的速度建模利用初至波层析反演完成的，射线密度决定了反演的精度，因此根据射线路径的密度选择一块反演精度较高的浅层区域，然后对融合界面进行选择，选择策略为:首先对初至波层析反演得到的速度场进行偏移成像得到一个浅层的地震剖面，之后在剖面上找到之前确定的高精度区域，在这个区域的最大深度附近找到一条较为平缓的层面，并以这个平面为基准面。然后将浅层速度模型与中深层速度体的横向上每隔一定的间距都抽取一个速度剖面，在速度剖面上对基准面范围附近进行分析，之后将两模型的剖面做差，找出速度值差异最小的层位，将这些层位信息记录下来，统计与上一步基准面位置最接近的层位，将这些层位采用三次函数插值法进行连接从而得到融合面，这种方法不仅确保了融合层位的正确性，并且使融合界面更加平滑。拾取的基准面如图1所示。
图1 融合面示意图
Fig.1 Schematic diagram of fusion surface
1.2 分区域加权融合
在三维速度模型融合时融合方法与过程对模型的形态与构造有至关重要的影响，因此需要确保融合的合理性。鉴于两模型的精确范围不同，在充分考虑这些因素之后，选择分区域加权融合。 在融合的过程中对浅层模型与中深模型度赋予不同的权重，融合的结果可以表示为:
Vf(x，y，z)=
Vα(x，y，z)·W+Vβ(x，y，z)·(1-W)(1)
式中，Vf(x，y)表示速度模型中的点(x，y)融合后的速度，Vα表示初至波层析反演出的浅层速度模型，Vβ表示中深层速度模型，W是权重，加权的方式有两种可以选择:
自定义权重范围的反距离加权:
W(i，j)=F-(Wmax-Wmin)
(2)
三角函数加权:
W(i，j)=
(3)
式中，F代表最大权重，Wmax代表最大权值，Wmin最小权值，Zmax表示最大深度，Zi表示融合区内某一深度。
自定义权重的反距离加权法是根据深度不同而划分权重的，权重的范围根据实际情况进行调节，由于该函数为权值可变的线性函数，因此适用于速度随深度的变化不剧烈的地层。而 sin函数加权是非线性的加权函数，通过改变参数能够适应上下速度变化较为剧烈的地层，可以使融合后的速度更接浅层速度或中深层速度，这个选择就是根据断层的位置与对应的层位，或是大倾角地层的层速度对应深度。因此这两种方法需要根据实际情况进行选择。
在融合区如果不做出合适的处理的方法时，会存在明显的区域间痕迹，这些痕迹即为融合造成的速度突变区，或是层位断裂区，其中后者我们可以通过调整融合界面来进行规避，前者可以通过下一节所提出的融合方法优化来有效解决，也就是改善速度值的突变问题。
1.3 融合方法优化
使用上述提出的方法进行速度融合时会面临着三个问题。
(1)两速度模型在融合处的速度值会在交点处出现扰动。如图2所示，在两速度模型的差剖面中标红的地方为扰动区域，指两速度在这个区域速度差值较大，从中色标可以发现二者相差数百米，融合后会出现突变的情况。
图2 融合处扰动区域
Fig.2 Disturbed region at fusion
(2)直接融合会改变速度模型原本的形态。这个问题与第一个问题较为相似，虽然在上述融合界面选择中进行了一定程度上的解决，但在两模型的速度值相差较大的地方还是无法合理解决，因此需要融合前的处理来解决上述两个问题:抽出两模型等间隔的速度剖面，并在剖面中标定出融合界面的深度范围，对比两模型的速度值是否在此区域存在较大差值，若差值较大则需对此区域的中深层模型进行局部平滑与调整，因为中深层模型的速度值在近地表区精度是较低的，经过调整与平滑处理可以缩小与浅层速度模型之间的差值，处理结果如图3所示。当融合区域附近的速度值都较为接近，有效地避免了融合后的异常值。由于地下介质在融合区域没有出现速度异常体，但是浅层速度模型速度变化较为剧烈，导致速度值扰动。这种情况下需要先在融合范围内选择出速度值变化相对最小的区域并在此区域内进行一定的平滑处理，从而使融合后的速度值不会明显偏离两模型的速度值，如图3所示。第二种情况就是当某一速度模型在融合界面处突然出现速度值的大幅度抖动，这种情况也会造成融合后的速度值异常，因此也选定这个区域进行平滑处理。这个问题解决之后可以有效地避免速度值突变带来的伪速度。
图3 扰动点处理后
Fig.3 Disturbance point after processing
(3)融合范围划分问题。初至波层析反演在浅层的射线密集区精度是最高的，如果这一部分参与融合就会降低浅层速度模型精度。同理，在深层部分，由于射线无法达到而导致精度不足，参与融合后同样会降低深层速度模型精度。因此我们需要选定融合区域，选定方法如图4所示，按照初至波层析反演的射线密度与射线最大深度，将模型分为融合区与不融合区，融合区按原有方法进行融合，不融合区在浅层直接使用初至波层析反演的速度值，深层部分直接使用常规速度分析的速度值。
图4 (a)融合范围;(b)初至波层析射线密度
Fig.4 (a)Fusion range;(b)Ray density of first-to-wave wave tomography
1.4 融合方法测试
对反距离加权法与sin函数加权法分别进行融合测试，选择三维模型进行测试，如图5所示。图(a)为原模型，模型大小为1041×512×917，(b)、(c)为中深层模型与浅层速度模型，(c)为融合界面。在融合时为了方便验证融合结果的准确性，我们将原模型进行平滑并乘以某一系数作为中深层初始模型来进行测试。
(1)反距离加权法。加权函数与融合模型如图6所示，根据图中速度值的分布与深度模型的形态来看，融合后的速度模型在近地表更接近初至波层析反演得到的浅层速度模型，深层部分更接近常规速度分析模型，且融合后模型与真实速度模型形态大致相同，初步判断融合结果正确。
(2)三角函数加权法。从图7(a)的函数图像中可以看到，浅层速度模型在近地表区域具有更大的权重，随着深度增加，中深层模型的权重逐渐增值最大。通过观测模型形态可以发现sin函数融合法的模型在浅层部分模型改变的较为平缓，且复杂构造区的形态也较为接近原模型。
这两种方法的正确性和准确度需要进一步判断，通过对比与原模型的差值体进行判断，由于在融合过程中模型经过平滑处理，因此将原模型进行平滑处理后分别与反距离加权法和三角函数加权法得到的融合模型分别做差，最终得到速度模型如图8所示:对比差值体可以发现sin加权法(a)在全局的差值不大且较为平均。反距离函数加权法(b)在融合时最大误差大于sin函数加权法，全局差值方面与反距离加权法相差不大。反距离加权法在此模型中误差较大的原因浅层模型与中深层模型在融合区域的差值较大且存在复杂构造，而反距离加权法无法适应构造与速度的剧烈变化。因此可以得出结论:当融合区域两模型差值较大时采用sin函数加权法得出的结果更为精确。
2 全局走时层析反演方法
2.1 全局层析优化
对于地下断层发育的复杂地区来说，常规网格层析在更新速度模型时无法考虑断层的影响，导致断层周围的速度不准确，而加入断层与层位约束的层析则可以解决这个问题，这种方法可以细致的刻画出速度模型的断层构造与层位信息，由于目前大部分探区都属于三维复杂介质勘探，传统层析反演方法得到的速度模型精度已无法满足偏移成像的需求，因此需要一种提高层析反演精度的方法，本文提出一种可以提高剩余曲率拾取精度的方法。
图5 (a)原模型;(b)浅层速度模型;(c)融合界面
Fig.5 (a)original model;(b)shallow velocity model;(c)fusion interface
图6 (a)反距离加权函数示意图;(b)融合后模型
Fig.6 (a)Schematic diagram of the inverse distance weighting function(b)Post-fusion model
图7 (a)sin加权函数示意图;(b)融合后模型
Fig.7 (a)Schematic diagram of the sin weighting function;(b)post-fusion model
图8 (a)sin加权融合差值速度体;(b)反距离加权融合差值速度体
Fig.8 (a)sin-weighted fusion difference velocity body;(b)Inverse distance weighting fusion difference velocity body
三维全局走时层析是之前二维的基础上推导出三维反演公式，并加入正则化约束，选择合适的模型参数，通过不断的迭代反演，得到最终的速度模型。三维层析反演方程的建立由三维Radon转换而来，三维Radon变换为:
(p，ε)=
f(x)δ(p-εx)dx
=
f(x，y，z)
δ(xsinθcosϕ+ysinθsinϕ+zcosθ-p)dxdydz (4)
式中:
p=xsinθcosϕ+ysinθcosϕ+zcosθ
ε=(sinθcosϕ，sinθsinϕ，cosθ)(5)
在上述的Radon正变换是对函数的投影。反过来从对函数的投影值即线积分值来求被积函数，称之为Radon反变换。这种由线积分值求被积函数的思想，早在20世纪数学家 N.H.Able 就开始研究了。他对某些规则的几何图形给出了结果。J.Radon的重要贡献是对任意形状的函数，在某些假设条件下求出被积函数。给出二维 Radon 反变换公式，同时给出在空间推导 Radon反变换公式，反变换如下:
f(x)=
sinθdϕdθ(6)
空间坐标如图9所示。
图9 三维Radon变换
Fig.9 3D radon transform
层析反演方程:
LΔs=Δt(7)
=
(8)
由于角度与共成像点道集是叠前深度偏移后的数据，在反应速度与深度的耦合关系上具有较高的精度，并且相较于其他共成像点道集来说减少了多路径造成的成像假象等干扰，由角道集的剩余曲率与时深转换公式来计算得到的旅行时差更为准确，求取角度域共成像点道集的方法有许多传统网格层析使用的是三点法但由于其精度达不到工业需求已经逐渐被其他方法取代，目前工业上广泛使用基于方向矢量方法并且这个方法是目前求取ADCIGs最高效的方法，图10所示为矢量法求取入射角信息，在逆时偏移的过程中提取角度域共成像点道集但这类方法在每个时间和网格点只能估算波场传播的一个主要方向，当涉及介质复杂时，这类方法可能不准确。
图10 矢量法计算入射角
Fig.10 Vector method to calculate the angle of incidence
由图中的几何关系我们可以得到:
θ=
arccos
(9)
式中，Ss，Sr代表由矢量计算得到的传播方向，通过几何关系可以得到入射角β。
(10)
在复杂介质中计算角道集的角度往往容易受到干扰，为了提高计算角度的稳定性，可以采取一段时间的叠加计算入射角度，使用最小平方求解角度。因为在一段时间内传播的角度α满足:
α(s，r，t)=arctan
(11)
式中，t代表当前时刻;Δt代表时间间隔数目;Ss，Sz表示传播方向。叠加后虽然可以提高计算效率与稳定性但在叠前深度偏移的过程中并不容易得到许多时间间隔的矢量，因此需要将公式改进为一定距离间隔的计算方法:
α(s，r，t)=
arctan
(12)
式中，x代表当横坐标;Δx代表位移间隔;m，n代表移动步长。可以通过改变步长来提高计算效率。综上所述入射角β也可由此方法计算得到:
β(a，b，s，r，t)=arccos
(13)
最后根据射线路径与深度之间的几何关系计算出其相应的走时残差代入式(8)。从而求解方程更新速度场。
2.2 断控约束
该方法具体的实现过程为:首先利用融合后的速度模型进行叠前深度偏移;然后在深度偏移剖面上进行层位和断层解释，利用解释的层位和断层信息作为层析反演的约束条件;层析划分网格时将断层两侧划分为不同的网格，以断层为分界面展开分块速度更新，避免了常规网格层析方式在断层附近的速度模糊特性;然后对更新后的速度再次进行叠前深度偏移，重新拾取层位与断层信息以及ADCIGs剩余残差，迭代更新速度模型。
断层拾取首先进行剖面的解释，然后再逐条测线查看断裂的延伸方向进行人工组合，做到平面—剖面—平面相结合，断裂系统应该有一定的规律性，尤其区域断裂和控制局部构造的主断裂走向分布的落实清楚，钻井、地质等资料相互匹配、相互验证。断层组合时要考虑区域、局部构造应力场背景，根据相邻、相交剖面断层特征一致性原则，确保组合的断层分布符合地质规律。
在断控约束的过程中以下几个问题需要解决:
(1)断面之间削截关系处理:复杂断层系统中往往存在着多条断层，必然会对断层网建模以及后续的地层模型产生不利影响，因此必须考虑多个断层面的削截关系。断层削截是指断层在纵向上发生相交情况时，就要削截断层。断层之间的接触关系采取分层次处理，即首先处理相交的断层面之间的接触关系，然后交互处理其他断层之间的关系，实现断层网的快速更新。
(2)断面穿过层位处理:区域内的断层发育并不一致，有些断层断开多个层位，有些断层只断开一个层位。在做断层网时，对于纵向上断穿其他层位的断层，要对其轮廓线进行重新绘制，使断层不穿透层位，防止骨架网格模型混乱。
(3)断点与断面匹配:受地震数据分辨率限制，可能出现断点与相应的断层面有一定距离，这时需要调整断层模型使其与断点吻合。如果断点与断层面距离太大，则可能是断点本身存在问题，这时需要重新进行地层对比，找出原因，使断层模型与实际断点吻合。
(4)层面模型建立:层面建模就是生成关键层面模型。通过一定的插值算法，将地震解释层位进行插值加密，生成一张覆盖全区的规则网格。具体实现时可采用普通克里金插值算法，插值时需考虑断面对原始数据点选取的影响，将断层两侧一定范围内的数据点忽略。
模型测试如下，图11为断控约束前后的速度模型，从图中对比可以发现加入断控后速度模型的形态发生了一定改变，更加接近地下断层的实际结构。
图12为断控约束前后的剖面对比，从图中对比可以发现断层部分更加清晰且连续，明显改善了成像质量。
图11 断控约束前(左)、后(右)速度模型对比图
Fig.11 Comparison of velocity models before(left)and after(right)the fault control constraint
图12 断控约束前(左)、后(右)剖面对比图
Fig.12 Comparison of profile before(left)and after(right)of the fault control constraint
2.3 实际资料测试
采用南川地区三维实际资料，为了观察融合细节，抽出其中一个构造清晰的二维剖面，其横向采样点1651，纵向采样点801，属于复杂山地构造区，测试结果如图13所示。
图13 (a)~(f)浅层速度模型剖面，中深层初始速度模型剖面，融合剖面，断层与层位约束，融合后模型，断控约束层析后三维速度模型
Fig.13 (a)~(f)Shallow velocity model profile，Initial velocity model profile in middle and deep layers，Post-fusion velocity profile，Fault and hierarchical constraints，Post-fusion model，3D velocity model after Fault control constraint
由图中融合后速度模型可以看到在浅层部分的速度形态已经和初至波层析速度反演的形态一致，且速度值在浅层部分与其保持一直。中深层部分的速度值与中深层初始模型一致，证明融合结果正确。三维走时层析反演后速度模型的形态发生一定改变，速度值的正确性由图14的道集与剖面验证如下。
图14 (a)~(f)中深层初始角道集，融合与全局层析后角道集，L1720线融合前后偏移剖面，L1480线融合前后偏移剖面
Fig.14 (a)~(f)Initial angle gathers in the middle and deep layers，angle gathers after fusion and global tomography，migration profile before and after fusion of line L1720，migration profile before and after fusion of line L1480
由图14中角道集的拉平程度可以看出融合后的全局层析得到的速度模型基本正确，(e)、(f)为三维工区中两条测线，通过观察可以明显发现剖面在浅层部分的同相轴质量明显提升，同相轴更加连续且聚焦，与此同时由于浅层速度的改善，深层部分的成像也有一定程度的提高，偏移结果符合预期。
3 结论
本文提出了一种将浅中深层速度融合后的三维全局层析反演方法，近地表的浅层速度建模使用初至波层析反演得到精度较高的浅层速度场，并通过优化拾取的初至、多次迭代等手段可以使反演出的速度模型有效深度增加至预期层位[16-19]，并根据浅、中、深层速度模型信息随深度的变化情况出发，将初至波层析反演出的浅层速度场与常规速度分析出的中深层速度场进行三角函数加权与反比例加权融合，通过融合层位来确定并修改速度模型以及选取不同深度融合的权值;最后运用优化后的角道集剩余曲率拾取方法来进行全局层析。此方法在复杂山地的速度建模有良好的应用前景。
利益冲突: 作者声明无利益冲突。
[①] *通讯作者 Corresponding author:王非翊，fei.y.wang@stu.pku.edu.cn
收稿日期:2023-03-25; 录用日期:2023-03-27; 发表日期:2023-03-28
资助项目:中国石化科技开发部项目“常压页岩气地球物理评价技术研究”(P21087-3)资助
参考文献(References)
[1] 秦宁, 李振春, 杨晓东. 基于角道集的井约束层析速度反演[J]. 石油地球物理勘探, 2011, 46(5): 725-731.
https: //doi. org/10. 13810/j. cnki. issn. 1000-7210. 2011. 05. 011
[2] 李振春, 王华忠, 马在田. 共中心点道集偏移速度分析[J]. 石油物探, 2000, 39(1): 20-26.
[3] 王延光, 尚新民, 赵胜天, 等. 山前带联合层析反演地震成像[J]. 石油地球物理勘探, 2021. 56(4): 782-791.
https: //doi. org/10. 13810/j. cnki. issn. 1000-7210. 2021. 04. 011
[4] 蔡杰雄, 王静波. 一种基于改进快速扫描法的多尺度近地表层析方法[J]. 石油物探, 2019, 58(6): 819-827.
https: //doi. org/10. 3969/j. issn. 1000-1441. 2019. 06. 004
[5] 王川, 李振春, 李文燕, 等. 复杂近地表三维初至波走时层析方法研究[J]. 地球物理学进展, 2018, 33(5): 1967-1973.
https: //doi. org/10. 6038/pg2018BB0216
[6] 徐嘉亮, 张冰. 基于层位约束的不规则网格层析速度建模方法[J]. 地球物理学进展, 2020, 35(2): 669-673.
https: //doi. org/10. 6038/pg2020BB0400
[7] 郭振波, 孙鹏远, 钱忠平, 等. 快速回转波近地表速度建模方法[J]. 石油地球物理勘探, 2019, 54(2): 261-267.
https: //doi. org/10. 13810/j. cnki. issn. 1000-7210. 2019. 02. 003
[8] 王华忠, 张兵, 刘少勇, 等. 山前带地震数据成像处理流程探讨[J]. 石油物探, 2012, 51(06): 574-583.
https: //doi. org/10. 3969/j. issn. 1000-1441. 2012. 06. 005
[9] 崔栋, 张研, 胡英, 等. 近地表速度建模方法综述[J]. 地球物理学进展, 2014, 29(6): 2635-2641.
https: //doi. org/10. 6038/pg20140625
[10] 张兵, 徐兆涛, 王华忠, 等. 山前带地震数据共成像点道集层析速度反演建模方法研究[J]. 石油物探, 2012, 51(6): 590-597.
https: //doi. org/10. 3969/j. issn. 1000-1441. 2012. 06. 007
[11] 张兵, 王华忠. 数据域初至波走时与成像域反射波走时联合层析速度建模方法[J]. 地球物理学报, 2019, 62(7): 2633-2644.
https: //doi. org/10. 6038/cjg2019M0536
[12] 王潇, 杨锴. 基于二维VTI拟声波程函方程的椭圆各向异性立体层析反演[J]. 地球物理学报, 2018, 61(12): 4954-4964.
https: //doi. org/10. 6038/cjg2018M0212
[13] 张建明, 董良国, 王建华. 波动方程初至波多信息联合反演方法[J]. 地球物理学报, 2021, 64(7): 2447-2460.
https: //doi. org/10. 6038/cjg2021O0278
[14] Billette F, Lambare G. velocity macro-model estimationfrom seismic reflection data by stereotomography [J]. Geophysical Journal International, 1988, 135(2): 671-690.
https: //doi. org/10. 1046/j. 1365-246X. 1998. 00632. x
[15] Siwei L, Vladimirsky A, Fomel S. First-break traveltime tomography with the double-square-root eikonalequation [J]. Geophysics, 2013, 78(6): U89-U101. doi. : 10. 1190/geo2013-0058. 1
[16] Taikun S, Jianzhong Z, Zhonglai H, et al. A layer-strippingmethod for 3D near-surface velocity model buildingusing seismic first-arrival times[J]. Journal ofEarth Science, 2015, 26(4): 502-507.
https: //doi. org/10. 1007/s12583-015-0569-0
[17] Bergman B. Seismic tomography studies of cover thicknessand near-surface bedrock velocities [J]. Geophysics, 2006, 71(6): 1ND-Z126.
https: //doi. org/10. 1190/1. 2345191
[18] Stork C. Reflection tomography in the postmigrated domain[J]. Geophysics, 1992, 57(5): 680-692.
https: //doi. org/10. 1190/1. 1443282
[19] Brandsberg-Dahl S, de Hoop V, and Ursin B. Velocity analysisin the common scattering-angle/azimuth domain [C]. SEG Expanded Abstracts, 1999: 1715-1718.
https: //doi. org/10. 1190/1. 1820865
[20] Wang B, Dirks V, Guillaume P, et al. A fast and lowcost alternative to subsalt wave equation migration perturbationscans [C]. SEG Expanded Abstracts, 2005: 2375-2378.
https: //doi. org/10. 1190/1. 2148166
[21] Woodward M, Nichols D, Zdraveva O, et al. A decade oftomography [J]. Geophysics, 2008, 73(5): VE5-VE11.
https: //doi. org/10. 1190/1. 2969907
[22] Biondi B, and Sava P. P Wave-equation migration velocityanalysis [C]. SEG Expanded Abstracts, 1999: 1723-1726.
https: //doi. org/10. 1190/1. 1820867
[23] Nurul Kabir, Uwe Albertin, Min Zhou, et al. Use of refraction, reflection, and wave-equation-based tomography forimaging beneath shallow gas: A Trinidad field data example[J]. Geophysics, 2008, 73(5): VE281-VE289.
https: //doi. org/10. 1190/1. 2953064
[24] Hongbu Zhou, Yu Zhang, Samuel H Gray, et al. Regularizationalgorithms for seismic inverse problems[C]. SEG Expanded Abstracts, 2002: 942–945.
https: //doi. org/10. 1190/1. 1817420
[25] Husen S, and Kissling E. LOCAL Earthquake tomographybetween rays and waves: Fatray tomography[J]. Physics of the Earh and Planetary Interiors, 2001, 125(1-4): 171-191.
https: //doi. org/10. 1016/S0031-9201(00)00206-5
[26] Weiwei G, Fengyuan X, Xiao W, et al. Anisotropicslope-tomography in VTI media preconditioned by animage guided smoothing filter[C]. SEG Expanded Abstracts, 2019: 5199-5203. doi. org/10. 1190/segam2019-3215621. 1
[27] Sergey Fomel. Time-migration velocity analysis by velocitycontinuation[J]. Geophysics, 2003, 68(5): 1446-1761.
https: //doi. org/10. 1190/1. 1620640
[28] Borhan F, Stéphane O, Alessandra R, et al. Matrix-freeanisotropic slope tomography: Theory and application[J]. Geophysics, 2019, 84(1): R21-R43.
https: //doi. org/10. 1190/geo2018-0112. 1
[29] 卢江波, 龙海滨, 张铭, 等. 初至波层析成像正演模拟中多发多收模式计算方法研究[J]. 地球物理学进展, 2020, 35(4): 1590-1594.
https: //doi. org/10. 6038/pg2020DD0291
[30] 潘奕铭, 桑运云, 姚雪峰, 等. 一种多信息约束的初至波走时层析反演优化方法[J]. 石油物探, 2023, 62(2): 271-280.
https: //doi. org/10. 3969/j. issn. 1000-1441. 2023. 02. 008
[31] 周衍, 饶莹. 黄土塬覆盖区的层析反演静校正方法研究[J]. 地球物理学报, 2019, 62(11): 4393-4400.
https: //doi. org/10. 6038/cjg2019N0304
A 3D Tomographic Inversion Method Based on Velocity Model Fusion with Fault Constraints
(1. Research Institute of Petroleum Exploration & Development, Sinopec East China Oil and Gas Company, Nanjing 210000, China
2. School of Earth and Space Sciences, Peking University, Beijing 100871, China)
Abstract: Conventional grid-tomographic inversion is based on the travel time of the reflected waves and the residual curvature of the co-imaged points to construct vector equations and thus update the velocities, however, in the complex surface region, factors such as low signal-to-noise ratio of a single shot and insufficient reflection information lead to insufficient accuracy of velocity modelling in the near-surface region, and the lack of layers and fault constraints, which can not accurately restore the actual structure of the deeper layers. The current method widely used in industry is to fuse the shallow velocities obtained from the tomographic inversion by first breaks with the surface elevation plane at a certain depth downward as the fusion surface and the mid-deep velocities obtained from the grid tomography, but this method will show obvious fusion marks in the fusion area and the velocity values at the fusion point will change abruptly. For this reason, this paper starts from the principle of tomography inversion by first breaks, picks up the fusion surface based on the first-to-wave ray density and the difference between the two models, and constructs a fusion function for velocity fusion that fits the geological background of the work area, so as to eliminate the problem of fusion trace and abrupt changes in velocity values at the fusion trace to the maximum extent, and add fault constraints during the subsequent stratigraphic inversion. This method uses velocity fusion, fault-control constraints and improved residual curvature pickup method to improve velocity modelling and imaging accuracy in complex mountain regions. The complex model results prove the correctness, and practical data test results prove the applicability.
Keywords: Velocity fusion, layer constraints, extraction of 3D ADCIGS, global tomography, migration imaging
DOI: 10.48014/cpngr.20230425001
Citation: XIONG Yingjie, WANG Feiyi, MENG Qingli, et al. A 3D tomographic inversion method based on velocity model fusion with fault constraints [J]. Chinese Petroleum and Natural Gas Research, 2023, 2(1): 1-14