CN108844616B - 一种次声事件检测方法及系统 - Google Patents

一种次声事件检测方法及系统 Download PDF

Info

Publication number
CN108844616B
CN108844616B CN201810365902.3A CN201810365902A CN108844616B CN 108844616 B CN108844616 B CN 108844616B CN 201810365902 A CN201810365902 A CN 201810365902A CN 108844616 B CN108844616 B CN 108844616B
Authority
CN
China
Prior art keywords
pixel
fisher
infrasound
signals
attribute
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.)
Active
Application number
CN201810365902.3A
Other languages
English (en)
Other versions
CN108844616A (zh
Inventor
李夕海
牛超
曾小牛
陈蛟
刘志刚
刘继昊
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Rocket Force University of Engineering of PLA
Original Assignee
Rocket Force University of Engineering of PLA
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Rocket Force University of Engineering of PLA filed Critical Rocket Force University of Engineering of PLA
Priority to CN201810365902.3A priority Critical patent/CN108844616B/zh
Publication of CN108844616A publication Critical patent/CN108844616A/zh
Application granted granted Critical
Publication of CN108844616B publication Critical patent/CN108844616B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01HMEASUREMENT OF MECHANICAL VIBRATIONS OR ULTRASONIC, SONIC OR INFRASONIC WAVES
    • G01H17/00Measuring mechanical vibrations or ultrasonic, sonic or infrasonic waves, not provided for in the preceding groups

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Image Analysis (AREA)
  • Investigating Or Analyzing Materials By The Use Of Ultrasonic Waves (AREA)

Abstract

本发明公开了一种次声事件检测方法及系统。所述方法及系统首先使用小波包分解的方法,将待检测次声信号分解之后,在不同频段内进行重构,然后对每个频段内的重构信号进行Fisher检测,并且对多阵元同一时窗内的每个Fisher值作为一个属性保存,计算结果生成一幅时间-频率平面上的二维Fisher属性分布图,同时,考虑信号检测的相关时间尺度,以一定的尺度阈值作为衡量信号的标准,降低了环境噪声对Fisher属性值计算的影响,从而降低了次声事件检测的漏警率和虚警率,提高了次声事件检测的效率和精度。

Description

一种次声事件检测方法及系统
技术领域
本发明涉及次声监测技术领域,特别是涉及一种次声事件检测方法及系统。
背景技术
对空中核爆炸(乃至地下核爆炸)事件准确监测的需求使得次声监测技术受到非常大的重视。在次声监测技术中,由于次声信号是一种低频声信号,极易受到背景风和其它大气噪声的影响,因此从众多背景噪声中准确地检测出次声事件,是次声监测面临的首要任务。目前,常用的次声检测方法主要是基于多通道信号的初至点检测,如Melton和Bailey提出的基于F-统计的Fisher(费歇尔)检测算法;Smart和Flinn给出了改进的Fisher检测方法,并将F-统计的原理引入频域;Cansi从阵列信号处理的角度将互相关算法(Progressive Multi-Channel Correlation,PMCC)引入次声检测;国际数据中心(International Data Center,IDC)根据短时平均和长时平均比(Short Term Average/Long TermAverage,STA/LTA)研究出了Libinfra软件,并将上述算法进行了实际应用。除此之外,David结合次声传播特性提出了基于Hough变换(霍夫变换)的次声事件检测方法。在这些检测方法中,基于F-统计的Fisher的检测及其改进算法以及IDC中应用的STA/LTA等方法均是基于待检测信号的瞬时特性而提出的,没有考虑信号波形的前后关系,必然会增加次声事件检测的虚警率(尤其在强噪声情况下)。虽然PMCC及基于Hough变换的次声检测等方法均利用了次声信号波形的前后关系,通过对信号波形的部分重叠增强了信号检测能力,但是目前的PMCC及基于Hough变换的次声检测受信噪比影响仍然较大,在低信噪比条件下,检测性能不佳,且检测虚检率高,增加了工作量,如PMCC法,从统计上说,80%以上检测结果为虚假或无用信号。因此,如何在保证次生事件信号检验效果的同时降低误检率,是本领域亟待解决的技术问题。
发明内容
本发明的目的是提供一种次声事件检测方法及系统,能够降低次声事件检测的漏警率和虚警率,提高次声事件检测效率和检测精度。
为实现上述目的,本发明提供了如下方案:
一种次声事件检测方法,所述检测方法包括:
获取次声台阵的多个阵元接收到的多路待检测次声信号;
对多路所述待检测次声信号进行小波包分解和信号重构,生成不同阵元不同频段的重构信号;
将所述不同阵元不同频段的重构信号中相同频段的重构信号进行归类,获得多个归类信号;
采用滑动窗对多个所述归类信号进行分割,获得多个不同频段的窗内信号;
计算各个所述窗内信号的Fisher属性值;
将所述Fisher属性值作为像素值,生成二维Fisher属性分布图;
对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇;每个所述像素簇中包括多个所述像素点;
判断所述像素簇中像素点的数量是否超过尺度阈值,获得第一判断结果;
当所述第一判断结果为所述像素簇中像素点的数量超过尺度阈值,确定发生了次声事件;
当所述第一判断结果为所述像素簇中像素点的数量未超过尺度阈值,确定未发生次声事件。
可选的,所述对多路所述待检测次声信号进行小波包分解和信号重构,生成不同阵元不同频段的重构信号,具体包括:
对多路所述待检测次声信号进行三层小波包分解,获得分解后的小波包分解系数;
对所述小波包分解系数实施信号重构,生成不同阵元不同频段的重构信号。
可选的,所述计算各个所述窗内信号的Fisher属性值,具体包括:
采用公式
Figure BDA0001637039830000021
计算各个所述窗内信号的Fisher属性值;其中F表示所述窗内信号的Fisher属性值;
Figure BDA0001637039830000031
N表示次声台阵的阵元数量,M为时间窗内采样点的数量,
Figure BDA0001637039830000032
xij表示第j个阵元在第i时刻接收到的次声信号,
Figure BDA0001637039830000033
可选的,所述将所述Fisher属性值作为像素值,生成二维Fisher属性分布图,具体包括:
将多个所述Fisher属性值作为所述二维Fisher属性分布图中多个像素点的像素值,生成二维Fisher属性分布图;其中所述二维Fisher属性分布图的横坐标为时间,纵坐标为频率。
可选的,所述对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇,具体包括:
判断所述二维Fisher属性分布图中像素点的像素值是否高于筛选阈值,获得第二判断结果;
当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值低于筛选阈值,将所述像素点的像素值归零;
当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值高于筛选阈值,判断所述像素点的像素值与其相邻像素点的像素值是否相同,获得第三判断结果;
当所述第三判断结果为所述像素点的像素值与其相邻像素点的像素值相同,将所述像素点与其相邻像素点合并为一个像素簇;
当所述第三判断结果为所述像素点的像素值与其相邻像素点的像素值不同,则不合并。
可选的,所述对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇,还包括:
判断所述二维Fisher属性分布图中像素点的像素值与其相邻像素点的像素值是否相同,获得第四判断结果;
当所述第四判断结果为所述像素点的像素值与其相邻像素点的像素值相同,将所述像素点与其相邻像素点合并为一个初选像素簇;
当所述第四判断结果为所述像素点的像素值与其相邻像素点的像素值不同,则不合并;
判断所述初选像素簇中像素点的像素值是否高于筛选阈值,获得第五判断结果;
当所述第五判断结果为所述初选像素簇中像素点的像素值低于筛选阈值,将所述初选像素簇归零;
当所述第五判断结果为所述初选像素簇中像素点的像素值高于筛选阈值,保留所述初选像素簇作为所述像素簇。
本发明还公开了一种次声事件检测系统,所述检测系统包括:
待检测次声信号获取模块,用于获取次声台阵的多个阵元接收到的多路待检测次声信号;
重构信号生成模块,用于对多路所述待检测次声信号进行小波包分解和信号重构,生成不同阵元不同频段的重构信号;
归类信号获取模块,用于将所述不同阵元不同频段的重构信号中相同频段的重构信号进行归类,获得多个归类信号;
窗内信号获取模块,用于采用滑动窗对多个所述归类信号进行分割,获得多个不同频段的窗内信号;
Fisher属性值计算模块,用于计算各个所述窗内信号的Fisher属性值;
Fisher属性分布图生成模块,用于将所述Fisher属性值作为像素值,生成二维Fisher属性分布图;
像素簇获取模块,用于对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇;每个所述像素簇中包括多个所述像素点;
第一判断结果获取模块,用于判断所述像素簇中像素点的数量是否超过尺度阈值,获得第一判断结果;
次声事件发生确定模块,用于当所述第一判断结果为所述像素簇中像素点的数量超过尺度阈值时,确定发生了次声事件;
次声事件未发生确定模块,用于当所述第一判断结果为所述像素簇中像素点的数量未超过尺度阈值时,确定未发生次声事件。
可选的,所述Fisher属性分布图生成模块具体包括:
Fisher属性分布图生成单元,用于将多个所述Fisher属性值作为所述二维Fisher属性分布图中多个像素点的像素值,生成二维Fisher属性分布图;其中所述二维Fisher属性分布图的横坐标为时间,纵坐标为频率。
可选的,所述像素簇获取模块具体包括:
第二判断结果获取单元,用于判断所述二维Fisher属性分布图中像素点的像素值是否高于筛选阈值,获得第二判断结果;
像素点归零单元,用于当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值低于筛选阈值时,将所述像素点的像素值归零;
第三判断结果获取单元,用于当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值高于筛选阈值时,判断所述像素点的像素值与其相邻像素点的像素值是否相同,获得第三判断结果;
像素簇合并单元,用于当所述第三判断结果为所述像素点的像素值与其相邻像素点的像素值相同时,将所述像素点与其相邻像素点合并为一个像素簇。
可选的,所述像素簇获取模块还包括:
第四判断结果获取单元,用于判断所述二维Fisher属性分布图中像素点的像素值与其相邻像素点的像素值是否相同,获得第四判断结果;
初选像素簇合并单元,用于当所述第四判断结果为所述像素点的像素值与其相邻像素点的像素值相同时,将所述像素点与其相邻像素点合并为一个初选像素簇;
第五判断结果获取单元,用于判断所述初选像素簇中像素点的像素值是否高于筛选阈值,获得第五判断结果;
像素簇归零单元,用于当所述第五判断结果为所述初选像素簇中像素点的像素值低于筛选阈值时,将所述初选像素簇归零;
像素簇生成单元,用于当所述第五判断结果为所述初选像素簇中像素点的像素值高于筛选阈值时,保留所述初选像素簇作为所述像素簇。
根据本发明提供的具体实施例,本发明公开了以下技术效果:
本发明提供一种次声事件检测方法及系统,所述方法及系统首先使用小波包分解的方法,将待检测次声信号分解之后,在不同频段内进行重构,由于小波包本身具有降噪特性,对信号噪声具有抑制作用,因此能够降低环境噪声对Fisher属性值计算的影响;然后对每个频段内的重构信号进行Fisher检测,并且对多阵元同一时窗内的每个Fisher值作为一个属性保存,计算结果生成一幅时间-频率平面上的二维Fisher属性分布图,同时,考虑信号检测的相关时间尺度,以一定的尺度阈值作为衡量信号的标准,降低了次声事件检测的漏警率和虚警率,提高了次声事件检测的效率和精度。
此外,本发明提供的次声事件检测方法及系统在获得像素簇之前,先对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,通过筛选阈值的设置能够有效排除毛刺信号的干扰,从而降低了次声事件检测的虚警率,进一步提高了次声事件的检测效率和检测精度。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动性的前提下,还可以根据这些附图获得其他的附图。
图1为本发明提供的一种次声事件检测方法的方法流程图;
图2为本发明提供的小波包三层分解树结构示意图;
图3为本发明提供的四元三角阵时延对齐后的实际次声信号示意图;
图4为本发明提供的实际次声信号叠加白噪声后信噪比为0db的仿真信号示意图;
图5为本发明提供的1号阵元进行三层小波包分解后得到的重构信号示意图;
图6为本发明提供的对像素点进行筛选合并处理后的二维分布图;
图7为本发明提供的像素簇经尺度阈值筛选后的二维分布图;
图8为本发明提供的信噪比为-7dB时采用Fisher方法获得的检测结果示意图;
图9为本发明提供的信噪比为-1dB时采用Fisher方法获得的的检测结果示意图;
图10为本发明提供的SGAR台阵接收的原始次声信号示意图;
图11为采用Fisher检测方法对实测信号进行检测的检测结果示意图;
图12为采用PMCC方法对实测信号进行检测的检测结果示意图;
图13为采用本发明提供的次声事件检测方法对实测信号进行检测的检测结果示意图;
图14为本发明提供的一种次声事件检测系统的结构示意图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
本发明的目的是提供一种次声事件检测方法及系统,能够降低次声事件检测的漏警率和虚警率,提高次声事件检测效率和检测精度。
为使本发明的上述目的、特征和优点能够更加明显易懂,下面结合附图和具体实施方式对本发明作进一步详细的说明。
首先对基于方差齐性检验(F-test,F检验)的Fisher检测(Fishertest,费歇尔检验)原理进行简单介绍:
次声信号检测分为单阵元次声信号检测和多阵元次声信号检测。多阵元次声信号检测方法中,次声台阵中设置的次声传感器记录了各种数据,当时间序列的幅值、形状或频率成分相对于背景噪声有明显变化时,就认为出现了信号。从统计学的角度分析,信号变化越明显,方差越大,因此,信号局部变化与信号总体变化构成的比值(Fisher值)变化也就越明显,从而达到检测的目的。
次声监测台站主要用于远程监测核爆和化爆等微弱事件,次声监测台站中往往设置多个次声监测台阵,其台阵结构对目标次声信号的分辨率及后续信号处理(如定位、识别等)的结果都有影响,常用的台阵结构有四元三角阵(包括四个阵元)、六元五边形阵(包括六个阵元)等,其中每个阵元接收一路次声信号。本文中出现的子台站、台阵子台、子台均表示次声监测台阵中的一个阵元。
定义N个子台站(阵元)的阵列记录为:
X=[x1(t),x2(t),x3(t),···,xN(t)]T(1)
其中,xN(t)表示第N个阵元(子台站)在t时刻记录的次声信号;X表示N个子台站(阵元)在t时刻记录的次声信号序列。
子台阵列样本(假定每个台阵的样本长度为M)的总体均值为:
Figure BDA0001637039830000081
其中xij表示第j个子台站在第i时刻记录的次声信号。
t时刻,N个子台的样本均值为:
Figure BDA0001637039830000082
假设噪声在时间和空间上是平稳的且不相关的,噪声方差为σ2
T时间窗内有M个采样点的样本方差和表示为:
Figure BDA0001637039830000083
经过计算可以将V分解如下:
Figure BDA0001637039830000084
设V=Vw+VB。其中:
Figure BDA0001637039830000085
由上式可知,第一部分Vw代表信号和噪声的总体能量,第二部分VB代表信号能量,根据统计理论可得:
E[Vw]=M(N-1)σ2 (6)
Figure BDA0001637039830000086
式中E[Vw]、E[VB]分别表示Vw和VB的均值。其中符号E[]指均值。
当没有信号时,各子台站的平均变化量与总体变化量保持一致,当信号到达时,变化比较明显时,时间窗内的平均变化量将显著增加。综上所述,可以构造F统计量检测信号:
Figure BDA0001637039830000091
F是定义为式(8)的一个统计量,也即信号的Fisher属性值。当有事件信号到来时,统计量F将显著增大。
从上文介绍的Fisher检测原理可以看出,Fisher检测是在噪声不相关的前提下提出的的,当存在相关噪声时,必然导致虚警率较高。此外,由于Fisher检测是基于待检测信号的瞬时特性而实施的,存在噪声情况下必然严重影响Fisher瞬时属性值的计算,进而影响次声信号的检测精度。为降低噪声对Fisher瞬时属性值计算的影响,如果将信号中的噪声进行滤波处理,然后再采用Fisher检测的方法进行信号事件检测,且考虑相关噪声的影响,以一定的时间尺度作为衡量信号的标准,其检测效果就可能有所改善。与此同时,由于次声波记录中的噪声不仅仅包含高频随机噪声,同时还有低频扰动,所以不可能用一个简单的滤波器就能起到去噪的效果。由此,本发明使用小波包分解的方法,将信号分解之后,在不同频段内进行重构,然后对每个频段内的重构信号进行Fisher检测,并且对多台站同一时窗内的每个Fisher值作为一个属性保存,计算结果形成一组时间-频率平面上的属性矩阵,且考虑信号检测的相关时间尺度。在噪声成分占很大比重的频段内,Fisher方法将检测不到事件信号,在事件信号占较大比重的频段内,Fisher方法则会发挥出该事件检测的有效性能。基于此构思,本发明提出了基于小波包分解和Fisher属性的次声事件检测方法。
图1为本发明提供的一种次声事件检测方法的方法流程图。参见图1,本发明提供的一种次声事件检测方法,具体包括:
步骤101:获取次声台阵的多个阵元接收到的多路待检测次声信号。
次声监测台站主要用于远程监测核爆和化爆等微弱事件,次声监测台站中往往设置多个次声监测台阵,其台阵结构对目标次声信号的分辨率及后续信号处理(如定位、识别等)的结果都有影响,常用的台阵结构有四元三角阵(包括四个阵元)、六元五边形阵(包括六个阵元)等,其中每个阵元接收一路次声信号。
步骤102:对多路所述待检测次声信号进行小波包分解和信号重构,生成不同阵元不同频段的重构信号,具体包括:
对多路所述待检测次声信号进行三层小波包分解,获得分解后的小波包分解系数。对每路所述待检测次声信号进行小波包分解的过程如下:
将所述待检测信号表示为s(t),对待检测信号s(t)进行三层小波包分解,分解结构如图2所示。图2为本发明提供的小波包三层分解树结构示意图,图2中(i,j)表示第i层小波包分解的第j个结点,每个结点都代表一定的信号特征。其中,(0,0)结点代表原始信号S,(1,0)代表小波包分解的第一层低频系数X10,(1,1)为小波包分解第一层的高频系数X11,(3,0)表示第三层第0个结点的系数,其它依此类推。
其中Xij表示次声信号s(t)进行小波包分解后,得到的第i层第j个节点(小波包)的系数。
对所述小波包分解系数Xij实施信号重构,生成不同阵元不同频段的重构信号,具体为:
对每个阵元检测到的每路次声信号s(t)分解后的8个小波包分解系数Xij实施信号重构,提取第三层各频带范围的重构信号。
以Sij表示小波包分解系数Xij的重构信号,例如S30表示X30的重构信号,S31表示X31的重构信号,其它依此类推。总信号可以表示为:
S=S30+S31+S32+S33+S34+S35+S36+S37 (9)
此处将分解后的小波包系数再重构信号相当于将信号分解成不同频率段信号的和。由于小波包本身具有降噪特性,对信号噪声具有抑制作用,因此能够降低环境噪声对后续Fisher属性值计算的影响,提高Fisher属性值计算精度,从而提高次声事件检测的效率和精度。
得到不同阵元的不同频段的重构信号Sij后,下一步要将不同阵元的相同频段的重构信号进行归纳。
步骤103:将所述不同阵元不同频段的重构信号中相同频段的重构信号进行归类,获得多个归类信号。
步骤104:采用滑动窗对多个所述归类信号进行分割,获得多个不同频段的窗内信号。
将相同频段的不同阵元(子台)信号进行归类,获得不同频段的多个归类信号,并以一滑动窗对所述归类信号进行分割,分割后得到不同频段的窗内信号。
步骤105:计算各个所述窗内信号的Fisher属性值。
采用公式(8)
Figure BDA0001637039830000111
计算各个所述窗内信号的Fisher属性值;其中F表示所述窗内信号的Fisher属性值;
Figure BDA0001637039830000112
N表示次声台阵的阵元数量,M为时间窗内采样点的数量,
Figure BDA0001637039830000113
xij表示第j个阵元在第i时刻接收到的次声信号,
Figure BDA0001637039830000114
由前文介绍的Fisher检测原理可知,当所述待检测次声信号是次声事件信号时,Fisher属性值比较大,因此可以准确检测出次声事件。
步骤106:将所述Fisher属性值作为像素值,生成二维Fisher属性分布图,具体包括:
以信号的时间长度为横轴,频率幅度为纵轴,窗长度范围内信号(窗内信号)的Fisher属性值为像素点,构成二维Fisher属性分布图,即以时间—频率构成二维Fisher属性分布图,所述二维Fisher属性分布图的横坐标为时间,纵坐标为频率。将得到的每一个所述Fisher属性值作为所述二维Fisher属性分布图中的一个像素点的像素值,一个Fisher属性值对应所述二维Fisher属性分布图中的一个像素点。
步骤107:对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇;每个所述像素簇中包括多个所述像素点。
在对所述二维Fisher属性分布图中所有像素点进行筛选合并处理时,可以先筛选后合并,也可以先合并后筛选,得到多个所述像素簇。
其中,对所述二维Fisher属性分布图中的像素点进行先筛选后合并处理的过程如下:
判断所述二维Fisher属性分布图中像素点的像素值是否高于筛选阈值,若否,将所述像素点的像素值归零;若是,判断所述像素点的像素值与其相邻像素点的像素值是否相同;若相同,将所述像素点与其相邻像素点合并为一个像素簇;若不相同则不合并。
其中,对所述二维Fisher属性分布图中的像素点进行先合并后筛选处理的过程如下:
判断所述二维Fisher属性分布图中像素点的像素值与其相邻像素点的像素值是否相同;若相同,将所述像素点与其相邻像素点合并为一个初选像素簇;若不同,则不合并;
判断所述初选像素簇中像素点的像素值是否高于筛选阈值,若否,将所述初选像素簇归零;若是,保留所述初选像素簇作为所述像素簇。
本发明提供的次声事件检测方法通过对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,有效排除了毛刺信号的干扰,从而降低了次声事件检测的虚警率,进一步提高了次声事件的检测效率和检测精度。
步骤108:判断所述像素簇中像素点的数量是否超过尺度阈值,获得第一判断结果。
步骤109:当所述第一判断结果为所述像素簇中像素点的数量超过尺度阈值,确定发生了次声事件。
步骤110:当所述第一判断结果为所述像素簇中像素点的数量未超过尺度阈值,确定未发生次声事件。
设定一像素簇连续长度标准,当像素簇中连续相邻的像素点数量达到该尺度阈值时,就判断所述待检测次声信号为次声事件信号,即表示检测到了次声事件。
次声波记录中的噪声不仅仅包含高频随机噪声,同时还有低频扰动,所以不可能用一个简单的滤波器就能起到去噪的效果效果,本发明基于小波包分解和Fisher联合检测提出的次声事件检测方法,由于小波包本身具有降噪特性,对信号噪声具有抑制作用,能够降低环境噪声对Fisher属性值计算准确度的影响;并且Fisher属性值比相关系数属性值变化梯度更加明显,能够更有利于次声信号的检测。因此发明提供的次声事件检测方法能够有效降低次声事件检测的漏警率和虚警率,提高次声事件检测的效率和准确度。
为验证本发明提出的一种次声事件检测方法的优劣,下面以不同信噪比条件下的次声仿真信号对该检测方法进行实验。
次声仿真信号检测实验的过程如下:
(1)构建一段次声仿真数据s(t)。图3为本发明实施例提供的四元三角阵时延对齐后的实际次声信号示意图。图4为本发明实施例提供的实际次声信号叠加白噪声后信噪比为0db的仿真信号示意图。图3和图4的横坐标为时间窗长度,纵坐标为信噪比。图3和图4中信号窗长度1024采样点,时间窗长度3000采样点,采样频率为1Hz。图3中的四路信号分别为四元三角阵中四个阵元时延对齐后的四个实际次声信号,图4中的四路信号分别为四路实际次声信号叠加白噪声后信噪比为0db的四路次声仿真信号。
(2)对每路次声仿真信号s(t)进行3层小波包分解,然后进行各频段的信号重构,得到8个重构信号。图5为1号阵元进行三层小波包分解后得到的重构信号示意图,横坐标为时间窗长度,纵坐标为信噪比。四个阵元接收的四路次声仿真信号共分解重构得到8×4个重构信号。
(3)将四个阵元的8×4个重构信号中相同频段的重构信号进行归类,得到多个归类信号,用一滑动窗(窗长50s,帧长10s)对归类信号进行分割,得到多个窗内信号。计算各个窗内信号的Fisher属性值,并以时间-频率构成二维Fisher属性分布图,其中每一个Fisher值作为一个像素点。
(4)设置筛选阈值α=8(经验阈值),低于所述筛选阈值的像素点进行归0处理;高于所述筛选阈值的像素点保留,并进行相同Fisher值的相邻像素点的合并,得到多个像素簇。图6为对像素点进行筛选合并处理后的二维分布图,图6横坐标为时间,纵坐标为频率。
(5)设置信号的尺度阈值β=5,当所述像素簇中连续相邻像素点的数量低于尺度阈值时进行滤除,最终得到如图7所示的像素簇经尺度阈值筛选后的二维分布图,图7横坐标为时间,纵坐标为频率。图7所示的像素簇中合并像素点的尺度高于所述尺度阈值,视为检测到了次声事件。
为了进一步阐述本发明基于小波包分解和Fisher联合检测提出的次声事件检测方法的有效性,我们设置相同的滑动窗长和检测步长,在不同的信噪比下对三种方法的检测结果进行对比分析,得到以下数据信息(其中原始信号的实际初至点为500,信号范围为500~1523;信号完整度是指接收信号的完整程度),如表1所示:
表1 Fisher、PMCC和本发明次声事件检测方法的基本参数及其实验结果
Figure BDA0001637039830000141
根据表1得到的实验结果可以得出以下结论:
1、根据Fisher理论,我们可以发现:当信噪比逐渐增大时,总体变化量也会相应增加,Fisher阈值作为一种衡量局部变化和总体变化的比值也会得到相应的增加。所以对于不同的信噪比,Fisher阈值是不断改变的。从表1中可以发现当我们选择一个固定Fisher阈值10进行检测时,在低信噪比条件下(信噪比为-7db、-5db、-3db时),信号完整度均低于90%,检测信号的完备性并不好,所以如果需要获得在低信噪比下的高精度检测必须实时地改变阈值。图8和图9分别为不同信噪比下,检测信号完整度达到93%的检测结果示意图。其中,图8为信噪比为-7dB时采用Fisher方法获得的检测结果。图8和图9的横坐标为时间,纵坐标为Fisher值。图8中的直线801表示所设Fisher阈值,图8中所设Fisher阈值为30。图9为信噪比为-1dB时采用Fisher方法获得的的检测结果,图9中的直线901表示所设Fisher阈值,图9中所设Fisher阈值为13。可见,Fisher检测方法在低信噪比(-7dB)下,如果想要获得和高信噪比(-1dB)下相同的检测精度,必须改变所设Fisher阈值。
2、PMCC方法虽然能够衡量阵列信号的一致性和相干性,但是实验表明相同检测精度下其所需时间窗长(80s)大于Fisher检测方法(50s),所以在相同的检测尺度(50s)下,PMCC方法的检测精度并不高。
3、本发明基于小波包分解和Fisher联合检测提出的次声事件检测方法,通过小波包分解,无论是高信噪比还是低信噪比,低频能量总是保持在第一个分解频段内,该频段信号经过重构,可以达到原始信号降噪的目的。所以在进行Fisher检测时,其噪声能量均值能够得到一定的抑制,从而可以保证Fisher检测阈值(8)即便不发生改变,也能获得较高的检测精度(采用本发明方法检测到的信号完整度均在90%以上)。
4、在相同的滑动窗长(50s)和步长(10s)下,本发明基于小波包分解和Fisher的次声信号检测方法相较于Fisher检测方法和PMCC方法来说,信号完整度最高,即次声事件信号的检测精度最高。
实际次声信号的检测与分析
下面以24组地面航天飞机再入大气层次声事件和10组噪声事件在6个台站的34个记录为基础数据库,应用Fisher、PMCC以及本文提出的次声事件检测方法分别进行事件检测。
其中34个地面真实次声事件的统计信息见表2。
表2地面真实事件目录
Figure BDA0001637039830000161
Figure BDA0001637039830000171
应用Fisher、PMCC以及本文提出的次声事件检测方法分别进行次声事件检测。图10、图11、图12、图13分别给出了2002年3月12日SGAR台阵记录的次声信号检测结果示意图。图10、图11、图12、图13的横坐标为时间窗长度,纵坐标为信噪比。其中图10为SGAR台阵接收的原始次声信号示意图。图11为采用Fisher检测方法对实测信号进行检测的检测结果示意图。图12为采用PMCC方法对实测信号进行检测的检测结果示意图。图13为采用本发明提供的次声事件检测方法对实测信号进行检测的检测结果示意图。
当对所有次声事件进行检测后,得到如下检测结果:
(1)Fisher方法共检测出了21个事件信号(其中真实事件18个,虚假事件3个),13组噪声事件。其中准确检测出18个事件,漏检6个事件,虚检3个事件信号。
(2)PMCC方法共检测出了23个事件信号,11组噪声事件。其中准确检测出19个事件,漏检5个事件,虚检4个事件信号。
(3)基于小波包和Fisher的联合检测方法共检测出了23个信号,11组噪声事件。其中准确检测出22个事件,漏检2个事件,虚检1个事件信号。
表3给出了三种方法的检测结果对比。
表3Fisher、PMCC、本发明方法的检测结果
检测方法 事件准确检测率 事件漏检率 事件虚检率
Fisher 75% 25% 12.5%
PMCC 79.2% 20.8% 16.7%
本文方法 91.7% 9.3% 4.2%
根据表3的检测结果可以看出,采用本发明提供的次声事件检测方法,有效降低了次声事件的漏检率和虚检率,大大提高了次声事件检测的准确率。
从仿真实验及实测次声信号的检测结果来看,本文提出的次声事件检测方法均具有较好的检测准确率,同时降低了虚警率和漏警率,适用于实际次声事件的高精度检测。
本发明还提供了一种次声事件检测系统。图14为本发明提供的一种次声事件检测系统的结构示意图。参见图14,所述次声事件检测系统包括:
待检测次声信号获取模块1401,用于获取次声台阵的多个阵元接收到的多路待检测次声信号;
重构信号生成模块1402,用于对多路所述待检测次声信号进行小波包分解和信号重构,生成不同阵元不同频段的重构信号;
归类信号获取模块1403,用于将所述不同阵元不同频段的重构信号中相同频段的重构信号进行归类,获得多个归类信号;
窗内信号获取模块1404,用于采用滑动窗对多个所述归类信号进行分割,获得多个不同频段的窗内信号;
Fisher属性值计算模块1405,用于计算各个所述窗内信号的Fisher属性值;
Fisher属性分布图生成模块1406,用于将所述Fisher属性值作为像素值,生成二维Fisher属性分布图;
像素簇获取模块1407,用于对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇;每个所述像素簇中包括多个所述像素点;
第一判断结果获取模块1408,用于判断所述像素簇中像素点的数量是否超过尺度阈值,获得第一判断结果;
次声事件发生确定模块1409,用于当所述第一判断结果为所述像素簇中像素点的数量超过尺度阈值时,确定发生了次声事件;
次声事件未发生确定模块1410,用于当所述第一判断结果为所述像素簇中像素点的数量未超过尺度阈值时,确定未发生次声事件。
其中,所述重构信号生成模块1402具体包括:
小波包分解单元,用于对多路所述待检测次声信号进行三层小波包分解,获得分解后的小波包分解系数;
信号重构单元,用于对所述小波包分解系数实施信号重构,生成不同阵元不同频段的重构信号。
所述Fisher属性值计算模块1405具体包括:
Fisher属性值计算单元,用于采用公式
Figure BDA0001637039830000191
计算各个所述窗内信号的Fisher属性值;其中F表示所述窗内信号的Fisher属性值;
Figure BDA0001637039830000192
N表示次声台阵的阵元数量,M为时间窗内采样点的数量,
Figure BDA0001637039830000193
xij表示第j个阵元在第i时刻接收到的次声信号,
Figure BDA0001637039830000194
Fisher属性分布图生成模块1406具体包括:
Fisher属性分布图生成单元,用于将多个所述Fisher属性值作为所述二维Fisher属性分布图中多个像素点的像素值,生成二维Fisher属性分布图;其中所述二维Fisher属性分布图的横坐标为时间,纵坐标为频率。
像素簇获取模块1407具体包括:
第二判断结果获取单元,用于判断所述二维Fisher属性分布图中像素点的像素值是否高于筛选阈值,获得第二判断结果;
像素点归零单元,用于当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值低于筛选阈值时,将所述像素点的像素值归零;
第三判断结果获取单元,用于当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值高于筛选阈值时,判断所述像素点的像素值与其相邻像素点的像素值是否相同,获得第三判断结果;
像素簇合并单元,用于当所述第三判断结果为所述像素点的像素值与其相邻像素点的像素值相同时,将所述像素点与其相邻像素点合并为一个像素簇。
第四判断结果获取单元,用于判断所述二维Fisher属性分布图中像素点的像素值与其相邻像素点的像素值是否相同,获得第四判断结果;
初选像素簇合并单元,用于当所述第四判断结果为所述像素点的像素值与其相邻像素点的像素值相同时,将所述像素点与其相邻像素点合并为一个初选像素簇;
第五判断结果获取单元,用于判断所述初选像素簇中像素点的像素值是否高于筛选阈值,获得第五判断结果;
像素簇归零单元,用于当所述第五判断结果为所述初选像素簇中像素点的像素值低于筛选阈值时,将所述初选像素簇归零;
像素簇生成单元,用于当所述第五判断结果为所述初选像素簇中像素点的像素值高于筛选阈值时,保留所述初选像素簇作为所述像素簇。
本发明提供的一种次声事件检测系统,采用小波包分解将待检测次声信号分解之后,在不同频段内进行重构,由于小波包本身具有降噪特性,对信号噪声具有抑制作用,因此能够降低环境噪声对Fisher属性值计算的影响;然后对每个频段内的重构信号进行Fisher检测,并且对多阵元同一时窗内的每个Fisher值作为一个属性保存,计算结果生成一幅时间-频率平面上的二维Fisher属性分布图,同时,考虑信号检测的相关时间尺度,以一定的尺度阈值作为衡量信号的标准,降低了次声事件检测的漏警率和虚警率,提高了次声事件检测的效率和精度。
本说明书中各个实施例采用递进的方式描述,每个实施例重点说明的都是与其他实施例的不同之处,各个实施例之间相同相似部分互相参见即可。对于实施例公开的系统而言,由于其与实施例公开的方法相对应,所以描述的比较简单,相关之处参见方法部分说明即可。
本文中应用了具体个例对本发明的原理及实施方式进行了阐述,以上实施例的说明只是用于帮助理解本发明的方法及其核心思想;同时,对于本领域的一般技术人员,依据本发明的思想,在具体实施方式及应用范围上均会有改变之处。综上所述,本说明书内容不应理解为对本发明的限制。

Claims (6)

1.一种次声事件检测方法,其特征在于,所述检测方法包括:
获取次声台阵的多个阵元接收到的多路待检测次声信号;
对多路所述待检测次声信号进行小波包分解和信号重构,生成不同阵元不同频段的重构信号;
将所述不同阵元不同频段的重构信号中相同频段的重构信号进行归类,获得多个归类信号;
采用滑动窗对多个所述归类信号进行分割,获得多个不同频段的窗内信号;
计算各个所述窗内信号的Fisher属性值;
将所述Fisher属性值作为像素值,生成二维Fisher属性分布图;
对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇;每个所述像素簇中包括多个所述像素点;
所述对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇,具体包括:
判断所述二维Fisher属性分布图中像素点的像素值是否高于筛选阈值,获得第二判断结果;
当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值低于筛选阈值,将所述像素点的像素值归零;
当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值高于筛选阈值,判断所述像素点的像素值与其相邻像素点的像素值是否相同,获得第三判断结果;
当所述第三判断结果为所述像素点的像素值与其相邻像素点的像素值相同,将所述像素点与其相邻像素点合并为一个像素簇;
当所述第三判断结果为所述像素点的像素值与其相邻像素点的像素值不同,则不合并;
判断所述像素簇中像素点的数量是否超过尺度阈值,获得第一判断结果;
当所述第一判断结果为所述像素簇中像素点的数量超过尺度阈值,确定发生了次声事件;
当所述第一判断结果为所述像素簇中像素点的数量未超过尺度阈值,确定未发生次声事件。
2.根据权利要求1所述的次声事件检测方法,其特征在于,所述对多路所述待检测次声信号进行小波包分解和信号重构,生成不同阵元不同频段的重构信号,具体包括:
对多路所述待检测次声信号进行三层小波包分解,获得分解后的小波包分解系数;
对所述小波包分解系数实施信号重构,生成不同阵元不同频段的重构信号。
3.根据权利要求2所述的次声事件检测方法,其特征在于,所述计算各个所述窗内信号的Fisher属性值,具体包括:
采用公式
Figure FDA0002466805830000021
计算各个所述窗内信号的Fisher属性值;其中F表示所述窗内信号的Fisher属性值;
Figure FDA0002466805830000022
N表示次声台阵的阵元数量,M为时间窗内采样点的数量,
Figure FDA0002466805830000023
xij表示第j个阵元在第i时刻接收到的次声信号,
Figure FDA0002466805830000024
4.根据权利要求3所述的次声事件检测方法,其特征在于,所述将所述Fisher属性值作为像素值,生成二维Fisher属性分布图,具体包括:
将多个所述Fisher属性值作为所述二维Fisher属性分布图中多个像素点的像素值,生成二维Fisher属性分布图;其中所述二维Fisher属性分布图的横坐标为时间,纵坐标为频率。
5.一种次声事件检测系统,其特征在于,所述检测系统包括:
待检测次声信号获取模块,用于获取次声台阵的多个阵元接收到的多路待检测次声信号;
重构信号生成模块,用于对多路所述待检测次声信号进行小波包分解和信号重构,生成不同阵元不同频段的重构信号;
归类信号获取模块,用于将所述不同阵元不同频段的重构信号中相同频段的重构信号进行归类,获得多个归类信号;
窗内信号获取模块,用于采用滑动窗对多个所述归类信号进行分割,获得多个不同频段的窗内信号;
Fisher属性值计算模块,用于计算各个所述窗内信号的Fisher属性值;
Fisher属性分布图生成模块,用于将所述Fisher属性值作为像素值,生成二维Fisher属性分布图;
像素簇获取模块,用于对所述二维Fisher属性分布图中所有像素点进行筛选合并处理,获得多个像素簇;每个所述像素簇中包括多个所述像素点;
所述像素簇获取模块具体包括:
第二判断结果获取单元,用于判断所述二维Fisher属性分布图中像素点的像素值是否高于筛选阈值,获得第二判断结果;
像素点归零单元,用于当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值低于筛选阈值时,将所述像素点的像素值归零;
第三判断结果获取单元,用于当所述第二判断结果为所述二维Fisher属性分布图中像素点的像素值高于筛选阈值时,判断所述像素点的像素值与其相邻像素点的像素值是否相同,获得第三判断结果;
像素簇合并单元,用于当所述第三判断结果为所述像素点的像素值与其相邻像素点的像素值相同时,将所述像素点与其相邻像素点合并为一个像素簇;
第一判断结果获取模块,用于判断所述像素簇中像素点的数量是否超过尺度阈值,获得第一判断结果;
次声事件发生确定模块,用于当所述第一判断结果为所述像素簇中像素点的数量超过尺度阈值时,确定发生了次声事件;
次声事件未发生确定模块,用于当所述第一判断结果为所述像素簇中像素点的数量未超过尺度阈值时,确定未发生次声事件。
6.根据权利要求5所述的次声事件检测系统,其特征在于,所述Fisher属性分布图生成模块具体包括:
Fisher属性分布图生成单元,用于将多个所述Fisher属性值作为所述二维Fisher属性分布图中多个像素点的像素值,生成二维Fisher属性分布图;其中所述二维Fisher属性分布图的横坐标为时间,纵坐标为频率。
CN201810365902.3A 2018-04-23 2018-04-23 一种次声事件检测方法及系统 Active CN108844616B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810365902.3A CN108844616B (zh) 2018-04-23 2018-04-23 一种次声事件检测方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810365902.3A CN108844616B (zh) 2018-04-23 2018-04-23 一种次声事件检测方法及系统

Publications (2)

Publication Number Publication Date
CN108844616A CN108844616A (zh) 2018-11-20
CN108844616B true CN108844616B (zh) 2020-06-30

Family

ID=64212252

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810365902.3A Active CN108844616B (zh) 2018-04-23 2018-04-23 一种次声事件检测方法及系统

Country Status (1)

Country Link
CN (1) CN108844616B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113657464B (zh) * 2021-07-29 2024-04-19 中国电子科技集团公司第三研究所 一种降低pmcc次声检测虚警事件数量的方法及存储介质

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102213611A (zh) * 2010-04-09 2011-10-12 中国科学院微电子研究所 一种探测次声波的方法
CN103076177A (zh) * 2013-01-16 2013-05-01 昆明理工大学 一种基于振动检测的滚动轴承故障检测方法
CN103369446A (zh) * 2012-03-28 2013-10-23 中国科学院声学研究所 一种基于pmcc算法的次声源定向方法
CN103941224B (zh) * 2013-01-21 2016-08-03 中国科学院声学研究所 一种次声源的定向定速方法及系统
CN107272061A (zh) * 2017-06-29 2017-10-20 禁核试北京国家数据中心 一种次声信号与地震事件的自动关联方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7079449B2 (en) * 2003-02-18 2006-07-18 Batelle Energy Alliance, Llc Methods and systems for low frequency seismic and infrasound detection of geo-pressure transition zones

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102213611A (zh) * 2010-04-09 2011-10-12 中国科学院微电子研究所 一种探测次声波的方法
CN103369446A (zh) * 2012-03-28 2013-10-23 中国科学院声学研究所 一种基于pmcc算法的次声源定向方法
CN103076177A (zh) * 2013-01-16 2013-05-01 昆明理工大学 一种基于振动检测的滚动轴承故障检测方法
CN103941224B (zh) * 2013-01-21 2016-08-03 中国科学院声学研究所 一种次声源的定向定速方法及系统
CN107272061A (zh) * 2017-06-29 2017-10-20 禁核试北京国家数据中心 一种次声信号与地震事件的自动关联方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Regional monitoring of infrasound events using multiple arrays: application to Utah and Washington State;Stephen J A , et al;《Geophysical Journal International》;20081002;第175卷(第1期);第292页左栏第28行-第295页左栏第22行 *
基于小波包分解和FK分析的次声信号参数估计;孟亮 等;《应用声学》;20160331;第35卷(第2期);第152页右栏-第155页右栏第6行 *

Also Published As

Publication number Publication date
CN108844616A (zh) 2018-11-20

Similar Documents

Publication Publication Date Title
Mousavi et al. Hybrid seismic denoising using higher‐order statistics and improved wavelet block thresholding
Wang et al. Multiwavelet denoising with improved neighboring coefficients for application on rolling bearing fault diagnosis
Küperkoch et al. Automated determination of P-phase arrival times at regional and local distances using higher order statistics
Wang et al. Maximum envelope-based autogram and symplectic geometry mode decomposition based gear fault diagnosis method
Kurzon et al. Real‐time automatic detectors of P and S waves using singular value decomposition
CN109340586B (zh) 一种供水管道泄露的检测方法及系统
Dodge et al. Initial global seismic cross‐correlation results: Implications for empirical signal detectors
CN110717472B (zh) 基于改进的小波阈值去噪的故障诊断方法及系统
Akhouayri et al. Automatic detection and picking of P-wave arrival in locally stationary noise using cross-correlation
Papanicolaou et al. Wavelet based estimation of local Kolmogorov turbulence
Zurbenko et al. Kolmogorov–Zurbenko filters in spatiotemporal analysis
CN108844616B (zh) 一种次声事件检测方法及系统
CN112539887A (zh) 一种基于wt-lcd-wd的管道泄漏信号去噪方法
Zou et al. Toward accurate extraction of bearing fault modulation characteristics with novel time–frequency modulation bispectrum and modulation Gini index analysis
CN117908120A (zh) 基于高分辨率时频分析的地震信号能量统计初至检测算法
Wu et al. Statistical significance test of intrinsic mode functions
CN104200815B (zh) 一种基于相关分析的音频噪声实时检测方法
CN116304735A (zh) 基于相依性度量的钢轨表面损伤诊断方法、系统及设备
CN114090949A (zh) 一种鲁棒的振动信号特征值计算方法
Yang et al. Fractal Slope‐Based Seismic Wave Detection Method
CN113050167B (zh) 地震数据的噪声压制方法及装置
Neocleous et al. Identification of possible Δ14C anomalies since 14 ka BP: A computational intelligence approach
Eldwaik et al. Microphone wind noise reduction using singular spectrum analysis techniques
CN112444877B (zh) 一种基于标准道的互相关信噪比的估算方法
El Mountassir et al. Detection of structural damage using novelty detection algorithm under variational environmental and operational conditions

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
GR01 Patent grant
GR01 Patent grant