一种充填型地层古地貌恢复方法、装置及存储介质

1.本发明涉及古地貌恢复技术领域，具体而言，涉及一种充填型地 层古地貌恢复方法、装置及存储介质。


背景技术：

2.古地貌的研究起步于20世纪50年代，主要是在古地貌的野外描 述和地貌年龄的恢复方面的研究。我国古地貌的研究工作起于20世 纪70年代，到目前为止取得了相当不错的成绩。对我国油气资源的 勘探和开采做出了巨大的贡献。
3.目前恢复古地貌的方法主要有：沉积学分析法、残余厚度法、印 模法、层拉平法和地震古地貌法、压实恢复法、地球物理法等以及上 述几种方法的综合。上述方法针对地址填充体地貌恢复，仅有少数能 够取得较好的结果，如压实法，印模法等，但是其工作量太大，效率 低，或者结果容易受一些客观因素或资料的影响，准确度较低。


技术实现要素：

4.本发明的目的在于提供一种充填型地层古地貌恢复方法、装置及 存储介质，其在地球物理法和层拉平法的基础上建立，以层拉平法为 核心思想，地球物理法的属性提取和聚类分析为主要工具进行古地貌 的恢复，既继承了层拉平法的直观、立体化的优点，又避免了压实校 正曲线难以实现的缺点，能够高效、准确的解决充填型古地貌的恢复 问题。
5.本发明的实施例通过以下技术方案实现：一种充填型地层古地貌 恢复方法，包括如下步骤：
6.根据二维地震数据，对目标地层顶底界面进行解释，获取顶层层 位解释数据h
top
和底层层位解释数据h
bot
，将顶层层位解释数据h
top
和底层层位解释数据h
bot
做差，得到目标地层的残余厚度h
res
；
7.对二维地震数据分别计算关于地质异常体的属性，基于对地质异 常体的属性剖面进行聚类分析处理，得到地质体厚度h
fil
；
8.将残余厚度h
res
和地质体厚度h
fil
做差，得到古地貌厚度h
anc
。
9.根据一种优选实施方式，所述地质异常体的属性包括相干体、灰 度共生矩阵差异性、灰度共生矩阵均一性、灰度共生矩阵能量以及灰 度共生矩阵熵。
10.根据一种优选实施方式，所述基于对地质异常体的属性剖面进行 聚类分析处理，得到地质体厚度h
fil
包括：
11.对五种属性剖面进行聚类处理，并对目标地层充填的地质体进行 聚类结果标定；
12.按道统计目标地层充填的地质体所属聚类类型所占的时间的采 样点数，再根据采样点数与采样频率做除，将对应的采样点数转化成 地质体厚度h
fil
。
13.根据一种优选实施方式，所述聚类分析处理采用som聚类或k
‑ꢀ
means聚类其中之一。
14.根据一种优选实施方式，在对地质异常体的属性剖面进行聚类分 析处理之前，还
包括：
15.对地质异常体的属性剖面进行滤波处理。
16.根据一种优选实施方式，所述滤波处理采用mlv滤波、kuwah ara滤波或mcv滤波其中之一。
17.本发明还提供一种充填型地层古地貌恢复装置，应用到如上述 所述的方法，所述装置包括：
18.残余厚度计算模块，用于根据二维地震数据，对目标地层顶底界 面进行解释，获取顶层层位解释数据h
top
和底层层位解释数据h
bot
， 将顶层层位解释数据h
top
和底层层位解释数据h
bot
做差，得到目标地 层的残余厚度h
res
；
19.属性提取模块，用于对二维地震数据分别计算关于地质异常体的 属性；
20.聚类分析模块，用于基于对地质异常体的属性剖面进行聚类分析 处理，得到地质体厚度h
fil
；
21.古地貌厚度计算模块，用于将残余厚度h
res
和地质体厚度h
fil
做差，得到古地貌厚度h
anc
。
22.本发明还提供一种计算机可读存储介质，所述计算机可读存储介 质上存储有计算机程序，所述计算机程序被处理器执行时实现如上述 所述的方法。
23.本发明实施例的技术方案至少具有如下优点和有益效果：本发明 在层拉平法和地球物理法的基础上综合起来建立，既继承了地球物理 法需要的数据资料少，简单、高效的优点，又有调整少量参数就能运 算出结果，使古地貌恢复实现简单，高效、容易操作的优点；同时也 继承了层拉平法的直观、立体化地反映古地貌的优点；也巧妙地避免 了这两种方法中地震相的划分等人为因素影响大的部分，使得工作量 相对减小，准确度高，能够更好地为后续的油气勘探开发做指引。
附图说明
24.图1为本发明实施例1提供的充填型地层古地貌恢复方法的流程 示意图；
25.图2为本发明实施例1提供的充填型地层古地貌恢复方法的原理 流程图。
具体实施方式
26.为使本发明实施例的目的、技术方案和优点更加清楚，下面将结 合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、 完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是 全部的实施例。通常在此处附图中描述和示出的本发明实施例的组件 可以以各种不同的配置来布置和设计。
27.实施例1
28.经申请人研究发现，目前恢复古地貌的方法主要有：沉积学分析 法、残余厚度法、印模法、层拉平法和地震古地貌法、压实恢复法、 地球物理法等以及上述几种方法的综合。上述方法针对地址填充体地 貌恢复，仅有少数能够取得较好的结果，如压实法，印模法等，但是 其工作量太大，效率低，或者结果容易受一些客观因素或资料的影响， 准确度较低。
29.基于此，本发明实施例提供一种充填型地层古地貌恢复方法，其 在地球物理法和层拉平法的基础上建立，以层拉平法为核心思想，地 球物理法的属性提取和聚类分析为主
要工具进行古地貌的恢复，既继 承了层拉平法的直观、立体化的优点，又避免了压实校正曲线难以实 现的缺点，能够高效、准确的解决充填型古地貌的恢复问题。
30.需要说明的是，区别于其他的充填型地层古地貌恢复方法，考虑 在对某时期古地貌恢复时，会存在古地貌低部位，在后期过程中，低 部位可能会被其他地质体充填，但此时的古地貌并不能代表沉积时期 真实的古地貌，需要将充填的地质体抠除之后才能体现真实的古地貌。
31.基于此，本发明实施例所采用的方案如下：
32.输入二维地震资料，一方面，根据二维地震数据，对目标地层顶 底界面进行解释，获取顶层层位解释数据h
top
和底层层位解释数据h bot
，将顶层层位解释数据h
top
和底层层位解释数据h
bot
做差h
res
＝h
bot
ꢀ‑htop
，得到目标地层的残余厚度h
res
。
33.另一方面，基于地球物理法进行属性提取，在本实施例的一种实 施方式中，提取五种关于地质异常体的属性，分别为相干体atv1、 灰度共生矩阵差异性atv2、灰度共生矩阵均一性atv3、灰度共生矩 阵能量atv4以及灰度共生矩阵熵atv5并分别进行计算，具体如下：
34.对于二维地震剖面s(t)，对应的矩阵大小记为a
×
b，从灰 度为i的参考点(x，y)出发，统计其偏移角度为θ、偏移量为δ、 灰度为j的点同时出现的概率p：
35.p(i，j，δ，θ)＝{(x，y)，(x dx y dy)∈ a
×
b|f(x，y)＝i，f(x dx y dy)＝j}
36.上式中，(x dx，y dy)为偏移点，f(x，y)＝i表示参考点的灰度 值为i，f(x dx，y dy)＝j表示偏移点的灰度值为j。
37.以下对相干体atv1、灰度共生矩阵差异性atv2、灰度共生矩阵 均一性atv3、灰度共生矩阵能量atv4以及灰度共生矩阵熵atv5的 计算流程分别进行举例说明：
38.相干体atv1计算流程如下：取二维地震剖面s(t)上一点m (x，y)，则沿视倾角(p，q)的相干值c1为：
[0039][0040]
上式中，视倾角(p，q)中p和q分别为x方向和y方向上的地 震道之间的时移量，c
ii
(i＝1，2)为第i道的自相关量，c
ij
(i，j＝1， 2)为第i道和第j道的互相关量。
[0041]
遍历地震资料每一个点就可以得到s(t)的相干属性：
[0042]
atv1＝c1(pm，qn)
[0043]
上式中，(pm，qn)中的pm和qn分别代表地震剖面上各个点x方 向和y方向上的地震道之间的时移量，m的取值范围为1～a，n的取 值范围为1～b。
[0044]
灰度共生矩阵差异性atv2、灰度共生矩阵均一性atv3、灰度共 生矩阵能量atv4以及灰度共生矩阵熵atv5的计算公式分别如下：
[0045][0046]
[0047][0048][0049]
上式中，l表示灰度值的级数。
[0050]
进一步地，对上述五种属性剖面相干体atv1、灰度共生矩阵差 异性atv2、灰度共生矩阵均一性atv3、灰度共生矩阵能量atv4以 及灰度共生矩阵熵atv5分别进行滤波处理，可选的滤波方法包括mlv滤波、kuwahara滤波或mcv滤波，其中，mlv滤波器效果较好、 kuwahara滤波器效果最差，本实施例以mlv滤波为例进行举例说明， 具体如下：
[0051]
需要说明的是，value-and-criterion structure结构的滤波器有两个 元素，一个是值(value)函数，另一个是评判标准(criterion)函数。
[0052]
首先，计算值函数v(x)的表达式如下：
[0053][0054]
上式中，(x，f(x))表示属性剖面上一数据点，n表示滤波窗口 中的数据点个数，n
x
表示窗口中的数据点。
[0055]
根据评判标准函数选择最后的值函数的值，其中，评判标准函数 的表达式如下：
[0056][0057]
点(x，f(x))滤波后输出得到：
[0058][0059]
上式中，表示n
x
的平移，使之等于以位置x为中心，x
′
和 mlv(.)分别对应原数据点滤波后的输出，v和c分别表示值函数和评 判标准函数，v(x)和c(x)表示其输出。
[0060]
遍历每一个数据点得到滤波后的数据：
[0061][0062]
上式中，matv1、matv2、matv3、matv4、matv5分别代表滤 波后的5个地震属性剖面，xa表示数据点(x，f(x))中的x坐标，x
′a表示取数据点上所有x值，a的取值范围为1～a(a指数据矩阵长度)。
[0063]
进一步地，对上述5个地震属性剖面进行聚类分析处理，可选聚 类分析方法包括som聚类或k-means聚类，但是k-means聚类受初 始化影响较大，可视化效果也比som聚类差；本实施例以som聚 类为例进行举例说明，具体如下：
[0064]
对五种属性剖面进行聚类处理，并对目标地层充填的地质体进行 聚类结果标定，假设充填的地质体在som聚类结果中类型值为i， 其余地震响应特征聚类结果标定为ii；在二维地震数据上，按道统计 目标地层充填的地质体所属聚类类型所占的时间的采样点数，如按道 统计i类值所占的时间的采样点数ni，并循环到整个二维地震数据， 其中此处i表示道坐标。
[0065]
再根据采样点数ni与采样频率fs做除将对应的采 样点数转化成地质体厚度h
fil
。
[0066]
将残余厚度h
res
和地质体厚度h
fil
做差h
anc
＝h
res-h
fil
，得到古地 貌厚度h
anc
。
[0067]
按照上述步骤可以对三维地震数据进行循环操作，即可最终得到 整个三维地震数据目标地层最终去除充填地质体的古地貌。
[0068]
本发明实施例还提供一种充填型地层古地貌恢复装置，应用到 如上述所述的方法，所述装置包括：
[0069]
残余厚度计算模块，用于根据二维地震数据，对目标地层顶底界 面进行解释，获取顶层层位解释数据h
top
和底层层位解释数据h
bot
， 将顶层层位解释数据h
top
和底层层位解释数据h
bot
做差，得到目标地 层的残余厚度h
res
；
[0070]
属性提取模块，用于对二维地震数据分别计算关于地质异常体的 属性；
[0071]
聚类分析模块，用于基于对地质异常体的属性剖面进行聚类分析 处理，得到地质体厚度h
fil
；
[0072]
古地貌厚度计算模块，用于将残余厚度h
res
和地质体厚度h
fil
做差，得到古地貌厚度h
anc
。
[0073]
此外，本发明实施例还提供一种计算机可读存储介质，所述计算 机可读存储介质上存储有计算机程序，所述计算机程序被处理器执行 时实现如上述所述的方法。
[0074]
综上所述，本发明实施例的技术方案至少具有如下优点和有益效 果：本发明在层拉平法和地球物理法的基础上综合起来建立，既继承 了地球物理法需要的数据资料少，简单、高效的优点，又有调整少量 参数就能运算出结果，使古地貌恢复实现简单，高效、容易操作的优 点；同时也继承了层拉平法的直观、立体化地反映古地貌的优点；也 巧妙地避免了这两种方法中地震相的划分等人为因素影响大的部分， 使得工作量相对减小，准确度高，能够更好地为后续的油气勘探开发 做指引。
[0075]
以上仅为本发明的优选实施例而已，并不用于限制本发明，对于 本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明 的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含 在本发明的保护范围之内。