基于小波变换模极值算法的轨道周期识别方法与流程

本发明涉及基于小波变换模极值算法的轨道周期识别方法，属于地球物理
技术领域：
。
背景技术：
：小波变换理论是近十几年发展起来的新的信号处理技术，因其在时间域和频率域都可达到高的分辨率，被称为“数学显微镜”。本发明所涉及一维连续小波变换能将一维测井信号变换成时间域和尺度域二维信息能谱，虽然能很直观的表达出信号的频谱结构，但由于能谱图存在数据冗余现象，很难精确指出其频率(周期)定量值。现有技术中，从地层测井数据中提取轨道周期信息的研究国内外学者均做过许多工作。但过去要么没有详细探讨如何从一维连续小波变换波谱定量求取周期信息，要么需结合快速傅里叶变换(FFT)来进行频谱分析，要么应用离散小波变换，但需要研究确定分解层数的问题。技术实现要素：本发明的目的在于克服现有识别方法存在的上述缺陷，提出了一种基于小波变换模极值算法的轨道周期识别方法，通过对小波变换系数矩阵做如下操作运算，可以得到优势旋回的均值定量信息。本发明是采用以下的技术方案实现的：一种基于小波变换模极值算法的轨道周期识别方法，包括如下步骤：步骤一：按照一维测井曲线模型f(x)∈L2(R)计算出连续小波变换系数矩阵：根据其中，为ψ(x)的共轭函数，若ψ(x)为实函数，则得出小波变换后系数矩阵：步骤二：对步骤一中的小波变换系数矩阵取绝对值：步骤三：对步骤二中绝对值矩阵的每一行求平均值，然后组成一个列向量：Ci‾=ΣC1j/nΣC2j/n···ΣCmj/n=C1‾C2‾···Cm‾,i=1,2,...,m;j=1,2,...,n;]]>步骤四：求出步骤三中的列向量序列极值点；步骤五：步骤四中极值点所对应的尺度即为优势周期的响应，通过小波分析中频率与尺度的对应关系式(1)换算得到优势周期的大小，Fa=Fca·Δ]]>式中：a为极值点所对应的尺度，Δ为采样周期，Fc为小波的中心频率，Fa为尺度a的准频率；步骤六：经过小波分析得到能谱图，进而得到能量环个数的比值。进一步地，步骤一中，一维测井曲线模型通过如下步骤构造：第一步：取轨道周期的主要周期比为20：40：100：400ka≈1：2：5：20，指定2π、π/2分别代表地球公转轨道的偏心率周期400、100ka，π/5代表黄赤交角的变化周期，π/10代表岁差变化周期20ka；第二步：通过函数y＝sinx，y＝sin4x，y＝sin10x，y＝sin20x来实现，把四个函数复合在一起便得到米氏周期的测井曲线模型，即y＝sinx+sin4x+sin10x+sin20x，其中，米氏周期的测井曲线模型的x为均匀取自0到10π之间的2000个点值，点间距均为2π/400；指定取样点间隔2π/400代表实际采样间距0.125m，换算后即为模拟信号对应的深度列。进一步地，步骤一的第二步中，米氏周期的测井曲线模型中所包含的时间周期分别为0.314、0.628、1.571和6.283，对应的厚度周期分别是2.5m、5m、12.5m和50m，周期个数为100、50、20和5个。进一步地，步骤二中，获得的时间周期分别是0.3093、0.6187、1.5660和6.3798，对应的厚度周期是2.46m、4.923m、12.46m和50.77m。进一步地，步骤五中，在一定的误差范围内，准频率Fa即为尺度a的真频率。进一步地，步骤五中，当Δ为实际采样厚度间距时，计算所得频率为单位厚度内的周期个数。进一步地，步骤五基于Morlet小波展开，Fc＝0.8125Hz，准频率Fa的倒数即为优势周期的大小。进一步地，步骤五中，小波系数尺度值的模平均值的极大值点是周期的响应，而且可以通过换算得到周期的大小。进一步地，步骤六中，加噪后的数据经分析以后在能谱图的体现就是能量环的颜色变化，利用能量环个数的比值，构建的模型在一定误差范围内具有可信性。本发明的有益效果是：本发明所述的基于小波变换模极值算法的轨道周期识别方法，直接借助连续小波变换能谱图进行模极值运算求取优势旋回，不仅做到直观定性与定量相结合，而且实现起来方便快捷。附图说明图1是轨道周期信号模型构造过程图。图2是米氏周期模型的小波变换能谱图。图3是小波系数模均值曲线图。图4是小波系数曲线图。图5是加噪前后信号对比图之一。图6是加噪前小波能谱图之一。图7是加噪后小波能谱图之二。图8是加噪前后信号对比图之二。图9是加噪前小波能谱图之三。图10是加噪后小波能谱图之四。图11是加噪前后信号对比图之三。图12是加噪前小波能谱图之五。图13是加噪后小波能谱图之六。图14是本发明的流程示意图。具体实施方式下面结合附图对本发明作进一步说明。如图14所示，通过对小波变换系数矩阵做如下操作运算，可以得到优势旋回的均值定量信息：(1)求取一维连续小波变换系数矩阵；(2)对小波系数矩阵取绝对值；(3)对每一列求平均值，然后组成一个行矩阵；(4)求出行矩阵的极值点；(5)极值点所对应的尺度即为优势周期的响应，通过小波分析中频率与尺度的对应关系式(1)换算得到优势周期的大小。Fa=Fca·Δ---(1)]]>式中：a为尺度，Δ为采样周期，Fc为小波的中心频率，Fa为尺度a的准频率。在一定的误差范围内，可以认为准频率Fa即为尺度a的真频率。Δ若为实际采样厚度间距时，则计算所得频率为单位厚度内的周期个数。由于是基于Morlet小波展开的，所以Fc＝0.8125，单位是Hz。准频率Fa的倒数即为优势周期的大小。此方法称为小波变换模极值法。为了验证上述方法的有效性，建立了利用正弦函数或者余弦函数的周期性构造的米氏周期测井曲线模型。由于轨道周期的主要周期比为20：40：100：400ka≈1：2：5：20，故指定2π、π/2分别代表地球公转轨道的偏心率周期400、100ka，π/5代表黄赤交角的变化周期，π/10代表岁差变化周期20ka，分别通过函数y＝sinx，y＝sin4x，y＝sin10x，y＝sin20x来实现，把四个函数复合在一起便得到米氏周期的测井反应，如图1，即y＝sinx+sin4x+sin10x+sin20x。其中x为均匀取自0到10π之间的2000个点值，点间距均为2π/400。为了便于和每米8个采样数据的常规测井信号作对比，指定取样点间隔2π/400代表实际采样间距0.125m，换算后即为模拟信号对应的深度列。这样，模型中所包含的时间周期分别为0.314，0.628，1.571，6.283，对应的厚度周期分别是2.5m、5m、12.5m、50m，周期个数为100、50、20、5个。应用小波变换模极值法对所建立的模拟曲线进行周期识别，如图2和图3，得到的时间周期分别是0.3093、0.6187、1.5660、6.3798，对应的厚度周期是2.46m、4.923m、12.46m、50.77m，而周期个数跟模型设计的个数是一样的，如图4。综上所述，米氏周期测井曲线模型经小波一维连续变换后可以准确地识别其中存在的周期，可以确定周期之比，小波系数尺度值的模平均值的极大值点是周期的响应，而且可以通过换算得到周期的大小。以上是理想信号条件下的识别情况，但是在实际实践中，由于地质情况错综复杂，地质构造复杂多变，影响因素多种多样，地层中实际的米氏周期组合以及呈现情况要复杂的多。因此，本文又构建了加入噪音信号的米氏周期组合的模型，在这里只列出四种周期组合的加入噪音信号的模型构建过程及其与未加入噪音信号的模型对比情况。图形的构造是在y值(GR)的基础上加一随机数，此时得到加噪后的图形与原图形的对比如图5。可以看出数据波形变化不是很大。按小波分析过程得到的能谱图为，如图6和图7。由此看出，由于噪声数据不明显，加噪后的数据的分析结果并没有多少明显的改变，将噪声数据的影响力加大得到的图形，如图8。将这组数据再加以分析，得到如下图形，如图9和图10。可以看出，图形形状总体没有改变，但是能量环的颜色发生了明显改变。再一次加大噪声数据影响，这时的原图和加噪后的图，如图11。经过程序分析后得到如下的能谱图，如图12和图13。由上述分析过程可以看出，加噪后的数据经分析以后在能谱图的体现就是能量环的颜色发生了改变，导致能谱图出现模糊的状况，但是能量环的个数比基本没有什么变化。而在实际研究中，需要研究以及利用的就是能量环个数的比值，由此说明，在此所构建的模型在一定误差范围内具有可信性。当然，上述内容仅为本发明的较佳实施例，不能被认为用于限定对本发明的实施例范围。本发明也并不仅限于上述举例，本
技术领域：
的普通技术人员在本发明的实质范围内所做出的均等变化与改进等，均应归属于本发明的专利涵盖范围内。当前第1页1 2 3