施俊杰，张建中

（1.厦门大学信息科学与技术学院，福建 厦门361005；2.中国海洋大学海洋地球科学学院，山东 青岛266100）

井间地震层析成像是估算井间介质速度最有力的技术之一。自McMechan［1］提出井间地震层析成像以来，人们提出了不少井间地震走时层析成像方法。但这类方法几乎都是使用直达波走时数据，往往只能给出地层速度，不能直接准确给出地层界面展布，且反演速度的分辨率不能满足勘探和开发的需要。

反射波走时比直达波走时对地层分界面的变化敏感，而直达波走时对地层速度的变化更敏感，因此联合两者同时反演分界面和速度，将会发挥各自的优势，并增加观测数据个数，减少层析反演问题的多解性。Calnan和Schuster［2］利用垂直井、水平地层二维模型，在已知水平反射界面和直射线路径的假设下对地层速度进行联合层析反演，结果表明联合走时层析成像比直达波走时层析成像具有更高的速度分辨率。Wu等［3］用二维矩形单元离散地层，用斜线表示反射界面，提出了1种二维井间曲射线反射和透射层析方法，提高了层析成像的分辨率，可以分辨出薄层。Van Schaack［4］用规则矩形网格离散地层，用一元深度函数表示反射界面，建立了地层与反射界面不耦合的井间二维离散模型，以及与该离散模型对应的直达波走时与层速度、反射波走时与界面深度的层析反演方程，提出了直达波和反射波走时联合层析成像方法。Spetzler［5］利用物理实验数据，研究了有限频率直达波和反射波走时联合层析反演问题。Bube and Langan［6］利用两层均匀介质的简单井间模型，从理论上研究了井间直达波和反射波走时联合层析反演的分辨率问题，认为联合层析反演可以很好地确定反射层深度，提高慢度分辨率。这些研究使用简单的或者地层与反射界面不耦合的二维井间模型，分别用直达波走时约束层速度，用反射波走时推断反射界面，还没有达到联合使用2种走时，同时反演层速度和反射界面的更符合实际情形的程度。

实际井孔轨迹一般不是过井口的铅垂线，而是1条曲线，地震波传播路径也是在1个三维空间内展布，因此在三维空间进行射线追踪和层析反演，更符合实际情况。另一方面，由于激发点和接收点只能放在井孔内，观测孔径有限，地震波射线难以覆盖对应的三维空间，也难以为三维反演问题提供足够的旅行时数据，从而增大了反演问题的不适定性。综合考虑这些情况，本文建立了基于不规则单元的地层与反射界面完全耦合的井间三维离散模型，并假设在垂直于连井截面方向上模型参数不变，即在正演时采用三维模型，以模拟三维空间的射线分布和走时数据；反演时相当于采用二维模型，以减小反演问题的不适定性。在此不规则单元离散模型基础上，建立了联合使用直达波和反射波走时同时反演地层速度和反射界面的层析反演方程，实现了井间直达波和反射波走时的联合层析反演。本方法已在文献［7］中得到了初步应用，本文将详细讨论本方法并在模型应用上展开分析。

1 方法原理

采用张建中等［8］的方法建立井间离散模型，图1所示是三维离散模型的xoz截面图。

模型中地层速度和分界面在连接左右两井轨迹中线的xoz平面内变化，而在垂直于该平面的y轴方向上保持不变。对于正演计算，模型在y轴方向上被等间距地离散成几个单元，对于反演问题，在y轴方向上相当于只有1个单元。由于在反射界面上速度发生突变，每个节点需设置地层上下2个速度值。这种离散模型的地层与界面完全耦合。

图1 井间离散模型截面示意图Fig.1 Sketched section of crosshole discretized model

考虑到走时数据存在测量和拾取误差，反问题的不适定性等因素，构造了下列联合层析反演问题的目标函数：

对上述非线性目标函数在初始模型邻域用泰勒级数展开至线性项，利用函数极小点满足的关系式，得到反演的迭代方程：

其中：k代表迭代次数；A表示Frechet矩阵；Δb表示观测走时与当前迭代模型计算走时之差；x表示模型参数；Δx为x的修正量，即线性方程中的未知量。各量的具体表达式可以写为：

式中：si是地层第i个节点上的慢度；zj是界面上第j个节点上的深度。

直达波走时和反射波走时对地层单元节点慢度si的导数计算式可以统一为：

其中：K为射线穿过的以节点i为顶点的单元个数；lk为第k个单元内的射线长度，由于网格慢度由4个顶点慢度平均所得，所以节点慢度的系数需除以4。

直达波和反射波走时对界面单元两端节点深度的导数式也可以统一为：

其中：xc和zc是直达波或反射波射线与界面单元交点的坐标；xl，zl，xr，zr表示该 界 面 单 元左、右端 点 的 坐标；sin和sout分别表示入射射线和出射射线所在地层单元的慢度；α和β分别是入射射线和出射射线与z轴的夹角

实际计算中，（4）式和（5）式中的方差均折算到权重系数。权重系数的设置应基于模型预判和经验值：ρtd可作为基准值设置为1，ρtr一般不超过1；慢度平滑系数λs，1分为x轴系数和z轴系数，其x轴系数在各种模型中变化不大，但z轴系数可因模型层次疏密程度进行宽幅变化；深度平滑系数λz，1一般设置为1左右即可得到良好效果；考虑到先验信息的可靠性问题，先验系数一般较小，在缺少先验信息的模型中甚至可以不设置先验方程，下文的模拟模型使用了测井先验信息，实例模型则未使用任何先验信息。在迭代反演过程中，应使后四项系数的值随着迭代次数逐次降低，以期在迭代初期获得稳定且有意义的收敛解，而在迭代后期获得与走时数据更匹配的高分辨率解。在求解线性方程（3）时，本文选择阻尼LSQR方法。

对于本文不规则单元的离散模型，使用基于波前扩展和走时插值的三维射线追踪方法［9］确定射线路径，该方法适应性强，计算精度和效率高。

2 模型实验

图2是含有起伏反射界面的层状介质井间模型试验。图2a是模型的三维尺寸，蓝色线为炮点阵列，绿色线为检波点阵列；图2b是理论模型的xoz截面。模型尺寸为200m×10m×960m；含有5个匀速层，4个反射界面；炮点和检波点各300个，分别均匀布置在33～930m深度范围内。层析反演时，选取了射线与水平面的夹角在15°～45°范围内共计15 544个直达波和10 972个反射波走时数据。初始界面取为左右井理论界面深度的连线，即平面或斜面，各层的初始速度均取为2.0km／s。地层离散单元起始尺寸为8m×5m×5m，且高度随着迭代过程中反射界面的变化而变化。图2c是直达波走时层析反演第六次迭代结果，图2d是联合层析反演的第六次迭代结果。可以看出，单直达波走时反演给出了较好的速度分布，但只能以速度过渡带的形式反映界面，特别是2个斜截面和下凹界面的中段反演精度不高；联合反演既给出了清晰的反射界面，也给出了精度更高的地层速度，整体结果与理论模型很接近。

3 实例

本文对胜利油田某区J41－108井的实际井间地震数据进行了层析成像。炮点329个，深度从869～1 853m，检波点300个，深度从866.1～1 763.1m，炮（检）深度间距均为3m；经坐标轴旋转后，炮点阵列在左井，且为直井，检波点阵列在右井，为斜井；地震有效道数43 568个，拾取了36 693道直达波走时数据，并参考其他地质资料，选择了6个反射界面及其6 102个反射波走时数据，6个反射面分别为1 412、1 480、1 570、1 614、1 697和1 726m。通 过 射 线 出 射 角 限 定 后，有11 210个直达波和1 203个反射波走时数据参与了反演计算。图3是J41－108第163炮的初始拾取结果，其中反射波皆为上行波。

图2 含非水平界面模型合成走时数据的层析反演结果Fig.2 Tomographic inversion models using synthetic traveltimes of a crosshole model with nonlevel interfaces

图3 胜利油田J41－108第163炮部分拾取结果Fig.3 Part of the picking result of 145th shot of 41－108wells in Shengli Oilfield

根据井斜资料，建立194m×15m×1 000m的模型，深度范围从860～1 860m；设置6个水平初始反射面，深度如上文所述，各地层初始速度均为2.0km／s；初始离散单元为5m×5m×3m。图4a为只用直达波走时反演的第6次迭代结果，图4b为联合层析反演第6次迭代结果，反演中没有施加任何慢度和深度的先验信息约束。为了直观表示反演结果的可靠程度，速度图两侧附加了声波测井曲线和边界反演速度的对比，其中，青色为声波测井速度，黑色为对应的边界反演速度。由于声波测井曲线的记录深度有限，为了与测井速度进行对比，图4只显示出深度范围为1 000～1 860m的反演速度结果。

由图3可见，该井间反射记录复杂，能够清晰拾取的反射走时数据很少，参与反演计算的反射波走时数据比直达波走时数据少得多，使得图4a和图4b的结果整体上比较接近，但图4b中联合层析反演模型的速度分布更精细，垂向上速度变化层次更多，层间变化更剧烈，此现象在反射数据走时最集中的1 480、1 570和1 697m的3个界面附近最为明显。两种方案沿井径的层析速度与测井声波速度的变化趋势也基本吻合，而联合层析的速度曲线与测井声波速度曲线的吻合度更好。实例结果说明，井间地震直达波和反射波联合层析反演方法是可靠的，在层析成像中加入反射波走时，不仅能够确定反射界面的起伏，而且能够提高模型速度的分辨率。

图4 41－108井实际资料层析反演结果Fig.4 Tomographic model section from field seismic data between wells 41and 108

4 结语

本文提出了1种基于不规则六面体（和五面体）单元的井间直达波和反射波走时联合层析反演方法。建立的地层与界面完全耦合的三维井间离散模型更合理，可以很好的模拟连井截面内反射界面起伏和每个地层内的速度变化；使用的三维射线追踪算法，能同时获得井间三维空间的直达波和反射波射线路径和走时数据；使用二维反演方程，大大减小了反演问题的不适定性。与井间直达波层析反演相比，井间直达波和反射波走时联合层析反演，既可以确定反射界面，也提高了井间层析成像的分辨率和可靠性。

［1］ McMechan G A.Seismic tomography in boreholes［J］.Geophys J Roy Astr Soc，1983，74：601－612.

［2］ Calnan C，Schuster G T.Reflection ＋ transmission crosswell tomography［C］.［s.l.］：59th Annual International Meeting，SEG，Expanded Abstracts，1989，8：908－911.

［3］ Wu L，Song W，Zhang M.Seismic crosshole curved ray reflection＆transmission tomography（CCRTT）［J］.CT Theory and applications，1994，3（1）：42－47.

［4］ Van Schaack M A.Velocity estimation for crosswell reflection imaging using combined direct and reflected arrival traveltime tomography［D］.Ph.D.dissertation.USA：Stanford University：1997.

［5］ Spetzler J.Finite－frequency wavefield theory for high－resolution velocity estimation using transmission and reflection data［C］.［s.l.］：73th Annual International Meeting，SEG，Expanded Abstracts，2003，22：2315－2318.

［6］ Bube K P，Langan R T.Resolution of slowness and reflectors in crosswell tomography with transmission and reflection traveltimes［J］.Geophysics，2008，73（5）：321－335.

［7］ 左建军，林松辉，孔庆丰，等.井间地震直达波和反射波联合层析成像及应用［J］.石油地球物理勘探，2011，46（2）：226－231.

［8］ 张建中，丁兴号.一种2.5维井间透射波层析成像方法［J］.声学学报，2007，32（1）：91－96.

［9］ Huang Y Q，Zhang J Z，Liu Q H.Three－dimensional GPR ray tracing based on wavefront expansion with irregular cells［J］.IEEE Transactions on Geoscience and Remote Sensing，2011，49（2）：676－686.

［10］ 张建中，陈世军，徐初伟.动态网络最短路径射线追踪［J］.地球物理学报，2004，47（5）：899－904.