CN104459772B - 一种地震数字信号的提频方法及装置 - Google Patents

Abstract

本发明公开了一种数字信号提频处理方法，其包括以下步骤：S101、输入一定时间里采集到的实信号道；S102、对所述实信号道进行希尔伯特变换以获取所述实信号道的瞬时振幅道；S103、基于所述瞬时振幅道对实信号道进行提频和极性变换处理，以获取提频信号道。由于消除了事件信号的极性，并且提高了频率，有效的弱事件信号与无效的干扰信号更容易区分，因此不需要大量强事件便能识别出微弱信号源，这在页岩气压裂微地震监测领域中彰显了环境保护和降低成本的优点。此外，本发明的零极性变换和提频处理的步骤简单、通用性强，一旦常数k₁和k₂给定，可以实现任何信号的倍频和零极性处理。

Description

一种地震数字信号的提频方法及装置

技术领域

本发明涉及数字信号处理技术，具体而言，涉及一种地震数字信号的提频方法及装置，更进一步地说，尤其涉及针对页岩气压裂开采过程中产生的微地震监测数据的分析处理。

背景技术

页岩气是一种重要的非常规天然气资源，其开采主要使用水力压裂法，即将化学物质和大量水、泥沙的混合物，用高压注入地下井，压裂附近的岩石构造，进而收集天然气。而岩石破裂会伴随产生强度较弱的地震波,称为“微地震”。

微地震监测技术就是通过观测、分析生产活动中所产生的微小地震事件来监测生产活动的影响、效果及地下状态的地球物理技术。其基本做法是：通过在井中或地面布置检波器排列接收生产活动所产生或诱导的微小地震事件，并通过对这些事件的反演求取微地震震源位置等参数。在页岩气压裂微地震监测领域，微地震资料的信噪比较低，导致弱事件的识别极其困难，目前无法进行微地震弱事件的震源成像定位。对此，现有技术还没有切实可行的方法。

为了微地震的震源成像定位，在现有技术条件下，可通过延长压裂施工时间，增加压裂液等方式试图获得更多的容易识别的强事件。但这样又将突出经济成本和环境保护问题。

因此，目前在微地震监测中，针对采集的数据资料具有较低的信噪比的情况，迫切需要一种能够准确提取有用弱事件的方法。

发明内容

本发明针对现有在页岩气压裂微地震监测领域中不能准确识别弱事件的不足，提出了一种新的数字信号提频处理方法，该方法中的提频和极性变换处理方法本发明称为零极性变换。

本发明方法包括以下步骤：

S101、输入一定时间里采集到的实信号道；

S102、对所述实信号道进行希尔伯特变换以获取所述实信号道的瞬时振幅道；

S103、基于所述瞬时振幅道对实信号道进行提频和极性变换处理，以获取提频信号道。

根据本发明的一个实施例，在所述步骤S103之后，继续按照以下步骤对所述提频信号道进一步优化处理：

S104、对所述提频信号道进行希尔伯特变换以获取对应于所述提频信号道的瞬时余弦相位函数道；

S105、将所述实信号道的瞬时振幅道和所述提频信号道的瞬时余弦相位函数道进行重构以优化所述提频信号道。

根据本发明的一个实施例，在步骤S105中，按照以下公式来进行重构：

z(t)＝cosξ(t)·a(t)

其中，z(t)表示优化的零极性的提频信号道，cosξ(t)表示提频信号道的瞬时余弦相位函数道，ɑ(t)表示所述实信号道的瞬时振幅道。

根据本发明的一个实施例，按照以下公式进行提频和极性变换处理：

y(t)＝k₁·|x(t)|-k₂·a(t)

其中，y(t)表示提频信号道，x(t)表示实信号道，a(t)表示所述实信号道的瞬时振幅道，k₁和k₂分别为常数。

根据本发明的一个实施例，常数k₁与k₂的比值在1.2至2.0范围内，处理后的提频信号道的频率与原实信号道频率是倍频关系。

根据本发明的一个实施例，常数k₁和k₂的取值优选分别为4和π。

根据本发明的一个方面，还提供了一种数字信号提频处理装置，该装置包括以下模块：

输入模块，其用于输入一定时间里采集到的实信号道；

第一变换模块，其用于对所述实信号道进行希尔伯特变换以获取所述实信号道的瞬时振幅道；

提频和极性变换处理模块，其用于基于所述瞬时振幅道对实信号道进行提频和极性变换处理，以获取提频信号道。

根据本发明的一个实施例，该装置还包括以下模块以对所述提频信号道进一步优化处理：

第二变换模块，其用于对所述提频信号道进行希尔伯特变换以获取对应于所述提频信号道的瞬时余弦相位函数道；

重构模块，其用于将所述实信号道的瞬时振幅道和所述提频信号道的瞬时余弦相位函数道进行重构以优化所述提频信号道。

根据本发明的一个实施例，在所述重构模块中，按照以下公式来进行重构：

z(t)＝cosξ(t)·a(t)

其中，z(t)表示优化的零极性的提频信号道，cosξ(t)表示提频信号道的瞬时余弦相位函数道，ɑ(t)表示所述实信号道的瞬时振幅道。

根据本发明的一个实施例，在所述提频和极性变换处理模块中，按照以下公式进行提频和极性变换处理：

y(t)＝k₁·|x(t)|-k₂·a(t)

其中，y(t)表示零极性的提频信号道，x(t)表示实信号道，a(t)表示所述实信号道的瞬时振幅道，k₁和k₂分别为常数。

本发明带来了以下有益效果：

(1)由于消除了事件信号的极性，并且提高了频率，有效的弱事件信号与无效的干扰信号更容易区分，因此不需要大量强事件便能识别出微弱信号源，这在页岩气压裂微地震监测领域中彰显了环境保护和降低成本的优点；

(2)本发明的零极性变换和提频处理的步骤简单、通用性强，一旦常数k₁和k₂给定，可以实现任何信号的倍频和零极性处理；

(3)由于本发明的公式简单，因此采用计算机实现的自动化程度高。

本发明的其它特征和优点将在随后的说明书中阐述，并且部分地从说明书中变得显而易见，或者通过实施本发明而了解。本发明的目的和其他优点可通过在说明书、权利要求书以及附图中所特别指出的结构来实现和获得。

附图说明

图1是模拟微地震资料的理论模型道集及其叠加道；

图2是在图1所示的模型上增加了随机噪声后的理论模型道集及其叠加道；

图3是根据本发明的一个实施例进行零极性变换提频处理的方法步骤流程图；

图4a－4f分别对应于按图3所示的步骤执行信号变换后的结果；

图5a显示了含有一个正极性子波和一个负极性子波的信号道；

图5b显示了对图5a所示的信号道作零极性变换后的无极性子波信号道；

图6a显示了含有一个主频为30Hz子波的信号道；

图6b显示了对应于图6a的子波的频谱图；

图7a显示了对图6a的信号道进行零极性变换后的信号道；

图7b显示了对应于图7a的信号道的频谱图；

图8a显示了含有不同频率子波的信号道；

图8b显示了对图8a的各子波进行零极性变换后的信号道；

图9为根据本发明的实施例对图1所示的模型进行零极性变换的结果及其叠加道；

图10为根据本发明的实施例对图2所示的模型进行零极性变换的结果及其叠加道；

图11a－11c显示了现有技术中页岩气压裂施工某工区实际的微地震强事件识别图；

图12a－12c显示了根据本发明对图11进行零极性变换后的微地震强事件识别图；

图13a－13c显示了现有技术中页岩气压裂施工某工区实际的微地震弱事件识别图；

图14a－14c显示了根据本发明对图13进行零极性变换后的微地震弱事件识别图；

图15是压裂事件的发震时间与震源垂直深度关系对应图；

图16是强压裂事件的震源定位3D图；

图17是页岩气压裂施工井的测井图。

具体实施方式

以下将结合附图来详细说明本发明的实施方式，借此对本发明如何应用技术手段来解决技术问题，并达成技术效果的实现过程能充分理解并据以实施。需要说明的是，只要不构成冲突，本发明各实施例以及各实施例中的各个特征可以相互结合，所形成的技术方案均在本发明的保护范围之内。

另外，附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。

以下就页岩气压裂微地震监测领域中的微地震事件来说明本发明的原理，但是，本发明并不限于此。因为，对于数字信号处理技术领域的普通技术人员而言，看到本发明之后可以想到其可应用于任何数字信号处理中对于有用信号的处理上。

实施例一

图1是模拟微地震资料的理论模型道集及其叠加道。在该图中，0.2秒处有一组扭曲的40Hz规则干扰同相轴；在0.6秒处有另一组具有相同扭曲度的20Hz规则干扰同相轴；在0.4秒处则有一组水平的极性不统一的30Hz事件同相轴。通过图中的叠加道显示可以看出，在没有随机噪声的背景下，两组规则干扰同相轴均能叠加成像，而事件同相轴却不能叠加成像。事件同相轴不能叠加成像的原因主要是极性不统一，在叠加的过程中出现了抵消的现象。

图2是在图1基础上增加随机噪声后的理论模型道集及其叠加道。通过该图中的叠加道显示可看出，仅在0.6秒处的20Hz低频规则干扰同相轴能够叠加成像。很显然，这将造成事件识别的误判。但是，注意到在0.2秒处的40Hz高频规则干扰同相轴并没有叠加成像，这是本发明人所希望的。该现象可根据菲涅耳带原理加以解释，具有相同扭曲度的两组规则干扰同相轴，低频同相轴更容易叠加成像。

根据上述发现，可以得到以下启示：如果将事件同相轴的极性统一，并且提高含噪信号道的频率，就能够获得事件同相轴的叠加成像，减少干扰同相轴的叠加成像，从而降低事件识别的误判率。

如图3所示，其中显示了根据本发明的一个实施例的方法流程图。在该图中，如步骤S101所示，输入一定时间里采集到的实信号道。

然后，在步骤S102中，对实信号道进行希尔伯特变换，以获取实信号道的瞬时振幅道。

希尔伯特变换(HT)是信号分析中的重要工具。假定一个连续的时间信号为x(t)，其希尔伯特变换为h(t)，则希尔伯特变换表达式为：

瞬时振幅表达式为：

瞬时相位表达式为：

瞬时余弦相位函数为：

故：x(t)＝cosθ(t)·a(t) (5)

可以看到，x(t)可分解为瞬时余弦相位函数cosθ(t)和瞬时振幅a(t)。

图4b显示了图4a所示的信号道希尔伯特变换后得到的瞬时振幅。

接下来，在步骤S103中，基于瞬时振幅道对实信号道x(t)进行提频和极性变换处理，以获取零极性的提频信号道。根据本发明的一个实施例，可以基于实信号道x(t)的瞬时振幅道a(t)对实信号道进行提频和极性变换处理。更具体地说，可以按照以下公式进行提频和极性变换处理：

y(t)＝k₁·|x(t)|-k₂·a(t) (6)

其中，y(t)表示零极性的提频信号道，x(t)表示实信号道，a(t)表示实信号道x(t)的瞬时振幅道，k₁和k₂分别为常数。如上处理后的提频信号道的频率是原实信号道的频率的二倍。如图4c和4d分别显示了经过上述变换后的信号道。在此实施例中，常数k₁和k₂的取值优选分别为4和π。这里，常数的取值并不限于此具体值。事实上，常数k₁和k₂之间的比值可以在1.2～2.0之间。

在本发明的一个实施例中，步骤S103之后，还可以继续对零极性的提频信号道进一步优化处理。例如，在步骤S104中，对零极性的提频信号道进行希尔伯特变换以获取对应于零极性的提频信号道的瞬时振幅道。在步骤S105中，将实信号道的瞬时振幅道和零极性的提频信号道的瞬时余弦相位函数道进行重构以优化零极性的提频信号道。

具体地说，在步骤S105中，按照以下公式来进行重构，

z(t)＝cosξ(t)·a(t) (7)

其中，z(t)表示优化的零极性的提频信号道，cosξ(t)表示提频信号道y(t)的瞬时余弦相位函数道，ɑ(t)表示实信号道x(t)的瞬时振幅道。变换的过程可参见图4e和4f所示。

通常，在步骤S102之前，还会对实信号道x(t)进行动校正等预处理以校平曲线。NMO动校正的目的是消除同一地震道上的各子波到达地面的时差，从而校平共深度点震源波时距曲线的轨迹。这样，可以达到利用水平叠加技术压制干扰的目的。

通过上述处理后，无论信号的极性是正极性还是负极性均变换为无极性，信号的频率提高了一倍，且其物理位置保持不变。通过理论模型及实际微地震数据试验，本方法提频效果明显，具有较强的针对性。

如图5a所示，其显示了含有一个正极性子波和一个负极性子波的信号道，而图5b显示了零极性变换后的无极性子波信号道。图中可见，零极性变换没有改变子波的物理位置。

图6a是含有一个主频为30Hz子波的信号道，图6b是其频谱图。图7a是对图6a进行零极性变换后的信号道，图7b是其对应的频谱图，频谱显示子波主频为60Hz。按照本发明的方法进行零极性变换后将信号频率提高了一倍。

正如图7b所示，零极性变换同样受采样定理的最高频率f_max限制，当原信号频率f＞0.5f_max时,零极性变换后将出现假频现象。但这种假频现象并不影响实质。

图8a是含有不同频率子波的信号道，其中在0.3秒处有一个20Hz子波，在0.7秒处有一个40Hz子波。图8b是对图8a的零极性变换处理的信号道，其中在0.3秒处的子波频率变换为40Hz，而在0.7秒处的子波频率变换为80Hz。说明零极性变换是对整个信号道进行提频，信号道中的各个子波在各自频率基础上实现频率翻倍。

图9是对图1中的理论模型道集进行零极性变换的结果及其叠加道。其中，在0.2秒处的一组扭曲的规则干扰同相轴的频率由40Hz变换为80Hz；在0.6秒处的另一组具有相同扭曲度的规则干扰同相轴的频率由20Hz变换为40Hz；在0.4秒处的一组水平的极性不统一的事件同相轴的频率由30Hz变换为60Hz，且统一了事件同相轴的极性。图9中的叠加道显示，在没有随机噪声的背景下，三组同相轴均能叠加成像。

图10是在图2中的含噪理论模型道集基础上进行零极性变换的结果及其叠加道。其中，叠加道显示仅在0.4秒处的一组水平的事件同相轴能够叠加成像。图10和图2的对比结果表明，在低信噪比条件下，高频率规则干扰同相轴是不容易叠加成像的，如果在图10的基础上进行自动识别或人工识别，将大大降低事件识别的误判率。

图11a－11c是页岩气压裂施工某工区实际的微地震强事件识别图。其中，图11a为微地震道集的扫描叠加能量团速度谱。图11b是未进行正常时差校正(NMO)的微地震超道集，由于事件极性的不统一造成了部分极性抵消现象。图11c为图11b进行NMO后的水平叠加道。由于该事件是强事件，在高信噪比条件下，对图11中的强事件识别是非常容易的。

图12a－12c是对图11进行零极性变换后的微地震强事件识别图。其中，图12a为零极性变换后的微地震道集的扫描叠加能量团速度谱。其中，强事件能量团变得更加清晰。图12b为零极性变换后未作NMO的微地震超道集，极性抵消现象已消失。图12c为图12b经NMO后的水平叠加道，叠加成像结果的物理位置没有发生变化。

图13a－13c是页岩气压裂施工某工区实际的微地震弱事件识别图。其中，图13a为微地震道集的扫描叠加能量团速度谱。图13b为经过NMO后的微地震超道集，由于事件极性的不统一造成极性抵消现象，从而导致事件水平同相轴的不可见。图13c为图13b的水平叠加道，同理，因极性抵消现象，使该弱事件不能叠加成像。很显然，在低信噪比条件下，对图13中的弱事件识别是非常困难的。

图14a－14c是对图13进行零极性变换后的微地震弱事件识别图。其中，图14a为零极性变换后的微地震道集的扫描叠加能量团速度谱。在26.1秒处出现了一个清晰的能量团。图14b为经NMO和零极性变换后的微地震超道集，在相同时间处可见事件水平同相轴。图14c为图14b的水平叠加道，叠加成像结果进一步证实了在26.1秒处发生了一个压裂弱事件。

理论模型及实际微地震数据试验结果表明：零极性变换统一了事件同相轴的极性，消除了极性抵消现象，使事件同相轴能够获得正常的叠加成像；零极性变换将含噪信号道的频率提高了一倍，在高频率状态下，随机噪声或干扰同相轴均难以成像；零极性变换保持了信号的物理位置不变，从而保证了震源反演结果的正确性。

实施例二

根据本发明的一个方面，上述方法可以实现于计算机设备中。与计算机设备及其它外围电路构成一种数字信号处理装置。该装置包括以下模块：

输入模块，其用于输入一定时间里采集到的实信号道；

第一变换模块，其用于对所述实信号道进行希尔伯特变换以获取所述实信号道的瞬时振幅道；

提频和极性变换处理模块，其用于基于所述瞬时振幅道对实信号道进行提频和极性变换处理，以获取零极性的提频信号道。

优选地，该装置还包括以下模块以对所述零极性的提频信号道进一步优化处理：

第二变换模块，其用于对所述零极性的提频信号道进行希尔伯特变换以获取对应于所述零极性的提频信号道的瞬时余弦相位函数道；

重构模块，其用于将所述实信号道的瞬时振幅道和所述零极性的提频信号道的瞬时余弦相位函数道进行重构以优化所述零极性的提频信号道。

在重构模块中，按照以下公式来进行重构：

z(t)＝cosξ(t)·a(t)

其中，z(t)表示优化的零极性的提频信号道，cosξ(t)表示提频信号道的瞬时余弦相位函数道，ɑ(t)表示所述实信号道的瞬时振幅道。

在提频和极性变换处理模块中，基于实信号道的瞬时振幅道对实信号道进行提频和极性变换处理。具体说，按照以下公式进行提频和极性变换处理：

y(t)＝k₁·|x(t)|-k₂·a(t) (6)

其中，y(t)表示零极性的提频信号道，x(t)表示实信号道，a(t)表示实信号道的瞬时振幅道，k₁和k₂分别为常数。

实施例三

使用本发明的方法对某工区页岩气井压裂监测数据进行批量处理，共获得879个压裂事件及震源定位，其中，强压裂事件127个，其余均为弱压裂事件。图15为全部879个压裂事件的发震时间与震源垂直深度关系对应图。图16为127个强压裂事件的震源定位3D图。

图15表明，在垂直深度2132m(测量深度2141m)附近有一个低应力水平薄层，弱压裂事件集中发生在该层位中，这一结论得到该井测井图(图17)的验证(图中的深度为测量深度)。

虽然本发明所揭露的实施方式如上，但所述的内容只是为了便于理解本发明而采用的实施方式，并非用以限定本发明。任何本发明所属技术领域内的技术人员，在不脱离本发明所揭露的精神和范围的前提下，可以在实施的形式上及细节上作任何的修改与变化，但本发明的专利保护范围，仍须以所附的权利要求书所界定的范围为准。

Claims (8)

1.一种地震数字信号提频处理方法，其特征在于，所述方法包括以下步骤：

S101、输入一定时间里采集到的实信号道；

S102、对所述实信号道进行希尔伯特变换以获取所述实信号道的瞬时振幅道；

S103、基于所述瞬时振幅道对实信号道进行提频和极性变换处理，以获取提频信号道；

S104、对所述提频信号道进行希尔伯特变换以获取对应于所述提频信号道的瞬时余弦相位函数道；

S105、将所述实信号道的瞬时振幅道和所述提频信号道的瞬时余弦相位函数道进行重构以优化所述提频信号道。

2.如权利要求1所述的方法，其特征在于，在步骤S105中，按照以下公式来进行重构：

z(t)＝cosξ(t)·a(t)

其中，z(t)表示优化的提频信号道，cosξ(t)表示提频信号道的瞬时余弦相位函数道，ɑ(t)表示所述实信号道的瞬时振幅道。

3.如权利要求1或2所述的方法，其特征在于，在所述步骤S103中，按照以下公式进行提频和极性变换处理：

y(t)＝k₁·|x(t)|-k₂·a(t)

其中，y(t)表示提频信号道，x(t)表示实信号道，a(t)表示所述实信号道的瞬时振幅道，k₁和k₂分别为常数。

4.如权利要求3所述的方法，其特征在于，常数k₁与k₂的比值在1.2至2.0范围内，处理后的提频信号道的频率与原实信号道频率是倍频关系。

5.如权利要求4所述的方法，其特征在于，常数k₁和k₂的取值分别为4和π。

6.一种地震数字信号提频处理装置，其特征在于，所述装置包括以下模块：

输入模块，其用于输入一定时间里采集到的实信号道；

第一变换模块，其用于对所述实信号道进行希尔伯特变换以获取所述实信号道的瞬时振幅道；

提频和极性变换处理模块，其用于基于所述瞬时振幅道对实信号道进行提频和极性变换处理，以获取提频信号道；

第二变换模块，其用于对所述提频信号道进行希尔伯特变换以获取对应于所述提频信号道的瞬时余弦相位函数道；

重构模块，其用于将所述实信号道的瞬时振幅道和所述提频信号道的瞬时余弦相位函数道进行重构以优化所述提频信号道。

7.如权利要求6所述的装置，其特征在于，在所述重构模块中，按照以下公式来进行重构：

z(t)＝cosξ(t)·a(t)

其中，z(t)表示优化的提频信号道，cosξ(t)表示提频信号道的瞬时余弦相位函数道，ɑ(t)表示所述实信号道的瞬时振幅道。

8.如权利要求6所述的装置，其特征在于，在所述提频和极性变换处理模块中，按照以下公式进行提频和极性变换处理：

y(t)＝k₁·|x(t)|-k₂·a(t)

其中，y(t)表示提频信号道，x(t)表示实信号道，a(t)表示所述实信号道的瞬时振幅道，k₁和k₂分别为常数。