查小惠，呂堅，鮑志誠

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

?

SEED轉(zhuǎn)SAC的連續(xù)波形截取程序

查小惠，呂堅，鮑志誠

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

摘要：基于LINUX系統(tǒng)下的MATLAB平臺，應(yīng)用rdseed程序和matTaup、saclab程序包，開發(fā)了一套SEED單臺連續(xù)波形截取轉(zhuǎn)存為SAC格式的程序。詳細(xì)介紹了波形截取過程中的理論到時計算、時間轉(zhuǎn)換運算和SAC頭文件處理等問題，實現(xiàn)了從SEED格式的連續(xù)波形數(shù)據(jù)中自動截取SAC格式的事件波形數(shù)據(jù)。

關(guān)鍵詞：SEED； 波形截??；SAC；MATLAB；LINUX

0引言

SEED作為一種數(shù)字地震交換的國際標(biāo)準(zhǔn)已成為地震數(shù)據(jù)記錄和使用的標(biāo)準(zhǔn)格式，這一數(shù)據(jù)格式在我國也被廣泛應(yīng)用[1-5]，并在2003年被規(guī)定為地震行業(yè)標(biāo)準(zhǔn)。近年來，中國地震局在各地震臺站推廣使用JOPENS系統(tǒng)，對于臺站的連續(xù)波形歸檔格式也有了進一步的規(guī)范要求。自2010年以來，南昌中心地震臺的連續(xù)波形數(shù)據(jù)歸檔均采用SEED格式，測震設(shè)備運行穩(wěn)定，波形資料的觀測質(zhì)量較高，積累了大量寶貴的連續(xù)波形資料。為了提高波形資料的使用效率，以南昌中心地震臺2010年以來積累的SEED格式的原始連續(xù)波形數(shù)據(jù)為基礎(chǔ)，嘗試性地編寫了固定單臺SEED連續(xù)波形數(shù)據(jù)截取程序[6]，自動截取出給定地震目錄、給定時間長度的事件波形數(shù)據(jù)，為后續(xù)開展相關(guān)的研究工作提供數(shù)據(jù)準(zhǔn)備。

1程序運行平臺和功能

本程序可以在個人計算機上安裝運行，基于LINUX系統(tǒng)下的MATLAB平臺開發(fā)，主要基于以下2點考慮:①LINUX系統(tǒng)下有成熟的可以截取SEED數(shù)據(jù)并同時轉(zhuǎn)換為SAC格式的RDSEED程序。②MATLAB平臺有很多高級函數(shù)，特別是相關(guān)的時間函數(shù)，可以很方便地解決波形截取過程中遇到的時間運算問題。基于MATLAB編寫的matTaup和saclab程序包，在理論走時計算和SAC文件處理等方面都非常方便。而且在MATLAB中可以利用UNIX命令方便的調(diào)用LINUX系統(tǒng)的程序。基于LINUX系統(tǒng)下的MATLAB平臺編寫程序，可以同時利用MATLAB高級函數(shù)和LINUX系統(tǒng)程序兩者的優(yōu)勢，高效的實現(xiàn)程序要完成的功能。

本程序的功能就是從固定單臺SEED連續(xù)波形數(shù)據(jù)中截取事件波形。連續(xù)波形數(shù)據(jù)是使用JOPENS軟件歸檔生成，以小時為單位記錄的SEED格式數(shù)據(jù)。如果連續(xù)波形的歸檔格式發(fā)生變化，如以天為單位記錄數(shù)據(jù)，則SEED連續(xù)波形截取程序需要修改。目前運行程序的LINUX系統(tǒng)為Ubuntu12.04，MATLAB版本為R2009b。程序的運行時間和地震目錄的個數(shù)及截取波形的長度有關(guān)，單個臺站截取500個地震目錄1 300 s長度的地震波形，用時大約5分鐘，認(rèn)為該程序的運行效率是可以接受的。

2程序運行流程

程序運行前，要進行相關(guān)的數(shù)據(jù)文件準(zhǔn)備。首先是原始連續(xù)波形數(shù)據(jù)，本程序?qū)蝹€臺站的所有連續(xù)波形數(shù)據(jù)存放到一個文件夾下，并將該文件夾按阿拉伯?dāng)?shù)字編號，對于單臺可以將該文件夾用阿拉伯?dāng)?shù)字命名為1，這樣處理是考慮方便文件夾的訪問和將來程序的擴展。文件夾下是以年為名稱的文件夾，如2010表示該文件夾下存放的是2010年的連續(xù)波形數(shù)據(jù)，年文件夾的子文件夾是月文件夾，如3表示該文件夾下存放的是3月份的連續(xù)波形數(shù)，月份文件夾下存放的才是SEED格式的連續(xù)波形數(shù)據(jù)。自2010年以來，南昌地震臺使用JOPENS軟件進行連續(xù)波形歸檔，連續(xù)波形數(shù)據(jù)每一個小時存為一個SEED文件。文件名格式為年+月+日+時.JX.00.seed，這里的文件名所用的時間為北京時間，如2011090116.JX.00.seed，表示的是北京時間2011年9月1日16時到17時的數(shù)據(jù)，但是SEED數(shù)據(jù)內(nèi)部時間是國際時間(記錄測震數(shù)據(jù)時使用GPS授時，故為國際時間)，所以2011090116.JX.00.seed文件記錄的是國際時間8點到9點的地震數(shù)據(jù)。

然后是地震目錄的整理，地震目錄不限行數(shù)，每行8列，前4列表示地震發(fā)震時刻，1到3列分別為年、月、日，第4列為時、分、秒，5到8列分別為震源的緯度、經(jīng)度、深度和震級，各列之間以空格分開。相關(guān)目錄可以在中國地震科學(xué)數(shù)據(jù)共享中心或USGS GLOBAL SEARCH下載。

最后是準(zhǔn)備臺站信息，共5列，分別為臺站編號、臺站經(jīng)度、臺站緯度、臺站高程和臺站代碼。臺站信息文件中的臺站編號和存放原始連續(xù)波形的文件夾命名編號一致，對于單臺通常編號為1，這樣方便數(shù)據(jù)的訪問。

程序的運行流程如下：

(1)程序的運行需要輸入4個參數(shù)，即截取波形相對于理論P波到時的提前時間，截取波形的總長度，地震事件相對于臺站的最小震中距和最大震中距。震中距的選擇可以為相關(guān)研究限定使用的地震目錄，如希望使用截取到的地震波形事件進行接收函數(shù)方面的研究，可將震中距范圍限制在30°到95°。

(2)讀入臺站信息，根據(jù)臺站的經(jīng)緯度信息和輸入的最大最小震中距信息，對地震目錄進行重新的篩選。

(3)挑選完地震目錄后，就是到時目錄的計算。主要就是計算出需要截取波形的前后2個時間點的時間。利用matTaup程序包中的taupTime計算出P波的理論走時，發(fā)震時刻加上P波理論走時就得到P波的理論到時，也就是波形截取的參考時間。然后根據(jù)調(diào)用函數(shù)時輸入的截取波形提前時間和截取波形總時間長度，計算出待截取波形的前后2個端點的時間。

(4)得到截取連續(xù)波形的2個時間端點后，就開始使用rdseed程序截取事件波形，輸出的數(shù)據(jù)格式為SAC。由于原始連續(xù)波形是1個小時1個記錄文件，所以需要判斷截取波形的2個時間端點是否跨小時，對于在1個小時內(nèi)的時間段，直接使用rdseed程序轉(zhuǎn)換。對于跨小時的時間段，需要再計算1個中間整點時間點，分2段截取，截取完成之后，再將2段SAC文件連接為單個的SAC文件。目前程序最多可以連接2段SAC文件，所以截取的波形長度最好不要超過3 600 s。

(5)對SAC格式數(shù)據(jù)的頭文件，如B、E、A、O等進行重寫。

(6)利用相關(guān)文件處理命令，對SAC文件進行存儲，循環(huán)截取下一個地震事件。最后獲得的SAC數(shù)據(jù)仍然以臺站為單位進行存放，也就是單臺截取的不同地震事件的SAC文件存放在同一個文件夾下。

3程序中的幾個關(guān)鍵問題

3.1理論走時的計算

在理論走時計算方面，本程序采用了matTaup程序包[7-8]中的taupTime子程序，調(diào)用語句為taupTime(′ak135′,depth,′p,P,PKP,Pdiff,PKIKP′,′deg′,dist)。在波形截取過程中，初至波的理論到時為發(fā)震時刻加上初至波的理論走時，以初至波的理論到時為參考時間，截取初至波到時前后特定時間長度的連續(xù)波形。在不同的震中距，初至波的震相并不一致，為了在不同震中距都可以計算得到初至波的理論走時，根據(jù)TAUP震相命名規(guī)則，在調(diào)用taupTime程序時，在震相列表中填入所有可能的初至震相，選取輸出結(jié)果的第一行震相走時就是初至波的理論走時。

3.2截取時間的計算

MATLAB中表示日期的時間有3種格式：日期字符串、連續(xù)的日期數(shù)值和日期向量[8]。日期字符串格式使用較多，也較為直觀，如2013-10-1 10:10:10，但是不方便進行時間的加減運算。連續(xù)的日期數(shù)值格式是以公元元年1月1日為起點的，用一個數(shù)值表示當(dāng)前時間到這個起點的時間距離，該數(shù)值為浮點數(shù)，所以可以很方便的進行時間的加減運算。datenum函數(shù)就可以把某種日期時間轉(zhuǎn)換成連續(xù)的日期數(shù)值輸出。datestr函數(shù)可以把某種日期格式轉(zhuǎn)換成日期字符串格式。波形截取程序中，需要進行儒略日的計算，利用datenum函數(shù)使用一條語句就可以加以實現(xiàn)， julianday=datenum(year,month,day)-datenum(year,1,1)+1。波形截取程序中需要進行時間的加減運算，以確定截取波形的前后時間點,日期字符串表示的時間進行加減運算非常困難，本程序采取的方法是，在matlab中先將日期字符串格式轉(zhuǎn)化為連續(xù)的日期數(shù)值格式，然后進行相關(guān)的時間運算，運算完以后再轉(zhuǎn)換回日期字符串的格式，這樣就使得時間的運算非常方便。

3.3波形數(shù)據(jù)的截取

本程序在MATLAB中通過打印字符串的方式輸出一個文本格式的調(diào)用rdseed的腳本文件，然后使用UNIX命令調(diào)用LINUX系統(tǒng)下的source命令，完成腳本的執(zhí)行，這樣就可以在MATLAB中間接調(diào)用rdseed程序。這里必須注意，JOPENS歸檔連續(xù)波形時對文件的命名使用的是北京時間，而SEED文件中數(shù)據(jù)的時間是國際時間，本程序使用的地震目錄時間以及截取波形的端點時間都采用國際時間，所以在搜索待截取的SEED文件時，需要將國際時轉(zhuǎn)換為北京時后再根據(jù)SEED文件名對SEED波形文件進行定位。

3.4SAC頭文件的處理

在波形截取和數(shù)據(jù)格式轉(zhuǎn)換過程中，要保持?jǐn)?shù)據(jù)時間的不變性，同一個數(shù)據(jù)點在SEED格式中的時間和SAC格式中的時間必須一致。在SAC文件中，數(shù)據(jù)點的時間計算涉及2個量，參考時間和時間偏移量。SAC頭文件中包含了1個參考時間，用6個整形數(shù)據(jù)來存儲，分別是 NZYEAR、NZJDAY、NZHOUR、NZMIN、NZSEC和NZMSEC。這個參考時間可以任意修改，一般是文件第一個數(shù)據(jù)點的時間，也可以是地震事件的時間或其他任意時間。SAC頭文件中的其他時間都是以該參考時間為起點的以秒計的偏移量，它們以浮點型數(shù)據(jù)存儲于頭文件中。

如果設(shè)置文件的第一個數(shù)據(jù)點時間為參考時間，則B值等于零，E值等于截取波形文件的長度。如果設(shè)置發(fā)震時刻T1為參考時間，假設(shè)P波理論走時為T2，截取理論P波前的波形時間為T3，截取波形總時間長度為T4，那么截取波形后的B值，即文件第一個數(shù)據(jù)點的時間相對T1的偏移量為T2-T3,E值為T2+T4-T3，O值為0，A值為T2。

JOPENS產(chǎn)生的SEED文件內(nèi)臺站控制頭中的方位角和傾角填寫有誤[2]，使用rdseed轉(zhuǎn)換得到SAC文件后，頭文件中的方位角和傾角參數(shù)需要修改，BHE、BHN和BHZ分量的CMPINC頭文件數(shù)值應(yīng)分別設(shè)置為90、90、0；CMPAZ頭文件數(shù)值應(yīng)分別設(shè)置為90、0、0。

4結(jié)束語

SEED連續(xù)波形截取程序，簡單來說就是一個MATLAB函數(shù)，它實現(xiàn)的功能就是從SEED格式連續(xù)波形數(shù)據(jù)中截取事件波形并轉(zhuǎn)存為SAC格式。使用該程序需要的SEED文件歸檔格式是以小時為單位的單臺記錄文件，最大截取的波形長度不能超過2個小時。對于省地震局臺網(wǎng)的數(shù)據(jù)歸檔格式，一般單個SEED文件包含多個臺站的波形文件，或者數(shù)據(jù)歸檔格式是以天為單位，則需要對程序做細(xì)微的改動以適應(yīng)新的數(shù)據(jù)歸檔格式(目前均已實現(xiàn))。該程序使用了rdseed這一核心子程序，并使用相關(guān)的MATLAB工具包，分別完成特定的功能，利用MATLAB的高端平臺特性，以簡單的編程高效的完成了波形截取的任務(wù)。在LINUX系統(tǒng)下利用MATLAB平臺連接各種處理程序的編程方式，可以減少編程時間，快速完成特定功能，但程序的移植性和易用性都較差，不過作為個人的科研使用還是有一定的價值。

參考文獻：

[1]王秀娟,和躍時,武利華.SEED格式地震數(shù)據(jù)的快速轉(zhuǎn)換軟件[J]．地震地磁觀測與研究，2001, 22(5):42-47.

[2]何加勇，李松林，陳會忠．地震波形歸檔格式分析和轉(zhuǎn)換[J]．震災(zāi)防御技術(shù)，2009, 4(4)：461-465.

[3]王秀文,姚立乎,賴德倫,等.地震數(shù)據(jù)交換標(biāo)準(zhǔn)[J].地震地磁觀測與研究, 1994, 15(2):1-3.

[4]張旸，滕云田，王喜珍，等.地震觀測儀器的MiniSEED數(shù)據(jù)格式的實現(xiàn)[J]. 地震地磁觀測與研究，2007,28(3):110-113.

[5]何加勇，吳建平，王未來.中國地震科學(xué)探測臺陣數(shù)據(jù)服務(wù)系統(tǒng)[J]．地震學(xué)報，2010, 32(4)：490-494．

[6]查小惠,張廣偉,楊雪超.REFTEK連續(xù)波形截取程序[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].北京：電子工業(yè)出版社, 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年本科畢業(yè)于中國海洋大學(xué)，主要從事前兆監(jiān)測和接收函數(shù)方面的研究工作。

中圖分類號：P315.39

文獻標(biāo)識碼：A

文章編號：1005-586X(2016)02-0054-04