SEED转SAC的连续波形截取程序

查小惠，吕坚，鲍志诚

(江西省地震局, 江西南昌　330039)



SEED转SAC的连续波形截取程序

查小惠，吕坚，鲍志诚

(江西省地震局, 江西南昌330039)

摘要：基于LINUX系统下的MATLAB平台，应用rdseed程序和matTaup、saclab程序包，开发了一套SEED单台连续波形截取转存为SAC格式的程序。详细介绍了波形截取过程中的理论到时计算、时间转换运算和SAC头文件处理等问题，实现了从SEED格式的连续波形数据中自动截取SAC格式的事件波形数据。

关键词：SEED； 波形截取；SAC；MATLAB；LINUX

0引言

SEED作为一种数字地震交换的国际标准已成为地震数据记录和使用的标准格式，这一数据格式在我国也被广泛应用[1-5]，并在2003年被规定为地震行业标准。近年来，中国地震局在各地震台站推广使用JOPENS系统，对于台站的连续波形归档格式也有了进一步的规范要求。自2010年以来，南昌中心地震台的连续波形数据归档均采用SEED格式，测震设备运行稳定，波形资料的观测质量较高，积累了大量宝贵的连续波形资料。为了提高波形资料的使用效率，以南昌中心地震台2010年以来积累的SEED格式的原始连续波形数据为基础，尝试性地编写了固定单台SEED连续波形数据截取程序[6]，自动截取出给定地震目录、给定时间长度的事件波形数据，为后续开展相关的研究工作提供数据准备。

1程序运行平台和功能

本程序可以在个人计算机上安装运行，基于LINUX系统下的MATLAB平台开发，主要基于以下2点考虑:①LINUX系统下有成熟的可以截取SEED数据并同时转换为SAC格式的RDSEED程序。②MATLAB平台有很多高级函数，特别是相关的时间函数，可以很方便地解决波形截取过程中遇到的时间运算问题。基于MATLAB编写的matTaup和saclab程序包，在理论走时计算和SAC文件处理等方面都非常方便。而且在MATLAB中可以利用UNIX命令方便的调用LINUX系统的程序。基于LINUX系统下的MATLAB平台编写程序，可以同时利用MATLAB高级函数和LINUX系统程序两者的优势，高效的实现程序要完成的功能。

本程序的功能就是从固定单台SEED连续波形数据中截取事件波形。连续波形数据是使用JOPENS软件归档生成，以小时为单位记录的SEED格式数据。如果连续波形的归档格式发生变化，如以天为单位记录数据，则SEED连续波形截取程序需要修改。目前运行程序的LINUX系统为Ubuntu12.04，MATLAB版本为R2009b。程序的运行时间和地震目录的个数及截取波形的长度有关，单个台站截取500个地震目录1 300 s长度的地震波形，用时大约5分钟，认为该程序的运行效率是可以接受的。

2程序运行流程

程序运行前，要进行相关的数据文件准备。首先是原始连续波形数据，本程序将单个台站的所有连续波形数据存放到一个文件夹下，并将该文件夹按阿拉伯数字编号，对于单台可以将该文件夹用阿拉伯数字命名为1，这样处理是考虑方便文件夹的访问和将来程序的扩展。文件夹下是以年为名称的文件夹，如2010表示该文件夹下存放的是2010年的连续波形数据，年文件夹的子文件夹是月文件夹，如3表示该文件夹下存放的是3月份的连续波形数，月份文件夹下存放的才是SEED格式的连续波形数据。自2010年以来，南昌地震台使用JOPENS软件进行连续波形归档，连续波形数据每一个小时存为一个SEED文件。文件名格式为年+月+日+时.JX.00.seed，这里的文件名所用的时间为北京时间，如2011090116.JX.00.seed，表示的是北京时间2011年9月1日16时到17时的数据，但是SEED数据内部时间是国际时间(记录测震数据时使用GPS授时，故为国际时间)，所以2011090116.JX.00.seed文件记录的是国际时间8点到9点的地震数据。

然后是地震目录的整理，地震目录不限行数，每行8列，前4列表示地震发震时刻，1到3列分别为年、月、日，第4列为时、分、秒，5到8列分别为震源的纬度、经度、深度和震级，各列之间以空格分开。相关目录可以在中国地震科学数据共享中心或USGS GLOBAL SEARCH下载。

最后是准备台站信息，共5列，分别为台站编号、台站经度、台站纬度、台站高程和台站代码。台站信息文件中的台站编号和存放原始连续波形的文件夹命名编号一致，对于单台通常编号为1，这样方便数据的访问。

程序的运行流程如下：

(1)程序的运行需要输入4个参数，即截取波形相对于理论P波到时的提前时间，截取波形的总长度，地震事件相对于台站的最小震中距和最大震中距。震中距的选择可以为相关研究限定使用的地震目录，如希望使用截取到的地震波形事件进行接收函数方面的研究，可将震中距范围限制在30°到95°。

(2)读入台站信息，根据台站的经纬度信息和输入的最大最小震中距信息，对地震目录进行重新的筛选。

(3)挑选完地震目录后，就是到时目录的计算。主要就是计算出需要截取波形的前后2个时间点的时间。利用matTaup程序包中的taupTime计算出P波的理论走时，发震时刻加上P波理论走时就得到P波的理论到时，也就是波形截取的参考时间。然后根据调用函数时输入的截取波形提前时间和截取波形总时间长度，计算出待截取波形的前后2个端点的时间。

(4)得到截取连续波形的2个时间端点后，就开始使用rdseed程序截取事件波形，输出的数据格式为SAC。由于原始连续波形是1个小时1个记录文件，所以需要判断截取波形的2个时间端点是否跨小时，对于在1个小时内的时间段，直接使用rdseed程序转换。对于跨小时的时间段，需要再计算1个中间整点时间点，分2段截取，截取完成之后，再将2段SAC文件连接为单个的SAC文件。目前程序最多可以连接2段SAC文件，所以截取的波形长度最好不要超过3 600 s。

(5)对SAC格式数据的头文件，如B、E、A、O等进行重写。

(6)利用相关文件处理命令，对SAC文件进行存储，循环截取下一个地震事件。最后获得的SAC数据仍然以台站为单位进行存放，也就是单台截取的不同地震事件的SAC文件存放在同一个文件夹下。

3程序中的几个关键问题

3.1理论走时的计算

在理论走时计算方面，本程序采用了matTaup程序包[7-8]中的taupTime子程序，调用语句为taupTime(′ak135′,depth,′p,P,PKP,Pdiff,PKIKP′,′deg′,dist)。在波形截取过程中，初至波的理论到时为发震时刻加上初至波的理论走时，以初至波的理论到时为参考时间，截取初至波到时前后特定时间长度的连续波形。在不同的震中距，初至波的震相并不一致，为了在不同震中距都可以计算得到初至波的理论走时，根据TAUP震相命名规则，在调用taupTime程序时，在震相列表中填入所有可能的初至震相，选取输出结果的第一行震相走时就是初至波的理论走时。

3.2截取时间的计算

MATLAB中表示日期的时间有3种格式：日期字符串、连续的日期数值和日期向量[8]。日期字符串格式使用较多，也较为直观，如2013-10-1 10:10:10，但是不方便进行时间的加减运算。连续的日期数值格式是以公元元年1月1日为起点的，用一个数值表示当前时间到这个起点的时间距离，该数值为浮点数，所以可以很方便的进行时间的加减运算。datenum函数就可以把某种日期时间转换成连续的日期数值输出。datestr函数可以把某种日期格式转换成日期字符串格式。波形截取程序中，需要进行儒略日的计算，利用datenum函数使用一条语句就可以加以实现， julianday=datenum(year,month,day)-datenum(year,1,1)+1。波形截取程序中需要进行时间的加减运算，以确定截取波形的前后时间点,日期字符串表示的时间进行加减运算非常困难，本程序采取的方法是，在matlab中先将日期字符串格式转化为连续的日期数值格式，然后进行相关的时间运算，运算完以后再转换回日期字符串的格式，这样就使得时间的运算非常方便。

3.3波形数据的截取

本程序在MATLAB中通过打印字符串的方式输出一个文本格式的调用rdseed的脚本文件，然后使用UNIX命令调用LINUX系统下的source命令，完成脚本的执行，这样就可以在MATLAB中间接调用rdseed程序。这里必须注意，JOPENS归档连续波形时对文件的命名使用的是北京时间，而SEED文件中数据的时间是国际时间，本程序使用的地震目录时间以及截取波形的端点时间都采用国际时间，所以在搜索待截取的SEED文件时，需要将国际时转换为北京时后再根据SEED文件名对SEED波形文件进行定位。

3.4SAC头文件的处理

在波形截取和数据格式转换过程中，要保持数据时间的不变性，同一个数据点在SEED格式中的时间和SAC格式中的时间必须一致。在SAC文件中，数据点的时间计算涉及2个量，参考时间和时间偏移量。SAC头文件中包含了1个参考时间，用6个整形数据来存储，分别是 NZYEAR、NZJDAY、NZHOUR、NZMIN、NZSEC和NZMSEC。这个参考时间可以任意修改，一般是文件第一个数据点的时间，也可以是地震事件的时间或其他任意时间。SAC头文件中的其他时间都是以该参考时间为起点的以秒计的偏移量，它们以浮点型数据存储于头文件中。

如果设置文件的第一个数据点时间为参考时间，则B值等于零，E值等于截取波形文件的长度。如果设置发震时刻T1为参考时间，假设P波理论走时为T2，截取理论P波前的波形时间为T3，截取波形总时间长度为T4，那么截取波形后的B值，即文件第一个数据点的时间相对T1的偏移量为T2-T3,E值为T2+T4-T3，O值为0，A值为T2。

JOPENS产生的SEED文件内台站控制头中的方位角和倾角填写有误[2]，使用rdseed转换得到SAC文件后，头文件中的方位角和倾角参数需要修改，BHE、BHN和BHZ分量的CMPINC头文件数值应分别设置为90、90、0；CMPAZ头文件数值应分别设置为90、0、0。

4结束语

SEED连续波形截取程序，简单来说就是一个MATLAB函数，它实现的功能就是从SEED格式连续波形数据中截取事件波形并转存为SAC格式。使用该程序需要的SEED文件归档格式是以小时为单位的单台记录文件，最大截取的波形长度不能超过2个小时。对于省地震局台网的数据归档格式，一般单个SEED文件包含多个台站的波形文件，或者数据归档格式是以天为单位，则需要对程序做细微的改动以适应新的数据归档格式(目前均已实现)。该程序使用了rdseed这一核心子程序，并使用相关的MATLAB工具包，分别完成特定的功能，利用MATLAB的高端平台特性，以简单的编程高效的完成了波形截取的任务。在LINUX系统下利用MATLAB平台连接各种处理程序的编程方式，可以减少编程时间，快速完成特定功能，但程序的移植性和易用性都较差，不过作为个人的科研使用还是有一定的价值。

参考文献：

[1]王秀娟,和跃时,武利华.SEED格式地震数据的快速转换软件[J]．地震地磁观测与研究，2001, 22(5):42-47.

[2]何加勇，李松林，陈会忠．地震波形归档格式分析和转换[J]．震灾防御技术，2009, 4(4)：461-465.

[3]王秀文,姚立乎,赖德伦,等.地震数据交换标准[J].地震地磁观测与研究, 1994, 15(2):1-3.

[4]张旸，滕云田，王喜珍，等.地震观测仪器的MiniSEED数据格式的实现[J]. 地震地磁观测与研究，2007,28(3):110-113.

[5]何加勇，吴建平，王未来.中国地震科学探测台阵数据服务系统[J]．地震学报，2010, 32(4)：490-494．

[6]查小惠,张广伟,杨雪超.REFTEK连续波形截取程序[J]. 地震地磁观测与研究，2015，36(1):145-148.

[7]Crotwell H P.Owens T.J and Ritsema J.The Taup ToolKit: Flexible Seismic Travel-Time and Raypath Utilities[J].Seismological Research Letters, 1999,70(2):154-160.

[8]王正林,刘明.精通MATLAB7[M].北京：电子工业出版社, 2007:10-150.

THE INTERCEPTION PROGRAM OF SEED CONTINUOUS WAVEFORM TO SAC

ZHA Xiaohui,LV Jian,BAO Zhicheng

(Earthquake Administration Of Jiangxi Province,Nangchang 330039,China)

Abstract:The interception program of SEED continuous waveform is developed based on the MATLAB platform under LINUX system corporate with rdseed program and matTaup、saclab MATLAB toolbox.The problems of theory time computation,time conversion and SAC head file processing are elaborated.

Key words:SEED;Waveform interception;Linux;Matlab;SAC

收稿日期：2015-12-19

作者简介：查小惠(1989—)，男，江西鄱阳县人，硕士，2010年本科毕业于中国海洋大学，主要从事前兆监测和接收函数方面的研究工作。

中图分类号：P315.39

文献标识码：A

文章编号：1005-586X(2016)02-0054-04