一种超低频次声异常信号判别方法与流程

本发明涉及一种异常信号判别方法，特别涉及一种超低频次声异常信号判别方法，属于试验测试技术领域。

背景技术：

在大地震前、海啸、流星、空气爆炸等自然与人工事件中会出现一些0.1Hz以下的超低频次声信号。目前普遍采用的监测方法多是通过多个站点全天无休的工作方式监测次声信号，这导致关注频段信号将淹没在大量次声数据中。又由于超低频次声波具有时间尺度大、传播距离远、信号失真严重、来源与背景信号复杂等特点，同一次声声源的发出的信号能被多个站点检测到，但不同站点监测到的次声信号相似且存在一些差别，常见的时域与频域分析方法无法有效的实现这些相似超低频次声异常是否来自同一事件的判别。

技术实现要素：

本发明的目的是为解决在次声监测中超低频异常次声信号的判断和多监测站点的超低频次声信号是否来自同一事件的判断问题，提出一种超低频次声异常信号判别方法。

本发明的思想是根据超低频次声信号时间尺度大、传播距离远、信号失真严重、来源与背景信号复杂等特点，通过能量曲线分析初步判断是否存在异常次声信号，若存在则通过广域多点长时数字离散相关分析，判断不同站点间的异常次声信号是否来自同一次声事件。

本发明的目的是通过以下技术方案实现的：

一种超低频次声异常信号判别方法，首先，通过对输入信号进行能量曲线分析初步判断是否存在异常次声信号，若存在则通过广域多点长时数字离散相关分析，判断不同站点间的异常次声信号是否来自同一次声事件，即实现异常一致性验证。

作为优选，所述能量曲线(STA/LTA)分析依据短时窗内次声信号能量在长时窗内次声信号能量的百分比来判断是否存在幅值较大的异常次声信号。

作为优选，所述能量曲线分析的具体过程为：

①对输入信号进行带通滤波，得到滤波后的信号X(n)；

②通过公式(1)计算滤波后信号X(n)的能量曲线S(n)，其中，L表示滤波后信号X(n)长度，N表示短时窗长度，n表示采样点序号；

③若能量曲线S(n)中存在大于预设的异常信号阈值TH的采样点，则判定这段能量曲线中存在异常次声信号，并取能量曲线S(n)大于TH部分的最大值为异常信号的粗略到达时间。

作为优选，所述广域多点长时数字离散相关分析通过对能量曲线离散化以及将来自不同站点的离散曲线进行相关性分析完成。

作为优选，所述对能量曲线离散化是对能量曲线分析中的能量曲线S(n)进行数字离散化，具体过程为：

①设定数字离散化的长时窗长度Ls,短时窗长度Ss；

②对能量曲线S(n)中长度为Ls的部分进行短时窗峰值提取，在长时窗中每Ss个点合成为1个数据点；将短时窗中具有较大能级的部分，转为一个能反应其能级大小的数，而对于不具有较大能级的部分，进行归零处理，从而得到离散化能量曲线Ss(n)。

作为优选，为使能量曲线Ss(n)中短时窗内能量曲线峰值点唯一化且将低能量的部分归零，所述Ss应为目标特征周期的2～100倍，Ls不小于信号特征周期的100倍。

作为优选，所述每Ss个点合成为1个数据点，合成方法为平均或加权平均。

作为优选，所述将来自不同站点的离散曲线进行相关性分析是计算两个不同站点的离散化能量曲线的相关系数，并判断其是否大于预设的相关阈值T，若是则判定两个离散化能量曲线对应的异常信号来自同一次声声源。

作为优选，所述相关系数为皮尔逊相关系数。

有益效果

传统人工选取异常次声信号存在工作繁重、判断主观的问题，本发明方法可以从次声监测数据中智能提取特定频段的高幅值信号，减少人工识别的工作量，为异常次声信号的判断提供了一个可以量化的依据，并对次声信号间的相关程度进行量化计算，提出一种判别超低频次声波是否来自同一事件的判定方法。

附图说明

图1是本发明实施例一种超低频次声异常信号判别方法流程示意图。

图2是本发明实施例经带通滤波后的输入信号示意图。

图3是通过基于能量曲线法对图2信号处理后的结果示意图。

图4是对图3中能量曲线离散化的处理结果示意图。

具体实施方式

下面结合附图与实施例对本发明加以详细说明。

实施例1

一种超低频次声异常信号判别方法，流程如图1所示，包括以下步骤：

1.选取任一监测站点一段监测信号记为Xs(n)，按照需要选取关注的信号频率范围F1～F2，本实施例中F1＝0.001Hz，F2＝0.01Hz,则周期大致在100s～1000s之间，取目标特征为1000s。

2.对原始信号Xs(n)进行频带范围为F1～F2的带通滤波，得到滤波后信号X(n)。

3.通过公式(1)计算滤波后信号X(n)的能量曲线Sx(n)，本实施例中，L取信号X(n)长度，短时窗长度N＝300。

4.将能量曲线Sx(n)与阈值TH比较，若则说明该段信号中不含有超低频异常次声信号，返回步骤1分析其他段信号；否则，通过取能量曲线Sx(n)最大值的出现时间为异常信号的粗略到达时间teX，能量曲线最大脉冲的持续时间(可看作目标特征周期)的一半为半周期持续时间TeX；本实施例中TH＝0.05。

5.选取其他任一监测站点一段监测信号记为Ys(n)，重复步骤2～4，得到滤波后信号Y(n)，能量曲线Sy(n)，异常信号的粗略到达时间teY，半周期持续时间TeY。

6.对能量曲线Sx(n)进行数字离散化，取长时窗长度Ls为能量曲线Sx(n)的信号长度,短时窗长度Ss为半周期持续时间TeX的5倍，得到离散化的能量曲线Ssx(n)。

7.对能量曲线Sy(n)进行数字离散化，取长时窗长度Lsy为能量曲线Sy(n)的信号长度,短时窗长度Ss为半周期持续时间TeY的5倍，得到离散化的能量曲线Ssy(n)。

8.计算离散化的能量曲线Ssx(n)和离散化的能量曲线Ssy(n)的相关系数，本实施例中采用皮尔逊相关系数，若系数大于阈值T,则信号Xs(n)和Ys(n)为同一声源发出的信号，反之，则不是；本实施例中T＝0.3。

试验结果

试验条件：

采用CASI-ISM系列传感器采集次声信号，在室温下工作。

如图2所示为位于重庆、厦门、南京、泰州、北京的次声监测站点监测到的次声信号，这些信号已经经过0.001Hz～0.01Hz的带通滤波，横坐标为时间，纵坐标为声波大小。

如图3所示为通过能量曲线分析(改进后的STA/LTA)方法对图2中的信号进行处理后得到的能量曲线与原始信号的对比图。图中横坐标为采样点序号，未标注的纵坐标表示短视窗内能量占整段能量的百分比。从图3中可以看出，在次声幅值较大的部分在能量曲线中也具有较大能量，因此，通过能量曲线可以很大程度上反应次声信号幅值大小，因此，可以通过阈值判断异常次声的到达时间。

如图4所示为在图3能量曲线的基础上通过长时数字离散相关分析得到的离散能量曲线。图中横坐标与图3一样，为采样点序号，纵坐标表示离散点附近的百分比均值。从图中可以看出，离散后的能量曲线能在很大程度上反应能量曲线的变化情况，并能较好的消除小能量的影响，有效提高各信号间的皮尔逊相关系数。

在试验中发现，计算得到次声信号间的皮尔逊相关系数并不理想，即使是南京、泰州这两个间距在100km左右的次声监测站点的信号，相关系数也常在0.3左右，相关程度为弱相关，考虑到次声波具有穿透能力强、传播损失低、传输距离远的特点，选择的南京、泰州的信号理论上应该具有同源性，却在距离如此近的情况下得到这么低的相关系数，判断可能是次声波在传输过程中，因为传播介质分布不均、空气流动等因素导致次声波的某些部分被拉伸或者压缩了，从而导致监测到的次声信号间的相关系数低的情况。

根据以上判断，我们计算了STA/LTA方法得到的各站点次声信号能量曲线间的相关系数，得到了较好的结果，南京、泰州的信号的能量曲线间的相关系数可以高达0.8以上，属于极强相关。能量曲线表示次声信号能量随时间的分布情况，能量曲线相关性高说明对各信号能量随时间分布情况类似，而我们之前判断的次声波在传播过程中因各种因素使得某些部分被压缩或者拉伸的情况对监测到的次声信号的能量分布影响较小，因此，这也从侧面证明之前的判断具有一定的正确性。

即使这样，北京、重庆这种远距离站点监测到的信号能量曲线间的相关系数依然不够理想，结合地震信号、次声信号处理中常用的指纹法的启发，对已有的能量曲线进行离散化，既保留能量曲线中能量的幅值信息，又保留了能量曲线中各能量峰值出现的时间信息。

以上所述的具体描述，对发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。