断控层析速度建模方法研究与应用

邓勇， 李芳， 彭海龙， 任婷， 刘国昌

1 中国石油大学(北京)油气资源与探测国家重点实验室， 北京 102249 2 中海石油(中国)有限公司湛江分公司， 湛江 524057 3 中海石油(中国)有限公司海南分公司，海口 570000

0 引言

随着勘探开发力度的逐步加大，准确获取地下地质构造形态及展布特征是当前油气藏开采工作的重要基础.尤其在复杂断块区域，真实的地震成像是推动区域油气藏勘探开发的重要保障(张昌榕等，2018).然而，常规的成像方法难以获取该区域准确的成像，使复杂断块区域的勘探开发陷入瓶颈.为更好的应对上述挑战，目前采用叠前深度偏移(Prestack Depth Migration，PSDM)方法解决上述问题(秦宁，2020).

在叠前深度偏移方法中，速度模型的精度制约着偏移成像的准确度，常用的速度建模方法包括速度分析法、层析反演方法、波动方程类方法等(王艳香等，2021).分析类方法存在条件假设，速度分析结果精度低，一般只用于初始速度建模(刘定进等，2019).波动方程方法的反演结果依赖于初始速度模型的精度，且需要多次迭代才能完成速度反演，效率较低(曹呈浩等，2020).而层析反演方法应用过程简单方便，具有较好的反演精度和反演效率，是当前应用最为广泛的速度建模方法(徐嘉亮等，2021).

传统的射线类走时层析反演方法一般采用规则网格单元的模型参数化方式表征空间坐标问题，在复杂地质构造速度场反演过程中，需要对网格进行精细剖分(郭振波等，2019)，但是在不规则阻抗界面区域中，该方法容易出现误差，给反演结果带来假象(徐嘉亮等，2020).目前常规解决策略是采用不规则网格单元进行模型参数化，可以较好地表征不规则异常体的波阻抗界面(何雷宇等，2017).但上述方案并没有完全解决断层引起的速度差异性问题，限制了断层区域速度建模精度的提升.除此之外，应用模型正则化方法，将地质构造特征引入层析进行约束，能够提升层析质量，但是在复杂断块区域中，该方法由于缺少断层的约束控制，反演结果中缺少高波数的速度成分，因此该方法有待进一步提升速度反演精度(李辉等，2015).

为解决断层区域的速度精度问题，Birdus(2007)以断层作为条件，在密集的叠前深度偏移道集网格上进行非双曲剩余曲率分析，然后约束速度更新过程.该方法需要高密度道集，且只能在与断层对应的有限区域内进行约束控制.Rodriguez等(2011)将断层上下盘分为两个区域分别进行速度反演，但该方法仅更新断层附近的速度，且在更新过程中难以实现有效迭代.Guo等(2015)提出用断层附近道集的剩余曲率来约束层析正则化过程，能够更新整个模型.但是准确拾取断层附近反射信息的剩余曲率较为困难，因此，该方法在应用中存在一定的挑战.

为更好地改进层析成像在复杂断块区域速度建模的精度问题，本文在模型正则化层析的基础上，提出一种基于断层控制的层析方法，该方法对成像过程中的断层进行识别判断，利用断层数据构建预条件矩阵，然后采用TV正则化方法改进目标泛函，从而避免层析过程中断层两侧速度的平滑过渡，提升断层两侧的速度场精度，得到分辨率较高的断块区域速度.经实际钻井数据验证，该方法得到的速度场和实际速度场匹配度高，具有较好的应用效果.

1 方法理论

在地震反演过程中，多个不同的解都适合作为同一个地球物理问题的解，使得反演问题具有不确定性.为解决上述问题，层析反演一般采用先验约束降低多解性，常用的方法是正则化(张力起等，2019).正则化是将逆问题的先验信息和逆解估计的期望加入到反演过程中，从而消除反演过程的病态问题，提高反问题求解的鲁棒性和适定性(代荣获等，2019).地球物理反问题中，通过增加正则约束来解决病态反演问题取得了很好的效果.尽管无法获得描述真实地下介质的完全精确解，但可以得到满足应用需求的近似最大似然解(苏扬等，2021).

正则化技术是获得有地质意义解的有效方法，按照应用条件一般可分为模型正则化和数据正则化.模型正则化利用地下构造信息约束模型参数，实现改善估计模型参数的目的(张兵和王华忠，2019)；数据正则化是在层析中加入空间介质参数的数据关联性(李辉等，2015)，但是收敛效率较低，收敛效果较差，目前应用较少.按照正则化函数的不同，正则化技术可以分为Tikhonov正则化和总变差(TV，Total Variation)正则化.相比较而言，Tikhonov正则化计算过程稳定，但是难以保持层析过程中的速度边界，TV正则化能够较好的保持速度的边缘信息，但是稳定性较差(王永奎等，2022).

在反演过程中引入模型的相关认识，即构成模型正则化.模型正则化的实现策略主要包括Tikhonov模型正则化和预条件模型正则化，预条件模型正则化方式收敛性又优于模型约束的Tikhonov模型正则化，且计算效率较高，因此本文采用预条件正则化进行层析矩阵构建，得到精度相对较高的初始速度.然后以此为基础，对成像数据中的断层进行识别判断，通过识别的断层数据建立预条件矩阵并采用TV正则方法改进目标泛函，实现速度场的更新，得到复杂断块区域的高精度速度场.

(1)预条件模型正则化

层析反演采用匹配射线走时的方法来得到速度值，应用射线追踪方法拟合地震数据中的反射信息走时，并通过走时差的反投影得到速度的更新量(范白涛等，2022).对于实际地震数据而言，拾取的地震走时可以表示为：

(1)

式中，τobs表示拾取的走时，s表示观测数据慢度，dl表示沿着射线路径ab的射线长度.如果已知背景慢度，则计算的走时可以表示为：

(2)

式中，τcal表示计算的走时，s0表示背景慢度，则走时扰动和慢度扰动具有如下关系：

(3)

将式(3)转为矩阵，则基于射线的线性走时层析反演表达式为：

Lδs=δτ，

(4)

式中，L表示线性化算子.由(4)式可知，每一条射线都存在一个对应的公式，所有射线所代表的公式则构成一个射线层析方程组，因此基于贝叶斯反演理论的最大后验概率密度思想，层析反演的目标泛函为：

(5)

层析反演在迭代过程中，采用线性化方式对目标函数进行优化，则线性迭代方程为：

(6)

式中，矩阵A表示A(m)的线性化算子，Δm表示当前迭代次数下的模型更新量，Δd表示当前迭代中的数据残差.

通过公式(5)和公式(6)中的数据协方差矩阵和正则化参数，能够对层析反演过程实现正则化.在不考虑正则化的情况下，即协方差矩阵为单位矩阵且阻尼因子为0，意味着仅用数据残差进行层析反演，则公式(6)退化为一般的层析方程：

Δd=AΔm.

(7)

ASu=Δd，

(8)

Δm=Su，

(9)

式中：S表示预条件算子，u为预条件的解.考虑阻尼因子ε，则层析方程的阻尼最小二乘方程可以表示为：

STATASu+εu=STATΔd.

(10)

将预条件下的解u=S-1A-1Δd代入方程(10)并适当变换得到：

SSTATAΔm+εΔm=SSTATΔd，

(11)

对于公式(11)而言，如果将地质构造信息引入预条件算子S，则公式(11)就是构造约束正则化的层析方程.所以，提取地质构造信息，并将其引入预条件算子所表示的光滑矩阵是上述步骤的关键点.

在同一个空间坐标系下，不同的散射点的位置和倾角信息具有较大的差异性，为构建一个空间上的光滑矩阵，可以通过坐标变换的方式将上述地质信息融入到变换后的光滑矩阵中(万弘等，2017)，空间坐标的转换示意图如图1所示.

图1 坐标变换示意图Fig.1 Schematic diagram of coordinate transformation

图1中，逆时针方向表示旋转的正方向，局部地质坐标系原点为散射点中心，记为T(u,v,w)，u轴的方向就是反射界面的走向，w轴方向是反射界面垂直方向；空间坐标系记为T′(x,y,z).然后参照李辉等人提出的坐标转换方法，空间坐标经过一系列的旋转、平移后变换为局部地质坐标系(李辉等，2015).记散射点的原始空间坐标为T′(x0,y0,z0)，则旋转后的局部地质坐标可以表示T(u0,v0,w0)，两个坐标系存在以下关系：

采访正值“第十九届国产高浓度磷复肥产销会”期间，杨三可来到瓮福集团展台，向国内外客商推介瓮福集团“DAP+”、磷酸二氢钾、聚磷酸铵等新型增效肥料。

(u0,v0,w0)=T(x-x0,y-y0,z-z0)，

(12)

式中，T为旋转矩阵，满足条件：

(13)

式中：φ和θ分别表示坐标变换中x和z坐标轴的旋转角度.两参数在地质意义中分别表示地层方位角和倾角.

光滑矩阵的某一行元素就是一个光滑函数，用于描述不同介质点之间的联系，令此光滑函数为高斯光滑函数，则光滑矩阵中第i行第j列的元素为：

(14)

式中，σui,σvi,σwi分别是高斯函数在局部地质坐标系中不同方向上的标准差.实际应用过程中，通常默认地下介质为均匀层状结构，认为反射界面方向模型参数基本一致，垂直于反射界面方向的模型参数变化剧烈，因此水平方向对应的光滑范围相对较大，即σui、σvi要大于σwi，则相应的标准差可以通过式(15)、(16)计算：

σui=σvi=abs(sin(φ))，

(15)

σwi=ξσui，

(16)

式中，ξ表示参数因子，通常取0.2～0.5的一个常数，至此将地质构造信息引入光滑矩阵.由于高斯函数的方差求取简单快速，可以作为一个已知量，所以构建光滑矩阵的关键点就是求取地质构造的倾角和方位角信息.

(2)地层倾角方位角计算

在成像剖面中，可以拾取主要反射同相轴作为某一套地层的顶底界面，以此为基线进行界面之间层位的插值，通过计算顶底界面的相关地层信息，然后结合地层接触关系则可以依次得到整套地层的倾角和方位角信息.开展上述工作时，首先要判断顶底界面之间地层的接触关系，通过地层关系实现整套地层信息的计算.根据地层沉积规律，顶底界面之间的地层应当存在连续沉积、侵入式沉积、剥蚀沉积等地层关系，具体如图2所示.

确定好顶底地层接触关系后，顶底地层之间的其他地层可以通过插值进行恢复，从而建立整个工区的地层架构.地层架构一旦确立，可以根据几何定理得到人工解释层位数据的倾角方位角，地层倾角方位角的几何关系如图3所示.图3中，θx表示x方向的视倾角，θy表示y方向的视倾角，θ表示地层的真倾角，ω表示地层的方位角.利用人工解释地层数据信息，分别在x方向和y方向建立网格，目标点位于中心，解释数据位于网格交点上，则数据位置分布如图4所示.

通过图3和图4，则目标点的视倾角、真倾角和方位角的公式为：

tanθx=Δx/Δz，

(17)

tanθy=Δy/Δz，

(18)

(19)

(20)

式中，Δx、Δy和Δz分别表示目标点和四周8个点不同方向上坐标的差值加权结果.完成顶底解释地层的方位角倾角计算后，则解释地层之间的地层倾角方位角可以采用图5的原理进行计算.

图2 不同沉积模式示意图(a) 等比例连续沉积； (b) 等厚底面侵入沉积； (c) 等厚顶面剥蚀沉积.Fig.2 Schematic diagram of different sedimentary modes(a) Equal proportion continuous deposition； (b) Equal thickness bottom intrusive sedimentation； (c)Equal thickness top denudation deposition.

图3 地层倾角方位角几何信息示意图Fig.3 Schematic diagram of geometric information of formation dip angle and azimuth angle

图4 层位网格示意图Fig.4 Schematic diagram of horizon grid

图5 不同沉积模式下地层倾角计算示意图(a) 等比例连续沉积地层信息关系； (b) 等厚底面侵入沉积地层信息关系； (c) 等厚顶面剥蚀沉积地层信息关系.Fig.5 Calculation diagram of formation dip angle under different sedimentary modes(a) Information relationship of equal proportion continuous sedimentary strata； (b) Information relationship of sedimentary strata invaded by equal thickness bottom； (c) Information relation between denudation sedimentary strata of equal thickness top surface.

在图5a表示等比例连续沉积，设定顶界面在平面坐标某位置处的角度为φ1，同一平面坐标位置处，底界面角度为φ3，则同一平面坐标位置处的中间某一深度的地层角度φ2为：

(21)

式中，H1表示顶界面的深度值，H3表示底界面的深度值，H2表示目标位置处的深度值，在图5b中，地层接触关系为等厚底面侵入关系，则中间地层的地质信息和顶界面地层一致.在图5c中，展示的是等厚顶面剥蚀地层的接触关系，中间地层的地质信息和底界面地层一致.结合图5所示的不同地层接触关系，运用公式(21)可以计算整个工区的构造信息.

(3)断层控制正则化

公式(11)引入构造信息后，可以得到较为准确的速度反演结果，但断层区域速度精度依旧有待提升，因此如何构建一个断层控制的目标泛函是关键.

预条件模型正则化方法会使得断层区域的速度边界模糊，难以保持断层两侧的速度差异性，因此有必要采用新方法提高断层两侧速度的精度.考虑TV正则化具有保持边缘信息的特征，可以在预条件模型正则化的基础上引入TV正则化.TV正则化在重构模型参数过程中，能够得到更为准确和丰富的参数信息，但是需要异常边界信息的支持，从而加快收敛速度，提高收敛精度，因此可以在重构算法中引入断层边界约束控制方法.根据当前迭代成像结果中断层边界位置信息，构建权重矩阵并将其用于下一次迭代，从而形成边界约束的TV正则化方法，则改进后正则化的目标泛函可以表示为：

(22)

式中，‖‖1表示L1范数，D表示三维梯度算子稀疏矩阵，其主要作用是对速度模型中的每一个速度值做梯度运算，在每次迭代过程中，使探测得到的参数边界避免正则化，从而保持该位置附近的参数边界；S表示预条件光滑矩阵.然而，在实际计算过程中，TV正则化存在着不稳定性，因此需要将公式(22)转为线性化方式进行求解，然后通过多次迭代并不断修正模型参数，使最终的反演结果不断逼近真实值.

对于模型参数的修正，需要采用迭代加权方式求解，通过对一系列权系数不断更新的加权最小二乘问题进行计算，得到目标函数的最优解.采用加权方式将L1范数松弛转为L2范数，则公式(22)转变为：

(23)

式中，ε和β表示正则化参数，取值为[0,1]，用于调节不同约束项的权重，w表示对角矩阵，其公式可以表示为：

(24)

对于一个Nx×Ny×Nz的层析网格，当i=1,2,…,Nx，j=1,2,…,Ny，k=1,2,…,Nz时，则梯度算子矩阵可以表示为：

(25)

由实际地质特征可知，断面左右的速度具有一定的差异性，要想保持断层两侧速度突变的真实性，就需要公式(23)在计算速度值梯度算子时不能跨越断层两侧，也就是位于断面网格位置处的速度值不能和断面两侧的速度值同时参与梯度计算，因此可以将上述断层代表的边界信息以权重矩阵的形式参与下一次迭代.即将此边界与公式(23)正则化过程结合，基于上述思想，则构建的对角权重数组表达式为：

(26)

式中P(xi,yj,zk)表示xi,yj,zk位置处网格坐标，F(Φ)表示解释的断层数据坐标集，W表示三维矩阵.

由公式(26)可知，矩阵的建立不受断层类型的影响，一旦获取到断层的位置坐标，就可以建立对角权重矩阵.定义断层控制的三维诱导矩阵Γ为：

(27)

将上述权重矩阵公式代入公式(23)，则公式变为：

(28)

随着速度场更新，断层在空间上会发生飘移，断层更新后的实际位置已经偏离最初解释的断层位置，因此需要对断层位置进行推测识别，采用新的断层信息对层析过程进行约束控制.

为提升断层的识别精度和提取效率，本文采用各向异性扩散方程中的梯度结构张量来提取地质构造特征，通过特征向量的情况来识别断层，则三维网格中梯度结构张量的计算公式为：

(29)

(30)

图6 特征值和空间地质结构关系图Fig.6 Relationship between eigenvalues and spatial geological structure

为提升流程效率，只需对目标区域的主断层区域进行断层控制模型更新，因此需要应用人工解释的断层数据F(Φ)对计算得到的断层位置集合Fn进行预估，其中F表示断层数据(三维)集合，Φ表示人工解释断层，n表示第n次迭代.速度模型的更新对于断层位置偏移量影响相对较小，因此可以定义偏移量来时预估新的断层位置集，如果推测出的断层位置集合Fn∈F(Φ+Δf)，则认为是目标断层区域.根据该处定义可知，即使人工解释断层刻画位置存在误差，只要在一定范围内，可以通过后续更新后的断层来约束断层位置.其中Δf表示定义的空间位置误差项，然后将计算得到的断层数据引入公式(24)和公式(27)，并结合公式(23)，则建立的断层控制矩阵的目标泛函可以表示为：

(31)

LW(mn)=(TΓn)Tdiag[max(ΓnTmn,γ)p-2]ΓnT，

(32)

+βLW(mn)mn，

(33)

(35)

(36)

通过上述方法，速度更新过程可以很好的保护断层区域速度的差异性，可以解决常规层析反演方法在该区域中遇到的难题，得到断层区域较为准确的速度值.

2 理论模型

为验证本文方法，建立一个含有复杂断层构造的地质模型，该模型含有多个断层，然后在此基础上开展正演，进行速度分析后执行叠前深度偏移.其中，真实速度模型剖面以及采用常规速度层析方法得到的偏移结果剖面如图8所示.

图7 断控层析速度建模流程算法示意图Fig.7 Schematic diagram of fault constrained tomography flow algorithm

对比图8可知，常规层析速度建模方法得到的成像剖面，其构造形态和原始速度模型大体一致，但是断层下方区域成像质量较差.因此该区域的成像质量有待进一步改善.在图8成像剖面的基础上，解释层位和断层数据，然后计算地层的倾角(所有地层方位角为同一数值)，则计算的地层倾角数据剖面如图9所示.在上述解释数据、倾角和方位角数据的基础上，建立断层约束算子，则速度更新前后的剖面对比如图10所示.最后将基于断控层析速度建模方法得到的速度应用于PSDM，则两种不同速度反演方法得到偏移成像结果如图11所示.

对比图10和11可知，采用断控层析反演方法得到速度，能够更好的描述断层区域的速度分布情况，断层两边的速度截断和实际情况也较为匹配.采用断控层析速度进行PSDM，断面成像更清晰，断层下部的地层成像更为连续，成像质量更好.

3 实际数据应用

选取南海某工区的实际地震数据进行测试，该区域断裂系统发育，断裂结构复杂，导致很难获取该工区准确的速度场，从而使得偏移成像结果不可靠，严重影响该区域的勘探开发进程.为解决上述问题，本文采用断控层析速度建模方法获取准确的速度场，从而有效提高PSDM成像准确度.

图8 速度模型及PSDM成像剖面Fig.8 Velocity model and PSDM imaging profile

图9 计算的地层倾角数据剖面Fig.9 Calculated dip data profile

(1)地层倾角方位角数据的计算

为准确计算地层的倾角和方位角，首先需要根据地质沉积规律来准确判断地层之间的接触关系，其次分析区域的实际地质情况，考虑地层和断层、地层与地层、断层与断层之间的接触关系，图12展示一条二维线上主要断层解释数据分布情况.

由图12可知，该区域地层和地层之间都是等比例沉积，所以只需要计算出区域主地层的倾角和方位角即可，然后利用公式(17)计算该区域所有地层倾角和方位角，计算的地层倾角和方位角如图13和图14所示.

图10 不同速度反演方法结果对比(a) 常规层析方法反演结果； (b) 断控层析反演结果.Fig.10 Comparison of results of different velocity inversion methods(a) Inversion results of conventional tomography； (b) Fault controlled tomography inversion results.

图11 不同速度建模方法的到的PSDM剖面(a) 常规层析反演方法得到PSDM剖面； (b) 断控层析反演方法得到PSDM剖面.Fig.11 PSDM profiles with different velocity modeling methods(a) PSDM profile obtained by conventional tomography inversion method； (b) PSDM profile obtained by fault controlled tomography inversion method.

图12 解释数据剖面Fig.12 Interpretation data profile

图13 计算的地层倾角剖面示意图Fig.13 Calculated dip profile

在速度更新过程中，需要采用梯度结构张量对断层的位置进行重新识别.随着速度的变化，部分断层的位置发生飘移，但是整体偏移位置不大.经过速度更新后，采用梯度结构张量来识别工区内的断层，则本研究工区中识别的断层位置和原始人工解释断层的叠合情况如图15所示.

由图15可知，速度更新后，梯度结构张量自动识别的断层位置和原始人工解释断层位置存在一些差异，但是位置的差异量不大.且根据识别结果可知，在部分区域，即使没有人工解释断层，很多小断层也被识别出来.为了尽可能提升计算效率，目前只需将解释断层附近被识别出断层数据引入层析过程.根据解释断层的位置，选取附近更新后并被识别出的断层数据，将其用于断控层析速度建模，从而实现速度的更新.对于识别出来其他非目标区域的断层，其断层数据则不用建立断层约束算子.根据上述方式，可以利用初始解释断层数据匹配自动识别的断层，采用新的断层数据建立断层约束算子，将其引入改进的目标函数用于速度场的迭代更新.

图14 计算的地层方位角剖面示意图Fig.14 Calculated stratigraphic azimuth profile

图15 层析反演前后的断层边界示意图Fig.15 Fault boundary before and after tomographic inversion profile

(2)断控层析速度建模以及偏移成像

利用提出的断层控制正则化层析速度反演方法对南海复杂断块区进行速度建模，在此基础上进行叠前深度偏移.图16和图17分别展示了常规层析速度、断控层析速度的剖面对比和水平切片.

断层区域由于断层两侧对接岩性的年代差异，会表现出速度突变现象.然而，在图16a和图17a中，断层两侧速度是平滑渐变，断层两侧速度差异不明显.断层约束层析反演速度建模结果中，断层两侧速度突变明显，断面等其他物理性质影响体现的更为准确.为进一步分析两种层析结果，对其成像结果再次对比，结果如图18所示.

对比图18可以发现，在常规方法的成像剖面中，断层区域(蓝色线框和黑色线框区域)的反射同相轴出现扭曲错段的现象，甚至出现一条垂直的虚假断层，上述现象的存在，严重干扰对断层区域地质状况的认识.经过分析可知，常规速度建模方法导致断层附近的速度不准确，使得成像剖面出现虚假的垂直断层和地层扭曲错断的现象.而在断控层析反演成像结果中，虚假的垂直断层消失，反射同相轴更加清楚和平整，表明断层附近的速度更为准确，剖面成像质量提升明显.

图16 不同速度建模方法的地震-速度叠合剖面(a) 常规层析速度； (b) 断控层析速度.Fig.16 Profile of different velocity inversion methods(a) Conventional chromatographic velocity；(b) Fault constrained tomography velocity.

图17 不同速度建模方法的地震-速度叠合水平切片(a) 常规层析速度； (b) 断控层析速度.Fig.17 Time slices of different velocity inversion methods(a) Conventional chromatographic velocity； (b) Fault constrained tomography velocity.

图18 不同速度建模方法的成像剖面对比示意图(a) 常规层析PSDM剖面； (b) 断控层析的PSDM剖面.Fig.18 Comparison diagram of imaging sections with different inversion velocities(a) Conventional tomographic PSDM profile；(b) PSDM profile of fault controlled tomography.

为进一步分析本文方法对于成像的影响，提取两种方法偏移成像结果的相干体属性进行对比，结果如图19所示.

从相干体切片中，断控层析反演速度建模方法的成像结果优势更加明显.在图19a的常规层析速度成像结果中，圆圈中断层的识别较为模糊，地层的反射也较为紊乱，难以识别断层的组合.然而，在图19b结果中，地层反射规律性强，断层位置识别清晰.

为检查速度反演质量，进一步量化反演速度结果，利用穿过断层的某实际钻井数据进行验证，井轨迹位置如图20所示.然后采用相同的地震子波和不同的地震速度模拟地震记录，将模拟结果和测井速度的合成地震记录进行相关，通过相关性来判断当前速度的精确度.图21展示的是测井数据合成记录同两种速度的正演结果对比，结果表明，断控层析速度的合成地震记录和测井地震记录的相关系数为0.75，常规层析速度的合成地震记录和测井地震记录的相关系数为0.58，证实断控层析速度建模方法的速度结果和实际钻井速度吻合度更高，最符合实际速度场分布情况.

4 结论

本文提出了一种基于断层控制的层析反演速度建模方法，该方法通过地震解释数据计算相关模型参数，建立断层约束的正则化约束算子，以此为基础建立断层控制的层析反演目标泛函，利用断层数据作为速度建模的先验信息，采用迭代最小二乘优化方法逼近L1范数的目标泛函，提高速度建模的精度.应用本文方法获取的速度模型，断层区域的速度分布更加符合地质认识，模型与真实情况吻合度更高.将反演得到的速度场用于PSDM，提高了复杂断块区域的成像质量，为后续目标评价提供了可靠的地震成像数据，从而有助于提高复杂区油气勘探成功率.由于相关研究条件限制，本方法的模型参数表征采用了矩形网格方式，因此在描述模型边界的精度方面依旧存在提升空间，下一步研究计划采用三角网格进行模型参数表征，开展层析反演的正则化研究工作.

图19 不同反演速度的成像剖面相干体提取结果对比(a) 常规层析成像结果提取的相干体； (b) 断控层析成像结果提取的相干体.Fig.19 Comparison of coherent volume extraction results of imaging profiles with different inversion velocities(a) Coherent volume extracted from conventional tomography results； (b) Coherent volume extracted from fault control tomography results.

图20 穿越断层的井轨迹示意图Fig.20 Schematic diagram of well trajectory crossing fault

图21 不同速度的合成记录对比示意图(a) 测井声波数据； (b) 测井密度数据； (c) 测井合成记录； (d) 常规层析反演速度合成记录； (e) 断控层析速度合成记录.Fig.21 Forward comparison diagram of different velocity(a) Logging acoustic data； (b) Logging density data； (c) Logging synthetic record； (d) Conventional tomographic inversion velocity synthesis record； (e) Fault controlled tomography velocity synthesis record.