CN106886044A - 一种基于剪切波与Akaike信息准则的微地震初至拾取方法 - Google Patents

Abstract

本发明涉及一种基于剪切波与Akaike信息准则的微地震初至拾取方法，属于实际地震勘探环境下微地震初至自动提取方法研究。针对微地震信号表现出的统计特性，把背景噪声和微地震信号看作两个统计特性不同的平稳过程，利用剪切波变换对微地震信号进行分解、重构，同时结合AIC模型对重构的各层信号进行初至拾取，最后综合各层的拾取结果确定微地震二维数据初至拾取结果。本发明所提出的方法无需对复杂背景噪声进行预处理，节约时间、计算成本，免除了去噪过程中微地震信号的能量损失，更有效的保留微地震二维记录中有效信号的幅频特性，为油气藏动态监测、地热动态监测、煤田开采和工程动态监测提供重要依据。

Description

一种基于剪切波与Akaike信息准则的微地震初至拾取方法

技术领域

本发明涉及的是实际地震勘探环境下微地震初至自动提取方法研究。

背景技术

在监测到的微地震记录中，除包含微弱的微地震信号外，还存在大量的随机背景噪声，人工拾取微地震有效事件需要很大的工作量，低信噪比情况下微地震记录初至的自动拾取十分困难，误差很大。目前，微地震初至拾取方法主要有长短时窗比法，分形分维法及自回归模型法等，但这些方法处理较高信噪比记录尚可，处理低信噪比记录时就极难做到准确拾取。

实际的微地震记录中，背景噪声的成分和性质十分复杂，几乎完全淹没微地震信号，使微地震初至的自动拾取极为困难。单一化的预处理(如：去噪)过程无法消除成分复杂的背景噪声，甚至还会在处理过程中对微地震记录的幅频特性造成损失。

近几年发展起来一种具有多尺度多方向特性的剪切波变换。剪切波变换不但能够保留传统小波变换的数学框架，还具有检测方向的能力。剪切波通过对一些生成函数进行抛物尺度缩放、剪切及位移而产生。因此，除位移参数以外，还是关于尺度参数a>0及剪切系数的变换。对于连续剪切波系统，抛物尺度矩阵A_a及剪切矩阵S_s分别定义如下：

A_a用于改变分辨率，S_s则用于改变方向，位移参数t用于改变位置。剪切波系统可定义为：

相应的剪切波变换为

剪切波变换的位置参数t可以让剪切波在二维图像上任意移动，即可以检测到二维图像上任何位置的信息。剪切系数s决定了剪切波方向，通过改变剪切系数可以获得任意方向的剪切波，即可以检测到任意方向的奇异性曲线。尺度参数a可以使剪切波伸缩，即可以检测到任意大小的线状奇异性形状。

剪切波变换作为一种新构造出来的多尺度几何分析函数，集合了许多函数的优良特性,同时还具有更简单的数学结构,为解决图像处理中的难点问题提供了新的研究思路。与传统小波变换不同,剪切波具有更好的局部化特性和很强的方向敏感性，使得剪切波对图像有很好的稀疏表示能力。此外，具有更好的光滑性和多分辨分率析能力,同时利用其良好的衰减性可以得到滤波性能更好的滤波器，有利于实现特征提取、图像消噪、水印嵌入与识别等任务。

发明内容

本发明提供一种基于剪切波与Akaike信息准则的微地震初至拾取方法，目的是研究复杂随机噪声环境下，既无需进行预处理，又节约时间、计算成本的微地震初至拾取方法。本发明利用剪切波变换对微地震信号进行分解、重构，同时结合AIC模型拾取微地震二维数据初至的算法，该算法针对微地震信号表现出的统计特性，把背景噪声和微地震信号看作两个统计特性不同的平稳过程，免除了去噪过程中微地震信号的能量损失，更有效的保留微地震二维记录中有效信号的幅频特性，能够准确的自动拾取微地震二维记录初至，为微地震震源定位及精确探明地质构造提供了重要依据。

本发明采取的技术方案是，包括下列步骤：

1)微地震二维记录获取：

微地震监测分为井中监测和地面监测两种方式。通过在井中(比如压裂井)或周围的地面布置检波器排列接收生产活动所产生或诱导的微小地震事件，注入作业期间引发的微地震事件在空间和时间上的分布是复杂的，但不是随机的，1Km范围内用适当的灵敏仪器就可以接收到。大多数微地震事件频率范围介于50～1500Hz之间，持续时间小于1s，在微地震记录上一般表现为清晰的脉冲；越弱的微地震事件，其频率越高，持续时间越短；能量越小，破裂的长度就越短。因此微地震信号很容易受其周围噪声的影响或遮蔽，对低信噪比微地震记录进行初至拾取更困难，误差也更大，

2)微地震二维记录的分频处理

(1)基于剪切波的多尺度分解

根据Shearlet变换的多尺度分解特性，将原始二维记录u(i,j)进行拉普帕斯金字塔(LP)分解，得到一次尺度分解结果的高频部分H₁及低频成分L₁如下式：

L₁(i,j)＝∑h(m,n)u(2i+m,2j+n) (1)

其中h为分析滤波器的冲击响应，g为综合滤波器的冲击响应，g与h具有对称性。对上述低频成分L₁继续进行第二次尺度分割，产生高频成分H₂及低频成分L₂；之后每一次尺度分解都针对上一次分解的低频结果而进行，可通过下式迭代得到各层的的高频部分H_k及低频部分L_k，k≥2。

L_k(i,j)＝∑h(m,n)L_k-1(2i+m,2j+n) (3)

经n次分解后，共产生n+1个尺度层，分别为H₁、H₂…、H_n及L_n；

(2)重构记录的获取

根据微地震信号及噪声的频率特性对上述多尺度分解结果进行重构，方案为：高频成分可以选择m个尺度层H₁、H₂、…、H_m叠加后整体经剪切波反变换重构为u_a，低频成分H_m+1…H_n及L_n可以分别重构u_b、u_c、…；

3)基于AR-AIC的初至拾取

(1)单道记录的初至拾取

若单道微地震信号x(t)是一个平稳过程，则可以用AR模型表示为：

其中M是模型阶数，a(m)是系数，ε(t)是具有零均值和方差为的平稳白噪声序列。某时刻t_k对应的Akaike信息准则值(AIC)应为：

其中和分别为t_k时刻前后的信号方差，N为整个序列的长度，C为常数；

假设t_p时刻为初至，即认为系统模型在t_p点前后分别对应随机噪声和含噪微地震信号，用AR模型表示分别为：

式中M₁和M₂分别为噪声和微地震信号的阶数，N为微地震记录的长度，t_p时刻前后的序列分属不同的平稳过程，即整个序列为非平稳过程，故利用(6)式求得的AR-AIC值会对应局部最小值，在(6)的求解中，往往需要不断的调试和修正去得到各段AR模型的阶数，而且在模型阶数过大的情况下需要对高阶矩阵进行求逆，使得运算量过大，故参照Maeda提出的简便算法，利用下式求取AR-AIC值：

根据(9)式的结果找到第一个局部极值点并将该点判定为初至时刻；

(2)微地震二维记录的初至判断

首先利用式(9)对二维记录的分频处理结果u_a、u_b、u_c、…的各道分别进行判断，将各道AR-AIC极小值对应的时刻标记到各分频处理结果中。其次根据不同地区数据的信噪比及信号的频率成分，从标记好的分频处理结果中选择差异较小的AR-AIC结果，经适当加权最终得出微地震二维记录的初至判断结果。

本发明旨在利用剪切波变换的多尺度多方向性结合AIC模型对微地震记录进行分解、重构，在低信噪比情况下精确拾取微震初至。利用剪切波变换将含噪图像分解，检测，可在变换域进行硬阈值收缩，分离出噪声分量，达到噪声消减的目的；利用剪切波对奇异曲线的敏感性，可进行图像的边缘检测；利用剪切波的稀疏性，可在压缩感知方法中选择剪切波字典对图像进行复原。

本发明利用剪切波变换对微地震信号进行分解、重构，同时结合AIC模型对重构的各层信号进行初至拾取，最后综合各层的拾取结果确定微地震二维数据初至。本发明可比现有方法更高效更准确的进行初至拾取，为微地震震源定位及精确探明地质构造提供了重要依据，有利于准确估计页岩气的储量及分布。

初至拾取是静校正处理的前提和基础，初至拾取的精度不但直接影响静校正的好坏，也影响着后续的噪声处理，甚至是微地震事件的定位精度。本发明在低信噪比情况下进行精确的初至拾取，有利于正确进行微地震静校正，为微地震资料的噪声处理提供依据，并为后续的定位工作奠定基础。

本发明的优点：微地震震源定位要求精确确定初至，人工检测和拾取初至效率较低，且在低信噪比条件下容易出现误判。本发明针对微地震信号表现出的统计特性，把背景噪声和微地震信号看作两个统计特性不同的平稳过程，利用剪切波变换对微地震信号进行分解、重构，同时结合AIC模型对重构信号进行初至拾取，无需对复杂背景噪声进行预处理，节约时间、计算成本，免除了去噪过程中微地震信号的能量损失，更有效的保留微地震二维记录中有效信号的幅频特性，自动拾取速度快，可比现有方法更高效更准确的拾取微地震初至，更有利于后续处理及解释工作，具有较强的实用性。

附图说明

图1是基于剪切波变换的多尺度多方向剖分示意图，对原始二维微地震记录采用非下采样拉普拉斯金字塔滤波器(NSP)经1级尺度分解(k＝1)后，得到高频成分H₁及低频成分L₁共2个部分，对上述低频成分L₁继续进行第2次尺度分解(k＝2)，继续产生高频成分H₂及低频成分L₂；以此类推，经n次分解(k＝n)后共产生n+1个成分，分别为H₁、H₂…、H_n及L_n。上述多尺度分解结果可进一步通过剪切滤波器(SF)来实现方向局部化；

图2是模拟的21道微地震合成记录，共有4个事件，每道采样点数为512点，人工拾取的初至点已经用红色的箭头标出；

图3是图2记录加入-10dB高斯白噪声的合成记录，人工拾取较困难，且由于随机噪声的影响，拾取结果会存在一定的误差；

图4是利用传统AIC方法对图3进行初至拾取的结果(红色箭头)，由于随机噪声能量远大于微地震信号能量，AIC无法对初至进行准确的拾取，约半数道存在较大的误差；

图5是本发明所述方法重构后第一层记录部分道的初至拾取结果，初至位置已用红色脉冲标出；

图6是本发明所述方法重构后第二层记录部分道的初至拾取结果，红色脉冲标出；

图7是本发明所述方法重构后第三层记录部分道的初至拾取结果，红色脉冲标出；

图8是综合各层结果得出对图3最终的初至拾取结果(红色箭头)；

图9是井中实际微地震数据，共15道，每道512采样点，利用本发明方法所拾取的初至如图中红色箭头。

具体实施方式

包括下列步骤：

1)微地震二维记录获取：

微地震监测分为井中监测和地面监测两种方式。通过在井中(比如压裂井)或周围的地面布置检波器排列接收生产活动所产生或诱导的微小地震事件，注入作业期间引发的微地震事件在空间和时间上的分布是复杂的，但不是随机的，1Km范围内用适当的灵敏仪器就可以接收到。大多数微地震事件频率范围介于50～1500Hz之间，持续时间小于1s，在微地震记录上一般表现为清晰的脉冲；越弱的微地震事件，其频率越高，持续时间越短；能量越小，破裂的长度就越短。因此微地震信号很容易受其周围噪声的影响或遮蔽，对低信噪比微地震记录进行初至拾取更困难，误差也更大，

2)微地震二维记录的分频处理

(1)基于剪切波的多尺度分解

根据Shearlet变换的多尺度分解特性，将原始二维记录u(i,j)进行拉普帕斯金字塔(LP)分解，得到一次尺度分解结果的高频部分H₁及低频成分L₁如下式：

L₁(i,j)＝∑h(m,n)u(2i+m,2j+n) (1)

其中h为分析滤波器的冲击响应，g为综合滤波器的冲击响应，g与h具有对称性。对上述低频成分L₁继续进行第二次尺度分割，产生高频成分H₂及低频成分L₂；之后每一次尺度分解都针对上一次分解的低频结果而进行，可通过下式迭代得到各层的的高频部分H_k及低频部分L_k，k≥2。

L_k(i,j)＝∑h(m,n)L_k-1(2i+m,2j+n) (3)

经n次分解后，共产生n+1个尺度层，分别为H₁、H₂…、H_n及L_n；

(2)重构记录的获取

根据微地震信号及噪声的频率特性对上述多尺度分解结果进行重构，方案为：高频成分可以选择m个尺度层H₁、H₂、…、H_m叠加后整体经剪切波反变换重构为u_a，低频成分H_m+1…H_n及L_n可以分别重构u_b、u_c、…；

3)基于AR-AIC的初至拾取

(1)单道记录的初至拾取

若单道微地震信号x(t)是一个平稳过程，则可以用AR模型表示为：

其中M是模型阶数，a(m)是系数，ε(t)是具有零均值和方差为的平稳白噪声序列。某时刻t_k对应的Akaike信息准则值(AIC)应为：

其中和分别为t_k时刻前后的信号方差，N为整个序列的长度，C为常数；

假设t_p时刻为初至，即认为系统模型在t_p点前后分别对应随机噪声和含噪微地震信号，用AR模型表示分别为：

式中M₁和M₂分别为噪声和微地震信号的阶数，N为微地震记录的长度，t_p时刻前后的序列分属不同的平稳过程，即整个序列为非平稳过程，故利用(6)式求得的AR-AIC值会对应局部最小值，在(6)的求解中，往往需要不断的调试和修正去得到各段AR模型的阶数，而且在模型阶数过大的情况下需要对高阶矩阵进行求逆，使得运算量过大，故参照Maeda提出的简便算法，利用下式求取AR-AIC值：

根据(9)式的结果找到第一个局部极值点并将该点判定为初至时刻；

(2)微地震二维记录的初至判断

首先利用式(9)对二维记录的分频处理结果u_a、u_b、u_c、…的各道分别进行判断，将各道AR-AIC极小值对应的时刻标记到各分频处理结果中。其次根据不同地区数据的信噪比及信号的频率成分，从标记好的分频处理结果中选择差异较小的AR-AIC结果，经适当加权最终得出微地震二维记录的初至判断结果。

初至拾取是静校正处理的前提和基础，初至拾取的精度不但直接影响静校正的好坏，也影响着后续的噪声处理，甚至是微地震事件的定位精度。本发明在低信噪比情况下进行精确的初至拾取，有利于正确进行微地震静校正，为微地震资料的噪声处理提供依据，并为后续的定位工作奠定基础。

应用举例：

模拟21道地面微地震记录，主频20Hz，如图2所示。该段记录共有4个事件，每道采样点数为512点。人工拾取的初至点已经用红色的箭头标出。在该记录中加入-10dB能量的高斯白噪声，合成记录如图3所示，由于信噪比较低，人工拾取初至困难且误差较大。利用传统AIC方法对含噪记录进行初至拾取是较常用的方法，图4给出了传统AIC方法的拾取结果，有半数道拾取较准确，其余偏差较大。本发明提出基于剪切波与AIC准则的微地震初至拾取方法，步骤如下：

1.基于剪切波变换的微地震记录分频处理

利用图1对含噪记录图3进行基于剪切波变换的多尺度分解。采用非下采样拉普拉斯金字塔滤波器组(NSPS)经4级尺度分解后得到5个频率成分，分别为H₁、H₂、H₃、H₄及L₄。根据信号的主频对上述成分进行重构，其中L₄重构为第一层记录，H₄重构为第二层记录，H₁、H₂、H₃合并重构为第三层记录。图5、图6及图7分别为第一、第二及第三层重构后的结果(只给出1-13道)。

2.基于AR-AIC的微地震初至拾取

对各层重构结果利用AR-AIC方法(式13)进行初至拾取，各层的拾取结果见图5、图6及图7中的红色脉冲。为去除极端值的影响，对每道数据进行处理时，将三层拾取结果中偏差较大的结果舍去，其余两个结果求取均值作为该道数据的最终初至拾取结果，如图8中红色箭头标识。为进一步说明结果，表1给出了传统AIC方法、小波变换-AIC方法和本发明方法的拾取结果对比。

表1不同方法对21道合成微地震数据初至拾取结果对比

计算上述方法与理论值的均方误差可知，AIC方法为245.7，小波-AIC方法为105.1，本发明方法仅为13.2。说明本发明基于剪切波变换的AIC初至拾取方法结果更准确。

对如图9所示的微地震井中实际数据进行处理，仍利用上述步骤先采用非下采样拉普拉斯金字塔滤波器组(NSPS)经2级尺度分解后得到3个频率成分，分别为H₁、H₂及L₂。根据信号的频带在重构时将L₂重构为第一层记录，H₂重构为第二层记录，H₁重构为第三层记录，再基于AR-AIC准则拾取初至。对每一层拾取初至后，将3层结果中偏差较大的结果舍去，再综合其余两结果得出各道的最终初至点，如图9红色箭头所示。可见本发明方法能够在较低信噪比记录中拾取微地震事件初至，尤其是第1道、第7道局部信噪比很低的情况下初至拾取较准确。

综上所述，基于剪切波与AR-AIC准则的微地震初至拾取方法是有效的，实用的，它无需预先对噪声进行消减处理，可以在信噪比较低的情况下准确、自动的识别微地震事件的初至，有利于正确进行微地震静校正，为微地震资料的噪声处理提供依据，并为后续的定位工作奠定基础。

Claims (4)

1.一种基于剪切波与Akaike信息准则的微地震初至拾取方法，其特征在于，包括下列步骤：

1)微地震二维记录获取：

微地震监测分为井中监测和地面监测两种方式，通过在井中或周围的地面布置检波器排列接收生产活动所产生或诱导的微小地震事件，注入作业期间引发的微地震事件在空间和时间上的分布是复杂的，但不是随机的，1Km范围内用适当的灵敏仪器就可以接收到，大多数微地震事件频率范围介于50～1500Hz之间，持续时间小于1s，在微地震记录上一般表现为清晰的脉冲；越弱的微地震事件，其频率越高，持续时间越短；能量越小，破裂的长度就越短；

2)微地震二维记录的分频处理：

(1)基于剪切波的多尺度分解：

(2)重构记录的获取：

根据微地震信号及噪声的频率特性对上述多尺度分解结果进行重构，方案为：高频成分可以选择m个尺度层H₁、H₂、…、H_m叠加后整体经剪切波反变换重构为u_a，低频成分H_m+1…H_n及L_n可以分别重构u_b、u_c、…；

3)基于AR-AIC的初至拾取：

(1)单道记录的初至拾取

(2)微地震二维记录的初至判断

首先利用式(9)对二维记录的分频处理结果u_a、u_b、u_c、…的各道分别进行判断，将各道AR-AIC极小值对应的时刻标记到各分频处理结果中，其次根据不同地区数据的信噪比及信号的频率成分，从标记好的分频处理结果中选择差异较小的AR-AIC结果，经适当加权最终得出微地震二维记录的初至判断结果。

2.根据权利要求1所述的一种基于剪切波与Akaike信息准则的微地震初至拾取方法，其特征在于，监测井采用压裂井。

3.根据权利要求1所述的一种基于剪切波与Akaike信息准则的微地震初至拾取方法，其特征在于，步骤2)微地震二维记录的分频处理中(1)基于剪切波的多尺度分解，具体包括：

根据Shearlet变换的多尺度分解特性，将原始二维记录u(i,j)进行拉普帕斯金字塔(LP)分解，得到一次尺度分解结果的高频部分H₁及低频成分L₁如下式：

L₁(i,j)＝∑h(m,n)u(2i+m,2j+n) (1)

$H_1(i,j)=u(i,j)-\mathrm{\Σ}g(m,n)L_1(\frac{i-m}{2},\frac{j-n}{2})---\left(2\right)$

其中h为分析滤波器的冲击响应，g为综合滤波器的冲击响应，g与h具有对称性，对上述低频成分L₁继续进行第二次尺度分解，产生高频成分H₂及低频成分L₂；之后每一次尺度分解都针对上一次分解的低频结果而进行，可通过下式迭代得到各层的的高频部分H_k及低频部分L_k，k≥2；

L_k(i,j)＝∑h(m,n)L_k-1(2i+m,2j+n) (3)

$H_k(i,j)=L_{k-1}(i,j)-\mathrm{\Σ}g(m,n)L_k(\frac{i-m}{2},\frac{j-n}{2})---\left(4\right)$

经n次分解后，共产生n+1个尺度层，分别为H₁、H₂…、H_n及L_n。

4.根据权利要求1所述的一种基于剪切波与Akaike信息准则的微地震初至拾取方法，其特征在于，步骤3)基于AR-AIC的初至拾取中(1)单道记录的初至拾取，包括：

若单道微地震信号x(t)是一个平稳过程，则可以用AR模型表示为：

$x\left(t\right)=\underset{m=1}{\overset{M}{\mathrm{\Σ}}}a\left(m\right)x(t-m)+\mathrm{\ϵ}\left(t\right)---\left(5\right)$

其中M是模型阶数，a(m)是系数，ε(t)是具有零均值和方差为的平稳白噪声序列；某时刻t_k对应的Akaike信息准则值(AIC)应为：

$A\left(k\right)=(k-M)\lg \left({\mathrm{\σ}}_1^2\right)+(N-M-k)\lg \left({\mathrm{\σ}}_2^2\right)+C---\left(6\right)$

其中和分别为t_k时刻前后的信号方差，N为整个序列的长度，C为常数；

假设t_p时刻为初至，即认为系统模型在t_p点前后分别对应随机噪声和含噪微地震信号，用AR模型表示分别为：

$x\left(t\right)=\underset{m=1}{\overset{M_1}{\mathrm{\Σ}}}{a}_1\left(m\right)x(t-m)+{\mathrm{\ϵ}}_1\left(t\right),t_1\≤t\≤t_p---\left(7\right)$

$x\left(t\right)=\underset{m=1}{\overset{M_2}{\mathrm{\&Sigma;}}}{a}_2\left(m\right)x(t-m)+{\mathrm{\&epsiv;}}_2\left(t\right),t_p<t\&le;t_N---\left(8\right)$

式中M₁和M₂分别为噪声和微地震信号的阶数，N为微地震记录的长度，t_p时刻前后的序列分属不同的平稳过程，即整个序列为非平稳过程，故利用(6)式求得的AR-AIC值会对应局部最小值，在(6)式的求解中，参照Maeda提出的简便算法，利用下式求取AR-AIC值：

$A\left(k\right)=k\lg \left({\mathrm{\&sigma;}}_1^2\right)+(N-k-1)\lg \left({\mathrm{\&sigma;}}_2^2\right)---\left(9\right)$

根据(9)式的结果找到第一个局部极值点并将该点判定为初至时刻。