单台数字测震仪地震实时准确报警方法研究*

李万金，陈 阳

(1.中国科学技术大学地球和空间科学学院，安徽合肥230026;2.云南省地震局个旧地震台，云南个旧661000)

0 引言

地震台站测震仪器系统在记录地震波形数据的同时，必须能够实时准确检测地震并报警，才能从根本上保障测震人员在地震发生后及时速报地震。目前单台仪器不能对所有记录到的地震进行准确检测，但至少应该保证一定强度地震(如ML≥3.0的地方震，ML≥5.0的近震，MS≥6.0的远震，MS≥7.0的极远震)的实时检测报警。但目前国家数字地震台使用的地震波形监视记录软件的地震报警功能很不理想，经常出现“有震不报警，无震乱报警”的情况，严重影响了台站日常地震速报工作。

由于台站观测仪器易受自然界背景噪声干扰(如机械振动、大风、爆破、人为活动等)(孟晓春，2005)和运行环境恶劣(如受潮、气流扰动、雷电、50 Hz交流电串入等)的影响，导致地震记录波形干扰严重，甚至因仪器运行不正常而使记录波形失真。如果不采取措施来改善波形记录质量以及不对干扰严重的波形数据进行除噪处理，直接使用原始记录波形进行报警，就无法清晰呈现出地震的记录特征，这样的波形有时人工都难以识别，地震记录软件就更不可能实现准确报警。

目前，我国已经完全实现从模拟测震到数字化测震的转变，数字测震仪具有记录频带宽、分辨率高、动态范围大并且便于计算机进行资料处理等优点(刘瑞丰等，1997)。例如，个旧地震台高灵敏度、宽频带、高性能的CTS－1E型甚宽频带地震计，其频带是50 Hz～120 s，对等速度输入响应平坦(蔡亚先等，2004)。可见它的频带基本涵盖了地震波的周期，即现在的一套数字测震仪记录的地震信息相当于过去短、中、长周期3套模拟测震仪记录信息的总和，因此只要把改善地震波记录质量、波形实时仿真滤波技术以及合适的地震检测算法相结合，就能有效地提高数字化测震仪地震检测报警的准确率。

测震台站承担着所在区域范围内有感及以上地震的速报工作。实现从干扰波形中自动检测地震的技术，是实现震相自动识别技术的前提和基础。因此，研究如何提高单台地震报警的准确率具有重要的现实意义和应用价值。

1 地震触发检测原理及改进措施

1.1 STA/LTA触发算法原理及参数设置

目前最常使用的地震触发检测算法是短时绝对均值与长时绝对均值之比(即STA/LTA)(彼得·鲍曼，2006)。STA测量信号的短时间变化，随时监视着地震;而LTA测量台站当前地震噪声的振幅水平。

(1)STA/LTA算法原理

首先分别计算移动的STA和LTA时间窗中振幅的绝对均值，再计算这两个均值的比值(STA/LTA)。程序将此比值与用户设定的触发阈值相比较。若大于触发阈值时，就认为检测到地震并触发报警。当STA/LTA比值低于解除触发阈值时，就认为地震信号结束并解除报警。STA/LTA解除触发阈值要小于STA/LTA触发阈值(图1)。

STA、LTA都作为滑动均值进行计算:

式中:xi是地震信号(经滤波或未经滤波);STA是短时均值;LTA是长时均值;NSTA和NLTA分别是STA窗和LTA窗中的点数。

图1 STA/LTA触发算法原理(a)一个连续的(经过滤波)的地震信号;(b)在STA、LTA两个窗口里平均绝对信号随时间的变化;(c)STA/LTA随时间的变化Fig.1 The principle of STA/LTA trigger algorithm(a)a continuous(filtered)seismic signal;(b)mean absolute signals in STA and LTA windows change with time;(c)STA/LTA changes with time

(2)STA/LTA基本触发参数选择

触发参数设置是否合适对能否成功检测地震事件很关键。STA/LTA的基本触发参数有:STA时间窗长、LTA时间窗长、STA/LTA触发阈值、STA/LTA解除触发阈值。对于仿真成短周期的信号，STA用来测量地震信号的“瞬时”振幅，STA窗长一般为0.5 s;LTA测量台站当前地震噪声的振幅水平，通常LTA窗长为50～500 s。安静的台站，STA/LTA触发阈值常设置为4;人为噪声强的台站，触发阈值可取8以上。STA/LTA解除触发阈值决定何时解除报警，能记录弱地震活动的安静台站，STA/LTA解除触发阈值的取值是2～3;对噪声大些的台站，该阈值应设置高一点。

(3)改进STA/LTA地震检测算法

为提高地震检测灵敏度，并有效减少地震的漏触发、虚触发和重复触发，本研究采用冻结LTA值和设置触发有效的最小持续时间、相临两次有效的地震检测的最小间隔时间来改进传统的STA/LTA地震检测算法。

冻结LTA:一旦检测到事件后就冻结LTA，从而使背景噪声参考水平不受地震事件信号的影响;当事件检测结束时，再解冻LTA来跟踪背景噪声变化。这样STA/LTA可以突出地震信号的变化。

检测触发有效的最小持续时间:为了有效减少尖锐脉冲或爆破等引发的虚触发，根据此类脉冲持续时间短的特点，采取一旦达到触发阈值，检测触发有效的最小持续时间，比如持续时间小于10 s可认定为干扰，持续时间超过10 s就认定为地震，从而能排除绝大部分因干扰引发的报警。

相临两次有效的地震检测最小间隔时间:为了避免同一次地震因波形起伏大引起的重复触发，可采取设置相临两次有效地震的最小间隔时间，即后一次地震的触发时刻与前一次地震解除触发时刻之间的时差(比如为240 s)，间隔小于最小间隔时间的触发就认定为同一次地震，间隔超过最小间隔时间的触发则认定为另一次地震。

1.2 提高地震波形记录质量

台站观测仪器工作环境好、仪器运行正常，才能产出正常、连续、可靠的地震波形数据，这是提高地震检测准确性的前提和基础。如果因仪器工作环境恶劣以及管理不善而导致仪器长期运行不正常、产出断续甚至畸形的波形，那么即便采用好的地震检测算法也无法显著提高地震报警准确率。

(1)台址选择

台址的选择，对地震波形记录质量尤为关键。必须按《地震台站观测规范》(国家地震局，1990)、《数字地震观测技术》(中国地震局监测预报司，2003)和《地震学与地震观测》(中国地震局监测预报司，2007)来选择。选址的注意事项主要有以下几条:测震台的选址，必须远离各种振动干扰源;台基应选择在岩性坚硬致密、完整的、并大面积出露的基岩上;台址的地势起伏要小，应避开风口;台站观测环境地噪声水平要符合要求;摆房要防潮保温。

(2)改善仪器运行环境

大多数测震仪器要在干燥的环境中才能正常工作。对潮湿的观测山洞，可以采取综合措施来降低潮湿程度，如采取工程措施对洞室内壁作防水防潮处理;用密闭的仪器罩罩住地震计，罩内放置能吸湿的硅胶干燥剂;不定期用除湿机对山洞除湿。对于室外摆线的铺设或架设，要采取加装外套或穿管埋地铺设等防护措施，防止鼠类动物啃咬线路。对于雨季雷患严重的地区，要制作良好的接地、电源防雷和GPS天线防雷，防止雷电串入干扰甚至烧毁设备。稳定的供电也是保证仪器正常运行的一个关键因素。

1.3 地震波形实时仿真

实践表明，如果将数字测震仪波形仿真成模拟短周期测震仪波形，就能滤除干扰的主要成分，突出每个地震的纵波，从而把各个地震区分出来。本研究采用无延时的实时仿真算法来仿真原始记录波形，原理如下:

结合宽频带数字测震仪的采样频率和模拟测震仪的传递函数，使用双线性变换法(王洪体等，2006;李万金，邓存华，2011)将模拟测震仪器的传递函数滤波器变换为不同采样率下的数字滤波器传递函数，从而构造出递归仿真滤波器;再用滤波的方法便可将数字地震波形仿真成模拟仪器记录波形，从而实现实时波形仿真。

由于模拟地震仪记录的是地动位移信号，而数字地震仪记录的是地动速度或地动加速度信号，因此需要将数字波形信号先转换成位移量。要将速度型数字地震仪波形转换成位移量，需要积分一次。

一般的做法是:先对原始数字地震记录数据进行积分，再用模拟仪器的传递函数来仿真。本研究的做法是:将速度型数字地震仪波形转换成位移量，需要积分一次，只需要将模拟地震仪的传递函数乘1/S后得到的传递函数来仿真;将加速度型数字地震仪波形转换成位移量，要积分两次，只需要将模拟地震仪的传递函数乘1/S2后得到的传递函数来仿真。这样可以大幅提高波形仿真的速度。

将减少了零点个数的模拟地震仪器传递函数变换为不同采样率下的数字滤波器传递函数，从而构造出递归仿真滤波器，可实现实时波形仿真。

使用MATLAB软件构造递归仿真滤波器(万永革，2007)的步骤:

假设模拟仪器的传递函数形式为H(s)=f(s)/g(s)，则

第一步，使用MATLAB函数f=poly(z)*k和g=poly(p)便可将仪器剩下的零极点转换成传递函数形式。

第二步，使用MATLAB双线性变换函数［a，b］=bilinear(f，g，Fs)将模拟仪器传递函数离散成数字滤波器。式中:f、g为模拟滤波器传递函数分子和分母多项式向量系数;Fs为采样频率，单位为Hz;a、b为数字滤波器传递函数分子和分母多项式向量系数。

如:将DD－1型模拟地震仪传递函数转换成采样率为100 Hz的速度信号的数字仿真滤波器:

a=0.2440 －0.8741 0.8128 0.5086－1.2196 0.4672 0.1628 －0.1，

b=1.0000 －5.3974 12.2810 －15.3370 11.4941 －5.3084 1.4621 －0.1944.

第三步，将构造出的递归仿真滤波器系数代入下式的滤波算法，即可对波形进行实时仿真:

式中:系数 ai(i=0，1，…，n)、bi(i=0，1，…，n)，n为递归滤波器阶数;x(k)、y(k)分别为输入的数字测震仪采集数据和输出的DD－1型模拟短周期测震仪仿真波形数据。

2 地震波形仿真和检测实例

2.1 波形仿真的作用

运用上述仿真算法将波形仿真成DD－1型模拟短周期测震仪波形，主要有以下作用:

(1)滤除仪器不正常的波形

测震仪器因放大电路板工作受潮、接地不良或雷电感应等原因引起波形大幅振荡或毛刺脉冲，仿真成短周期波形后可以显著的滤除这些干扰，从而减少地震虚报。

(2)消除仪器因零漂产生的偏移

利用零漂周期长的特点在仿真成短周期波形过程中就能有效地消除数据中的仪器零漂，突出地震信号。

(3)滤除长周期大幅度干扰

每年台风季节，宽频带地震仪的波形上常记录到卓越周期为3～8 s的大幅度台风干扰。若不滤除，STA/LTA经常检测不到地震，导致地震漏报。

(4)从强震重叠波形中识别出余震

强震波形会湮没其中的余震波形，仿真成短周期波形后可以显露出这些余震，从而减少地震漏报。

地震仿真实例:以个旧地震台记录的2011年3月11日日本本州东海岸附近海域9.0级大地震(震中距为36.6°)为例，采用MATLAB编程对地震信号进行仿真以及地震检测判别，从仿真后的波形中可以发现主震之后又发生了多次余震(图2)。

图2 2011年3月11日日本9.0级大地震波形(a)CTS－1仪器记录的原始地震波形;(b)仿真成DD－1测震仪的波形Fig.2 Japan MS9.0 earthquake on Mar.11，2011(a)original seismic waveform recorded by CTS－1 seismograph;(b)simulatied to waveforms of DD－1 seismograph

2.2 地震波形检测

对近年来个旧地震台记录到的一些地震波形和干扰波形采用上述算法进行重新检测，既要有效排除大部分干扰波形，又要能保证一定级别的地震报警。经过比较，个旧地震台的STA/LTA报警检测参数合理设置为:STA窗长为0.5 s，LTA窗长为50 s，触发阈值是4，解除触发阈值是2，检测有效最小持续时间为8 s，相临两次有效地震检测的最小间隔时间为180 s。以图2b中的垂直向波形为例进行检测演示，如图3所示。

图3 2011年3月11日日本9.0级大地震的检测过程(a)仿真后的DD－1波形;(b)STA;(c)LTA;(d)STA/LTA;(e)地震事件触发及解除触发Fig.3 The detection process of Japan MS9.0 earthquake on Mar.11，2011(a)simulation waveforms of DD－1 seismograph;(b)STA;(c)LTA;(d)STA/LTA;(e)seismic event triggering and releasing triggering

与采用传统算法(图1)对比，从STA/LTA曲线上可看出改进后的STA/LTA算法，当STA/LTA达到触发阈值后，冻结LTA能有效抑制地震信号对LTA的影响，从而提高STA/LTA值;检测触发有效的最小持续时间参数能够有效地排除短时干扰;相临两次有效的地震检测最小间隔时间参数能显著抑制因地震波形起伏较大造成对同一次地震的重复检测。

2.3 地震波形检测效果统计

为了充分验证本算法的可行性和准确性，笔者从个旧地震台数字地震仪近年来记录到的地震中选取5个地震事件进行检验(表1)，其中涵盖了近震、远震和极远震，波形中或有尖锐的爆破波形干扰、或有大周期台风振幅干扰、或余震完全被主震波形湮没，对台站需要检测的地震和波形所受的干扰具有显著的代表性。从表1的地震检测效果统计结果可看出，在仿真成短周期的地震波形上，使用改进后的STA/LTA触发检测算法来检测地震效果显著，只要地震间的间隔时间不要过短，完全可以实现一定级别地震的准确检测，并有效减少地震的虚报、漏报和重复触发，只是有时触发持续时间较短。

表1 地震检测效果统计Tab.1 Statistic results of seismic detection effect

3 讨论及结论

(1)本文首先分析了采用传统STA/LTA地震检测算法在原始记录波形上检测地震时，地震的判别受各种振动干扰或主震波形的影响常常造成地震的漏报或虚报，地震检测效果很不理想，而且对于选择偏重检测近震还是偏重检测远震，参数设置差别大，严重影响了台站的地震速报工作。然后提出了综合解决措施:采取工程技术措施改善仪器的运行环境，提高波形产出质量;将数字测震仪波形仿真成模拟短周期测震仪波形，就能滤除干扰的主要成分，突出每个地震的纵波，从而把各个地震区分出来，然后对仿真后的波形运用改进后的STA/LTA算法来检测地震，只要一套参数，就能对一定强度级别的近震、远震和极远震敏感，从而实现较为准确的地震检测报警，地震漏报、虚报大幅降低。

(2)从检测效果统计来看，有时地震检测的持续时间与地震波形的完整持续时间相比较短，对波形间隔太近(如本文采用的180 s)甚至重叠的两次地震在触发检测时会作为一次地震来处理会造成地震漏报，但考虑到从地震触发报警到横波甚至面波出现后才能准确分析需要一段时间，且分析地震也需要时间，分析人员能从波形浏览中注意到软件漏报的地震，因而不会造成真正的漏报。

(3)只要在台站的地震数据实时记录计算机上设计一个软件程序来从地震数据采集器输出的地震波形数据流中实时接收垂直向的数据，并对接收到的每个数据点应用本文提出的实时仿真算法和改进的STA/LTA检测算法，就能将此地震检测报警算法应用到实际观测中。

(4)要严格区分天然地震波形和爆破等干扰波形需要从波形的形态、P波初动方向、频谱等多个特征综合辨别，检测辨别算法复杂以及处理时间长;本方法只根据STA/LTA比值超过触发阈值的持续时间(如本文采用的8 s)来区分地震和干扰，判别算法简单，并不能完全排除干扰，但能排除绝大部分短时干扰且判别时间短，应该说这种取舍是合适的。

(5)提高单台地震报警准确率，除了采用好的地震检测判别算法外，最重要的是要实时产出质量好的地震波形数据，这需要台站地震观测仪器系统长期正常运转。在影响仪器正常运行因素中，有些问题需要采用工程技术措施来解决，但归根结底还要依靠对台站的精心管理、观测人员高度的责任心以及观测技能的提升。

彼得·鲍曼.2006.新地震观测实践手册［M］.中国监测预报司，译.北京:地震出版社，804－820.

蔡亚先，吕永清，周云耀，等.2004.CTS－1甚宽频带地震计［J］.大地测量与地球动力学，24(3)109－114.

国家地震局.1990.地震台站观测规范［M］.北京:地震出版社，18－20.

李万金，邓存华.2011.将宽频带数字地震仪波形实时仿真成模拟地震仪波的研究［J］.地震通讯，(4):9－15.

刘瑞丰，陈培善，党京平，等.1997.宽频带数字地震记录仿真的应用［J］.地震地磁观测与研究，18(3):7 －12.

孟晓春.2005.地震信息分析技术［M］.北京:地震出版社，12－14.

万永革.2007.数字信号处理的MATLAB实现［M］.北京:科学出版社，154－157.

王洪体，陈阳，庄灿涛.2006.使用双线性变换构造实时地震波形仿真滤波器的方法研究［J］.地震地磁观测与研究，27(1):68－73.

中国地震局监测预报司.2003.数字地震观测技术［M］.北京:地震出版社，563－565.

中国地震局监测预报司.2007.地震学与地震观测［M］.北京:地震出版社，205－226.