CN101872024A - 一种利用时移地震进行井位部署的方法 - Google Patents

Abstract

本发明公开了一种利用时移地震进行井位部署的方法。该方法包括如下步骤：1)按照步骤a和b的方法分别获得两个不同时间点的地震数据：a、测定待测区域的地层基本参数；制作待测区域的油藏剖面，依据油藏剖面建立地层剖面模块；根据地层基本参数和地层剖面模块得到纵波速度、横波速度和密度的三个数据体；b、进行AVO正演模拟，得到各个入射角的叠加的地震数据；2)处理两次地震数据，得地震数据差异；3)将地震数据差异进行反演，得纵波阻抗差、横波阻抗差、密度差和纵横波速度比差；4)预测油藏参数变化；5)进行井位的部署。本发明形成了基于时移地震差异反演实现油藏参数变化预测的过程并确定井位的部署，这是运用常规方法不能实现的。

Description

一种利用时移地震进行井位部署的方法

技术领域

本发明涉及时移地震领域，特别涉及一种利用时移地震进行油井部署的方法。

背景技术

时移地震是目前进行油气藏监测，合理调整开发方案、提高油气采收率的有效手段，其主要目的在于分析研究油气藏开采过程中造成的储层流体运动、流体成分变化、流体饱和度变化、压力变化、孔隙度和温度变化等油气藏特性的变化所引起的地震响应的变化，并由地震响应的变化反演油气藏特性的变化。而地震响应的变化本质是地下储层弹性参数(速度和密度)的变化，建立弹性参数与油藏参数之间的统计关系，便可以实现由地震数据预测油藏参数(见图1)。其中时移地震反演是目前常用的方法，但通常做法是基于两次独立的地震数据体进行反演，分别得到弹性参数，对弹性参数求差得到弹性参数差异数据，即先反演后求差的过程，由此弹性参数差异再进行油藏参数的预测(见图2，即图1中的方案一实现过程)，能否直接利用时移地震数据的差异进行反演并完成油藏参数预测的过程(见图3，即图1中方案二实现过程)没有资料可查。

发明内容

本发明的目的在于提供一种利用时移地震进行井位部署的方法。

本发明提供的利用时移地震进行井位部署的方法，包括如下步骤：

1)按照下述步骤a和b的方法分别获得两个不同时间点的地震数据：

a、测定目标区域的地层基本参数；制作所述目标区域的油藏剖面，依据所述油藏剖面建立所述目标区域的地层剖面模块；根据所述地层基本参数和所述地层剖面模块得到所述目标区域的纵波速度、横波速度和密度的三个数据体；

b、根据步骤a)得到的三个数据体进行AVO正演模拟，得到入射角分别为0-15、15-30、30-45度的叠加的地震数据；

2)处理步骤1)得到的所述两个不同时间点的两次地震数据，求得所述两个不同时间点的地震数据差异；

3)将步骤2)得到的地震数据差异进行测井约束反演，得到所述目标区域地下介质的如下四个弹性参数差：纵波阻抗差、横波阻抗差、密度差和纵横波速度比差；

4)根据步骤3)得到的四个弹性参数差，预测油藏参数变化；

5)依据步骤4)获得的油藏参数变化进行井位的部署。

步骤1)和步骤2)中，上述两个不同时间点的间隔没有限制，时间间隔不等，以观测到油藏的变化为依据，很多地区实施时移地震的动态监测，3个月，6个月就可能采集一次数据，本发明中的两个时间点一个是2003年、另一个是2008年。

步骤1)中，上述地层基本参数是孔隙度、压力、含水饱和度、含气饱和度和含油饱和度。

步骤1)中，上述油藏剖面的生成方法是在ECLIPSE数值模拟软件中，针对目标地质模型，经过数值计算，可以输出任意开发阶段、任意剖面的网格属性参数(孔隙度、渗透率、饱和度、压力、油藏深度等)，带有不同属性值的剖面图即为油藏剖面图。

步骤3)中，上述测井约束反演是将两次的地震数据差异(δS)代入下述公式，计算得到弹性参数差(δL)：

(B^TB+αD^TD+βM^TM)δL＝B^TδS+βM^TMδL^*    (公式I)

公式I中，

δS地震数据差异；

δL是弹性参数差，又称为差异波阻抗对数矩阵；

δL^＊是差异波阻抗对数逆矩阵的转置；

B是地震子波矩阵；

α是光滑因子，决定光滑的强弱；

D是光滑矩阵；

β是低频约束因子，决定约束的强弱；

M是低频矩阵；

B^T、D^T、M^T分别是B、D、M矩阵的转置矩阵。

步骤5)中，依据油藏参数变化进行井位的部署是通过含油饱和度，压力等油藏参数变化可以得到剩余油气的分布规律，判断油气富集的区域。含油气饱和度高或者气藏压力高的区域就是井位部署的有力区域。

任一上述的方法是在油藏监测中的应用也属于本发明的保护范围之内。

上述公式I的推导过程如下：

地震记录与对数波阻抗之间的关系可写为：

S＝BL

展开写为：

$\left(\begin{array}{c}{s}_0\\ s_1\\ s_2\\ s_3\\ s_4\\ s_5\\ s_6\\ s_7\\ s_8\\ s_9\\ s_{10}\end{array}\right)=\left(\begin{array}{ccccccccccc}0& & & & & & & & & & \\ -b_1& b_1& & & & & & & & & \\ -b_2& b_2-b_1& b_1& & & & & & & & \\ -b_3& b_3-b_2& b_2-b_1& b_1& & & & & & & \\ -b_4& b_4-b_3& b_3-b_2& b_2-b_1& b_1& & & & & & \\ -b_5& b_5-b_4& b_4-b_3& b_3-b_2& b_2-b_1& b_1& & & & & \\ & -b_5& b_5-b_4& b_4-b_3& b_3-b_2& b_2-b_1& b_1& & & & \\ & & & & & & b_2-b_1& b_1& & & \\ & & & & & & & b_2-b_1& b_1& & \\ & & & & & & & & b_2-b_1& b_1& \\ & & & & & -b_5& b_5-b_4& b_4-b_3& b_3-b_2& b_2-b_1& b_1\end{array}\right)\left(\begin{array}{c}{\mathrm{li}}_0\\ {\mathrm{li}}_1\\ {\mathrm{li}}_2\\ {\mathrm{li}}_3\\ {\mathrm{li}}_4\\ {\mathrm{li}}_5\\ {\mathrm{li}}_6\\ {\mathrm{li}}_7\\ {\mathrm{li}}_8\\ {\mathrm{li}}_9\\ {\mathrm{li}}_{10}\end{array}\right)$

假设两次勘探得到的地震资料分别记为S₁，S₂，对应的对数弹性阻抗分别为L₁，L₂，则有如下关系成立：

S₁＝BL₁

S₂＝BL₂

S₂-S₁＝B(L₂-L₁)

δS＝BδL

采用优化的办法求解上述问题。求如下泛函的极小：

||δS-BδL||

该极小问题一般是不适定的，为了增加反演过程的稳定性，在上述泛函中加入稳定项得：

||δS-BδL||+α||DδL||

这里D是光滑矩阵，α是光滑因子，决定光滑的强弱，该目标泛函的极小可由如下方程的解得到：

(B^TB+αD^TD)δL＝B^TδS

加入模型约束后，目标泛函变为

||δS-BδL||α|dδL||+β||MδL-MδL^*||

这里M为取低频矩阵，该矩阵乘上模型向量的结果是取模型的低频分量，β为低频约束因子，决定约束的强弱，该目标泛函的极小可由如下方程的解得到：

(B^TB+αD^TD+βM^TM)δL＝B^TδS+βM^TMδL^*    (公式I)

利用此公式就可以实现差异数据体的弹性波阻抗反演。

有益效果：常规的基于时移地震叠前弹性参数反演进行油藏参数预测的技术，是基于两次地震数据体分别计算的过程。时移地震资料通过匹配处理后可以求得两次地震资料的差，如果能够直接通过两次地震数据体的差异，完成上述的反演和预测过程，则可以极大地减少计算量，同时提高预测过程的精度。

本发明形成了基于时移地震差异反演实现油藏参数变化预测的过程并确定井位部署，这是运用常规方法不能实现的。通过中国海上某气田模型数据的计算，以及同常规方法计算结果的对比分析，验证了该方法的有效性和实用性。

附图说明

图1说明基于地震叠前反演进行油藏参数预测的过程，其中，方案1是基于两次独立地震数据体的反演，方案2是基于两次地震数据体差异的反演。

图2表示图1中方案1的反演过程，是传统的方法：即基于两次独立的地震数据体进行弹性参数反演，由此分别预测油藏参数，从而求取油藏参数变化的过程。

图3表示图1中方案2的反演过程，是本文发明的方法，即基于两次地震数据的差异，直接进行弹性参数变化的反演，从而预测得到油藏参数的变化。

图4是油藏剖面与地层剖面模块图，A是油藏剖面，也是2003年和2008年的孔隙度剖面；B是地层剖面模块；C是压力剖面(左为2003年，右为2008年)；D是含气饱和度剖面(左为2003年，右为2008年)。

图5是纵波速度、横波速度和密度的模型，A是纵波速度，B是横波速度、C是密度，左边均是2003年的结果，右边均是2008年的结果。

图6角度道集叠加剖面，左边均是2003年的结果，右边均是2008年的结果，从上至下依次是0～15、15～30、30～45度叠加结果。

图7常规流程反演时移前后纵横波阻抗差，A为纵波阻抗差，B为横波阻抗差。

图8本发明提供的方法(差异反演)时移前后纵横波阻抗差，A为纵波阻抗差，B为横波阻抗差。

图9模型压力差异剖面。

图10本发明提供的方法与常规流程预测压力变化，A是本发明提供的方法，B是常规方法。

具体实施方式

下面结合具体实施例对本发明作进一步说明，但本发明并不限于以下实施例。

下述实施例中，如无特殊说明，均为常规方法。

时移地震反演(波阻抗反演和弹性阻抗反演)原理

波阻抗反演

反演本身是一个数学上的概念，在地震资料处理的许多环节都会提到反演，在广义上讲，整个地震资料的处理过程本身就是一个反演过程。地震资料的反演，包括叠后地震资料波阻抗反演、叠前地震资料AVO反演、角道集叠加地震资料弹性阻抗反演。

叠后波阻抗反演是，已知叠加偏移后的地震资料，求取波阻抗剖面。对应的正问题是已知波阻抗，求叠后地震记录。长期以来，叠后地震资料主要用于构造解释，研究地下反射界面的几何形态，在构造油气勘探中发挥了巨大的作用。随着油气田勘探程度的深入，勘探开发成本的急剧上升，构造勘探已不能满足油气田勘探开发的需要。地震勘探技术已由构造勘探研究深入到岩性勘探研究。所谓构造勘探，是指通过在地震剖面上研究地质构造，进而圈闭油气藏。而岩性勘探是指通过地震资料研究地下介质的性质(如：速度，密度，波阻抗，孔隙度，渗透率等)，进一步圈闭油气藏。波阻抗反演，就是利用叠后地震资料求取地下介质的波阻抗，进而对地下介质的速度、孔隙度、岩性等参数进行预测，实现岩性勘探。

叠后地震资料正演就是，假设已知地下介质的波阻抗，来合成叠后地震资料。假设地下介质是水平层状介质，地下介质的波阻抗为Z(t)，这里t是时间深度，用双程旅行时表示的深度。离散后变为Z(t_i)，i＝1，2，...，n，t_i代表离散后第i网格点的时间。设相应的速度为v(t_i)，密度为ρ(t_i)，则有

Z(t_i)＝ρ(t_i)v(t_i)，i＝1，2，...，n                                          (111)

根据弹性波理论，在平面波垂直入射的情况下，界面的反射系数由界面两侧的波阻抗决定，即

$R\left(t_i\right)=\frac{Z\left(t_{i+1}\right)-Z\left(t_i\right)}{Z\left(t_{i+1}\right)+Z\left(t_i\right)}---\left(\mathrm{1.1.2}\right)$

假设地震子波b(t)已知，叠后地震记录是反射系数与地震子波的褶积：

$s\left(t\right)=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}R\left(t_i\right)b(t-t_i)=R\left(t\right)*b\left(t\right)---\left(\mathrm{1.1.3}\right)$

其中，s(t)是合成的地震记录，b(t)代表地震子波。这样，利用公式(1.1.2)、(1.1.3)就能实现由波阻抗Z(t)到叠后地震记录s(t)的正演过程。

波阻抗反演是已知地震记录s(t)，求波阻抗Z_i，将(1.1.3)离散后可用如下的方程来表示：

$\left(\begin{array}{c}{s}_0\\ s_1\\ s_2\\ s_3\\ s_4\\ s_5\\ s_6\\ s_7\\ s_8\\ s_9\end{array}\right)=\left(\begin{array}{cccccccccc}{b}_0& & & & & & & & & \\ b_1& b_0& & & & & & & & \\ b_2& b_1& b_0& & & & & & & \\ b_3& b_2& & b_0& & & & & & \\ b_4& b_3& & & b_0& & & & & \\ & b_4& & & & b_0& & & & \\ & & & & & & b_0& & & \\ & & & & & & & b_0& & \\ & & & & & & & & b_0& \\ & & & & & b_4& b_3& b_2& b_1& b_0\end{array}\right)\left(\begin{array}{c}{r}_0\\ r_1\\ r_2\\ r_3\\ r_4\\ r_5\\ r_6\\ r_7\\ r_8\\ r_9\end{array}\right)---\left(\mathrm{1.1.4}\right)$

假设地下介质的波阻抗是时间深度的连续函数时，在波阻抗与反射系数关系公式(1.1.2)中，可以认为

2Z(t_i)≈Z(t_i+1)+Z(t_i)                                              (1.1.5)

$Z\left(t_{i+1}\right)-Z\left(t_i\right)\≈\mathrm{\Δt}\frac{\mathrm{dZ}\left(t_i\right)}{\mathrm{dt}}---\left(\mathrm{1.1.6}\right)$

其中Δt代表采样间隔。对于常规地震资料来说，Δt是常数，在讨论相对波阻抗时，在(1.1.6)式中，不乘Δt不会影响波阻抗的相对关系，(1.1.6)也可以写为：

$Z\left(t_{i+1}\right)-Z\left(t_i\right)\≈\frac{\mathrm{dZ}\left(t_i\right)}{\mathrm{dt}}---\left(\mathrm{1.1.7}\right)$

这样，波阻抗与反射系数关系公式(1.1.2)变为：

$r\left(t\right)\≈\frac{\mathrm{dZ}\left(t\right)}{\mathrm{dt}}\frac{1}{2Z\left(t\right)}---\left(\mathrm{1.1.8}\right)$

$r\left(t\right)\≈\frac{1}{2}\frac{d\ln Z\left(t\right)}{\mathrm{dt}}---\left(\mathrm{1.1.9}\right)$

对(1.1.9)离散后，记li_i＝lnZ_i，则反射系数与波阻抗之间的关系可写为：

$\left(\begin{array}{c}{r}_0\\ r_1\\ r_2\\ r_3\\ r_4\\ r_5\\ r_6\\ r_7\\ r_8\\ r_9\end{array}\right)=\left(\begin{array}{ccccccccccc}-1& 1& & & & & & & & & \\ & -1& 1& & & & & & & & \\ & & -1& 1& & & & & & & \\ & & & -1& 1& & & & & & \\ & & & & -1& 1& & & & & \\ & & & & & -1& 1& & & & \\ & & & & & & -1& 1& & & \\ & & & & & & & -1& 1& & \\ & & & & & & & & -1& 1& \\ & & & & & & & & & -1& 1\end{array}\right)\left(\begin{array}{c}{\mathrm{li}}_0\\ {\mathrm{li}}_1\\ {\mathrm{li}}_2\\ {\mathrm{li}}_3\\ {\mathrm{li}}_4\\ {\mathrm{li}}_5\\ {\mathrm{li}}_6\\ {\mathrm{li}}_7\\ {\mathrm{li}}_8\\ {\mathrm{li}}_9\\ {\mathrm{li}}_{10}\end{array}\right)---\left(\mathrm{1.1.10}\right)$

这样地震记录与对数波阻抗之间的关系可写为：

$\left(\begin{array}{c}{s}_0\\ s_1\\ s_2\\ s_3\\ s_4\\ s_5\\ s_6\\ s_7\\ s_8\\ s_9\end{array}\right)=\left(\begin{array}{cccccccccc}{b}_0& & & & & & & & & \\ b_1& b_0& & & & & & & & \\ b_2& b_1& b_0& & & & & & & \\ b_3& b_2& & b_0& & & & & & \\ b_4& b_3& & & b_0& & & & & \\ & b_4& & & & b_0& & & & \\ & & & & & & b_0& & & \\ & & & & & & & b_0& & \\ & & & & & & & & b_0& \\ & & & & & b_4& b_3& b_2& b_1& b_0\end{array}\right)\left(\begin{array}{ccccccccccc}-1& 1& & & & & & & & & \\ & -1& 1& & & & & & & & \\ & & -1& 1& & & & & & & \\ & & & -1& 1& & & & & & \\ & & & & -1& 1& & & & & \\ & & & & & -1& 1& & & & \\ & & & & & & -1& 1& & & \\ & & & & & & & -1& 1& & \\ & & & & & & & & -1& 1& \\ & & & & & & & & & -1& 1\end{array}\right)\left(\begin{array}{c}{\mathrm{li}}_0\\ {\mathrm{li}}_1\\ {\mathrm{li}}_2\\ {\mathrm{li}}_3\\ {\mathrm{li}}_4\\ {\mathrm{li}}_5\\ {\mathrm{li}}_6\\ {\mathrm{li}}_7\\ {\mathrm{li}}_8\\ {\mathrm{li}}_9\\ {\mathrm{li}}_{10}\end{array}\right)$

$\left(\begin{array}{c}{s}_0\\ s_1\\ s_2\\ s_3\\ s_4\\ s_5\\ s_6\\ s_7\\ s_8\\ s_9\\ s_{10}\end{array}\right)=\left(\begin{array}{ccccccccccc}0& & & & & & & & & & \\ -b_1& b_1& & & & & & & & & \\ -b_2& b_2-b_1& b_1& & & & & & & & \\ -b_3& b_3-b_2& b_2-b_1& b_1& & & & & & & \\ -b_4& b_4-b_3& b_3-b_2& b_2-b_1& b_1& & & & & & \\ -b_5& b_5-b_4& b_4-b_3& b_3-b_2& b_2-b_1& b_1& & & & & \\ & -b_5& b_5-b_4& b_4-b_3& b_3-b_2& b_2-b_1& b_1& & & & \\ & & & & & & b_2-b_1& b_1& & & \\ & & & & & & & b_2-b_1& b_1& & \\ & & & & & & & & b_2-b_1& b_1& \\ & & & & & -b_5& b_5-b_4& b_4-b_3& b_3-b_2& b_2-b_1& b_1\end{array}\right)\left(\begin{array}{c}{\mathrm{li}}_0\\ {\mathrm{li}}_1\\ {\mathrm{li}}_2\\ {\mathrm{li}}_3\\ {\mathrm{li}}_4\\ {\mathrm{li}}_5\\ {\mathrm{li}}_6\\ {\mathrm{li}}_7\\ {\mathrm{li}}_8\\ {\mathrm{li}}_9\\ {\mathrm{li}}_{10}\end{array}\right)---(1.1.11)$

(1.1.11)可以简写为矩阵的形式：

S＝BL                                                       (1.1.12)

采用优化的办法求解上述问题。求如下泛函的极小：

||S-BL||

该极小问题一般是不适定的，为了增加反演过程的稳定性，在上述泛函中加入稳定项得

||S-BL||+α||DL||                                             (1.1.13)

这里D是光滑矩阵，α是光滑因子，决定光滑的强弱，该目标泛函的极小可由如下方程的解得到：

(B^TB+αD^TD)L＝B^TS                                             (1.1.14)

加入模型约束后，目标泛函变为

||S-BL||+α||DL||+β||ML-ML^*||                                (1.1.15)

这里M为取低频矩阵，该矩阵乘上模型向量的结果是取模型的低频分量，β为低频约束因子，决定约束的强弱，该目标泛函的极小可由如下方程的解得到：

(B^TB+αD^TD+βM^TM)L＝B^TS+βM^TML^*                               (1.1.16)

利用公式就可以实现波阻抗反演。

弹性阻抗反演

和AVO反演相类似，弹性阻抗反演也是由叠前共中心点道集计算地下地层的物性参数，实际上是AVO反演的另一种实现方式，它的反演过程和叠后波阻抗反演过程完全一致，可以认为是在不同入射角的部分叠加剖面上做波阻抗反演。弹性阻抗EI(elastic impedance)是纵波速度、横波速度、密度的函数和入射角的函数，通过弹性阻抗反演可以求得不同入射角的弹性阻抗值，进一步分析弹性阻抗值随入射角的变化规律，相对于波阻抗，更有利于岩性分析。

振幅随入射角的变化规律可用如下Zoeppritz方程来描述：

r(θ)＝A+Bsin²θ+Csin²θtan²θ                                 (1.1.17)

上式中，r(θ)是随入射角变化的反射系数，A，B，C是参数。要想用叠后波阻抗反演的方法来反演物性参数，就要把上式写成波阻抗与反射系数之间关系的形式，即要把r(θ)写成如下的形式：

$r\left(\mathrm{\θ}\right)=\frac{f\left(t_{i+1}\right)-f\left(t_i\right)}{f\left(t_{i+1}\right)+f\left(t_i\right)}---\left(\mathrm{1.1.18}\right)$

这里函数f(t)称为弹性阻抗(EI)，和(1.1.9)式类似，我们利用对数差分的形式表示上式

$r\left(\mathrm{\θ}\right)\≈\frac{1}{2}\frac{\mathrm{\ΔEI}}{\stackrel{\‾}{\mathrm{EI}}}\≈\frac{1}{2}\mathrm{\Δ}\ln \left(\mathrm{EI}\right)---\left(\mathrm{1.1.19}\right)$

其中ΔEI＝f(t_i)-f(t_i-1)、

，Δ代表求差，以下公式中的Δ也代表同样的含义。将Zoeppritz方程中A、B、C的具体形式带入(1.1.17)式，结合(1.1.19)式得

$\frac{1}{2}\mathrm{\Δ}\ln \left(\mathrm{EI}\right)=\frac{1}{2}(\frac{\mathrm{\Δ}{V}_p}{\stackrel{\‾}{V_p}}+\frac{\mathrm{\Δ\ρ}}{\stackrel{\‾}{\mathrm{\ρ}}})+$

$(\frac{\mathrm{\Δ}{V}_p}{2\stackrel{\‾}{V_p}}-4\frac{V_s^2}{V_p^2}\frac{\mathrm{\Δ}{V}_s}{\stackrel{\‾}{V_s}}-2\frac{V_s^2}{V_p^2}\frac{\mathrm{\Δ\ρ}}{\stackrel{\‾}{\mathrm{\ρ}}}){\sin }^2\mathrm{\θ}+\frac{1}{2}\frac{\mathrm{\Δ}{V}_p}{\stackrel{\‾}{V_p}}{\sin }^2\mathrm{\θ}{\tan }^2\mathrm{\θ}---\left(\mathrm{1.1.20}\right)$

式中，

分别是纵波速度、横波速度和密度的变化率，令

有

$\frac{1}{2}\mathrm{\Δ}\ln \left(\mathrm{EI}\right)=\frac{1}{2}(\frac{\mathrm{\Δ}{V}_p}{\stackrel{\‾}{V_p}}(1+{\sin }^2\mathrm{\θ})+\frac{\mathrm{\Δ\ρ}}{\stackrel{\‾}{\mathrm{\ρ}}}(1-4K{\sin }^2\mathrm{\θ})$

$-\frac{\mathrm{\Δ}{V}_s}{\stackrel{\‾}{V_s}}8K{\sin }^2\mathrm{\θ}+\frac{\mathrm{\Δ}{V}_p}{\stackrel{\‾}{V_p}}{\sin }^2\mathrm{\θ}{\tan }^2\mathrm{\θ})---\left(\mathrm{1.1.21}\right)$

因为sin²θtan²θ＝tan²θ-sin²θ，所以有

Δln(EI)＝(1+tan²θ)Δln(V_p)-8K sin²θΔln(V_s)+(1-4K sin²θ)Δln ρ

$\mathrm{\Δ}\ln \left(\mathrm{EI}\right)=\mathrm{\Δ}\ln \left(V_p^{(1+{\tan }^2\mathrm{\θ})}\right)-\mathrm{\Δ}\ln \left(V_s^{8K{\sin }^2\mathrm{\θ}}\right)+\mathrm{\Δ}\ln {\mathrm{\ρ}}^{(1-4K{\sin }^2\mathrm{\θ})}$

$\mathrm{\Δ}\ln \left(\mathrm{EI}\right)=\mathrm{\Δ}\ln \left(V_p^{(1+{\tan }^2\mathrm{\θ})}{V}_s^{(-8K{\sin }^2\mathrm{\θ})}{\mathrm{\ρ}}^{(1-4K{\sin }^2\mathrm{\θ})}\right)$

$\mathrm{EI}=V_p^{(1+{\tan }^2\mathrm{\θ})}{V}_s^{(-8K{\sin }^2\mathrm{\θ})}{\mathrm{\ρ}}^{(1-4K{\sin }^2\mathrm{\θ})}---\left(\mathrm{1.1.22}\right)$

$\mathrm{EI}=V_P\left(V_p^{\left({\tan }^2\mathrm{\θ}\right)}{V}_s^{(-8K{\sin }^2\mathrm{\θ})}{\mathrm{\ρ}}^{(1-4K{\sin }^2\mathrm{\θ})}\right)---\left(\mathrm{1.1.23}\right)$

这样我们就得到了弹性阻抗的表达式，由此式可以看出，弹性阻抗是纵波速度、横波速度、密度和入射角的函数，随入射角变化的反射系数和弹性阻抗的关系可写为：

$r(\mathrm{\θ},t_i)=\frac{\mathrm{EI}\left(t_i\right)-\mathrm{EI}\left(t_{i-1}\right)}{\mathrm{EI}\left(t_i\right)+\mathrm{EI}\left(t_{i-1}\right)}---\left(\mathrm{1.1.24}\right)$

某一入射角的地震记录d(θ，t)和弹性阻抗的关系，和叠后地震记录类似，有下面的公式：

$d(\mathrm{\θ},t)=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}r(\mathrm{\θ},t_i)w(t-t_i)=r(\mathrm{\θ},t)*w\left(t\right)---\left(\mathrm{1.1.25}\right)$

此公式和叠后地震记录的正演公式完全一致，如果已知某一入射角的地震记录d(θ，t)就可以完全按叠后波阻抗反演的方法求得弹性阻抗。

某一入射角的地震记录d(θ，t)可以通过角道集部分叠加得到。共中心点道集通过射线追踪可以变换为角道集，即把时间-偏移距域的地震记录变换到时间-角度的地震记录，然后进行部分角道集叠加，如将15°-25°范围内的地震记录动校正后进行叠加作为角度20°的地震记录，在不同角度范围进行叠加就可得到不同角度的角道集叠加记录。有了角道集叠加记录就可以按波阻抗反演的方法进行反演得到不同角度的弹性阻抗。

弹性阻抗是对波阻抗的推广，它是入射角的函数，波阻抗是入射角为零时的弹性阻抗的特例，弹性阻抗反演使得波阻抗反演从叠后发展到叠前，角道集叠加剖面可保留地震波的许多AV0特征，弥补了从传统叠加资料里无法得到岩性参数这一缺点，结合弹性阻抗和波阻抗可以更好的解释地下介质的岩性及其含油气性。

由于弹性阻抗和地震记录的关系和波阻抗反演是一致的，弹性阻抗反演方法和波阻抗反演方法也是一致的。模型约束下弹性阻抗反演的目标泛函为

||S-BL||+α||DL||+β||ML-ML^*||(1.1.26)

求解上述泛函极小的方程组为：

(B^TB+αD^TD+βM^TM)L＝B^TS+βM^TML^*(1.1.27)

利用这个公式就可以实现弹性阻抗反演。

实施例1、预测油藏参数

一、地震数据的获得

1、地质模型的建立

1)建立地层剖面模块

用femlab软件依据研究区块的油藏剖面(图4A)建立研究区块地层剖面模块(图4B)。其中油藏剖面是在ECLIPSE数值模拟软件中，针对目标地质模型，经过数值计算，可以输出任意开发阶段、任意剖面的网格属性参数(孔隙度、渗透率、饱和度、压力、油藏深度等)，带有不同属性值的剖面图即为油藏剖面图。

2)输入基本参数

测定待测地层的基本参数(孔隙度Por、压力Pp、含水饱和度Sw、含气饱和度Sg、含油饱和度So)，依次为各个小区块(如图4B中依据图4A中不同颜色建立的块体)输入基本参数，从而输入到步骤1)建立好的地层剖面模块(图4B)中。

其中孔隙度剖面如图4A(2003年和2008年这个两次时间的孔隙度剖面一样)，图4C是压力剖面(左为2003年，右为2008年)，图4D是含气饱和度剖面(左为2003年，右为2008年)

3)测井统计

(1)根据测井数据，区分砂层和泥层，分别统计出砂层和泥层的泥质含量SH和孔隙度Por的关系。

(2)对于砂层，根据其流体情况，应用流体替换技术，把砂岩的纵波速度Vp和密度ρ数据转化为饱水砂岩的纵波速度Vp_Sat和密度ρ_sat，统计出纵波速度Vp_Sat和密度ρ_sat与孔隙度Por的关系；对于泥层，直接统计出其纵波速度及密度与孔隙度间的关系。

(3)对于砂层，应用泥岩公式得到饱水砂岩的横波数据Vs_Sat与纵波速度Vp_Sat的关系；对于泥层，应用泥岩公式得到岩石横波数据和纵波速度的关系。

4)模型计算

对于砂岩，将步骤3)统计出的泥质含量、饱水砂岩的纵横波速度及密度(Vp_Sat，Vs_Sat，ρ_sat)公式输入到建立好的地层剖面模块中，应用流体替换插值计算得到实际地层的纵横波速度和密度数据(Vp，Vs，ρ)；对于泥岩，步骤3)统计出的纵横波速度及密度即为实际地层数据。(纵波速度如图5A，横波速度如图5B、密度如图5C)。(图5的左边是2003年的结果)

重复上述步骤同时考虑压力变化的影响，可得到第二次(2008年)的实际地层数据(纵波速度、横波速度和密度)(图5的右边是2008年的结果)。

2、通过AVO正演模拟和射线追踪得到两次角度部分叠加数据

利用步骤1获得的两次模型(纵波速度、横波速度和密度)进行AVO正演，从而可以得到两次角度部分叠加数据。

时移地震叠前AVO属性变化能表征油藏流体变化，识别流体前缘，而时移地震AVO正演是分析油藏不同参数变化时地震波AVO曲线变化规律，从而进行可行性分析和有效提高时移地震数据处理、解释精度的重要工具。时移地震叠前正演模拟的核心是Zoeppritz方程理论和射线追踪技术。以Zoeppritz方程计算各个界面的反射系数和透射系数，以射线追踪技术计算各个界面的入射角和透射角，由此与子波褶积得到时间域的地震记录，如图6所示，左图为模拟2003年采集的角度部分剖面，从上至下分别为0～15、15～30、30～45度叠加，右图为模拟2008年采集的角度部分剖面，从上至下分别为0～15、15～30、30～45度叠加。

二、弹性参数的差异的获得

本发明的实验组：操作步骤如图3所示，利用地震数据的差异和公式I，直接反演弹性参数差异，具体步骤如下：

1)首先对两次地震数据(即角度部分叠加数据)做匹配处理，求差，得到差异地震数据；

2)利用公式I对差异数据做测井约束反演得到差异弹性参数，包括纵波阻抗差、横波阻抗差、密度差、纵横波速度比差。

其中，测井约束反演是将两次地震数据差异δS代入下述公式，计算得到差异弹性参数δL：

(B^TB+αD^TD+βM^TM )δL＝B^TδS+βM^TMδL^*

反演结果如图8(A为纵波阻抗差，B为横波阻抗差)所示，差异反演所得到的纵横波阻抗差与原始模型(图9)(原始剖面从油藏数值模拟剖面直接可以得到)有较好的对应关系。

本发明的对照组：以传统的方法作为对照。操作步骤如图2所示，利用两次独立的地震数据和公式II，分别进行反演，得到两次独立的弹性参数，由此预测得到两次油藏参数，进而得到油藏参数的变化。

(B^TB+αD^TD+βM^TM)L＝B^TS+βM^TML^*                                   (公式II)

从图7(A为纵波阻抗差，B为横波阻抗差)反演结果看，常规反演流程所得到的纵横波阻抗差与原始模型(图9)并和本发明的反演结果有较好的对应关系。

从理论模型试算知，时移地震差异反演方法与常规流程一样有效可行，且降低了计算量，另外时移地震差异反演不用做时移测井预测，从而具有更广泛的适用性。

三、预测油藏压力的变化

基于步骤二的反演得到的弹性阻抗参数的差异，运用统计方法可以进行油藏参数变化的预测。

实验组油藏压力的变化的预测：

1)根据油藏岩石物理模型分析，油藏有效压力变化时，纵波速度和横波速度的变化有固定的比例，而油藏有效压力变化对密度变化可以忽略。

2)根据Gassmann方程油藏含油饱和度的变化对纵波阻抗变化影响明显，而对剪切模量没有影响，对横波阻抗变化影响也很小。

3)构建方程

进行饱和度变化影响预测，式中

通过

与进行曲线拟合得到，

分别为含水饱和度，纵横波阻抗变化，a压力变化时，纵横波阻抗变化比例因子，由岩石物理分析得到。

4)利用反演得到的纵横波阻抗，基于以上方法即可以得到油藏饱和度变化。

5)根据岩石物理关系，油藏含油饱和度变化与纵横波阻抗通过曲线拟合都可以得到相应的曲线关系，从而计算饱和度变化引起的纵横波阻抗。

6)从纵横波阻抗变化中减去饱和度变化引起的纵横波阻抗，就为压力变化引起的纵横波阻抗。通过曲线拟合建立压力变化与纵横波阻抗变化的关系，利用纵横波阻抗差异计算有效压力变化。

由于油藏参数变化与时移地震纵波阻抗和横波阻抗差异之间不是简单的显示函数关系，课题研究基于BP神经网络的时移地震反演数据解释方法，利用实验室岩石物理测量数据和测井数据获取油藏参数变化与地震纵波阻抗和横波阻抗差异关系，并作为样本进行神经网络训练，然后利用实际时移地震反演的纵波阻抗和横波阻抗差异预测不同油藏参数的变化。

利用两种方案反演得到的开发前后纵、横波阻抗差异数据，结合岩石物理分析，通过多参数统计油气藏参数预测技术，最终预测出开发前后气藏的主要参数压力的变化如图10所示(图中A是本发明提供的方法，B是常规方法)，与模型实际压力变化(图9)剖面有较好的对应关系。

综合原始油藏模型(图9)可知，本发明提供的方法与传统的方法均可以预测得到和原始油藏模型相吻合的油藏参数变化剖面。

四、油井位置的部署

根据步骤三获得的油藏参数变化进行井位的部署。通过含油饱和度，压力等油藏参数变化可以得到剩余油气的分布规律，判断油气富集的区域。含油气饱和度高或者气藏压力高的区域就是井位部署的有力区域。

Claims (4)

1.一种利用时移地震进行井位部署的方法，包括如下步骤：

1)按照下述步骤a和b的方法分别获得两个不同时间点的地震数据：

a、测定目标区域的地层基本参数；制作所述目标区域的油藏剖面，依据所述油藏剖面建立所述目标区域的地层剖面模块；根据所述地层基本参数和所述地层剖面模块得到所述目标区域的纵波速度、横波速度和密度的三个数据体；

b、根据步骤a)得到的三个数据体进行AVO正演模拟，得到入射角分别为0-15、15-30、30-45度的叠加的地震数据；

2)处理步骤1)得到的所述两个不同时间点的两次地震数据，求得所述两个不同时间点的地震数据差异；

3)将步骤2)得到的地震数据差异进行测井约束反演，得到所述目标区域地下介质的如下四个弹性参数差：纵波阻抗差、横波阻抗差、密度差和纵横波速度比差；

4)根据步骤3)得到的四个弹性参数差，预测油藏参数变化；

5)依据步骤4)获得的油藏参数变化进行油井位置的部署。

2.根据权利要求1所述的方法，其特征在于：步骤1)中，所述地层基本参数是孔隙度、压力、含水饱和度、含气饱和度和含油饱和度。

3.根据权利要求1或2所述的方法，其特征在于：步骤3)中，所述测井约束反演是将地震数据差异代入下述公式，计算得到弹性参数差：

(B^TB+αD^TD+βM^TM)δL＝B^TδS+βM^TMδL^*

其中，

δS是地震数据差异；

δL是弹性参数差；

δL^＊是差异波阻抗对数逆矩阵的转置；

B是地震子波矩阵；

α是光滑因子，决定光滑的强弱；

D是光滑矩阵；

β是低频约束因子，决定约束的强弱；

M是低频矩阵；

B^T、D^T、M^T分别是B、D、M矩阵的转置矩阵。

4.权利要求1-3中任一所述的方法在油藏监测中的应用。

Images