大地電場(chǎng)數(shù)據(jù)頻譜分析及幅值提取軟件的研發(fā)

王利兵，李細(xì)順，王 靜，胡秀娟，羅 娜，暢國(guó)平，賈 華，于春頌

（1.河北省地震局紅山基準(zhǔn)地震臺(tái)，河北 邢臺(tái) 054000；2.河北省地震局秦皇島中心地震臺(tái)，河北 秦皇島 066000）

0 引言

地電場(chǎng)是地球物理基本場(chǎng)之一，分布于地表的地電場(chǎng)可分為大地電場(chǎng)和自然電場(chǎng)［1］。大地電場(chǎng)是由電離層（高約100km）中的電流體系在地下感應(yīng)產(chǎn)生的電場(chǎng)，強(qiáng)度在地面上隨時(shí)間、地點(diǎn)而異。它的變化周期有11年、1年、1月、24小時(shí)、12小時(shí)、8小時(shí)、6小時(shí)以及更短的周期，幅度一般以毫伏／千米計(jì)算。近十幾年來(lái)，電極技術(shù)的改善，為地電場(chǎng)觀測(cè)的發(fā)展創(chuàng)造了有利條件。在此基礎(chǔ)上，實(shí)施“九五”、首都圈示范工程及“十五”計(jì)劃，已建成數(shù)字化地電場(chǎng)觀測(cè)網(wǎng)。

紅山基準(zhǔn)地震臺(tái)ZD9A－II地電場(chǎng)儀是2007年1月開(kāi)始試運(yùn)行，2008年正式運(yùn)行并參加評(píng)比。ZD9AII地電場(chǎng)儀是用于觀測(cè)地表地電場(chǎng)強(qiáng)度的專用儀器，具有高靈敏度、大動(dòng)態(tài)范圍及良好的抗干擾能力和較高的可靠性等特點(diǎn)。主要用來(lái)觀測(cè)地電場(chǎng)強(qiáng)度及其時(shí)空變化，為地球物理場(chǎng)研究和地震監(jiān)測(cè)預(yù)測(cè)服務(wù)。

國(guó)內(nèi)外已開(kāi)展大地電場(chǎng)觀測(cè)用于地震預(yù)測(cè)的研究多年，并記錄到部分地震電信號(hào)。從理論研究和室內(nèi)巖石力學(xué)實(shí)驗(yàn)的結(jié)果出發(fā)，不少地震學(xué)家預(yù)言地震前有異常電信號(hào)存在，但這些信號(hào)混雜于正常大地電場(chǎng)變化的背景中不易識(shí)別。顯然，尋找并分離正常變化與異常信號(hào)的方法是大地電場(chǎng)法預(yù)報(bào)地震的一項(xiàng)重要研究?jī)?nèi)容［2］。充分利用數(shù)字化觀測(cè)資料，開(kāi)展地電場(chǎng)預(yù)測(cè)地震的研究，對(duì)于推動(dòng)我國(guó)地電場(chǎng)觀測(cè)與研究，具有重要的意義［3］。對(duì)觀測(cè)數(shù)據(jù)進(jìn)行傅立葉變換可以提取電信號(hào)中各頻率成分的幅值，可對(duì)大地電場(chǎng)的強(qiáng)度、頻率分布進(jìn)行研究。在此基礎(chǔ)上，對(duì)數(shù)據(jù)進(jìn)行分時(shí)間段處理，提取不同時(shí)間段同一頻率的幅值，可對(duì)大地電場(chǎng)的某一頻率信號(hào)強(qiáng)度隨時(shí)間的變化進(jìn)行分析，研究其震前異常信息。為此，研發(fā)“大地電場(chǎng)數(shù)據(jù)頻譜分析及幅值提取軟件”。

1 觀測(cè)數(shù)據(jù)文件格式

ZD9A－II地電場(chǎng)儀觀測(cè)數(shù)據(jù)為文本格式，分兩列，為時(shí)間和其對(duì)應(yīng)的觀測(cè)值（見(jiàn)表1）。時(shí)間序列數(shù)據(jù)前四位數(shù)字代表年，之后每?jī)晌粩?shù)字分別代表月、日、時(shí)、分。ZD9A－II地電場(chǎng)儀采樣率為1次／（分鐘·通道），測(cè)量準(zhǔn)確度優(yōu)于±（0.1%讀數(shù)＋0.02%滿度），測(cè)量分辨力優(yōu)于10μV，頻帶范圍為DC～0.005Hz，測(cè)量范圍為±1 000.000mV，動(dòng)態(tài)范圍不小于100dB。

表1 ZD9A－Ⅱ地電場(chǎng)儀觀測(cè)數(shù)據(jù)格式Table.1 Observational data format of ZD9A－Ⅱgeoelectric meter

2 軟件設(shè)計(jì)與開(kāi)發(fā)

2.1 軟件設(shè)計(jì)

軟件需要有一個(gè)良好的用戶界面，通過(guò)一些選項(xiàng)來(lái)進(jìn)行數(shù)據(jù)路徑及必要功能的設(shè)置，可以方便、大量地進(jìn)行數(shù)據(jù)運(yùn)算及處理。

軟件的功能主要包括：

（1）觀測(cè)數(shù)據(jù)分段頻譜分析。

（2）各段數(shù)據(jù)頻譜同頻幅值提取。

2.2 軟件開(kāi)發(fā)

Matlab是一門(mén)計(jì)算機(jī)編程語(yǔ)言，取名來(lái)源于Matrix Laboratory，本意是專門(mén)以矩陣的方式來(lái)處理計(jì)算機(jī)數(shù)據(jù)，它把數(shù)值計(jì)算和可視化環(huán)境集成到一起，非常直觀，而且提供大量的函數(shù)，工具箱越來(lái)越多，應(yīng)用范圍也越來(lái)越廣。

利用Matlab語(yǔ)言可方便地建立軟件圖形用戶界面（見(jiàn)圖1）。圖中各個(gè)按鈕功能如下：

（1）“數(shù)據(jù)文件名”按鈕用于指定數(shù)據(jù)文件的名稱。

（2）“數(shù)據(jù)載入”按鈕用于加載觀測(cè)數(shù)據(jù)，主要由“textread”函數(shù)完成。

（3）“傅氏變換”按鈕用于數(shù)據(jù)分段頻譜分析并存儲(chǔ)，主要由“fft”函數(shù)和“fprintf”函數(shù)完成。

（4）“數(shù)據(jù)提取”按鈕用于從各段頻譜數(shù)據(jù)文件中提取某一頻率對(duì)應(yīng)的幅值并存儲(chǔ)，主要由“dir”函數(shù)、“textread”函數(shù)、“fprintf”函數(shù)完成。

（5）“退出”按鈕用于退出軟件。

圖1 軟件圖形用戶界面Fig.1 Graphical user interface of the software

2.3 軟件運(yùn)行環(huán)境

此軟件由 Matlab7.0語(yǔ)言研發(fā)，其運(yùn)行需要有Matlab7.0或以上版本語(yǔ)言平臺(tái)；程序運(yùn)行處理大量的觀測(cè)數(shù)據(jù)，需要有足夠的內(nèi)存，至少1G以上。

2.4 軟件操作步驟及注意事項(xiàng)

（1）在程序目錄下雙擊“waveFormAnalysis.fig”文件或在Matlab命令行中輸入字符“waveFormAnalysis”回車(chē)，軟件運(yùn)行顯示用戶界面。軟件運(yùn)行時(shí)，圖形界面文件（＊.fig）和各功能文件（＊.m）應(yīng)在同一目錄中。

（2）單擊“數(shù)據(jù)文件名”按鈕，打開(kāi)數(shù)據(jù)文件選擇窗口，指定要處理的原始觀測(cè)數(shù)據(jù)文件。

（3）單擊“數(shù)據(jù)載入”按鈕，加載觀測(cè)數(shù)據(jù)，并在圖形顯示區(qū)顯示數(shù)據(jù)波形。

（4）單擊“傅氏變換”按鈕，對(duì)加載的數(shù)據(jù)進(jìn)行分段傅里葉變換，在圖形顯示區(qū)顯示變換結(jié)果，存儲(chǔ)結(jié)果數(shù)據(jù)。變換結(jié)果數(shù)據(jù)文件名前3個(gè)字符為“fft”，其后為數(shù)據(jù)段的起始日期，年為4位數(shù)，月、日均為2位數(shù)，如，“fft20121221.txt”。文件數(shù)據(jù)分兩列，分別為頻率值及其對(duì)應(yīng)的幅值（見(jiàn)表2）。

表2 傅里葉變換結(jié)果文件數(shù)據(jù)格式Table.2 File data format of Fourier transform results

（5）單擊“數(shù)據(jù)提取”按鈕，將各段數(shù)據(jù)的傅氏變換結(jié)果文件中的某一指定頻率對(duì)應(yīng)的幅值提取出來(lái)，與其時(shí)間一起存入數(shù)據(jù)文件中（見(jiàn)表3）。數(shù)據(jù)文件名以頻率代碼字符開(kāi)頭，其后為觀測(cè)數(shù)據(jù)文件名，如，“FHW20110101－20120111.txt”。“DW”為日波數(shù)據(jù)，“HDW”為半日波數(shù)據(jù)，“EHW”為8小時(shí)波數(shù)據(jù)，“SHW”為6小時(shí)波數(shù)據(jù)，“FHW”為4小時(shí)波數(shù)據(jù)。

（6）單擊“退出”按鈕，退出軟件。

3 軟件數(shù)據(jù)應(yīng)用

利用EXCEL辦公軟件，打開(kāi)此軟件提取的各段數(shù)據(jù)傅氏變換結(jié)果文件中的某一指定頻率對(duì)應(yīng)的幅值，所形成的數(shù)據(jù)文件（以“固定寬度”方式打開(kāi)），則第1列為日期值，第2列為幅值。以日期值為橫坐標(biāo)，以某一頻率幅值或同一測(cè)道不同頻率幅值的比值為縱坐標(biāo)畫(huà)圖（見(jiàn)圖2）。圖中，圖2a為高郵臺(tái)原始數(shù)據(jù)波形圖，圖2b為同一測(cè)道的半日波與日波幅值的比值曲線圖。從圖中可以清楚地看到，比值隨時(shí)間的變化曲線，從曲線形態(tài)上可直觀地分析研究大地電場(chǎng)背景信號(hào)與異常信號(hào)，總結(jié)異常信號(hào)的特點(diǎn)及規(guī)律。同時(shí)，也可對(duì)同一頻率不同測(cè)道的結(jié)果數(shù)據(jù)進(jìn)行對(duì)比分析，從同步形態(tài)或其差分?jǐn)?shù)值上找出異常信號(hào)。

表3 同一頻率幅值數(shù)據(jù)文件格式Table.3 Data format of the same frequency amplitude

4 結(jié)語(yǔ)

該軟件作為地電場(chǎng)儀觀測(cè)數(shù)據(jù)的一個(gè)處理工具，在大地電場(chǎng)數(shù)據(jù)分析、研究過(guò)程中，可起到一定的輔助作用。大地電場(chǎng)信號(hào)成分復(fù)雜，只有某些電信號(hào)成分對(duì)地震孕育過(guò)程反應(yīng)敏感，利用此軟件可以完成信號(hào)頻率的抽取，方便工作人員和分析研究人員針對(duì)性地對(duì)電信號(hào)進(jìn)行辨別、分析、研究，有利于大地電場(chǎng)前兆異常信號(hào)的研究工作，提高數(shù)據(jù)的利用價(jià)值。由于數(shù)據(jù)分析是相通的，利用此軟件還可以對(duì)具有相同格式的其他觀測(cè)儀器數(shù)據(jù)進(jìn)行分析處理，得到相似的研究成果。如有需要，對(duì)軟件稍做修改，可以處理多道數(shù)據(jù)，并進(jìn)行對(duì)比分析，得到推廣應(yīng)用。

圖2 同測(cè)道半日波與日波幅值比變化曲線Fig.2 Amplitude ratio variation of semidiurnal wave and daily wave at the same trace

［1］高玉芬，錢(qián)家棟.地震及前兆數(shù)字觀測(cè)技術(shù)規(guī)范（電磁觀測(cè)）［M］.北京：地震出版社，2002：20－70.

［2］田 山，王建國(guó)，徐學(xué)恭，等.大地電場(chǎng)觀測(cè)地震前兆異常提取技術(shù)研究［J］.地震學(xué)報(bào)，2009，31（4）：424－431.

［3］馬君釗，張 磊，關(guān)華平，等.大地電場(chǎng)觀測(cè)各類干擾源的調(diào)研與分析［J］.地震地磁觀測(cè)與研究，2010，31（5）：65－72.