CN114355441A - 一种微震震动波到时拾取方法 - Google Patents
一种微震震动波到时拾取方法 Download PDFInfo
- Publication number
- CN114355441A CN114355441A CN202210011651.5A CN202210011651A CN114355441A CN 114355441 A CN114355441 A CN 114355441A CN 202210011651 A CN202210011651 A CN 202210011651A CN 114355441 A CN114355441 A CN 114355441A
- Authority
- CN
- China
- Prior art keywords
- time window
- microseismic
- length
- time
- picking
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
Images
Landscapes
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明属于微震监测领域,公开了一种微震震动波到时自动拾取方法,该方法主要为一种时窗法和AIC法联合使用的综合法,我们将此综合方法取名为W‑AIC。本发明提供一种微震震动波到时自动拾取系统,使用简便且精确,适用于微震定位。该方法通过微震数据采集系统,将采集到的波形文件传输到计算机内,为了有效去除噪声的影响,采用改进后的FIR滤波进行处理,将处理后的微震信号用时窗能量特征值法进行事件的初至拾取,对拾取并划分范围内的数据使用改进后的W‑AIC算法精确拾取到时,再进行峰值判断,精确得到初至点。该方法提高了震动波到时自动拾取的精度和科学性,大大提高了工作效率。
Description
技术领域
本发明涉及微震监测领域,具体涉及一种微震震动波到时拾取方法。
背景技术
随着信息化、数字化、智能化的发展,微震监测技术在高速公路边坡、地下空间、地铁开发、水电大坝、高层建筑、岩爆、矿山透水、非常规油气压裂开采、煤与瓦斯突出和冲击地压等领域的应用越来越广泛。
如何准确地拾取微震信号,科学准确的判断该信号的到时,从而提高震源定位精度和预测预报岩体破裂的准确性已成为亟待解决的问题。自2000年以来,国内一些大型矿山引进了国外的微震监测产品,并与相关高校和科研机构的合作获得了很多研究成果,为微震监测技术的发展做出了重要贡献,但同时也面临着产品成本较高和售后技术服务短缺的问题,现有的微震监测系统中的到时识别与拾取技术还处在发展阶段,现阶段到时拾取比较显著的问题是精度较低、操作复杂,对操作人员的经验和分析能力要求较高,无法实现微震监测系统的精准性、实时性和科学性。
发明内容
本发明的目的是提供一种可以快速、精准和科学地拾取到时的微震震动波到时自动拾取方法。
为了实现上述目的,本申请实施例采用如下技术方案:
一种微震震动波到时拾取方法,包括:S1:输入及处理数据;利用微震数据采集装置,采集一段震动波,将所采集到的微震事件波形文件传输到计算机内的到时拾取系统中,为了有效压制底噪的影响,确保震动波到时拾取的精准度,采用FIR滤波器进行滤波;S2:自动选取时窗长度L;对经FIR滤波器滤波后的微震信号,自动选取一个合适时窗长度L,并截取在该时窗的长度范围内的微震数据;S3:时窗能量特征值法精确拾取;
在选取的一个合适时窗长度L内,采用改进后时窗能量特征值算法精确拾取到时点;S4:划分精确拾取范围;通过步骤S3时窗能量特征值法拾取到的到时点,为更加精确拾取到时点,在所拾取的到时点周围划分得到一个合适的拾取范围M,该M远远小于本文中的时窗长度L,即M<<L;S5:AIC精确拾取;对划分精确拾取范围M,使用改进后W-AIC算法计算AIC值,然后进行峰值判断,在此范围M的数据内更加精确的拾取出到时点;S6:峰值判断;经峰值判断为是,可得到精确到时点;若拾取到事件结束的点,峰值判断为否,则进行步骤S4、S5,重新划分精确的拾取范围Mi,并在重新划分精确的拾取范围Mi进行AIC精确拾取,直到峰值判断为是。
进一步的,若微震震动波为周期性的波,所述时窗长度L为微震震动波的两个周期长度。
进一步的,若微震震动波为非周期性的波,所述时窗长度L为:当最大振幅位于波谷时,找到其前、后相邻的波谷位置,并把这两个位置间的跨度长度作为时窗长度L;所述时窗长度L计算如下:
a.求取|x(t0)|max
b.如果x(t0)·x(t0-1)<0或者x(t0)=0则a1=a1+1
c.当a1=2时,取x(t1-n)first,min,n=1,2,3,4…
d.当a1=2时,取x(t2-n)first,max,n=1,2,3,4…
式中,|x(t0)|max为最大幅值;a1为幅值异号点次数,t0为最大幅值点;x(t1-n)first,min为次波谷;x(t2-n)first,max为次波峰;L为选取的时窗长度。
进一步的,若微震震动波为非周期性的波,所述时窗长度L为:当最大振幅位于波峰时,找到其前、后相邻的波峰位置,并把这两个位置间的跨度长度作为时窗长度L;所述时窗长度L计算如下:
a.求取|x(t0)|max
b.如果x(t0)·x(t0-1)β<0或者x(t0)=0则a1=a1+1
c.当a1=2时,取x(t1-n)first,min,n=1,2,3,4…
d.当a1=2时,取x(t2-n)first,max,n=1,2,3,4…
式中,|x(t0)|max为最大幅值;a1为幅值异号点次数,t0为最大幅值点;x(t1-n)first,min为次波谷;x(t2-n)first,max为次波峰;L为选取的时窗长度。
具体的,所述用改进后时窗能量特征值算法为:
式中,取一个稳定因子λ,其中稳定因子λ取2或3,A为能量特征值,CF(t)为微震记录幅值;T1为时窗起点;T0为时窗中点;T2为时窗终点;n为单道采样点总数。
具体的,改进后W-AIC算法为:
式中:g为AR过程阶数,i为2个局部统计时段的分界点,var(x1,i)和var(x[i+1,N])为2个局部统计时段的拟合误差,RC为稳定常数。
具体的,所述FIR滤波器加入了不同窗函数。
本发明的有益效果在于,
1、本发明系统安全性高,操作简单,整体成本低;
2、本发明提高了微震信号到时拾取的精度及效率,比传统拾取方法的精确度更高,提高微震监测效果;
3、本发明在整体上提高了工作效率,并提高其经济效益。
附图说明
图1是本发明是一种微震震动波到时自动拾取方法的流程图;
图2是本发明到时拾取后的软件运行的系统界面图;
图3是本发明拾取到时对比。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚明了,下面结合具体实施方式并参照附图,对本发明实施例中的技术方案进行清楚、完整的描述,显然,所描述的实施例是本发明一部分实施例,而不是全部实施例。因此,以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围,而是仅仅表示本发明的选定实施例。基于本发明中的实施例,本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
此外,在以下说明中,省略了对公知结构和技术的描述,以避免不必要地混淆本发明的概念。
实施例1
本发明主要是一种微震震动波到时自动拾取方法,本发明方法在具体的实例中通过各个设备之间信息的相互传递进行收集数据,即通过专业微震数据采集设备采集一段微震震动波,然后将采集的微震事件传递到系统中的控制器内,传递方式有多种,可以用硬盘或U盘及储存卡等多种储存设备进行数据的传递,将传递后的数据经系统控制器进行数据的处理,该系统控制器可用CPU控制器,为确保系统运行效率,建议震动波到时自动拾取系统装配在一台工业计算机中。
结合附图对本发明方法进行再一步的说明,包括以下步骤:
步骤S1:输入及处理数据;
利用微震数据采集装置,采集一段震动波,将所采集到的微震事件波形文件传输到计算机内的到时拾取系统中,为了有效压制底噪的影响,确保震动波到时拾取的精准度,采用FIR滤波器进行滤波;其中,FIR滤波器广泛应用于数字信号处理中,主要功能就是将不感兴趣的信号滤除,留下微震信号。FIR滤波器是全零点结构,系统永远稳定;并且具有线性相位的特征,在有效频率范围内所有信号相位上不失真,同时本发明中的FIR滤波器,为了避免频率泄露,加入了不同的窗函数,可以针对不同的微震波,使用不同的窗函数来减少频率泄露。
步骤S2:自动选取时窗长度L;对经FIR滤波器滤波后的微震信号,自动选取一个合适时窗长度L,并截取在该时窗的长度范围内的微震数据。
具体的,对于周期性的波,时窗长度L恰好等于两倍周期时最为合理,因为时窗前、后能量比特征反映的恰好是波相邻两个完整周期长度的情况。
对于非周期性的微震波,要完全实现时窗长度L的自适应性是很困难的,这是一个复杂的模式识别过程,本发明中,采用了当最大振幅位于波谷时,找到其前、后相邻的波谷位置,并把这两个位置间的跨度L作为时窗的长度;当最大振幅位于波峰时,方法同样。时窗后移范围控制在最大峰值点前面相邻的波谷或波峰位置,这样既避免了能量特征值的双峰现象,又实现了时窗长度L与不同主频微震波的近似匹配。自动选取时窗长度L算法如式:
a.求取|x(t0)|max
b.如果x(t0)·x(t0-1)<0或者x(t0)=0则a1=a1+1
c.当a1=2时,取x(t1-n)first,min,n=1,2,3,4…
d.当a1=2时,取x(t2-n)first,max,n=1,2,3,4…
式中,|x(t0)|max为最大幅值;a1为幅值异号点次数,t0为最大幅值点;x(t1-n)first,min为次波谷;x(t2-n)first,max为次波峰;L为选取的时窗长度。
这样当微震数据庞大时,因短时间内无法对全部数据进行拾取到时,为了提高效率,自动选取微震数据存在的合适时窗L,很大程度上提高了拾取效率。
步骤S3:时窗能量特征值法精确拾取;
在选取的一个合适时窗长度L内,采用改进后时窗能量特征值算法精确拾取到时点;
具体的,时窗的能量特征对地震波初始到时拾取是最有帮助的,其中时窗内后、前能量比的最大值是检测初始到时时刻的一个很好标志,它优于其他能量特征和分维特征,为避免特殊情况,改进后的时窗内后、前能量比特征值算法如式:
式中,取一个稳定因子λ,其中稳定因子取2或者3比较合适,A为能量特征值,CF(t)为微震记录幅值;T1为时窗起点;T0为时窗中点;T2为时窗终点;n为单道采样点总数。
步骤S4:划分精确拾取范围;
通过步骤S3时窗能量特征值法拾取到的到时点,为更加精确拾取到时点,在所拾取的到时点周围划分得到一个合适的拾取范围M,该M远远小于本文中的时窗长度L,即M<<L。
步骤S5:AIC精确拾取;
对划分精确拾取范围M,使用改进后W-AIC算法计算AIC值,然后进行峰值判断,在此范围M的数据内更加精确的拾取出到时点。
具体的,AIC法,拾取利用时间序列分析,试图把信号数据建立成AR模型(自回归模型),写出AR方程(自回归方程)(类似曲线拟合),而AIC最初的目的在于确定AR方程的阶数和系数(定阶方法),来描述AR模型的拟合程度,AIC具体表现为描述拟合程度的矩和参数个数两部分,是一种惩罚方程(当拟合程度过高或者过低,同时参数个数又不合理时,AIC值就会过高,是对模型的一种惩罚),当一段微震波形中存在较高噪声信号时,AIC可以写成两段不同信号的矩取对数乘以信号的长度再加上一个稳定常数的形式。当AIC最小时,说明噪声信号和微震信号得以区分,而此时分界点即到时点,为符合需求,改进后的W-AIC算法可描述为:
式中:g为AR过程阶数,i为2个局部统计时段的分界点,var(x1,i)和var(x[i+1,N])为2个局部统计时段的拟合误差,Rc为稳定常数。
步骤S6:峰值判断;
经峰值判断为是,可得到精确到时点;若拾取到事件结束的点,峰值判断为否,则进行步骤S4、S5,重新划分精确的拾取范围Mi,并在重新划分精确的拾取范围Mi进行AIC精确拾取,直到峰值判断为是。
进行峰值判断时,为避免拾取到事件结束的点,规定峰值前的点为事件开始的点。
本发明以科学的方法将震动波到时得到精准地拾取,大大减少了人工识别的工作量及误差,所述震动波到时自动拾取系统装配在一台工业计算机中,以确保系统运行效率。
现有技术中对微震事件震动波的数据进行到时拾取一般采用的方法为AIC法和STA/LTA法,由图3可知,使用同一个微震事件震动波动的数据,本发明方法在精度上面有明显的优势,本发明方法在计算速度方面也要高于这2种常规方法,在处理采样点极多的微震震动波数据,会增加计算效率,节省时间。
为了验证本发明方法的明显优势,利用这2种常规方法和本发明方法对模拟数据进行对比。
传统AIC拾取到时方法,进行进一步说明,步骤如下:
步骤S1:输入微震震动波数据
利用微震数据采集装置,采集一段震动波,将所采集到的微震事件波形文件传输到计算机内的已编写好的AIC到时拾取系统中。
步骤S2:AIC拾取到时
AR-AIC算法是一种基于AIC信息准则的到时拾取方法,该方法基于微震信号的非平稳特性,将信号划分为多个固定长度的波形段,然后进行自回归AR处理,求取AR模型的阶数和系数,这是整个到时拾取算法的前提。
该方法的计算过程实际是求取微震信号AIC函数局部最小值的过程,AR-AIC算法可描述为
STA/LTA法,进行进一步说明,步骤如下:
步骤S1:输入微震震动波数据
利用微震数据采集装置,采集一段震动波,将所采集到的微震事件波形文件传输到计算机内的已编写好的AIC到时拾取系统中。
步骤S2:计算STA/LTA长度
选波峰或者波谷:x(i1)=|x(i)|max
选左侧次波峰或波谷:
x(iguodu1)=(|x(i1)|-|x(i-a)|)max a=1,2,3…
x(i2)=(|x(iguodu1)|-|x(iguodu-a)|)min a=1,2,3…
选右侧次波峰或波谷:
x(iguodu2)=(|x(i1)|-|x(i+a)|)max a=1,2,3…
x(i3)=(|x(iguodu2)|-|x(iguodu+a)|)min a=1,2,3…
STA长度为:L1=a×|i3-i2| a取2~3之间
LTA长度为:L2=b×L1 a取5~10之间
步骤S3:选取CF(j),这里仅仅以绝对值举例
CF(j)=|x(i)|
步骤S4:选取大干λ的范围
步骤S5:到时点和微震事件区间长度
微震事件长度:L=ib-ia
到时点:ia
应当理解的是,本发明的上述具体实施方式仅仅用于示例性说明或解释本发明的原理,而不构成对本发明的限制。因此,在不偏离本发明的精神和范围的情况下所做的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。此外,本发明所附权利要求旨在涵盖落入所附权利要求范围和边界、或者这种范围和边界的等同形式内的全部变化和修改例。
Claims (6)
1.一种微震震动波到时拾取方法,其特征在于,包括:
S1:输入及处理数据;
利用微震数据采集装置,采集一段震动波,将所采集到的微震事件波形文件传输到计算机内的到时拾取系统中,为了有效压制底噪的影响,确保震动波到时拾取的精准度,采用FIR滤波器进行滤波;
S2:自动选取时窗长度L;
对经FIR滤波器滤波后的微震信号,自动选取一个合适时窗长度L,并截取在该时窗的长度范围内的微震数据;
S3:时窗能量特征值法精确拾取;
在选取的一个合适时窗长度L内,采用改进后时窗能量特征值算法精确拾取到时点;
S4:划分精确拾取范围;
通过步骤S3时窗能量特征值法拾取到的到时点,为更加精确拾取到时点,在所拾取的到时点周围划分得到一个合适的拾取范围M,该M远远小于本文中的时窗长度L,即M<<L;
S5:AIC精确拾取;
对划分精确拾取范围M,使用改进后W-AIC算法计算AIC值,然后进行峰值判断,在此范围M的数据内更加精确的拾取出到时点;
S6:峰值判断;
经峰值判断为是,可得到精确到时点;
若拾取到事件结束的点,峰值判断为否,则进行步骤S4、S5,重新划分精确的拾取范围Mi,并在重新划分精确的拾取范围Mi进行AIC精确拾取,直到峰值判断为是。
2.根据权利要求1所述的微震震动波到时拾取方法,其特征在于,若微震震动波为周期性的波,所述时窗长度L为微震震动波的两个周期长度。
3.根据权利要求1所述的微震震动波到时拾取方法,其特征在于,若微震震动波为非周期性的波,所述时窗长度L为:
当最大振幅位于波谷时,找到其前、后相邻的波谷位置,并把这两个位置间的跨度长度作为时窗长度L;所述时窗长度L计算如下:
a.求取|x(t0)|max
b.如果x(t0)·x(t0-1)<0或者x(t0)=0则a1=a1+1
c.当a1=2时,取x(t1-n)first,min,n=1,2,3,4…
d.当a1=2时,取x(t2-n)first,max,n=1,2,3,4…
式中,|x(t0)|max为最大幅值;a1为幅值异号点次数,t0为最大幅值点;x(t1-n)first,min为次波谷;x(t2-n)first,max为次波峰;L为选取的时窗长度。
4.根据权利要求1所述的微震震动波到时拾取方法,其特征在于,若微震震动波为非周期性的波,所述时窗长度L为:
当最大振幅位于波峰时,找到其前、后相邻的波峰位置,并把这两个位置间的跨度长度作为时窗长度L;所述时窗长度L计算如下:
a.求取|x(t0)|max
b.如果x(t0)·x(t0-1)<0或者x(t0)=0则a1=a1+1
c.当a1=2时,取x(t1-n)first,min,n=1,2,3,4…
d.当a1=2时,取x(t2-n)first,max,n=1,2,3,4…
式中,|x(t0)|max为最大幅值;a1为幅值异号点次数,t0为最大幅值点;x(t1-n)first,min为次波谷;x(t2-n)first,max为次波峰;L为选取的时窗长度。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210011651.5A CN114355441A (zh) | 2022-01-06 | 2022-01-06 | 一种微震震动波到时拾取方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210011651.5A CN114355441A (zh) | 2022-01-06 | 2022-01-06 | 一种微震震动波到时拾取方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN114355441A true CN114355441A (zh) | 2022-04-15 |
Family
ID=81106932
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210011651.5A Pending CN114355441A (zh) | 2022-01-06 | 2022-01-06 | 一种微震震动波到时拾取方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114355441A (zh) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN114994791A (zh) * | 2022-05-27 | 2022-09-02 | 中国矿业大学 | 一种用于评估井地一体微震监测系统监测能力的方法 |
-
2022
- 2022-01-06 CN CN202210011651.5A patent/CN114355441A/zh active Pending
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN114994791A (zh) * | 2022-05-27 | 2022-09-02 | 中国矿业大学 | 一种用于评估井地一体微震监测系统监测能力的方法 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN106353792B (zh) | 一种适用于水力压裂微震震源定位的方法 | |
CN104266894B (zh) | 一种基于相关性分析的矿山微震信号初至波时刻提取方法 | |
CN102928873B (zh) | 基于四维能量聚焦的地面微地震定位方法 | |
CN104216008A (zh) | 一种井中压裂微地震事件识别方法 | |
CN114063153B (zh) | 一种自动反演震源机制解的方法与装置 | |
CN103399346B (zh) | 一种井震联合初始波阻抗建模方法 | |
CN103365916A (zh) | 地震事件参数估计获取方法和系统,地震事件搜索引擎 | |
CN104570076A (zh) | 一种基于二分法的地震波初至自动拾取方法 | |
CN110910613A (zh) | 一种岩体微震无线监测接收预警系统 | |
CN208872879U (zh) | 一种微地震数据采集系统 | |
CN103742131A (zh) | 随钻声波井下信号采集与处理系统的时差实时提取方法 | |
CN114355441A (zh) | 一种微震震动波到时拾取方法 | |
CN111929728A (zh) | 一种三维三分量超前精细化地质预报方法 | |
CN102879813B (zh) | 一种微地震信号到时自动拾取的方法及装置 | |
CN102721523B (zh) | 一种区间风谱模型的建立方法 | |
CN104182651B (zh) | 用于三分量检波器接收的微地震事件方位角自动质控方法 | |
CN114002733B (zh) | 微震波信号初至到时自动拾取方法及微震监测装置 | |
CN111308548A (zh) | 一种高精度微震数据初至拾取装置、系统及方法 | |
CN109298451B (zh) | 一种改进偏斜度的自动拾取s波震相方法 | |
CN110146920B (zh) | 基于幅值相对变化的微震事件检测方法和系统 | |
CN106646592B (zh) | 单炮记录实时叠加方法 | |
Takanami et al. | The SIL seismological data acquisition system—As operated in Iceland and in Sweden— | |
CN105572733B (zh) | 一种地震速度谱自动拾取方法 | |
JIA et al. | Joint arrival-time picking method of microseismic P-wave and S-wave based on time-frequency analysis | |
CN105549066A (zh) | 生命信息探测方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination |