基于三分量偏振梯度的微地震初至识别方法及装置

摘要

公开了基于三分量偏振梯度的微地震初至识别方法及装置。该方法可以包括：基于原始三分量数据构建瞬时三分量协方差矩阵CXYZ(t)；基于所述瞬时三分量协方差矩阵CXYZ(t)构建瞬时偏振度函数PF(t)；基于所述瞬时偏振度函数PF(t)构建目标梯度函数Δ(t)；以及通过将所述目标梯度函数Δ(t)与给定门槛值K进行比较来识别微地震初至时间。本发明直接从原始微地震三分量数据出发，通过构建瞬时偏振函数的目标梯度函数，实现了微地震初至的快速、自动识别。

权利要求

1.一种基于三分量偏振梯度的微地震初至识别方法，所述方法包括：基于原始三分量数据X0(t)、Y0(t)、Z0(t)构建瞬时三分量协方差矩阵CXYZ(t)；

基于所述瞬时三分量协方差矩阵CXYZ(t)构建瞬时偏振度函数PF(t)；

基于所述瞬时偏振度函数PF(t)构建目标梯度函数Δ(t)；以及通过将所述目标梯度函数Δ(t)与给定门槛值K进行比较来识别微地震初至时间；

其中，基于所述瞬时三分量协方差矩阵CXYZ(t)构建瞬时偏振度函数PF(t)包括：对所述瞬时三分量协方差矩阵CXYZ(t)进行奇异值分解以计算特征值：其中，λ1、λ2、λ3为瞬时三分量协方差矩阵的特征值且λ1＞λ2＞λ3，对应的特征向量为V1、V2、V3，利用计算出的特征值构建瞬时偏振度函数PF(t)：其中，γ为设定的指数，γ∈(0,1)；

基于所述瞬时偏振度函数PF(t)构建目标梯度函数Δ(t)包括：对瞬时偏振函数PF(t)求偏导计算标梯度函数Δ(t)：

2.根据权利要求1所述的基于三分量偏振梯度的微地震初至识别方法，其中，基于原始三分量数据X0(t)、Y0(t)、Z0(t)构建瞬时三分量协方差矩阵CXYZ(t)包括：以当前时间t0为中心，利用hanning窗截取原始三分量数据X0(t)、Y0(t)、Z0(t)，获得新的三分量数据X(t)、Y(t)、Z(t)，利用所获得的新的三分量数据构建协方差矩阵CXYZ(t)：其中，

N为hanning窗半长度。

3.根据权利要求1所述的基于三分量偏振梯度的微地震初至识别方法，其中，通过将所述目标梯度函数Δ(t)与给定门槛值K进行比较来识别微地震初至时间包括：如果出现所述目标梯度函数Δ(t)大于所述给定门槛值K，则寻找所述目标梯度函数Δ(t)的区域峰值，将所述区域峰值对应的时间位置确定为微地震初至时间。

4.一种基于三分量偏振梯度的微地震初至识别装置，所述装置包括：协方差矩阵构建单元，用于基于原始三分量数据X0(t)、Y0(t)、Z0(t)构建瞬时三分量协方差矩阵CXYZ(t)；

偏振度函数构建单元，用于基于所述瞬时三分量协方差矩阵CXYZ(t)构建瞬时偏振度函数PF(t)；

目标梯度函数计算单元，用于基于所述瞬时偏振度函数PF(t)构建目标梯度函数Δ(t)；以及微地震初至时间识别单元，用于通过将所述目标梯度函数Δ(t)与给定门槛值K进行比较来识别微地震初至时间；

其中，基于所述瞬时偏振度函数构建目标梯度函数包括：对所述瞬时三分量协方差矩阵CXYZ(t)进行奇异值分解以计算特征值：其中，λ1、λ2、λ3为瞬时三分量协方差矩阵的特征值且λ1＞λ2＞λ3，对应的特征向量为V1、V2、V3，利用计算出的特征值构建瞬时偏振度函数PF(t)：其中，γ为设定的指数，γ∈(0,1)；

其中，基于所述瞬时偏振度函数构建目标梯度函数包括：对瞬时偏振函数PF(t)求偏导计算标梯度函数Δ(t)：

5.根据权利要求4所述的基于三分量偏振梯度的微地震初至识别装置，其中，构建瞬时三分量协方差矩阵包括：以当前时间t0为中心，利用hanning窗截取原始三分量数据X0(t)、Y0(t)、Z0(t)，获得新的三分量数据X(t)、Y(t)、Z(t)，利用所获得的新的三分量数据构建协方差矩阵CXYZ(t)：其中，

N为hanning窗半长度。

6.根据权利要求4所述的基于三分量偏振梯度的微地震初至识别装置，其中，通过将所述目标梯度函数与给定门槛值进行比较来识别微地震初至时间包括：如果出现所述目标梯度函数Δ(t)大于所述给定门槛值K，则寻找所述目标梯度函数Δ(t)的区域峰值，将所述区域峰值对应的时间位置确定为微地震初至时间。

说明书

基于三分量偏振梯度的微地震初至识别方法及装置

技术领域

[0001] 本发明涉及微地震三分量信号处理领域，更具体地，涉及一种基于三分量偏振梯度的微地震初至识别方法和一种基于三分量偏振梯度的微地震初至识别装置。

背景技术

[0002] 精确拾取有效的事件初至时间，对于微地震事件的震源定位以及解释震源机制具有重要意义。人工查看震动图拾取初至是非常耗时的，并且由于人工操作带有主观性，即使是专业人员，在不同时间拾取的初至也会有所不同。因此，研究微地震三分量初至自动识别并且提高拾取时间准确性，成为微地震处理关键技术之一。

[0003] 由于微地震事件信号和背景噪声在振幅特征和频率组成上有很大的不同，这就为在三分量震动图上识别事件初至时间提供了可能。微地震事件的特征是明显的初至点、高频以及频率随着时间快速衰减，而背景噪音信号可以根据它们的低振幅、低频来进行区分。

[0004] 发明人发现，由于微地震三分量数据振幅和极性具有杂乱无序的特点，而在现有技术中，只有在微地震三分量数据前期大量处理的基础上才能提高微地震事件初至识别能力。因此，有必要开发一种能够直接从原始三分量数据快速、准确地识别微地震事件初至的方法。

[0005] 公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域技术人员所公知的现有技术。

发明内容

[0006] 本发明提出了一种能够直接从原始三分量数据快速、准确地识别微地震事件初至的方法，本发明还公开了相应的装置。

[0007] 根据本发明的一方面，提出了一种基于三分量偏振梯度的微地震初至识别方法，该方法包括：基于原始三分量数据构建瞬时三分量协方差矩阵；基于所述瞬时三分量协方差矩阵构建瞬时偏振度函数；基于所述瞬时偏振度函数构建目标梯度函数；以及通过将所述目标梯度函数与给定门槛值进行比较来识别微地震初至时间。

[0008] 根据本发明的另一方面，提出了基于三分量偏振梯度的微地震初至识别装置，该装置包括：协方差矩阵构建单元，用于基于原始三分量数据构建瞬时三分量协方差矩阵；偏振度函数构建单元，用于基于所述瞬时三分量协方差矩阵构建瞬时偏振度函数；目标梯度函数计算单元，用于基于所述瞬时偏振度函数构建目标梯度函数；以及微地震初至时间识别单元，用于通过将所述目标梯度函数与给定门槛值进行比较来识别微地震初至时间。

[0009] 现有技术通常是基于对地震三分量数据前期大量处理来提高微地震事件初至识别能力。本发明直接从原始微地震三分量数据出发，通过构建瞬时偏振函数的目标梯度函数，实现了微地震初至的快速、自动识别。

附图说明

[0010] 通过结合附图对本发明示例性实施方式进行更详细的描述，本发明的上述以及其它目的、特征和优势将变得更加明显，其中，在本发明示例性实施方式中，相同的参考标号通常代表相同部件。

[0011] 图1示出了根据本发明的基于三分量偏振梯度的微地震初至识别方法的流程图；

[0012] 图2A-2C分别示出了微地震模型三分量单道数据的X分量、Y分量、Z分量；

[0013] 图3示出了微地震模型三分量偏振度；

[0014] 图4示出了对微地震模型数据应用本发明得到的偏振度梯度；

[0015] 图5A-5C分别示出了微地震实际资料三分量单道数据的X分量、Y分量、Z分量；

[0016] 图6示出了微地震实际资料三分量偏振度；

[0017] 图7示出了对微地震实际资料数据应用本发明得到的偏振度梯度。

具体实施方式

[0018] 下面将参照附图更详细地描述本发明的优选实施方式。虽然附图中显示了本发明的优选实施方式，然而应该理解，可以以各种形式实现本发明而不应被这里阐述的实施方式所限制。相反，提供这些实施方式是为了使本发明更加透彻和完整，并且能够将本发明的范围完整地传达给本领域的技术人员。

[0019] 实施例1

[0020] 图1示出了根据本发明的一个实施例的基于三分量偏振梯度的微地震初至识别方法的流程图，该方法包括：

[0021] 步骤101，基于原始三分量数据X0(t)、Y0(t)、Z0(t)构建瞬时三分量协方差矩阵CXYZ(t)；

[0022] 步骤102，基于所述瞬时三分量协方差矩阵CXYZ(t)构建瞬时偏振度函数PF(t)；

[0023] 步骤103，基于所述瞬时偏振度函数PF(t)构建目标梯度函数Δ(t)；以及[0024] 步骤104，通过将所述目标梯度函数Δ(t)与给定门槛值K进行比较来识别微地震初至时间。

[0025] 微地震信号在三分量X、Y、Z矢量数据上存在两两耦合关系，也即偏振关系，而有效信号之外的随机信号无此特征。本发明利用这一点，直接从原始微地震三分量数据出发，通过构建瞬时偏振函数的目标梯度函数实现了微地震初至的快速、自动识别。

[0026] 在一个示例中，在步骤101中，基于原始三分量数据X0(t)、Y0(t)、Z0(t)构 建瞬时三分量协方差矩阵CXYZ(t)可以包括：以当前时间t0为中心，利用hanning窗截取原始三分量数据X0(t)、Y0(t)、Z0(t)，获得新的三分量数据X(t)、Y(t)、Z(t)，利用所获得的新的三分量数据构建协方差矩阵CXYZ(t)。

[0027] 具体地，hanning窗可以是：

[0028]

[0029] 其中，N为hanning窗半长度，t0为当前时间。

[0030] 将hanning窗与原始三分量数据X0(t)、Y0(t)、Z0(t)相乘，获得新的三分量数据X(t)、Y(t)、Z(t)：

[0031] X(t)＝X0(t)h(t)、Y(t)＝Y0(t)h(t)、Z(t)＝Z0(t)h(t)  (2)。

[0032] 利用所获得的新的三分量数据构建协方差矩阵CXYZ(t)：

[0033]

[0034] 其中，

[0035]

[0036]

[0037]

[0038]

[0039]

[0040]

[0041]

[0042]

[0043] 在一个示例中，在步骤102中，基于所述瞬时三分量协方差矩阵CXYZ(t)构建瞬时偏振度函数PF(t)可以包括：对所述瞬时三分量协方差矩阵CXYZ(t)进行奇异值分解以计算特征值，利用计算出的特征值构建瞬时偏振度函数PF(t)。

[0044] 具体地，对以当前时间t0为中心的三分量协方差矩阵CXYZ(t)进行奇异值分解：

[0045]

[0046] 其中，λ1、λ2、λ3为瞬时三分量协方差矩阵的特征值且λ1＞λ2＞λ3，对应的特征向量为V1、V2、V3。

[0047] 利用计算出的特征值构建瞬时偏振度函数PF(t)：

[0048]

[0049] 其中，指数γ∈(0,1)，其具体数值根据实际数据测试而得。

[0050] 在一个示例中，在步骤103中，基于所述瞬时偏振度函数PF(t)构建目标梯度函数Δ(t)可以包括：对瞬时偏振度函数PF(t)求偏导以计算目标梯度函数Δ(t)。

[0051] 具体地，通过直接对瞬时偏振函数PF(t)求偏导计算标梯度函数Δ(t)：

[0052]

[0053] 根据该步骤对所有样点时间计算对应的梯度函数Δ(t)，若前后数据不足则补零值。

[0054] 在一个示例中，在步骤104中，通过将所述目标梯度函数Δ(t)与给定门槛值K进行比较来识别微地震初至时间可以包括：如果出现所述目标梯度函数Δ(t)大于所述给定门槛值K，则寻找所述目标梯度函数Δ(t)的区域峰值，将所述区域峰值对应的时间位置确定为微地震初至时间。

[0055] 具体地，如果出现目标函数Δ(t)大于给定门槛值K，说明存在微地震事件，则在大于K的范围内寻找目标函数Δ(t)的区域峰值对应的时间，实现微地震事件初至时间识别。如果没有出现目标函数Δ(t)大于给定门槛值K，说明不存在微地震事件。

[0056] 实施例2

[0057] 本发明还公开了基于三分量偏振梯度的微地震初至识别装置。在本实施例中，该装置包括：协方差矩阵构建单元，用于基于原始三分量数据X0(t)、Y0(t)、Z0(t)构建瞬时三分量协方差矩阵CXYZ(t)；偏振度函数构建单元，用于基于所述瞬时三分量协方差矩阵CXYZ(t)构建瞬时偏振度函数PF(t)；目标梯度函数计算单元，用于基于所述瞬时偏振度函数PF(t)构建目标梯度函数Δ(t)；以及微地震初至时间识别单元，用于通过将所述目标梯度函数Δ(t)与给定门槛值K进行比较来识别微地震初至时间。

[0058] 在一个示例中，基于原始三分量数据X0(t)、Y0(t)、Z0(t)构建瞬时三分量协方差矩阵CXYZ(t)可以包括：以当前时间t0为中心，利用hanning窗截取原始三分量数据X0(t)、Y0(t)、Z0(t)，获得新的三分量数据X(t)、Y(t)、Z(t)，利用所获得的新的三分量数据构建协方差矩阵CXYZ(t)。

[0059] 在一个示例中，基于所述瞬时三分量协方差矩阵CXYZ(t)构建瞬时偏振度函数PF(t)可以包括：对所述瞬时三分量协方差矩阵CXYZ(t)进行奇异值分解以计算特征值，利用计算出的特征值构建瞬时偏振度函数PF(t)。

[0060] 在一个示例中，基于所述瞬时偏振度函数PF(t)构建目标梯度函数Δ(t)可以包括：对瞬时偏振度函数PF(t)求偏导以计算目标梯度函数Δ(t)。

[0061] 在一个示例中，通过将所述目标梯度函数Δ(t)与给定门槛值K进行比较来识别微地震初至时间可以包括：如果出现所述目标梯度函数Δ(t)大于所述给定门槛值K，则寻找所述目标梯度函数Δ(t)的区域峰值，将所述区域峰值对应的时间位 置确定为微地震初至时间。

[0062] 应用示例

[0063] 为便于理解本发明实施例的方案及其效果，以下给出一个具体应用示例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体细节并非意在以任何方式限制本发明。

[0064] 下面通过实例来验证如图1所示的本发明的微地震初至识别方法的效果和准确性。

[0065] 首先对微地震模型数据应用本发明进行测试。图2A-2C分别为含一定噪音的微地震模型三分量单道数据的X分量、Y分量、Z分量。已知P波初至时间为100ms，S波初至时间为250ms。

[0066] 根据如图1所示的本发明的微地震初至识别方法，设置hanning窗长度为21，根据公式(3)对模型三分量数据构建瞬时协方差矩阵CXYZ(t)，并根据公式(4)进行奇异值分解；然后根据公式(5)构建瞬时三分量偏振度函数PF(t)，结果如图3所示；最后，通过公式(6)获得最终目标梯度函数Δ(t)，如图4所示。对于给定门槛值0.3，搜索出两个微地震事件，寻找出区域峰值对应的时间位置，分别为99ms、249ms，即为微地震事件初至识别的结果。

[0067] 应用根据本发明的基于三分量偏振梯度的微地震初至识别方法识别的初至时间为99ms、249ms，与模型的理论初至时间100ms、250ms相比，可知本发明的初至识别准确性比较高，有一定抗噪音能力。

[0068] 下面对微地震实际资料数据应用本发明进行测试。图5A-5C为微地震实际资料三分量单道数据的X分量、Y分量、Z分量。其中，Z分量振幅较大、信噪比最高，其次是X分量，数据最差的是Y分量。这些差异跟三分量观测系统与野外复杂环境有关。

[0069] 根据如图1所示的本发明的微地震初至识别方法，设置hanning窗长度为81，根据公式(3)对微地震实际资料三分量数据构建瞬时协方差矩阵CXYZ(t)， 并根据公式(4)进行奇异值分解；然后根据公式(5)构建瞬时三分量偏振度函数PF(t)，结果如图6所示；最后，通过公式(6)获得实际三分量事件识别目标函数，即目标梯度函数Δ(t)，如图4所示。对于给定门槛值0.15，搜索出一个微地震事件，并识别出初至时间位置为697ms。

[0070] 通过以上两个实例可以看出，根据本发明的基于三分量偏振梯度的微地震初至识别方法能够准确地识别事件及其初至时间。该方法可以直接在原始数据上操作，过程简单，计算和搜索过程自动化处理，并且识别出的事件初至时间准确性较高。

[0071] 本发明可以是系统、方法和/或计算机程序产品。计算机程序产品可以包括计算机可读存储介质，其上载有用于使处理器实现本发明的各个方面的计算机可读程序指令。

[0072] 计算机可读存储介质可以是可以保持和存储由指令执行设备使用的指令的有形设备。计算机可读存储介质例如可以是——但不限于——电存储设备、磁存储设备、光存储设备、电磁存储设备、半导体存储设备或者上述的任意合适的组合。计算机可读存储介质的更具体的例子(非穷举的列表)包括：便携式计算机盘、硬盘、随机存取存储器(RAM)、只读存储器(ROM)、可擦式可编程只读存储器(EPROM或闪存)、静态随机存取存储器(SRAM)、便携式压缩盘只读存储器(CD-ROM)、数字多功能盘(DVD)、记忆棒、软盘、机械编码设备、例如其上存储有指令的打孔卡或凹槽内凸起结构、以及上述的任意合适的组合。这里所使用的计算机可读存储介质不被解释为瞬时信号本身，诸如无线电波或者其他自由传播的电磁波、通过波导或其他传输媒介传播的电磁波(例如，通过光纤电缆的光脉冲)、或者通过电线传输的电信号。

[0073] 这里所描述的计算机可读程序指令可以从计算机可读存储介质下载到各个计算/处理设备，或者通过网络、例如因特网、局域网、广域网和/或无线网下载到外部计算机或外部存储设备。网络可以包括铜传输电缆、光纤传输、无线传 输、路由器、防火墙、交换机、网关计算机和/或边缘服务器。每个计算/处理设备中的网络适配卡或者网络接口从网络接收计算机可读程序指令，并转发该计算机可读程序指令，以供存储在各个计算/处理设备中的计算机可读存储介质中。

[0074] 用于执行本发明操作的计算机程序指令可以是汇编指令、指令集架构(ISA)指令、机器指令、机器相关指令、微代码、固件指令、状态设置数据、或者以一种或多种编程语言的任意组合编写的源代码或目标代码，所述编程语言包括面向对象的编程语言—诸如Smalltalk、C++等，以及常规的过程式编程语言—诸如“C”语言或类似的编程语言。计算机可读程序指令可以完全地在用户计算机上执行、部分地在用户计算机上执行、作为一个独立的软件包执行、部分在用户计算机上部分在远程计算机上执行、或者完全在远程计算机或服务器上执行。在涉及远程计算机的情形中，远程计算机可以通过任意种类的网络—包括局域网(LAN)或广域网(WAN)—连接到用户计算机，或者，可以连接到外部计算机(例如利用因特网服务提供商来通过因特网连接)。在一些实施例中，通过利用计算机可读程序指令的状态信息来个性化定制电子电路，例如可编程逻辑电路、现场可编程门阵列(FPGA)或可编程逻辑阵列(PLA)，该电子电路可以执行计算机可读程序指令，从而实现本发明的各个方面。

[0075] 这里参照根据本发明实施例的方法、装置(系统)和计算机程序产品的流程图和/或框图描述了本发明的各个方面。应当理解，流程图和/或框图的每个方框以及流程图和/或框图中各方框的组合，都可以由计算机可读程序指令实现。

[0076] 这些计算机可读程序指令可以提供给通用计算机、专用计算机或其它可编程数据处理装置的处理器，从而生产出一种机器，使得这些指令在通过计算机或其它可编程数据处理装置的处理器执行时，产生了实现流程图和/或框图中的一个或多个方框中规定的功能/动作的装置。也可以把这些计算机可读程序指令存储在计算机可读存储介质中，这些指令使得计算机、可编程数据处理装置和/或其他设备以特定方式工作，从而，存储有指令的计算机可读介质则包括一个制造品，其包括实现流程图和/或框图中的一个或多个方框中规定的功能/动作的 各个方面的指令。

[0077] 也可以把计算机可读程序指令加载到计算机、其它可编程数据处理装置、或其它设备上，使得在计算机、其它可编程数据处理装置或其它设备上执行一系列操作步骤，以产生计算机实现的过程，从而使得在计算机、其它可编程数据处理装置、或其它设备上执行的指令实现流程图和/或框图中的一个或多个方框中规定的功能/动作。

[0078] 附图中的流程图和框图显示了根据本发明的多个实施例的系统、方法和计算机程序产品的可能实现的体系架构、功能和操作。在这点上，流程图或框图中的每个方框可以代表一个模块、程序段或指令的一部分，所述模块、程序段或指令的一部分包含一个或多个用于实现规定的逻辑功能的可执行指令。在有些作为替换的实现中，方框中所标注的功能也可以以不同于附图中所标注的顺序发生。例如，两个连续的方框实际上可以基本并行地执行，它们有时也可以按相反的顺序执行，这依所涉及的功能而定。也要注意的是，框图和/或流程图中的每个方框、以及框图和/或流程图中的方框的组合，可以用执行规定的功能或动作的专用的基于硬件的系统来实现，或者可以用专用硬件与计算机指令的组合来实现。

[0079] 以上已经描述了本发明的各实施例，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施例。在不偏离所说明的各实施例的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。本文中所用术语的选择，旨在最好地解释各实施例的原理、实际应用或对市场中的技术的改进，或者使本技术领域的其它普通技术人员能理解本文披露的各实施例。