CN105093279A

CN105093279A - 针对山前带的三维地震初至波菲涅尔体层析反演方法

Info

Publication number: CN105093279A
Application number: CN201410206859.8A
Authority: CN
Inventors: 李凌云; 张传强; 梁鸿贤; 王常波; 杜继修; 秦宁; 傅金荣; 石林光; 李建明; 马季; 孙治国; 孙朋朋; 王磊; 王敬阁; 高丽; 王静轩; 关键; 孙淑琴; 陈会芬; 武天祥
Original assignee: China Petroleum and Chemical Corp; Geophysical Research Institute of Sinopec Shengli Oilfield Co
Current assignee: China Petroleum and Chemical Corp; Geophysical Research Institute of Sinopec Shengli Oilfield Co
Priority date: 2014-05-16
Filing date: 2014-05-16
Publication date: 2015-11-25

Abstract

本发明公开了一种针对山前带的三维地震初至波菲涅尔体层析反演方法，包括：获取地震资料的初至波时间；建立离散模型；计算初至波前走时；确定射线路径；计算菲涅尔体；建立层析反演方程；解该反演方程；将第（2）步的速度模型换成第（7）步的速度模型，重复第（2）至（7）步，即可得到最终的速度模型。该方法考虑了射线传播是沿着菲涅尔带传播的特征，相比单独射线方法更符合实际，同时也更适合山前带这种特殊的复杂地表条件，所得到的结果效果明显。该方法流程及参数设置简单，运算速度快，适合应用于三维地震资料处理。

Description

针对山前带的三维地震初至波菲涅尔体层析反演方法

技术领域

本发明属于油气勘探地震资料处理领域，是一种地震资料表层速度建模的有效方法。

现有技术

目前，地震层析反演方法主要有走时层析和波形层析。

走时层析是通过解释地震波到时来反演地下速度分布，它是基于射线理论，地震波走时是对慢度沿激发点和接收点之间的射线路径的线积分。假设地震波频率足够高，地震能量沿零体积的射线传播，这种理想化的射线又被称为数学射线，在模型空间分布很不均匀，常常遇到低速体发生偏离，而在高速区聚集，使层析反演问题具有严重的不适定性。而山前带速度变化比较大，浅层资料覆盖次数低、信噪比低，许多模型单元上甚至没有一条射线通过，从而要求层析反演的离散模型单元往往要取得较大。当炮点和检波点布设稀疏，特别在连续废炮、坏道之处，层析反演结果往往较差。

波形层析是基于地震波的动力学特征，考虑了地震波的频率因素以及整个模型空间所有介质的效应，从理论上比走时层析具有更高的分辨率，但它需要计算地震波的正向传播及剩余波场的反向传播，还要计算二者的相关沿时间的积分，计算效率非常低，其反演的目标函数和速度摄动之间表现为强烈的非线性关系，对初始模型要求很高，再加上地震子波反演困难、地震信号的信噪比比较低，针对山前带速度变化大，浅层资料覆盖次数低、信噪比低等特点很难得到应有的效果。

山前带发育于洋-陆或陆-陆碰撞背景下，包括前陆褶皱冲断带和前陆盆地靠近造山带的部分。紧挨山前带的造山带通常由基底卷入高角度逆冲断层引起的变质结晶基底隆升和复杂、强烈的褶皱变形组成。构造非常复杂，山前带地震资料品质普遍很差，表现为噪音干扰强烈、速度变化大，浅层资料覆盖次数低、信噪比低等，传统反射波层析反演很难得到准确的浅层速度模型。

地震信号具有一定的频带宽度,不仅高频射线路径上的介质性质影响地震波的传播,射线附近的介质特征同样会对地震波的传播产生影响，只有定量的描述出空间中每一点对接收信息的影响程度才有可能得到更加精确的反演结果。

发明内容

本发明的目的是针对山前带速度变化大，浅层资料覆盖次数低、信噪比低等特点，提出了一种针对山前带的三维地震初至波菲涅尔体层析反演方法。该方法介于波动方程和射线理论之间的三维菲涅尔体层析反演方法，充分考虑了实际地震频率是有限的，地震能量主要在包括理想射线在内的第一菲涅耳带内传播；走时与菲涅耳带内的所有介质有关。菲涅尔体路径相当于物理路径，更合乎地震波实际传播规律；提高了反演稳定性和分辨率，可以提高对山前带复杂地形和地质情况的适应性。

本发明的一种针对山前带的三维地震初至波菲涅尔体层析反演方法包括

(1)获取地震资料的初至波时间；

(2)建立离散模型；

(3)计算初至波前走时；

(4)确定射线路径；

(5)计算菲涅尔体；

(6)建立层析反演方程；

(7)解该反演方程；

(8)将第(2)步的速度模型换成第(7)步的速度模型，重复第(2)-(7)步，即可得到最终的速度模型。

上述方案进一步细化为：

(1)获取地震资料的初至波时间；

(2)建立离散模型，采用规则网格把表层模型离散成若干立方体单元，用网格节点处的速度值表示整个模型速度的变化；

(3)计算初至波前走时，采用程函方程二阶差分和波前群推算法GMM；

(4)确定射线路径，采用波前扩展和双线性插值射线追踪方法，根据互换原理，提供射线追踪精度；

(5)计算菲涅尔体，先求出菲涅尔体半径，再以射线为中心轴线，用菲涅尔体半径确定菲涅尔体空间分布；

(6)建立层析反演方程，利用菲涅尔体空间及其能流密度的分布，建立初至走时残差与网格节点慢度修正量之间满足的线性方程；

(7)解该反演方程，对所有菲涅尔体建立单元慢度增量与走时残差之间的方程，或者表示成节点慢度增量与走时残差之间的方程，解该反演方程；

(8)将第(2)步的速度模型换成第(7)步的速度模型，重复第(2)-(7)步5-15次，即可得到最终的速度模型。

本发明具体技术方案：

(1)获取地震资料的初至波时间。

(2)建立离散模型，采用规则网格把表层模型离散成若干立方体单元，用网格节点处的速度值表示整个模型速度的变化。对于正演问题，网格节点处的速度值为已知；对于反演问题，网格节点处的速度值为未知量。

(3)计算初至波前走时，采用程函方程二阶差分和波前群推算法GMM(GroupMarchingMethod)。该方法不需要对波前窄带中的点进行排序来寻找具有最小波前时间的节点，而是在保证传播时间因果性的条件下，找出多个次级震源同时向外扩展波前，从而大大提高了运算效率。

(4)确定射线路径，创新性采用波前扩展和双线性插值射线追踪方法，这是根据互换原理，提供射线追踪精度。双线性插值函数表示为：

t_{S^{'}} = \frac{(t_{2} - t_{1}) - (t_{4} - t_{3})}{ΔyΔz} gh + \frac{t_{2} - t_{1} + t_{4} - t_{3}}{2 Δy} g + \frac{t_{1} - t_{3} + t_{2} - t_{4}}{2 Δz} h + \frac{t_{1} + t_{2} + t_{3} + t_{4}}{4}

①其中，g＝y_S′-y₀，h＝z_S′-z₀，y₀和z₀flR(l)ik分别为界面的中点坐标，

即y₀＝(y₁+y₂)/2，z₀＝(z₁+z₃)/2。

(5)计算菲涅尔体，先求出菲涅尔体半径，再以射线为中心轴线，用菲涅尔体半径确定菲涅尔体空间分布。下式是计算菲涅尔体半径的近似式：

R (l) \approx \sqrt{\frac{l (L - l)}{fLs (l)}}

②

式中，L是激发点与接收点之间的射线长度；f为地震波主频；l为从激发点到待计算菲涅尔体半径处的射线长度，s(l)为该处的平均慢度，R(l)为该处的菲涅尔体半径。

(6)建立层析反演方程，利用菲涅尔体空间及其能流密度的分布，建立初至走时残差与网格节点慢度修正量之间满足的线性方程。下式是每对激发点和接受点可获得的反演方程：

Σ_{j = 1}^{J} Σ_{k = 1}^{K} ω_{k} \frac{v_{k}}{V} l_{ij} {Δs}_{k} = Δ t_{i}

③

其中，Δt_i是第i条射线对应的走时残差；Δs_k是第k个单元的慢度增量，可以用该单元节点上的慢度增量表示，从而把该式可以整理成节点慢度增量与走时残差的关系式；ω_k是反映各单元能流密度大小的权系数，且ω_k可通过对一条菲涅耳体内各处的影响因子归一化得到。

(7)解该反演方程，对所有菲涅尔体用③式建立单元慢度增量与走时残差之间的方程，或者表示成节点慢度增量与走时残差之间的方程，解该反演方程，采用基于反投影方法和SIRT原理的方法。

(8)将第(2)步的速度模型换成第(7)步的速度模型，重复第(2)-(7)步10次左右，即可得到最终的速度模型。

发明效果

本发明的方法能够较好地获取浅层地表的速度模型，有着其他技术不具备的优势，其具体优势和特点表现在以下几个方面：

第一、技术效果的可靠性。该方法考虑了射线传播是沿着菲涅尔带传播的特征，相比单独射线方法更符合实际，同时也更适合山前带这种特殊的复杂地表条件，所得到的结果效果明显。

第二、操作简单易实现。该方法流程及参数设置简单，运算速度快，适合应用于三维地震资料处理。

附图说明

图1针对山前带的三维地震初至波菲涅尔体层析反演方法的处理流程图

图2为实际工区高程图。

图3为层析反演速度模型。

图4为纵线101线菲涅尔层析静校正前叠加剖面。

图5为纵线101线菲涅尔层析静校正后叠加剖面。

具体实施方式

下面结合附图对本发明的技术方案做进一步说明。

本实施例对XX油田XX区块山前带三维地震资料为靶区进行了应用，用该方法反演了近地表速度模型，以验证本方法的效果，具体流程见图1。该工区满次覆盖面积116平方公里，一次覆盖面积233平方公里，共4362炮。工区地形复杂，高程在400m与1300m之间变化，如图3所示。

1)依据步骤1，采用商业软件拾取初至波时间。

2)其依据步骤2，根据炮点和接收点的分布确定速度模型的范围，在接收线方向20560m，在炮线方向6960m，在垂直方向2040m。用40×40×40的网格离散化该模型，共计4686500个单元。地震波主频取为30Hz，初至模型取为速度为6000m/s的匀速模型。

3)依据步骤3，采用程函方程二阶差分和波前群推算法计算初至波前走时。

4)依据步骤4，采用波前扩展和双线性插值射线追踪方法确定射线路径。计算出离散介质所有网格节点上的波前传播时间后，根据互换原理，从接收点开始，反向确定满足Fermat原理的射线在相应单元界面上的位置，直至到震源所在单元为止，从而获得相应的射线路径。

5)依据步骤5，计算菲涅尔体。在正演计算得到射线路径后，

沿射线路径逐点确定射线的菲涅尔体半径，到射线的垂直距离小于该半径的点组成菲涅尔体。

6)依据步骤6，建立层析反演方程。

7)依据步骤7，解反演方程，得到速度模型。至此，第一次循环结束。

8)将步骤2的速度模型换成步骤7的速度模型，频率提高为40Hz，循环步骤2至步骤7，如此循环10次，每次循环将频率提高10Hz，就可以得到较为准确的速度模型。图3是第10次的反演出来的速度模型。

为了检验层析反演的近地表速度模型的有效性，我们利用该速度模型求取了静校正量，并将静校正应用于实际数据，用来验证速度模型的准确性。图4是菲涅尔体层析静校正前纵线101线叠加结果，图5菲涅尔体层析静校正后纵线101线叠加结果，从图中可以看出，静校正效果明显，说明采用三维菲涅尔体层析反演方法建立的近地表速度模型较准确。

Claims

1.针对山前带的三维地震初至波菲涅尔体层析反演方法，包括

(1)获取地震资料的初至波时间；

(2)建立离散模型；

(3)计算初至波前走时；

(4)确定射线路径；

(5)计算菲涅尔体；

(6)建立层析反演方程；

(7)解该反演方程；

2.根据权利要求1所述的针对山前带的三维地震初至波菲涅尔体层析反演方法，包括：

(1)获取地震资料的初至波时间；

3.根据权利要求2所述的针对山前带的三维地震初至波菲涅尔体层析反演方法，包括：

(1)获取地震资料的初至波时间；

(3)计算初至波前走时，采用程函方程二阶差分和波前群推算法GMM，在保证传播时间因果性的条件下，找出多个次级震源同时向外扩展波前；

(4)确定射线路径，创新性采用波前扩展和双线性插值射线追踪方法，根据互换原理，提供射线追踪精度；其中，双线性插值函数表示为：

t_{S^{'}} = \frac{(t_{2} - t_{1}) - (t_{4} - t_{3})}{ΔyΔz} gh + \frac{t_{2} - t_{1} + t_{4} - t_{3}}{2 Δy} g + \frac{t_{1} - t_{3} + t_{2} - t_{4}}{2 Δz} h + \frac{t_{1} + t_{2} + t_{3} + t_{4}}{4}

即y₀＝(y₁+y₂)/2，z₀＝(z₁+z₃)/2；

(5)计算菲涅尔体，先求出菲涅尔体半径，再以射线为中心轴线，用菲涅尔体半径确定菲涅尔体空间分布，下式是计算菲涅尔体半径的近似式：

R (l) \approx \sqrt{\frac{l (L - l)}{fLs (l)}}

②

式中，L是激发点与接收点之间的射线长度；f为地震波主频；l为从激发点到待计算菲涅尔体半径处的射线长度，s(l)为该处的平均慢度，R(l)为该处的菲涅尔体半径；

(6)建立层析反演方程，利用菲涅尔体空间及其能流密度的分布，建立初至走时残差与网格节点慢度修正量之间满足的线性方程，下式是每对激发点和接受点可获得的反演方程：

Σ_{j = 1}^{J} Σ_{k = 1}^{K} ω_{k} \frac{v_{k}}{V} l_{ij} {Δs}_{k} = Δ t_{i}

③

其中，Δt_i是第i条射线对应的走时残差；Δs_k是第k个单元的慢度增量，用该单元节点上的慢度增量表示，从而把该式整理成节点慢度增量与走时残差的关系式；ω_k是反映各单元能流密度大小的权系数，且ω_k通过对一条菲涅耳体内各处的影响因子归一化得到；

(7)解该反演方程，对所有菲涅尔体用③式建立单元慢度增量与走时残差之间的方程，或者表示成节点慢度增量与走时残差之间的方程，采用基于反投影方法和SIRT原理的方法解该反演方程；

(8)将第(2)步的速度模型换成第(7)步的速度模型，重复第(2)-(7)步8-12次，即可得到最终的速度模型。