一种基于DWT_STA/LTA的含噪信号P波初至峰度拾取方法技术

本发明专利技术公开了一种基于DWT_STA/LTA的含噪信号P波初至峰度拾取方法，包括如下步骤：提取矿山微震信号；离散小波变换(DWT)得到微震细节信号；计算各细节信号的STA/LTA值；判定STA/LTA最大值是否大于阀值a，若是则取首个触发点对应的采样点序号为该细节信号的大致到时，否则记该细节信号到时为0；再对到时非0的细节信号进行局部峰态拾取；最后采用权重分析确定P波最终到时。该方法有效排除了低信噪比、刺突、尾部震荡及纯噪音信号的干扰，极大地降低了错误拾取率，增强了P波拾取的准确性。

【技术实现步骤摘要】
一种基于DWT_STA/LTA的含噪信号P波初至峰度拾取方法
本专利技术属于信号处理
，尤其是涉及一种矿山微震信号P波初至时刻联合拾取方法。
技术介绍
P波初至拾取是矿山微震识别、震源定位及震源机制解释中最为基础和重要的一步。目前主要采用人工拾取P波到时，但其易受个人经验、情绪等因素影响，且拾取数据量大时非常耗时。鉴于上述原因，国内外提出了一系列P波初至自动拾取方法，但到目前为止还没有任何一种算法能够在不同的震源环境、震中距范围、噪声背景下进行微震初至的一致性检测，更不用说对后续震相初至的有效检测和拾取。常见的P波拾取方法包括长短时窗平均值比法、分形维数法、自回归模型、人工神经网络法、高阶统计量法等。长短时窗平均值比法利用短时窗平均值与长时窗平均值之比来反映信号幅值的变化，并取比值大于某一设定阈值时为P波到时。长短时窗平均值比法具有算法简单、计算速度快、平均值特征函数多样等特点，但其拾取低信噪比信号困难，且阀值较大时，则可能拾取不到P波初至或增大拾取误差；阀值较小时，则可能过早的被噪音触发。分形维数法认为噪音与信号叠加时分形维数发生变化，并以此作为P波拾取的依据。该算法抗噪性能好，拾取精度较高，但其对插值的准确性、时窗和步长的依赖性很高，稍有不慎就会严重影响拾取结果，且其计算速度较慢。自回归模型将微震信号分为两个局部统计时段，并取自相关最小值点作为P波到时。该算法拾取精度较高，但其对低信噪比和尾部震荡信号拾取不稳定，且其采用最小值点作为P波到时，因此对纯噪音信号仍会拾取到时(HaijiangZhang等，2003)。人工神经网络法采用峰值振幅、时窗内均方根振幅比、峰值与其前后峰值的包络斜率及噪声与信号的比值等作为神经网络输入因素，人工神经网络法是一种较为综合的P波初至拾取方法，对微震波形适应性强，但其参数计算工作量大，学习速度慢，实现复杂。高阶统计量法(PAI-S/K法)由Saragiotis等(2002,2004)提出，取峰度(或偏度)最大值点作为P波初至时刻，该算法借助于信号只包含噪音时，峰度值和偏度值趋于0，而当P波到达时峰度和偏度值开始增加，并当滑动窗口内包含部分微震信号时达到最大值。该算法拾取精度高，但未考虑尖刺和尾部震荡信号对拾取的影响，由此可能产生较大的拾取误差，且其采用最大值拾取P波到时，对纯噪音信号仍会拾取到时。可见现有的矿山微震信号P波初至拾取方法存在很大的局限，需要研究一种拾取精度高、稳定性好的自动拾取方法。
技术实现思路
本专利技术所要解决的技术问题是提供一种改进的含噪信号P波初至峰度拾取方法，解决峰度法(PAI-K法)拾取低信噪比、刺突、尾部震荡及纯噪音信号不稳定的技术问题，该矿山微震信号P波初至拾取方法适用性强、准确性高。针对矿山微震信号P波初至PAI-K法拾取精度较高，但对低信噪比、刺突、尾部震荡及纯噪音信号拾取不稳定的技术问题，本专利技术提出了W-S/L-K拾取法，该方法借助于小波分解高频细节信号能够较好的保留P波初至信息以及STA/LTA法能很好地克服刺突和尾部震荡的影响并确定P波大致到时，进而采用权重分析排除上述干扰，极大地降低了错误拾取率，增强了P波拾取的准确性。一种基于DWT_STA/LTA的含噪信号P波初至峰度拾取方法，包括以下几个步骤：步骤1：提取微震信号x(n)从矿山微震信号中提取待拾取信号x(n)，n＝1,2,…,N，其中，N为所述微震信号的采样点个数，取N＝400～600，微震信号采样频率600Hz；步骤2：对微震信号采用离散小波变换DWT得到细节信号Di，i＝1,2,3,4；步骤3：依据STA/LTA算法按照以下公式确定细节信号Di的P波大致到时ki；其中，STA(k,Di)和LTA(k,Di)分别为短、长时窗内微震信号的振幅平均值；x(n,Di)为细节信号Di的第n个采样信号；k为采样的第k个点，k＝WLTA,WLTA+1,…,N；WSTA和WLTA分别为短、长时窗的长度，分别取10和60个采样点的长度；若细节信号Di对应的短时窗和长时窗的振幅平均值的比值λ(k,Di)的最大值大于设定阀值a，则取首个触发点对应的采样点序号为该细节信号Di的P波大致到时ki，否则，记该细节信号Di的P波大致到时ki为0。其中，i＝1,2,3,4；步骤4：计算细节信号的P波大致到时ki＝0的个数L，i＝1,2,3,4；步骤5：利用权重分析计算P波初至准确到时；若步骤4中L>2，则认为该信号为纯噪音信号，反之则采用以下权重公式计算P波准确到时k'：其中，1/8、3/8、3/8和1/8分别为D4、D3、D2和D1的权重系数；Ai和ki'分别为Ri的最大振幅和局部峰度拾取到时，Ri为Di在区间[ki-b,ki+b]的信号，b为确定局部峰度拾取区间的常量，取b＝10～15；由公式计算x(n,Di)的滑动峰度值K(k,Di)，k为采样的第k个点，k＝ki-b,ki-b+1,…,ki+b，M为滑动时窗的长度，取M＝20～30，并取K(k,Di)最大值点作为局部峰度拾取到时ki'；步骤6：依据步骤5输出拾取结果，该信号为纯噪音信号或者P波初至准确到时为k'。所述设定阈值a＝2～3。利用Mallat算法快速实现离散小波变换。并采用低通滤波器h和高通滤波器g获得原始信号c0(等同于本专利技术中的x(n))的小波系数{d1,d2,…,dJ}，从而将微震信号x(n)分解得到细节信号Di(取i＝1,2,3,4)。式中：cj+1,k为第k个采样点的近似信号值cj+1，dj+1,k为第k个采样点的细节信号值dj+1，k为采样的第k个点，l为滤波器确定指标，j为尺度参数(j＝0,1,…,J-1)，J为最大分解层数。为简化说明，记Cj+1＝{cj+1,k}和Dj+1＝{dj+1,k}。设定ki＝0时，Ai＝0，ki'＝0，i＝1,2,3,4，简化计算过程，提高计算速度。有益效果本专利技术提供了一种基于DWT_STA/LTA的含噪信号P波初至峰度拾取方法，包括如下步骤：提取矿山微震信号；离散小波变换(DWT)得到微震细节信号；计算各细节信号的STA/LTA值；判定STA/LTA最大值是否大于阀值a，若是则取首个触发点对应的采样点序号为该细节信号的大致到时，否则记该细节信号到时为0；再对到时非0的细节信号进行局部峰态拾取；最后采用权重分析确定P波最终到时。本专利技术拾取纯噪音信号时，细节信号STA/LTA触发数通常小于2，设定STA/LTA触发数大于等于2即可排除噪音信号；拾取低信噪比信号时，STA/LTA对含噪音较大的细节信号无拾取，且细节信号无拾取的数目通常不大于2，因此可用其余细节信号拾取初至，从而排除低信噪比的影响；通常刺突信号对STA/LTA值影响较小，设定合适的STA/LTA触发值即可排除刺突的影响；STA/LTA通常在P波初至时已触发，不会受尾部震荡的影响。由此，通过一种基于DWT_STA/LTA的含噪信号P波初至峰度拾取方法有效排除了低信噪比、刺突、尾部震荡及纯噪音信号的干扰，极大地降低了错误拾取率，增强了P波拾取的准确性。附图说明图1是本专利技术所述方法流程图；图2是典型波形W-S/L-K法实现过程图，其中，(a)为矿山微震原始信号x(n)(n＝1,2,…,400)，虚线为人工拾取P波初至点本文档来自技高网...

【技术保护点】
一种基于DWT_STA/LTA的含噪信号P波初至峰度拾取方法，其特征在于，包括以下几个步骤：步骤1：提取微震信号x(n)从矿山微震信号中提取待拾取信号x(n)，n＝1,2,…,N，其中，N为所述微震信号的采样点个数，取N＝400～600，微震信号采样频率600Hz；步骤2：对微震信号采用离散小波变换DWT得到细节信号Di，i＝1,2,3,4；步骤3：依据STA/LTA算法按照以下公式确定细节信号Di的P波大致到时ki；λ(k,Di)=STA(k,Di)LTA(k,Di)=1WSTAΣn=k-WSTAk|x(n,Di)|1WLTAΣn=k-WLTAk|x(n,Di)|]]>其中，STA(k,Di)和LTA(k,Di)分别为短、长时窗的振幅平均值；x(n,Di)为细节信号Di的第n个采样信号；k为采样的第k个点，k＝WLTA,WLTA+1,…,N；WSTA和WLTA分别为短、长时窗的长度，分别取10和60个采样点；若细节信号Di对应的短时窗和长时窗的振幅平均值的比值λ(k,Di)的最大值大于设定阀值a，则取首个触发点对应的采样点序号为该细节信号Di的P波大致到时ki，否则，记该细节信号Di的P波大致到时ki为0。其中，i＝1,2,3,4；步骤4：计算细节信号的P波大致到时ki＝0的个数L，i＝1,2,3,4；步骤5：利用权重分析计算P波初至准确到时；若步骤4中L>2，则认为该信号为纯噪音信号，反之则采用以下权重公式计算P波准确到时k'：k′=18A4·k4′+38A3·k3′+38A2·k2′+18A1·k1′18A4+38A3+38A2+18A1]]>其中，1/8、3/8、3/8和1/8分别为D4、D3、D2和D1的权重系数；Ai和k′i分别为Ri的最大振幅和局部峰度拾取到时，Ri为Di在区间[ki‑b,ki+b]的信号，b为确定局部峰度拾取区间的常量，取b＝10～15；由公式K(k,Di)=(M-1)Σn=k-M+1k(x(n,Di)-m^(k,Di))4Σn=k-M+1k(x(n,Di)-m^(k,Di))2-3]]>计算x(n,Di)的滑动峰度值K(k,Di)，k为采样的第k个点，k＝ki‑b,ki‑b+1,…,ki+b，M为滑动时窗的长度，取M＝20～30，并取K(k,Di)最大值点作为局部峰度拾取到时k′i；步骤6：依据步骤5输出拾取结果，该信号为纯噪音信号或者P波初至准确到时为k'。...

【技术特征摘要】
1.一种基于DWT_STA/LTA的含噪信号P波初至峰度拾取方法，其特征在于，包括以下几个步骤：步骤1：提取微震信号x(n)从矿山微震信号中提取待拾取信号x(n)，n＝1,2,…,N，其中，N为所述微震信号的采样点个数，取N＝400～600，微震信号采样频率600Hz；步骤2：对微震信号采用离散小波变换DWT得到细节信号Di，i＝1,2,3,4；步骤3：依据STA/LTA算法按照以下公式确定细节信号Di的P波大致到时ki；其中，STA(k,Di)和LTA(k,Di)分别为短、长时窗内微震信号的振幅平均值；x(n,Di)为细节信号Di的第n个采样信号；k为采样的第k个点，k＝WLTA,WLTA+1,…,N；WSTA和WLTA分别为短、长时窗的长度，分别取10和60个采样点的长度；若细节信号Di对应的短时窗和长时窗的振幅平均值的比值λ(k,Di)的最大值大于设定阀值a，则取首个触发点对应的采样点序号为该细节信号Di的P波大致到时ki，否则，记该细节信号Di的P波大致到时ki为0，其中，i＝1,2,3,4；步骤4：计算细节信号的P波大致到时ki＝0的个数L，i＝1,2,3,...

【专利技术属性】
技术研发人员：李夕兵，尚雪义，董陇军，王泽伟，马举，张楚旋，
申请(专利权)人：中南大学，
类型：发明
国别省市：湖南;43