本发明涉及地球物理勘探技术领域，具体涉及与盆地波动分析有关的一种在沉积地层中识别与提取周期性波动的方法及系统。
背景技术：
沉积盆地波动分析是从一些周期的、似周期的、非周期的地质现象中分解出有规律的波动过程并建立波动方程，实现对地质过程定量研究的一种方法。波动周期识别作为盆地波动分析研究中的关键技术，目前常用的方法是滑动窗口滤波法，即：利用尺度变化的滑动窗口对沉积速率曲线进行滑动滤波分析，根据分析结果识别和提取具有周期性的波动曲线。
由于盆地沉积演化过程受多层次周期与非周期性因素控制和影响，盆地内充填的沉积往往表现为周期性兼随机性、突变性和“多时间尺度”等特征，因此要识别和提取出盆地沉积演化过程中的多层次周期性波动因素控制下的周期性波动过程，需要有更精细的沉积速率曲线。
使用滑动窗口滤波法对分辨率更高的精细沉积速率曲线进行周期性波动识别与提取时，由于更精细的沉积速率曲线的分辨率高，变化频率快，地质滤波滑动窗口尺寸变化需要以小步长变化，为得到更丰富的周期信息，还需要有更长的钻井深度，这增加了识别与提取周期性波动的难度和复杂性。因此针对更精细的盆地沉积速率曲线开展波动分析时，需要有新的、简便的方法对记录在沉积地层中的周期性波动进行识别与提取。
技术实现要素：
针对现有周期性波动识别与提取技术存在操作复杂、难度大的问题，本发明的目的是：应用小波分析技术，提出一种简便易行，且操作难度小、效率高的周期性波动识别与提取方法及系统。
本发明提供一种在沉积地层中识别与提取周期性波动的方法，包括以下步骤：
s10：建立地层沉积速率曲线；
s20：对步骤s10中的地层沉积速率曲线作离散小波分析，得到分析结果；
s30：在步骤s20得到的分析结果中寻找有稳定周期的曲线，并归纳出相应的数学方程；
s40：对步骤s30得到的所述曲线分别进行信号重构，此时得到的曲线即为相应周期性波动影响下的综合沉积波动曲线；
s50：根据步骤s40得到的波动曲线建立相应的的波动方程；
s60：根据步骤s50建立的波动方程做出钻井未钻达的地层的沉积波动曲线预测。
在本发明的一个实施例中，步骤s20中所述的离散小波分析是以meyer小波作为母小波进行离散小波变换分析。
在本发明的一个实施例中，步骤s20中所述分析结果为一系列分解得到的近似分量曲线an和细节分量曲线dn，n取值为一系列正整数，所述正整数取自正整数集{1，2，……，k-1，k}，其中k为分解的次数。
在本发明的一个实施例中，步骤s30中所述的曲线为周期性曲线或近周期性曲线。
在本发明的一个实施例中，步骤s40中所述的信号重构是根据离散小波变换mallat算法的信号重构过程进行的。
在本发明的一个实施例中，步骤s50中所述的波动方程依据离散小波变换mallat算法的信号重构公式得到，所述信号重构公式如下：
an-1＝an+dn，
其中，为a近似分量曲线，d为细节分量曲线，n为分解的次数。
在本发明的一个实施例中，步骤s60中所述的沉积波动曲线预测，通过沉积速率的变化对未钻达的地层的地层情况进行预测。
在本发明的一个实施例中，所述预测的内容为未钻达地层的岩性、沉积环境变化、埋藏史和/或构造变动。
本发明的另一个方面，还提供一种在沉积地层中识别与提取周期性波动的系统，包括：
沉积速率曲线单元，用于建立地层沉积速率曲线；
离散小波分析单元，用于对沉积速率曲线单元中的地层沉积速率曲线作离散小波分析，得到分析结果；
周期搜索单元，用于从离散小波分析单元中得到的分析结果中寻找有稳定周期的曲线，并归纳出相应的数学方程；
波动曲线恢复单元：用于对周期搜索单元中得到的所述曲线分别进行信号重构，此时得到的曲线即为相应周期性波动影响下的综合沉积波动曲线；
波动方程建立单元，用于根据波动曲线恢复单元中得到的波动曲线建立相应的波动方程；以及
波动曲线预测单元，用于波动方程建立单元中建立的波动方程做出钻井未钻达的地层的沉积波动曲线预测。
在本发明的一个实施例中，波动曲线预测单元通过沉积速率的变化对未钻达的地层的地层情况进行预测。
本发明中的在沉积地层中识别与提取周期性波动的方法及系统，扩充了波动分析的技术方法，可以从沉积速率曲线中快速分解和提取出周期波，操作简便易行，有利于波动分析研究在实际工作中的推广应用，解决生产中的实际问题。
附图说明
在下文中将基于实施例并参考附图来对本发明进行更详细的描述。其中：
图1是本发明实施例一在沉积地层中识别与提取周期性波动的方法的流程图；
图2是本发明实施例二中步骤s10选取的沉积速率曲线；
图3是本发明实施例二中步骤s20中分解得到的近似分量曲线和细节分量曲线；
图4是本发明实施例二中在沉积地层中识别与提取周期性波动的方法应用结果与采用滑动滤波法分析出的结果的对比图；
图5是本发明实施例三中的沉积地层中识别与提取周期性波动的方法应用结果与采用滑动滤波法分析出的结果的对比图；
图6是本发明实施例四中的沉积地层中识别与提取周期性波动的结构示意图。
在附图中，相同的部件使用相同的附图标记。附图并未按照实际的比例。
具体实施方式
为了使本发明所揭示的技术内容更加详尽和完备，下面参照附图和具体实施例详细地说明本发明的目的、技术方案和技术效果。需要特别说明的是，虽然是针对附图进行说明，但是本领域的技术人员应当理解，附图及实施例并非用来限制本发明所涵盖的范围。
本申请中涉及到的时间单位ma(millionanniversary)为百万年；mallat算法是由s.mallat在用于图像分解的金字塔算法启发下，结合多分辨率分析理论，提出的信号的塔式多分辨率分析的分解与重构算法，简称mallat快速算法，是一种快速正交小波分解和重构算法，它在小波分析中的地位等同于快速傅里叶变换在傅里叶分析中的地位，该算法甚至不需要尺度函数和小波函数的具体表达式，只需要二尺度方程中的滤波器系数，就可以进行小波分解与重构。
meyer小波是由法国数学家y.meyer于1986年构造的第一个具有任意阶光滑性的规范正交小波，具有任意阶连续性、对称性和指数衰减特性，它无时域表达式，是通过一对共轭正交镜像滤波器组的频谱来定义，其小波函数在频率域的定义为式(a)：
式(a)中，v(x)为构造meyer小波的辅助函数，满足式(b)
即：
且
ν(x)+ν(1-x)＝1(d)。
dmeyer小波是离散meyer小波，由于meyer小波用于离散小波变换时不能使用快速离散小波变换算法，因此在使用mallat算法时，需使用离散meyer小波，即dmeyer小波。
实施例一
如图1所示，为本实施例中的在沉积地层中识别与提取周期性波动的方法的流程图，本实施中的方法包括以下步骤：
s10：建立地层沉积速率曲线；
s20：对步骤s10中的地层沉积速率曲线作离散小波分析，得到分析结果；
s30：在步骤s20得到的分析结果中寻找有稳定周期的曲线，并归纳出相应的数学方程；
s40：对步骤s30得到的所述曲线分别进行信号重构，此时得到的曲线即为相应周期性波动影响下的综合沉积波动曲线；
s50：根据步骤s40得到的波动曲线建立相应的的波动方程；
s60：根据步骤s50建立的波动方程做出钻井未钻达的地层的沉积波动曲线预测。
实施例二
本实施例中，就实施例一中的方法的具体的应用，具体为对文章《柴达木盆地新生代波动过程及与油气关系》(金之钧，李京昌，汤良杰等.柴达木盆地新生代波动过程及与油气关系[j].地质学报，2006，80(3)：359-365。)中使用的沉积速率直方图通过离散采样得到的沉积速率曲线，使用本发明提出的方法进行周期波识别与提取，可根据使用本发明得到的识别结果与文章中应用滑动滤波法得到的结果进行对比验证。具体过程以及结果的对比在本实施例中进行详述。
步骤s10：建立地层沉积速率曲线；
层沉积速率曲线具体的建立方法现有技术中已有广泛的记载，本领域技术人员根据现有技术可以方便完成，在此不在赘述。在本实施例中，具体的选取的沉积速率曲线如图2所示，为文章《柴达木盆地新生代波动过程及与油气关系》(金之钧，李京昌，汤良杰等.柴达木盆地新生代波动过程及与油气关系[j].地质学报，2006，80(3)：359-365。)中使用的沉积速率直方图通过离散采样得到的沉积速率曲线。
步骤s20：对步骤s10中的地层沉积速率曲线作离散小波分析，得到分析结果；本实施例中离散小波变换使用mallat算法，母小波meyer需使用离散meyer小波，即dmeyer小波。其中的沉积速率曲线经本步骤离散小波变换后得到一系列分解的近似分量曲线an和细节分量曲线dn，n取值为一系列正整数，所述正整数取自正整数集{1，2，……，k-1，k}，其中k为分解的次数。
具体的曲线如图3所示，其中a为近似分量曲线和d为细节分量曲线：为沉积速率曲线经步骤s20作离散小波变换分解出的近似分量曲线an和细节分量曲线dn，在本实施例中，k为6。
步骤s30：在步骤s20得到的分析结果中寻找有稳定周期的曲线，并归纳出相应的数学方程；在本实施例中，所述的曲线为周期性曲线或近周期性曲线，该曲线对应的波动即为以该曲线周期值为周期的周期性波动；
具体的，由于地质过程具有复杂性、非平稳性和在某一趋势下演化的非线性等特点，导致周期性波动中的每一波动不是数学意义上的周期严格不变的绝对周期波，它们的周期会围绕某一固定值作微小变化，表现为基本或大致围绕某一周期作周期性旋回演化，这一周期在某一地质历史时期可能长点，在另一地质历史时期可能又变得短点，因此在细节分量曲线中找到周期值围绕某一数值在误差范围内基本稳定的曲线即可。
具体的，在本实施例中，本实施例在步骤s20得到的细节分量中找到周期分别为20ma(图3中曲线d4)和9.8ma(图3中曲线d3)的两条有稳定周期的曲线。从图3中相应曲线上读出初始相位值，并根据振幅变化情况拟合出振幅函数，然后写出这两条周期曲线的近似方程，如式(1)和式(2)所示：
步骤s40：对步骤s30得到的所述曲线分别进行信号重构，此时得到的曲线即为相应周期性波动影响下的综合沉积波动曲线；在本实施例中，信号重构是根据离散小波变换mallat算法的信号重构过程进行的。
根据离散小波变换mallat算法的信号重构公式进行信号重构，此时得到的曲线即为对应周期性波动影响下的综合沉积速率波动曲线。依据该重构过程，本实施例中20ma周期波动(图3中的曲线d4)影响下的综合沉积速率波动曲线为近似分量a3(图3中曲线a3)，9.8ma周期波动(图3中曲线d3)影响下的综合沉积波动曲线为近似分量a2(图3中曲线a2)。
步骤s50：根据步骤s40得到的波动曲线建立相应的的波动方程；本实施例中，波动方程依据离散小波变换mallat算法的信号重构公式得到，所述信号重构公式为：an-1＝an+dn，(其中，为a近似分量曲线，d为细节分量曲线，n为分解的次数。)
在本实施例中，根据步骤s40中的分解曲线，可以建立综合沉积速率波动曲线a3和a2的波动方程，分别如式(3)和(4)所示：
步骤s60：根据步骤s50建立的波动方程做出钻井未钻达的地层的沉积波动曲线预测。本实施例中，步骤s60中所述的沉积波动曲线预测，通过沉积速率的变化对未钻达的地层的地层情况进行预测。预测的内容为未钻达地层的岩性、沉积环境变化、埋藏史和/或构造变动等。
图4为本实施例中的在沉积地层中识别与提取周期性波动的方法应用结果与采用滑动滤波法分析出的结果的对比图：图4中纵坐标为时间，单位ma(millionanniversary，百万年)，横坐标为沉积速率，单位m/ma(米/百万年，米每百万年)。图4中间两栏为文章《柴达木盆地新生代波动过程及与油气关系》中使用滑动滤波法分析出的结果，即为原图；图4中左端栏和右端栏曲线为应用本发明实施例分析出的结果。在文章《柴达木盆地新生代波动过程及与油气关系》应用滑动滤波法分析出的结果，即原图中，曲线m为30ma滑动窗得到的沉积速率滑动平均曲线，曲线n为15ma滑动窗得到的沉积速率滑动平均曲线，曲线y为5ma滑动窗得到的沉积速率滑动平均曲线，曲线n1为周期为20ma的波动曲线，曲线n2为周期为9.6ma的波动曲线；在应用本发明实施例分析出的结果中，左端栏中的曲线m′与原图中的曲线m对应，曲线n′与原图中的曲线n对应，曲线y′与原图中的曲线y对应，右端栏中的周期波n1′与原图分解的周期波n1对应，周期波n2′与周期波n2对应，从对比结果看，应用本发明分析出的综合沉积波动曲线和周期波都与文章《柴达木盆地新生代波动过程及与油气关系》应用滑动滤波法分析出的结果基本一致。
实施例三
为了更进一步的验证本发明中的方法的有效性，本实施例将本发明中的方法应用于文章《塔里木盆地高频波识别及其意义》(范国章，金之钧，刘国臣，等.塔里木盆地高频波识别及其意义[j].沉积学报，2001，19(2)：245-248。)中使用的沉积速率曲线进行分析。具体方法与实施例二相同，在此不在赘述。本实施例中，如图5所示为应用本发明中的方法进行分析的结果与文章中应用滑动窗口滤波法得到的分析结果的对比图。如图5所示，曲线a′和曲线b′为应用本发明中的方法得到的波动曲线；曲线a和曲线b为文章中应用滑动窗口滤波法得到的波动曲线。其中：曲线a′与曲线a对应，曲线b′与曲线b对应。由图中的对比可知，二者基本一致。进一步证明本发明中的方法的有效性。
实施例四
如图6所示，为本发明中的在沉积地层中识别与提取周期性波动的系统的结构示意图，该系统包括：
沉积速率曲线单元1，用于建立地层沉积速率曲线；
离散小波分析单元2，用于对沉积速率曲线单元中的地层沉积速率曲线作离散小波分析，得到分析结果；
周期搜索单元3，用于从离散小波分析单元中得到的分析结果中寻找有稳定周期的曲线，并归纳出相应的数学方程；
波动曲线恢复单元4：用于对周期搜索单元中得到的所述曲线分别进行信号重构，此时得到的曲线即为相应周期性波动影响下的综合沉积波动曲线；
波动方程建立单元5，用于根据波动曲线恢复单元中得到的波动曲线建立相应的的波动方程；以及
波动曲线预测单元6，用于波动方程建立单元中建立的波动方程做出钻井未钻达的地层的沉积波动曲线预测。
优选的，离散小波分析单元2以dmeyer小波为母小波，采用mallat算法对沉积速率曲线进行离散小波变换分析。
优选的，离散小波分析单元2的分析结果为一系列分解得到的近似分量曲线an和细节分量曲线dn，其中n为分解的次数，n取值为一系列正整数，所述正整数取自正整数集{1，2，……，k-1，k}，其中k为分解的次数。
优选的，周期搜索单元3寻找的有稳定周期的曲线为周期性曲线或近周期性曲线。
优选的，波动曲线恢复单元4根据离散小波变换mallat算法的信号重构过程对周期搜索单元中得到的所述曲线分别进行信号重构，此时得到的曲线即为相应周期性波动影响下的综合沉积波动曲线；
优选的，波动方程建立单元5根据波动曲线恢复单元中得到的波动曲线建立相应的的波动方程，具体的依据离散小波变换mallat算法的信号重构公式得到，所述信号重构公式如下：
an-1＝an+dn，
其中，为a近似分量曲线，d为细节分量曲线，n为分解的次数。
优选的，波动曲线预测单元6通过沉积速率的变化对未钻达的地层的地层情况进行预测。预测的内容为未钻达地层的岩性、沉积环境变化、埋藏史和/或构造变动。
周期性波动识别作为盆地波动分析研究中的关键点，在实际波动分析研究中有较大难度。本发明中的在沉积地层中识别与提取周期性波动的方法及系统可以有效降低这一难度，使波动周期识别与提取简便易行且高效。
虽然已经参考优选实施例对本发明进行了描述，但在不脱离本发明的范围的情况下，可以对其进行各种改进并且可以用等效物替换其中的部件。尤其是，只要不存在结构冲突，各个实施例中所提到的各项技术特征均可以任意方式组合起来。本发明并不局限于文中公开的特定实施例，而是包括落入权利要求的范围内的所有技术方案。