CN113607268B - 一种区域次声事件自动关联方法 - Google Patents

一种区域次声事件自动关联方法 Download PDF

Info

Publication number
CN113607268B
CN113607268B CN202110102318.0A CN202110102318A CN113607268B CN 113607268 B CN113607268 B CN 113607268B CN 202110102318 A CN202110102318 A CN 202110102318A CN 113607268 B CN113607268 B CN 113607268B
Authority
CN
China
Prior art keywords
signal
infrasound
time
grid
signals
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
CN202110102318.0A
Other languages
English (en)
Other versions
CN113607268A (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.)
Ctbt Beijing National Data Centre
Original Assignee
Ctbt Beijing National Data Centre
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 Ctbt Beijing National Data Centre filed Critical Ctbt Beijing National Data Centre
Priority to CN202110102318.0A priority Critical patent/CN113607268B/zh
Publication of CN113607268A publication Critical patent/CN113607268A/zh
Application granted granted Critical
Publication of CN113607268B publication Critical patent/CN113607268B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

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
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Alarm Systems (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

本发明涉及一种区域次声事件自动关联方法,根据设定的次声监测区域,组建专用次声监测台网,通过建立符合设定区域事件特征的信号特征模型和信号走时模型,设计并创建包含信号走时、方位角、震中距、单位距离走时偏差等参数的区域次声事件关联格点,采用区域细化格点搜索方法,准确检索并匹配出与事件相关的全部次声信号。实现了利用国际监测系统次声台站监测数据及其信号检测结果,对设定区域次声事件进行自动关联处理,并形成与事件相匹配的次声信号集,对于观测大气层及地面次声事件具有重要价值,实现了区域次声事件的准确自动匹配。

Description

一种区域次声事件自动关联方法
技术领域
本发明属于次声监测领域,具体涉及一种区域次声事件自动关联方法。
背景技术
次声监测技术是有效监测大气层和地面次声事件的有效技术手段,全面禁止核试验条约国际监测系统计划建立包括60个次声台站的全球监测网络,目前已经完成52个次声台站核证,并向国际数据中心实时传输监测数据。次声台站关联事件已经成为全面禁止核试验条约组织审核事件公报的重要组成部分,公报收录次声台网关联事件3200余个,占次声台站参与关联事件的15%。次声台网关联事件在南纬60゜以北的广大地区几乎均有分布,但主要集中欧亚地区和太平洋东西海岸,具有明显的地域分布特征。
全面禁止核试验条约组织建立的国际次声监测台网外,美国建立了一定规模的次声监测网络,其他国家或组织建设和运行的次声台站数量较少,很少能构成多台共同监测的次声台网。因此国内外对于大规模次声台网关联方法研究较少,主要关联方法为全球均匀格点搜索法,该方法的针对性不强、关联形成的虚假事件较多。
发明内容
本发明的目的是提供一种区域次声事件自动关联方法,实现自动关联结果准确率高、漏检率和误检率较低、次声信号自动筛选匹配快速,不仅可以具备次声台网的自动关联与定位能力,还能提高区域重复事件的精准识别能力。
本发明的技术解决方案是:一种区域次声事件自动关联方法,其特征在于,包括以下步骤:
S1、选定监测区域,选择对该区域具有监测能力的次声台站组成专用次声监测台网;
S2、根据选定监测区域已有的历史事件,分别确定上述次声台站历次事件所关联的次声信号,提取包括信号的频率特征、能量强度特征、相关性特征在内的信号特征参数,并建立该区域次声信号特征模型;根据历史事件关联的次声信号,确定上述信号特征参数的最低阈值;
S3、针对监测区域创建格点并计算格点参数文件;
S4、针对给定的时间段,获取专用次声监测台网所属次声台站的检测信号及其特征参数,利用步骤2中根据历史事件次声信号建立的特征模型,逐个台站对检测信号进行识别筛选,通过一次筛选,剔除不符合目标监测区域次声信号特征的检测信号;对于一次筛选剩余的信号利用格点参数文件中定义的目标监测区域次声信号的方位角特征,进行二次筛选,从而得到符合该区域次声事件信号特征的最小检测信号集合;
S5、从专用次声监测台网中选择其中一个台站作为触发台站,其记录的符合区域信号特征的筛选信号作为触发信号S(i),先将触发信号与格点进行方位角匹配,筛选出与触发信号方位相匹配的格点,并作标识;
S6、从除触发台站以外的台站中,选择一个待关联台站,将其记录的符合区域信号特征的筛选信号作为待关联信号D(j),根据信号的实际到达时刻,计算触发信号S(i)与待关联信号D(j)的到时差,并按照触发台站和待关联台站至格点中心的实际距离、设定的次声信号传播速度,计算预期速度条件下,次声信号由格点传输至两个不同次声台站的信号走时差;
S7、对于触发台站和待关联台站,分别根据其震中距格点参数、设定的单位震中距最小传输走时time_delta_min和最大传输走时time_delta_max,分别计算触发台站和待关联台站同时记录到次声信号的情况下,信号走时差的最大值tt_res_max和最小值tt_res_min;
S8、对触发信号和待关联信号进行信号走时匹配分析,当满足以下条件时认为触发信号与待关联信号匹配:
其中,time_th为设定的信号走时偏差阈值,az_gridsec为待关联台站的方位角格点;
S9、对第一个待关联台站的所有待关联信号分别与触发信号进行匹配,筛选出所有满足公式(7)的匹配信号;而后再对下一待关联台站进行上述相同处理,直至所有待关联台站均与触发信号进行完匹配。
S1中,专用次声监测台网通常由3至4个次声台站组成,监测区域原则上小于10°×10°。
S3中,设定格点步长为grid_range,计算出整个区域的格点数N×M,格点文件结构为{STA台站,az_residual_min最小方位角偏差,az_residual_max最大方位角偏差,az_residual_th方位角偏差阈值,time_residual_min最小走时偏差,time_residual_max最大走时偏差,time_residual_th走时偏差阈值,time_delta_min最小单位距离走时偏差,time_delta_max最大单位距离走时偏差,time_delta_th单位距离走时偏差阈值,grid_azimuth方位角格点参数,grid_dist震中距格点参数},
其中,az_residual_th=max(abs(az_residual_max),abs(az_residual_min)),
time_residual_th=max(abs(time_residual_max),abs(time_residual_min)),time_delta_th=max(abs(time_delta_max),abs(time_delta_min)),grid_azimuth与grid_dist均为N×M参数文件,分别记录了台站至每个格点的方位角与震中距。
S4中,次声信号的方位角特征应满足:
min(grid_azimuth)-az_th≤Sig_az≤max(grid_azimuth)+az_th (1)
其中,grid_azimuth为格点方位角参数,az_th为设定的方位角偏差阈值,Sig_az为提取的检测信号的方位角参数。
S4中,当所有次声台站均含有符合区域次声事件信号特征要求的检测信号时,继续进行后续关联,否则结束整个关联进程。
S5中,对于触发信号S(i),其检测方位角为S_az(i),信号到时为S_time(i),对于N×M的特定区域方位角格点文件grid_azmiuth,将触发信号的检测方位角依次与格点文件中定义的方位角值进行匹配,当满足:
|S_az(i)-grid_azimuth(n,m)|≤az_th (2)
则认为触发信号S(i)能与格点匹配,其中i为时序信号序列,信号总数为I,1≤i≤I,1≤n≤N、1≤m≤M,az_th为设定的方位角偏差阈值,通过公式(2)从格点文件中筛选出与触发信号S(i)相匹配的格点grid_sel,筛选的格点分布在S_az(i)方位两侧,筛选出的格点总数为K。
S6中,对于待关联信号D(j),其信号总数为J,j为时序信号序列,1≤j≤J,检测方位角为D_az(j),信号到时为D_time(j);触发信号S(i)与待关联信号D(j)的到时差为:
time_arr_diff=S_time(i)-D_time(j) (3)
其中,grid_distfirst触发台站对应的区域格点震中距参数,grid_distsec为待关联台站对应的区域格点震中距参数,speed为设定的次声信号传播速度,time_travel_diff表示次声信号由格点传输至两个不同次声台站的信号走时差。
S7中,
其中,time_delta_minsec为待关联台站单位距离信号传输最小单位距离走时偏差,time_delta_maxsec为待关联台站单位距离信号传输最大走时偏差,time_delta_minfirst为触发台站单位距离信号传输最小走时偏差,time_delta_maxfirst为触发台站单位距离信号最大走时偏差。
S9中,当且仅当次声网中所有次声台站均有与区域内格点相匹配的次声信号时,形成临时关联事件,若含有匹配信号的台站数少于台网中次声台站总数,则认为无足够匹配台站,针对给定的监测区域无法形成有效的信号关联。
还包括S10:当有多个自动关联形成的临时关联事件时,需对临时关联事件进行重复关联事件筛选,挑选出触发信号不同、其余台站关联信号基本一致的事件,并进行信号合并,具体判别过程如下:假设临时关联事件共关联A个信号,除触发信号外关联的信号为A-1个,已关联事件共关联B个信号,除触发信号外关联的信号为B-1个,当满足(A-1)/(B-1)≥N_th时,将临时关联事件与已有的事件进行关联信号合并,并剔除重复信号;否则将临时关联事件作为一个独立事件而保存,其中N_th为预先设定的百分比阈值。
有益效果:本发明可根据设定的次声监测区域,组建专用次声监测台网,通过建立符合设定区域事件特征的信号特征模型和信号走时模型,设计并创建包含信号走时、方位角、震中距、单位距离走时偏差等参数的区域次声事件关联格点,采用区域细化格点搜索方法,准确检索并匹配出与事件相关的全部次声信号。实现了利用国际监测系统次声台站监测数据及其信号检测结果,对设定区域次声事件进行自动关联处理,并形成与事件相匹配的次声信号集,对于观测大气层及地面次声事件具有重要价值,具体优点如下:
1、本发明采用历史次声事件建立信号特征模型的方法,将实测信号参数提炼成具有可操作性的特征模型,有效解决了次声检测信号筛选,从而筛选出最符合特定区域信号特征的检测信号,提高了关联过程的针对性。同时设计了一种区域次声自动关联的格点,将关键参数在格点中预设定义,并根据不同台站创建特有的关联参数,提高了信号匹配搜索的效率。通过触发信号筛选出与触发信号方位相匹配的格点,而后将这些筛选出的格点与待关联信号进行方位角匹配、触发信号与待关联信号走时匹配,从而确定属于某次声事件的全部次声检测信号,实现了区域次声事件的准确自动匹配。
2、本发明实现了利用国家监测系统次声台站监测数据对区域次声事件进行精准自动监测,解决了区域次声信号特征模型建立,匹配次声信号筛选等技术难题,对于开展大气层次声事件的自动、快速、准确监测具有重要价值。测试结果表明,所述方法对区域次声事件的整体自动识别率达到90%以上,且保持了较低的误识别率和漏识别率。
附图说明
图1区域次声事件自动关联流程
具体实施方式
以下结合附图表对本发明做详细说明:
一种区域次声事件自动关联方法,其特征在于,包括以下步骤:
步骤1、首先,根据检测需求设定目标监测区域和次声监测台网。针对某次声事件设定监测区域,范围为北纬53゜~58゜,东经40゜~45゜,选定的次声监测台站分别为I31KZ、I43RU、I46RU,收集上述3个次声台站记录的选定区域内的历史事件。
步骤2、根据选定监测区域已有的历史事件,确定各次声台站记录历史事件次声信号的特征模型,主要包括信号的频率特征(fcfreq中心频率、fband信号频带宽度)、能量强度特征(SNR信噪比、Amp信号幅值)和相关性特征(Corr互相关系数),根据历史事件关联的次声信号,确定上述信号特征参数的最低阈值;具体如下:
图1为区域次声事件自动关联流程,根据历史事件分别,分别确定上述次声台站历次事件所关联的次声信号,提取包括频率特征、能量特征、相关性特征在内的信号特征参数,并建立该区域次声信号特征模型和信号走时模型;
步骤3、针对监测区域创建格点并计算格点参数文件:设定格点步长为grid_range,计算出整个区域的格点数N×M,格点文件结构为{STA台站,az_residual_min最小方位角偏差,az_residual_max最大方位角偏差,az_residual_th方位角偏差阈值,time_residual_min最小走时偏差,time_residual_max最大走时偏差,time_residual_th走时偏差阈值,time_delta_min最小单位距离走时偏差,time_delta_max最大单位距离走时偏差,time_delta_th单位距离走时偏差阈值,grid_azimuth方位角格点参数,grid_dist震中距格点参数},
其中,az_residual_th=max(abs(az_residual_max),abs(az_residual_min)),
time_residual_th=max(abs(time_residual_max),abs(time_residual_min)),
time_delta_th=max(abs(time_delta_max),abs(time_delta_min)),grid_azimuth与grid_dist均为N×M参数文件,分别记录了台站至每个格点的方位角与震中距;
具体格式如表1所示,其包含的方位角残差阈值、走时残差阈值、单位距离走时残差阈值、方位角格点参数、震中距格点参数等对于信号的自动关联尤为重要,方位角格点参数格式如表2所示,震中距格点参数格式如表3所示。
步骤4、根据监测时间段,提取出指定监测台网内全部次声台站的检测信号,并确保每个台站均有信号,否则结束关联进程。利用已建立的次声信号特征模型对提取的次声检测信号进行一次筛选,剔除不符合特征要求的检测结果。而后根据方位角格点参数文件,对剩余次声检测信号进行二次筛选,从而得到符合该区域次声事件信号特征的最小检测信号集合,并用于后续关联程序。具体如下:
次声信号的方位角特征应满足:
min(grid_azimuth)-az_th≤Sig_az≤max(grid_azimuth)+az_th (1)
其中,grid_azimuth为格点方位角参数,az_th为设定的方位角偏差阈值,Sig_az为提取的检测信号的方位角参数。
当所有次声台站均含有符合区域次声事件信号特征要求的检测信号时,继续进行后续关联,否则结束整个关联进程;
步骤5至步骤8为关联的主要过程:
步骤5、首先选择其中一个次声台站作为关联过程的触发台站,其对应的次声信号为触发信号,先将触发信号与格点进行方位角匹配,参考公式2,筛选出与触发信号方位相匹配的格点,并作标识。
对于触发信号S(i),其检测方位角为S_az(i),信号到时为S_time(i),对于N×M的特定区域方位角格点文件grid_azmiuth,将触发信号的检测方位角依次与格点文件中定义的方位角值进行匹配,当满足:
|S_az(i)-grid_azimuth(n,m)|≤az_th (2)
则认为触发信号S(i)能与格点匹配,其中i为时序信号序列,信号总数为I,1≤i≤I,1≤n≤N、1≤m≤M,az_th为设定的方位角偏差阈值,通过公式(2)从格点文件中筛选出与触发信号S(i)相匹配的格点grid_sel,筛选的格点分布在S_az(i)方位两侧,筛选出的格点总数为K。
步骤6、将除触发台站以外的次声台站作为待关联台站,选择一个待关联台站,其对应的次声信号为待关联信号,根据信号的实际到达时刻,计算触发信号S(i)与待关联信号D(j)的到时差time_arr_diff,如公式3所示,同时根据公式4按照触发台站和待关联台站至格点中心的实际距离、设定的次声信号传播速度,计算预期速度条件下,次声信号由格点传输至两个不同次声台站的信号走时差time_travel_diff。
对于待关联信号D(j),其信号总数为J,j为时序信号序列,1≤j≤J,检测方位角为D_az(j),信号到时为D_time(j);触发信号S(i)与待关联信号D(j)的到时差为:
time_arr_diff=S_time(i)-D_time(j) (3)
其中,grid_distfirst触发台站对应的区域格点震中距参数,grid_distsec为待关联台站对应的区域格点震中距参数,speed为设定的次声信号传播速度,time_travel_diff表示次声信号由格点传输至两个不同次声台站的信号走时差;
步骤7、对于触发台站和待关联台站,分别根据其震中距格点参数、设定的单位震中距最小传输走时time_delta_min和最大传输走时time_delta_max,分别计算触发台站和待关联台站同时记录到次声信号的情况下,信号走时差的最大值tt_res_max和最小值tt_res_min。
tt_res_min=time_delta_minsec*grid_distsec(grid_sel(k,1),grid_sel(k,2))-time_delta_maxfirst*grid_distfirst(grid_sel(k,1),grid_sel(k,2)) (5)
tt_res_max=time_delta_maxsec*grid_distsec(grid_sel(k,1),grid_sel(k,2))-time_delta_minfirst*grid_distfirst(grid_sel(k,1),grid_sel(k,2)) (6)
其中,time_delta_minsec为待关联台站单位距离信号传输最小单位距离走时偏差,time_delta_maxsec为待关联台站单位距离信号传输最大走时偏差,time_delta_minfirst为触发台站单位距离信号传输最小走时偏差,time_delta_maxfirst为触发台站单位距离信号最大走时偏差。
步骤8、对触发信号和待关联信号进行信号走时匹配分析,当满足公式7时认为触发信号与待关联信号匹配。
其中,time_th为设定的信号走时偏差阈值,az_gridsec为待关联台站的方位角格点;
步骤7至步骤8是重点利用信号走时关系进行信号关联匹配。
步骤9、对第一个待关联台站的所有待关联信号分别与触发信号进行匹配,筛选出所有满足公式(7)的匹配信号;而后再对下一待关联台站进行上述相同处理,直至所有待关联台站均与触发信号进行完匹配;
若所有待关联台站均有1个以上信号能与触发信号匹配,则认为包括触发信号在内的所有匹配信号均来自某一匹配格点所在区域的事件,从而形成临时关联事件。若所有待关联台站均有1个以上信号能与触发信号匹配,则认为包括触发信号在内的所有匹配信号均来自某一匹配格点所在区域的事件,从而形成临时关联事件。
若含有匹配信号的台站数少于台网中次声台站总数,则认为无足够匹配台站,针对给定的监测区域无法形成有效的信号关联;对于形成的临时关联事件参考步骤10进行筛选合并。
步骤10、当有多个自动关联形成的临时关联事件时,需对临时关联事件进行重复关联事件筛选,挑选出触发信号不同、其余台站关联信号基本一致的事件,并进行信号合并,具体判别过程如下:假设临时关联事件共关联A个信号,除触发信号外关联的信号为A-1个,已关联事件共关联B个信号,除触发信号外关联的信号为B-1个,当满足(A-1)/(B-1)≥N_th时,将临时关联事件与已有的事件进行关联信号合并,并剔除重复信号;否则将临时关联事件作为一个独立事件而保存,其中N_th为预先设定的百分比阈值。
实施实例:
利用全面禁止核试验条约国际监测系统记录的2020年10月7日I31KZ、I43RU、I46RU三个次声台站的数据,利用本方法对当日某次声事件进行自动关联测试,自动关联结果与国际数据中心公报记录如表1所示,对于国际数据中心记录的7次事件,本方法全部实现了自动检测,所有7次事件中,除发生于11:26:09的事件,自动处理未能实现I43RU台站公报信号匹配外,其余所有公报中记录的信号在自动处理结果中均有匹配,本方法对2020年10月7日某次声事件国际数据中心公报结果的100%自动关联。本专利方法对10月7日当日数据测试,除与国际数据中心匹配的7个事件外,另外形成了3个位于图1设定区域的事件,经验证,额外关联形成的3个事件均属于真实事件。
表1格点文件参数示例
台站 I43RU I31KZ I46RU
最小方位角残差 -1.66 -2.75 -5.42
最大方位角残差 3.39 9.65 6.24
方位角残差阈值 3.39 9.65 6.24
最小走时残差 -94.18 -492.58 -287.18
最大走时残差 98.47 349.61 474.70
走时残差阈值 98.47 492.58 474.70
最小单位距离走时残差 -30.00 -46.32 -11.84
最大单位距离走时残差 24.55 34.82 20.27
单位距离走时残差阈值 30.00 46.32 20.27
方位角格点参数 N*M数组 N*M数组 N*M数组
震中距格点参数 N*M数组 N*M数组 N*M数组
区域方位角最小值 48.40 290.06 284.03
区域方位角最大值 155.57 320.04 297.59
区域震中距最小值 1.53 8.47 22.31
区域震中距最大值 5.82 12.93 26.24
表2 I43RU次声台站选定区域内部分格点的方位角参数
表3特定次声台站选定区域内部分格点的震中距参数
表4某次声事件自动关联结果比较
*注:斜体标记信号为与REB审核公报结果相比,自动处理多关联的信号;粗体标记信号为与REB审核公报匹配的信号。
所述次声台站为全面禁止核试验条约组织所建国际监测系统次声台站,泛指各类安装微气压计、微麦克风等用于记录大气压力扰动的监测站点;次声数据是指次声传感器记录的数据;所述大气层或地表次声事件是指被多个次声台站同时记录的、源项位于大气层或地表的爆炸;信号走时是指次声信号由源项传输至监测台站所需时间;信号走时偏差是指次声信号由源项传输至监测台站的实际时长与理论时长的差值;信号方位角是指台站至源项连线与正北方向的夹角;信号方位角偏差是指信号实际方位角与理论方位角的差值;格点步长是指相邻两个格点的间隔距离;信号关联是指与某事件相关的所有次声台站记录的所有次声信号自动搜素匹配过程;次声审核事件是指由波形分析专家审核提交的真实次声事件;慢度是指信号沿地球表面传输1°所需时间。
所述次声台网关联格点是指本发明所定义的格点。该格点包含了次声台站、最小方位角偏差、最大方位角偏差、方位角偏差阈值、最小走时偏差、最大走时偏差、走时偏差阈值、最小单位距离走时偏差、最大单位距离走时偏差、单位距离走时偏差阈值、位角格点参数、震中距格点参数等信息;
区域细化格点是指将设定区域按照一定格点步长进行细化划分而形成的一系列小区域;
单位距离走时偏差是指信号传输1゜距离时,次声信号由源项传输至监测台站的实际时长与理论时长的差值;
最小单位距离走时偏差是指所有格点中单位距离走时偏差的最小值;
最大单位距离走时偏差是指所有格点中单位距离走时偏差的最大值;
单位距离走时偏差阈值是指信号传输1゜距离时,次声信号由源项传输至监测台站的实际时长与理论时长差值的最大允许值;
最小走时偏差是指所有格点中信号走时偏差的最小值;
最大走时偏差是指所有格点中信号走时偏差的最大值;;
信号走时偏差阈值是次声信号由源项传输至监测台站的实际时长与理论时长的差值的最大允许值;
最小方位角偏差是指所有格点中方位角偏差的最小值;
最大方位角偏差是指所有格点中方位角偏差的最大值;
方位角偏差阈值是指信号实际方位角与理论方位角差异的最大允许值。

Claims (6)

1.一种区域次声事件自动关联方法,其特征在于,包括以下步骤:
S1、选定监测区域,选择对该区域具有监测能力的次声台站组成专用次声监测台网;
S2、根据选定监测区域已有的历史事件,分别确定上述次声台站历次事件所关联的次声信号,提取包括信号的频率特征、能量强度特征、相关性特征在内的信号特征参数,并建立该区域次声信号特征模型;根据历史事件关联的次声信号,确定上述信号特征参数的最低阈值;
S3、针对监测区域创建格点并计算格点参数文件;设定格点步长为grid_range,计算出整个区域的格点数N×M,格点文件结构为{STA台站,az_residual_min最小方位角偏差,az_residual_max最大方位角偏差,az_residual_th方位角偏差阈值,time_residual_min最小走时偏差,time_residual_max最大走时偏差,time_residual_th走时偏差阈值,time_delta_min最小单位距离走时偏差,time_delta_max最大单位距离走时偏差,time_delta_th单位距离走时偏差阈值,grid_azimuth方位角格点参数,grid_dist震中距格点参数},
其中,az_residual_th=max(abs(az_residual_max),abs(az_residual_min)),
time_residual_th=max(abs(time_residual_max),abs(time_residual_min)),
time_delta_th=max(abs(time_delta_max),abs(time_delta_min)),grid_azimuth与grid_dist均为N×M参数文件,分别记录了台站至每个格点的方位角与震中距;
S4、针对给定的时间段,获取专用次声监测台网所属次声台站的检测信号及其特征参数,利用步骤2中根据历史事件次声信号建立的特征模型,逐个台站对检测信号进行识别筛选,通过一次筛选,剔除不符合目标监测区域次声信号特征的检测信号;对于一次筛选剩余的信号利用格点参数文件中定义的目标监测区域次声信号的方位角特征,进行二次筛选,从而得到符合该区域次声事件信号特征的最小检测信号集合;次声信号的方位角特征应满足:
min(grid_azimuth)-az_th≤Sig_az≤max(grid_azimuth)+az_th (1)
其中,grid_azimuth为格点方位角参数,az_th为设定的方位角偏差阈值,Sig_az为提取的检测信号的方位角参数;
S5、从专用次声监测台网中选择其中一个台站作为触发台站,其记录的符合区域信号特征的筛选信号作为触发信号S(i),先将触发信号与格点进行方位角匹配,筛选出与触发信号方位相匹配的格点,并作标识;对于触发信号S(i),其检测方位角为S_az(i),信号到时为S_time(i),对于N×M的特定区域方位角格点文件grid_azmiuth,将触发信号的检测方位角依次与格点文件中定义的方位角值进行匹配,当满足:
|S_az(i)-grid_azimuth(n,m)|≤az_th
则认为触发信号S(i)能与格点匹配,其中i为时序信号序列,信号总数为I,1≤i≤I,1≤n≤N、1≤m≤M,az_th为设定的方位角偏差阈值,通过上述公式从格点文件中筛选出与触发信号S(i)相匹配的格点grid_sel,筛选的格点分布在S_az(i)方位两侧,筛选出的格点总数为K;
S6、从除触发台站以外的台站中,选择一个待关联台站,将其记录的符合区域信号特征的筛选信号作为待关联信号D(j),根据信号的实际到达时刻,计算触发信号S(i)与待关联信号D(j)的到时差time_arr_diff,并按照触发台站和待关联台站至格点中心的实际距离、设定的次声信号传播速度,计算预期速度条件下,次声信号由格点传输至两个不同次声台站的信号走时差time_travel_diff;对于待关联信号D(j),其中j为待关联信号序列,1≤j≤J,J为待关联信号总数,检测方位角为D_az(j),信号到时为D_time(j);触发信号S(i)与待关联信号D(j)的到时差为:
time_arr_diff=S_time(i)-D_time(j)
其中,grid_distfirst触发台站对应的区域格点震中距参数,grid_distsec为待关联台站对应的区域格点震中距参数,speed为设定的次声信号传播速度,time_travel_diff表示次声信号由格点传输至两个不同次声台站的信号走时差;
S7、对于触发台站和待关联台站,分别根据其震中距格点参数、设定的单位震中距最小传输走时time_delta_min和最大传输走时time_delta_max,分别计算触发台站和待关联台站同时记录到次声信号的情况下,信号走时差的最大值tt_res_max和最小值tt_res_min;
S8、对触发信号和待关联信号进行信号走时匹配分析,当满足以下条件时认为触发信号与待关联信号匹配:
其中,time_th为设定的信号走时偏差阈值,az_gridsec为待关联台站的方位角格点;
S9、对第一个待关联台站的所有待关联信号分别与触发信号进行匹配,筛选出所有满足S8的匹配信号;而后再对下一待关联台站进行上述相同处理,直至所有待关联台站均与触发信号进行完匹配。
2.根据权利要求1所述的一种区域次声事件自动关联方法,其特征在于,S1中,专用次声监测台网通常由3至4个次声台站组成,监测区域原则上小于10°×10°。
3.根据权利要求1所述的一种区域次声事件自动关联方法,其特征在于,S4中,当所有次声台站均含有符合区域次声事件信号特征要求的检测信号时,继续进行后续关联,否则结束整个关联进程。
4.根据权利要求1所述的一种区域次声事件自动关联方法,其特征在于,S7中,
tt_res_min=time_delta_minsec*grid_distsec(grid_sel(k,1),grid_sel(k,2))-time_delta_maxfirst*grid_distfirst(grid_sel(k,1),grid_sel(k,2))
tt_res_max=time_delta_maxsec*grid_distsec(grid_sel(k,1),grid_sel(k,2))-time_delta_minfirst*grid_distfirst(grid_sel(k,1),grid_sel(k,2))
其中,time_delta_minsec为待关联台站单位距离信号传输最小单位距离走时偏差,time_delta_maxsec为待关联台站单位距离信号传输最大走时偏差,time_delta_minfirst为触发台站单位距离信号传输最小走时偏差,time_delta_maxfirst为触发台站单位距离信号最大走时偏差。
5.根据权利要求1所述的一种区域次声事件自动关联方法,其特征在于,S9中,当且仅当次声网中所有次声台站均有与区域内格点相匹配的次声信号时,形成临时关联事件,若含有匹配信号的台站数少于台网中次声台站总数,则认为无足够匹配台站,针对给定的监测区域无法形成有效的信号关联。
6.根据权利要求5所述的一种区域次声事件自动关联方法,其特征在于,还包括S10:当有多个自动关联形成的临时关联事件时,需对临时关联事件进行重复关联事件筛选,挑选出触发信号不同、其余台站关联信号基本一致的事件,并进行信号合并,具体判别过程如下:假设临时关联事件共关联A个信号,除触发信号外关联的信号为A-1个,已关联事件共关联B个信号,除触发信号外关联的信号为B-1个,当满足(A-1)/(B-1)≥N_th时,将临时关联事件与已有的事件进行关联信号合并,并剔除重复信号;否则将临时关联事件作为一个独立事件而保存,其中N_th为预先设定的百分比阈值。
CN202110102318.0A 2021-01-26 2021-01-26 一种区域次声事件自动关联方法 Active CN113607268B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110102318.0A CN113607268B (zh) 2021-01-26 2021-01-26 一种区域次声事件自动关联方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110102318.0A CN113607268B (zh) 2021-01-26 2021-01-26 一种区域次声事件自动关联方法

Publications (2)

Publication Number Publication Date
CN113607268A CN113607268A (zh) 2021-11-05
CN113607268B true CN113607268B (zh) 2024-01-09

Family

ID=78303267

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110102318.0A Active CN113607268B (zh) 2021-01-26 2021-01-26 一种区域次声事件自动关联方法

Country Status (1)

Country Link
CN (1) CN113607268B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116068494B (zh) * 2023-03-31 2023-07-18 中国人民解放军96901部队 一种基于次声传播模型的广域次声源定位方法

Citations (17)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
RU24890U1 (ru) * 2001-06-13 2002-08-27 Федеральное государственное унитарное предприятие "Центральный научно-исследовательский институт "Морфизприбор" Автономная станция гидроакустического наблюдения
CN101738611A (zh) * 2009-12-15 2010-06-16 中国科学院声学研究所 一种水声目标信号检测和识别方法
JP2010266392A (ja) * 2009-05-18 2010-11-25 Kochi Univ Of Technology インフラサウンド測定装置
CN105676287A (zh) * 2016-01-29 2016-06-15 禁核试北京国家数据中心 一种检测特定地区核爆炸地震事件的方法
CN105716707A (zh) * 2015-12-10 2016-06-29 中国航空工业集团公司北京长城计量测试技术研究所 一种超低频次声异常信号判别方法
CN106094021A (zh) * 2016-06-01 2016-11-09 北京科技大学 一种基于到时差数据库的微地震震源快速定位方法
CN108401562B (zh) * 2013-08-02 2017-01-11 中国人民解放军军事科学院防化研究院 一种爆炸次声波的三站定位算法
CN106448233A (zh) * 2016-08-19 2017-02-22 大连理工大学 基于大数据的公交线路时刻表协同优化方法
RU2623837C1 (ru) * 2016-03-25 2017-06-29 Федеральное государственное бюджетное учреждение науки Специальное конструкторское бюро средств автоматизации морских исследований Дальневосточного отделения Российской академии наук Способ экологического мониторинга и охраны районов нефтегазодобычи
CN107272061A (zh) * 2017-06-29 2017-10-20 禁核试北京国家数据中心 一种次声信号与地震事件的自动关联方法
CN107290787A (zh) * 2017-06-29 2017-10-24 禁核试北京国家数据中心 一种地震次声同址台站的监测信号关联方法
US9817925B1 (en) * 2014-02-04 2017-11-14 The United States Of America As Represented By The Secretary Of The Navy Probit method of cumulative distribution function determination of energetic sensitivity
CN109633548A (zh) * 2018-12-10 2019-04-16 禁核试北京国家数据中心 一种水声台网关联方法
CN109669185A (zh) * 2018-12-10 2019-04-23 禁核试北京国家数据中心 一种次声台网定向搜索关联方法
CN111060965A (zh) * 2019-12-05 2020-04-24 禁核试北京国家数据中心 一种基于卷积神经网络的地震震相拾取及事件检测方法
CN111208556A (zh) * 2020-01-14 2020-05-29 禁核试北京国家数据中心 一种水声台站监测信号到时区间估算方法
WO2020148547A1 (en) * 2019-01-18 2020-07-23 Gaiacode Ltd Infrasound detector

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20060248954A1 (en) * 2005-04-26 2006-11-09 Snieder Roelof K System for and method of monitoring structural integrity of a structure

Patent Citations (17)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
RU24890U1 (ru) * 2001-06-13 2002-08-27 Федеральное государственное унитарное предприятие "Центральный научно-исследовательский институт "Морфизприбор" Автономная станция гидроакустического наблюдения
JP2010266392A (ja) * 2009-05-18 2010-11-25 Kochi Univ Of Technology インフラサウンド測定装置
CN101738611A (zh) * 2009-12-15 2010-06-16 中国科学院声学研究所 一种水声目标信号检测和识别方法
CN108401562B (zh) * 2013-08-02 2017-01-11 中国人民解放军军事科学院防化研究院 一种爆炸次声波的三站定位算法
US9817925B1 (en) * 2014-02-04 2017-11-14 The United States Of America As Represented By The Secretary Of The Navy Probit method of cumulative distribution function determination of energetic sensitivity
CN105716707A (zh) * 2015-12-10 2016-06-29 中国航空工业集团公司北京长城计量测试技术研究所 一种超低频次声异常信号判别方法
CN105676287A (zh) * 2016-01-29 2016-06-15 禁核试北京国家数据中心 一种检测特定地区核爆炸地震事件的方法
RU2623837C1 (ru) * 2016-03-25 2017-06-29 Федеральное государственное бюджетное учреждение науки Специальное конструкторское бюро средств автоматизации морских исследований Дальневосточного отделения Российской академии наук Способ экологического мониторинга и охраны районов нефтегазодобычи
CN106094021A (zh) * 2016-06-01 2016-11-09 北京科技大学 一种基于到时差数据库的微地震震源快速定位方法
CN106448233A (zh) * 2016-08-19 2017-02-22 大连理工大学 基于大数据的公交线路时刻表协同优化方法
CN107272061A (zh) * 2017-06-29 2017-10-20 禁核试北京国家数据中心 一种次声信号与地震事件的自动关联方法
CN107290787A (zh) * 2017-06-29 2017-10-24 禁核试北京国家数据中心 一种地震次声同址台站的监测信号关联方法
CN109633548A (zh) * 2018-12-10 2019-04-16 禁核试北京国家数据中心 一种水声台网关联方法
CN109669185A (zh) * 2018-12-10 2019-04-23 禁核试北京国家数据中心 一种次声台网定向搜索关联方法
WO2020148547A1 (en) * 2019-01-18 2020-07-23 Gaiacode Ltd Infrasound detector
CN111060965A (zh) * 2019-12-05 2020-04-24 禁核试北京国家数据中心 一种基于卷积神经网络的地震震相拾取及事件检测方法
CN111208556A (zh) * 2020-01-14 2020-05-29 禁核试北京国家数据中心 一种水声台站监测信号到时区间估算方法

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
Ultra-Wideband Infrasonic Signals Generated by Series of Chemical Explosions;Leonid Chernogor等;2018 9th International Conference on Ultrawideband and Ultrashort Impulse Signals (UWBUSIS);第318-321页 *
基于广域次声传感器网络的地震本地次声波监测;郭泉;杨亦春;吕君;滕鹏晓;;地球科学(中国地质大学学报)(第12期);第1807-1817页 *
基于慢度估计的次声台阵信号自动检测算法及应用;陈虎虎等;声学技术;第34卷(第1期);第85-89页 *
天津"8·12"爆炸次声波分析;任文涛等;第十二届国家安全地球物理专题研讨会;第288-293页 *
朝鲜4.13"光明星3 号"卫星发射次声信号分析;唐伟等;环境工程;第31卷(第1期);第81-84页 *

Also Published As

Publication number Publication date
CN113607268A (zh) 2021-11-05

Similar Documents

Publication Publication Date Title
US7211994B1 (en) Lightning and electro-magnetic pulse location and detection for the discovery of land line location
CN109061317B (zh) 融合甚高频和雷声探测的雷电全过程监测方法及系统
CN103649778A (zh) 用于自动探测海洋动物的方法和装置
CN109633548B (zh) 一种水声台网关联方法
CN113607268B (zh) 一种区域次声事件自动关联方法
Fuchs et al. Seismic detection of rockslides at regional scale: examples from the Eastern Alps and feasibility of kurtosis-based event location
CN107290787B (zh) 一种地震次声同址台站的监测信号关联方法
CN103336298A (zh) 一种大地震断裂区域前兆数据的采集和分析方法
CN106382981A (zh) 一种单站次声波信号识别提取方法
CN107272061A (zh) 一种次声信号与地震事件的自动关联方法
CN114280541A (zh) 一种基于深海分布式垂直线列阵的目标被动定位方法
CN101893661A (zh) 一种光电信号同步监测的闪电数据处理方法
Simon et al. Recording earthquakes for tomographic imaging of the mantle beneath the South Pacific by autonomous MERMAID floats
CN109669185B (zh) 一种次声台网定向搜索关联方法
CN108196269A (zh) 卫星导航抗干扰天线系统内部弱谐波干扰信号检测方法
CN115828187B (zh) 一种星基和地基闪电数据融合方法
Visser et al. A comprehensive earthquake catalogue for the Fort St. John–Dawson Creek region, British Columbia
Kitov et al. Detection, estimation of magnitude, and relative location of weak aftershocks using waveform cross-correlation: The earthquake of August 7, 2016, in the town of Mariupol
Visser et al. A comprehensive earthquake catalogue for northeastern British Columbia: The northern Montney trend from 2017 to 2020 and the Kiskatinaw Seismic Monitoring and Mitigation Area from 2019 to 2020
CN109738050B (zh) 一种水声台网关联格点设计方法
Cataldi et al. Registration of Pre-Seismic Radio Signals Related To The Russian And Jamaican Earthquakes With The RDF System Developed By The Radio Emissions Project
CN105388497A (zh) 基于聚类的卫星导航欺骗攻击防御方法及系统
Ringdal Teleseismic event detection using the NORESS array, with special reference to low-yield Semipalatinsk explosions
Stevens et al. Infrasound scaling and attenuation relations from Soviet explosion data and instrument design criteria from experiments and simulations
Le Bras et al. Integration of seismic, hydroacoustic, infrasound and radionuclide processing at the prototype international data centre

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