一种深度域测井和地震数据的直接标定方法

1.本发明属于油气地震勘探领域，涉及一种深度域测井和地震数据的直接高效标定方法。


背景技术：

2.地球物理测井是把依据电、磁、声、热、核等物理原理制造的各种测井仪器放入钻井内，沿井连续记录随深度变化的各种物理参数，然后通过这些物理参数的变化曲线识别地下的岩性，以及油、气、水、煤、金属等矿藏。在测井点处，由于有测井时所取得的岩芯及测量的物理参数，人们对该点处地下地质情况有清晰的认识。若要在远离测井点处推测地下的地质情况，目前主要是用人工地震勘探所获得的地震数据来进行。根据地震数据的特征(如反射特征等)，赋予其明确的地质意义和概念模型以明确地下地质情况就是地震资料的解释工作。在地震资料解释中，测井和地震数据标定是将地下地质情况与观测到的地震数据关联起来的关键步骤，也是进行可靠地层解释和地震反演的基石。测井和地震数据标定的基本原理是，首先用测井所获得的速度和密度参数制作合成地震记录，然后将合成地震记录中明显的波形特征与测井处附近地震数据的波形特征进行对齐(对齐操作需要重复多次)，从而将测井物理参数曲线的特征与地震数据的特征关联起来，以实现测井点处地质认知的外推，即可以对远离测井点处的地下的地质情况进行推测。
3.由于以往用于地震资料解释的地震数据大多是时间域的地震数据，故测井和地震数据标定需要将随深度变化的测井物理参数根据测井速度信息数转换为随时间变换的测井物理参数，然后在时间域基于褶积模型合成时间域地震记录，进而与时间域地震数据进行标定。时间域的测井和地震数据标定技术已经很成熟。然而，随着对地下构造成像精度更高的深度域地震数据在油气地震勘探中的应用越来越普遍，地震资料解释也越来越多地基于深度域地震数据进行。目前，与深度域地震资料解释配套的深度域测井和地震数据标定技术还不成熟，主要做法是将深度域测井数据和深度域地震数据转换至时间域完成标定，然后再将标定结果转换回深度域；这样的标定流程复杂，存在多次域转换，会在测井数据和地震数据中引入转换误差，从而使得标定精度不高，甚至在错误的深度位置将测井数据与地震数据进行匹配。


技术实现要素：

4.本发明的目的在于提供一种深度域测井和地震数据的直接标定方法，其特点是在标定过程中采用了直接在深度域合成地震记录和提取地震子波的方法，即地震记录合成和地震子波提取过程中不存在域转换操作。因此，本发明方法在充分利用测井信息的同时，能够避免间接方法在不同域之间转换时产生的误差积累及由于重采样造成的数据信息丢失。本发明方法在第一次测井和地震数据标定后就能够在测井数据与地震数据之间建立起较为准确的测井深度-地震深度对应关系，从而可以减少后续在局部调整测井数据与地震数据之间的深度对应关系的次数。本发明方法包括以下主要步骤：
5.(1)剔除测井速度数据和测井密度数据中的空值和异常值，得到有效测井速度数据和有效测井密度数据。
6.(2)用有效测井速度数据和有效测井密度数据计算不同深度位置处的测井反射系数r＝[r(1),
…
,r(l),
…
,r(l)]，其中，l为有效测井速度数据和有效测井密度数据的采样点数量。
[0007]
(3)从目标深度段内且经过井位置的剖面地震数据s中提取零相位初始深度域地震子波，其过程按照以下步骤进行：
[0008]
①
计算剖面地震数据s的平均谱y：
[0009][0010]
其中，fft[
·
]表示向量的快速傅里叶变换，剖面地震数据s共包含n道地震记录，s(n)为剖面地震数据s中的第n道地震记录，且其包含m个采样点；
[0011]
②
构建目标函数并使用最小二乘法对目标函数求解，获得拟合参数α和β；其中，||
·
||2表示向量的l2范数，exp[
·
]表示向量的指数，地震波数向量k＝[k(1),
…
,k(m),
…
,k(m)]，单位是千米分之一，δd为剖面地震数据s的深度采样间隔，单位是米，拟合参数α和β的单位均为千米分之一。
[0012]
③
使用拟合参数β和测井波数向量k构建零相位初始深度域地震子波w0：
[0013][0014]
其中，k＝[k(1),
…
,k(l),
…
,k(l)]，δd为测井数据的深度采样间隔，单位是米。
[0015]
(4)用零相位初始深度域地震子波w0和测井反射系数r制作深度域合成地震记录α：
[0016][0017]
其中，real{
·
}表示取向量实部的运算，
⊙
表示向量的哈达玛积，矩阵p为一个l阶矩阵，矩阵为矩阵p的共轭矩阵，矩阵p中的元素p(l,l)具有如下表达式：
[0018][0019]
其中，为虚数单位，v
max
为有效测井速度数据v中的最大值，单位是米/秒。此时，深度域合成地震记录α对应的深度坐标为d
α
＝[d
α
(1),
…
,d
α
(l),
…
,d
α
(l)]，且d
α
(1)为与测井对应的深度域井旁地震记录的起始深度。
[0020]
(5)将测井对应的深度域井旁地震记录按照测井数据的深度采样间隔采用三次样
条插值方法进行插值。
[0021]
(6)将深度域合成地震记录α与插值后的深度域井旁地震记录进行互相关运算，并用互相关结果中最大值对应的深度位置d0将深度域合成地震记录α对应的深度坐标更新为d
α
＝d
α
+[d
0-d
α
(1)]；同时，从d0起，从插值后的深度域井旁地震记录中截取l个采样点得到深度域井旁地震记录片段β，对应的深度坐标为d
β
＝[d
β
(1),
…
,d
β
(l),
…
,d
β
(l)]。
[0022]
(7)基于深度域井旁地震记录片段β和测井反射系数r，使用直接深度域地震子波提取方法提取深度域地震子波；基于提取的深度域地震子波和测井反射系数r，使用深度域地震记录的直接合成方法更新深度域合成地震记录α；
[0023]
(8)选取深度域井旁地震记录片段β中与深度域合成地震记录α中波形特征相似且振幅绝对值较大的波峰或波谷位置，通过压缩或是拉伸深度域井旁地震记录片段β，将所选取的波峰或波谷位置与深度域合成地震记录α进行对齐，并产生相应的深度偏差点集合e：
[0024]
e＝{e(1),
…
,e(j),
…
,e(j)},j＜l，
[0025]
其中，e(j)＝d
β
(j)-d
α
(j)表示深度域井旁地震记录片段β与深度域合成地震记录α在所选取的第j点位置处的深度坐标之差。对深度偏差点集合中的元素使用三次样条插值方法进行插值，得到深度校正曲线c。
[0026]
(9)将深度校正曲线c分别与当前深度域合成地震记录α和深度域井旁地震记录片段β对应的深度坐标d
α
和d
β
进行叠加，得到各自的新深度坐标，即d
α
＝d
α
+c和d
β
＝d
β
+c，并根据新的深度坐标从插值后的深度域井旁地震记录中截取对应的深度域井旁地震记录片段β，然后重复步骤(7)；
[0027]
(10)如果深度域井旁地震记录片段β与深度域合成地震记录α在波形特征相似，振幅绝对值较大的波峰或波谷位置已经对齐，且二者的相关系数高于之前的标定过程，就结束标定操作，并生成测井深度-地震深度对应关系和深度域地震子波，否则重复步骤(8)至(9)。
附图说明
[0028]
图1是本发明实施例的经过处理后的有效测井数据(图1a-1c)及之对应的深度域测井反射系数(图1d)和深度域井旁地震记录(图1e)。有效测井数据包括伽玛数据(无量纲单位)、速度数据(单位是米/秒)和密度数据(单位是克/立方厘米)。有效测井数据及深度域测井反射系数(无量纲单位)的深度范围为2991米至3199.2米，深度采样间隔为0.15米；深度域井旁地震记录(无量纲单位)的深度范围为2400米至3895米，深度采样间隔为5米。
[0029]
图2是采用本发明方法对有效测井数据和对应的深度域井旁地震记录(灰色虚线)进行标定的过程示意图。图2a至图2e中的黑色实线均为深度域合成地震记录。图2a是由图2f的初始深度域地震子波和测井反射系数制作的深度域合成地震记录与深度域井旁地震记录的对比，此时，深度域合成地震记录的起始深度与深度域井旁地震记录的起始深度相同；图2b中，将深度域合成地震记录与整个深度域井旁地震记录进行互相关后，将互相关结果中最大值对应的深度位置(2941.2米)作为深度域合成地震记录的起始深度，同时，从该深度位置起截取深度域井旁地震记录片段(图2b中灰色区域内的灰色虚线)，基于此深度域井旁地震记录片段和测井反射系数提取的深度域地震子波如图2g所示；图2c是由图2g的深
度域地震子波和测井反射系数制作的深度域合成地震记录与深度域井旁地震记录的对比，选取深度域井旁地震记录片段中与深度域合成地震记录中波形相似且振幅值较大的波峰或波谷位置(如图3a所示)，通过压缩或是拉伸深度域井旁地震记录片段，将所选取的波峰或波谷位置与深度域合成地震记录进行对齐，进而产生相应的深度偏差点(如图3b中黑色标记点所示)，再由这些偏差点生成深度校正曲线(如图3b中黑色曲线所示)对深度域合成地震记录和深度域井旁地震记录片段的深度坐标进行更新，并截取新的深度域井旁地震记录片段，基于更新后的深度域井旁地震记录片段和测井反射系数提取的深度域地震子波如图2h所示；图2d是由图2h的深度域地震子波和测井反射系数制作的深度域合成地震记录与深度域井旁地震记录的对比，再次经过与图3类似的一次标定后，基于更新后的深度域井旁地震记录片段和测井反射系数提取的深度域地震子波如图2i所示；图2e是由图2i的深度域地震子波和测井反射系数制作的深度域合成地震记录与深度域井旁地震记录的对比，此时，深度域井旁地震记录片段中与深度域合成地震记录中波形特征相似且振幅值较大的波峰或波谷位置已经对齐，且二者的相关系数高于之前的标定，可以结束标定操作。
[0030]
图3是与图2c对应的一次标定示意图。其中，图3a是图2c灰色区域部分中深度域合成地震记录(黑色实线)与深度域井旁地震记录片段(灰色虚线)的对比，两条曲线上的黑色标记点是要进行对齐的波峰或波谷位置；图3b中的黑色标记点是通过压缩或拉伸图3a中的深度域井旁地震记录片段，将所选取的波峰或波谷位置与深度域合成地震记录进行对齐后产生的深度偏差点，对这些深度偏差点采用三次样条插值方法进行插值就得到了图3b中的黑色实线，也即深度校正曲线；图3c是将所选取的波峰或波谷位置对齐后，深度域井旁地震记录片段与深度域合成地震记录的对比，此时，二者的波形特征比对其前更加匹配。
[0031]
图4是与图3对应的测井深度与地震深度对应关系更新示意图，其纵坐标为测井深度，单位是米，其横坐标为地震深度，单位是米。其中，灰色虚线是未叠加深度校正曲线前的测井深度与地震深度对应关系，此时，测井深度的2991米至3199.2米对应地震深度的2941.2米至3149.4米，黑色实线是叠加深度校正曲线后的测井深度与地震深度对应关系，此时，测井深度的2991米至3199.2米对应地震深度的2920.3米至3154.2米。
[0032]
图5是本发明实施例的位于目标深度段内且过测井的inline方向深度域地震剖面，其纵坐标为深度，单位是米，其横坐标为道号，无量纲，色标表示地震数据的振幅，无量纲，深色越深，振幅值越小，反之，振幅值越大；该深度域地震剖面的深度范围为2800米至3300米，深度采样间隔为5米；测井位于剖面道号5449处，剖面中所展示的测井曲线为伽玛曲线；此外，在剖面上还标注了根据该工区内的已有地质认知及已有测井信息预测的两套砂体。
具体实施方式
[0033]
图2展示了将测井数据与深度域井旁地震记录(图1e)进行标定的过程示意图，其对应的具体实施方式包括以下步骤：
[0034]
(1)剔除测井速度数据和测井密度数据中的空值和异常值，得到有效测井速度数据v和有效测井密度数据ρ；
[0035]
(2)用有效测井速度数据v和有效测井密度数据ρ按照下式依次计算不同深度位置
处的测井反射系数r＝[r(1),
…
,r(l),
…
,r(l)]，其中，l＝1389为有效测井速度数据v和有效测井密度数据ρ的采样点数量。
[0036]
(3)从目标深度段2800米至3400米内且经过测井位置的inline方向剖面地震数据s中提取一个零相位初始深度域地震子波，其包括如下步骤：
[0037]
①
计算剖面地震数据s的平均谱其中，fft[
·
]表示向量的快速傅里叶变换，剖面地震数据s共包含n＝965道地震记录，s(n)为地震数据s中的第n道地震记录，且其包含m＝121个采样点；
[0038]
②
构建目标函数并使用最小二乘法对目标函数求解，获得拟合参数α和β；其中，||
·
||2表示向量的l2范数，exp[
·
]表示向量的指数，地震波数向量k＝[k(1),
…
,k(m),
…
,k(m)]，剖面地震数据s的深度采样间隔δd为5米，拟合参数α和β的单位均为千米分之一。
[0039]
③
使用拟合参数β和测井波数向量k构建零相位初始深度域地震子波w0：
[0040][0041]
其中，k＝[k(1),
…
,k(l),
…
,k(l)]，测井数据的深度采样间隔δd为0.15米。
[0042]
(4)用零相位初始深度域地震子波w0和测井反射系数r制作深度域合成地震记录α：
[0043][0044]
式中，real{
·
}表示取向量实部的运算，
⊙
表示向量的哈达玛积，矩阵p为一个l阶矩阵，矩阵为矩阵p的共轭矩阵，矩阵p中的元素p(l,l)具有如下表达式：
[0045][0046]
式中，为虚数单位，有效测井速度数据v中的最大值v
max
＝3954.2米/秒。此时，深度域合成地震记录α对应的深度坐标为d
α
＝[d
α
(1),
…
,d
α
(l),
…
,d
α
(l)]，其中，测井对应的深度域井旁地震记录的起始深度d
α
(1)＝2400米。
[0047]
(5)将测井对应的深度域井旁地震记录按照测井数据的深度采样间隔采用三次样条插值方法进行插值。
[0048]
(6)将深度域合成地震记录α与插值后的深度域井旁地震记录进行互相关运算，并用互相关结果中最大值对应的深度位置，即d0＝2941.2米，将深度域合成地震记录α对
应的深度坐标更新为d
α
＝d
α
+[d
0-d
α
(1)]；同时，从d0起，从插值后的深度域井旁地震记录中截取l个采样点得到深度域井旁地震记录片段β，对应的深度坐标为d
β
＝[d
β
(1),
…
,d
β
(l),
…
,d
β
(l)]。
[0049]
(7)基于深度域井旁地震记录片段β和测井反射系数r，使用基于模型的深度域地震子波提取方法(可参照申请号为2020105042059中国发明专利中所披露的方法)提取深度域地震子波；基于提取的深度域地震子波和测井反射系数r，使用深度域地震记录的直接合成方法(同样可参照申请号为2020105042059中国发明专利中所披露的方法)更新深度域合成地震记录α；
[0050]
(8)选取深度域井旁地震记录片段β中与深度域合成地震记录α中波形特征相似且振幅绝对值较大的波峰或波谷位置(如图3a所示)，通过压缩或是拉伸深度域井旁地震记录片段β，将所选取的波峰或波谷位置与深度域合成地震记录α进行对齐，并产生相应的深度偏差点集合e(如图3b的黑色标记点所示)：
[0051]
e＝{e(1),
…
,e(j),
…
,e(j)},j＜l，
[0052]
式中，e(j)＝d
β
(j)-d
α
(j)表示深度域井旁地震记录片段β与深度域合成地震记录α在所选取的第j点位置处的深度坐标之差。对深度偏差点集合中的元素使用三次样条插值方法进行插值，得到深度校正曲线c(如图3b的黑色实线所示)；
[0053]
(9)将深度校正曲线c分别与当前深度域合成地震记录α和深度域井旁地震记录片段β对应的深度坐标d
α
和d
β
进行叠加，得到各自的新深度坐标(如图3c中的横坐标)，即d
α
＝d
α
+c和d
β
＝d
β
+c，并根据新的深度坐标从插值后的深度域井旁地震记录中截取对应的深度域井旁地震记录片段β，然后重复步骤(7)；
[0054]
(10)如果深度域井旁地震记录片段β与深度域合成地震记录α在波形特征相似，振幅绝对值较大的波峰或波谷位置已经对齐，且二者的相关系数高于之前的标定过程(如图2e所示)，就结束标定操作，并生成测井深度-地震深度对应关系和地震子波(如图2i所示)，否则重复步骤(8)至(9)。
[0055]
从图2可以看出，使用本发明方法将测井数据与深度域井旁地震记录进行2次标定后(图2c和图2d所示)就得到了合理的测井深度-地震深度对应关系和地震子波(图2i所示)。
[0056]
测井伽玛曲线可以用于岩性划分，对于砂泥岩剖面，纯砂岩的伽玛值最低，伽玛值随着泥质含量的增加而升高。在图5所示的过井深度域地震剖面中，测井伽玛曲线低值区域均与预测的两套砂体吻合，说明由本发明方法是有效可靠的。
[0057]
本发明方法的优点在于：在深度域测井和地震数据标定过程中，采用了深度域地震记录的直接合成方法和地震子波提取方法，能够在充分利用测井信息的同时，避免间接方法在不同域(如深度域、时间域)之间转换时在测井数据和地震数据中引入误差积累及由于重采样造成的数据信息丢失；本发明方法在第一次测井和地震数据标定后，就能够在测井数据与地震数据之间建立起较为准确的测井深度-地震深度对应关系，从而可以减少后续在局部调整测井数据与地震数据之间的深度对应关系的次数。基于本发明方法产生的测井深度-地震深度对应关系能够进一步被用于建立更加准确的低频地质背景模型，进而得到更加可靠的地震反演结果。
[0058]
上述各实施例仅用于说明本发明，其中方法的各实施步骤等都是可以有所变化
的，凡是在本发明技术方案的基础上进行的等同变换和改进，均不应排除在本发明的保护范围之外。