洪旭瑜,林苗祿,關(guān)玉梅,李 強

(福建省地震局,福州 350003)

0 引言

隨著信息化建設(shè)日益成熟,依托計算機軟件編程、WebGIS、數(shù)據(jù)庫以及實時融合等技術(shù)實現(xiàn)定制需求包括自動產(chǎn)出在很大程度上都能夠?qū)崿F(xiàn)[1-3]。林巖釗等[4]基于地震預(yù)警與烈度速報系統(tǒng)實現(xiàn)自動產(chǎn)出更接近現(xiàn)實的地震專題圖;張博[5]利用 MATLAB自動產(chǎn)出地震活動性資料;劉堅等[6]采用自編程序,實現(xiàn)地震應(yīng)急與科學產(chǎn)品后臺自動加工處理與發(fā)布。因此,地震應(yīng)急及科技產(chǎn)品能夠由軟件程序?qū)崿F(xiàn)的,讓人工操作向自動化轉(zhuǎn)變是未來發(fā)展趨勢。

福建已建成前兆臺網(wǎng)技術(shù)觀測系統(tǒng),圖1顯示福建省前兆臺網(wǎng)各學科觀測臺站,到2018年底共有270套左右前兆觀測儀器,包括形變、重力、流體、電磁等學科的儀器,儀器觀測測項數(shù)達500個左右。隨著前兆儀器數(shù)據(jù)量及信息日益增多,這使得日常和應(yīng)急會商需花費較多人工精力對前兆曲線繪制及質(zhì)量分析,且圖件繪制無統(tǒng)一規(guī)格標準,而地震應(yīng)急會商時要在規(guī)定時間內(nèi)快速給出初步的震情專報及會商意見。因此,除震后趨勢判定外,通過統(tǒng)計并計算震中一定范圍內(nèi)相關(guān)歷史地震活動性資料,按照相對固定格式自動產(chǎn)出初步震情專報,同自動繪制產(chǎn)出的規(guī)范化前兆觀測數(shù)據(jù)曲線及異常干擾情況說明相結(jié)合,為震后趨勢判定提供參考,并可應(yīng)用于日常會商資料產(chǎn)出和各類圖件報告的制作。

圖1 福建省前兆臺網(wǎng)各學科觀測臺站Fig.1 Observation stations of all discipline of Fujian precursory network

Java發(fā)展到現(xiàn)在已經(jīng)很成熟,它主要由編程語言 (Java語言)、運行環(huán)境 (JVM:Java虛擬機)、框架 (Java API群)組成[7],具有面向?qū)ο?、可移植性、高性能、多線程、動態(tài)性等特點,相比MATLAB,JAVA在可視化界面、兼容性和跨平臺性具有較大的優(yōu)勢。MATLAB在算法開發(fā)、數(shù)據(jù)可視化、數(shù)值計算等高級技術(shù)計算語言和交互式環(huán)境擁有強大的功能。因此,自動產(chǎn)出程序設(shè)計在面向?qū)ο罂梢暬故炯八惴ǔ绦蛟O(shè)計語言分別選擇JAVA技術(shù)和MATLAB編程語言,充分利用二者之間的優(yōu)勢。

考慮用戶使用的便利性,采用JAVA調(diào)用MATLAB程序,JAVA作為可視化界面進行設(shè)計,通過JVM實現(xiàn)平臺無關(guān)性（可同時兼容32位及64位 Window7-10系統(tǒng)及各版本 Linux操作系統(tǒng)）,無需擔心用戶使用的操作系統(tǒng)版本,便于系統(tǒng)的安裝使用。操作系統(tǒng)在沒有安裝 MATLAB相關(guān)版本是無法使用內(nèi)置于JAVA的MATLAB程序,但是讓每個用戶單獨安裝較大文件 MATLAB程序顯然不太現(xiàn)實,采用將MATLAB程序代碼打包成 jar文件內(nèi)置到Java作為一個包,安裝 MATLAB軟件提供的精簡版MCR（MATLAB Compiler Runtime）軟件即可實現(xiàn)JAVA與MATLAB之間的代碼共享兼容。通過JAVA可視化界面調(diào)用MATLAB程序,便捷、高效產(chǎn)出美觀、規(guī)范化、標準化的會商資料。

1 測震信息產(chǎn)出及其計算原理

1.1 地震活動性參數(shù)

地震發(fā)生后,震后趨勢判定需要結(jié)合地震的震中位置,震中一定范圍內(nèi)歷史地震分布情況,歷史地震序列類型,序列類型時空分布情況等信息進行綜合分析。通過地震速報信息系統(tǒng)（EQIM）正式報與省級地震臺網(wǎng)速報內(nèi)容相結(jié)合,選擇合適的地震目錄,算出震區(qū)周邊歷史中強震活動性情況。

1.2 相關(guān)參數(shù)計算

震中附近歷史地震分布情況：需要產(chǎn)出震中附近50 km和100 km范圍內(nèi)歷史地震分布情況：計算震級在3.0～3.9級、4.0～4.9級、5.0～5.9級、6.0～6.9級以及7.0級以上地震個數(shù),要得出50 km和100 km范圍內(nèi)發(fā)生的最大地震的時間、地點以及震級。因此需要計算歷史地震與當前震中的震中距。震中距有多種算法,不同算法精度不一樣,其中精度較高、采用較廣的算法主要有Robbins法、Gauss法以及地心緯度法,Robbins法、Gauss法對震中緯度（高緯度或是低緯度）以及方位角的要求較高,而計算800 km以內(nèi)的震中距,受方位角影響較小的一般采用地心緯度法,獲取兩個點的經(jīng)緯度就可得出震中距,公式如下[8]：

式中：lat2為歷史地震震中緯度,lat1為當前地震震中緯度;long2為歷史地震經(jīng)度,long1為當前地震經(jīng)度,D為震中距。

歷史地震序列類型和序列類型時空分布情況計算原理同歷史地震分布情況,這里不再敘述。

1.3 自動觸發(fā)

自動觸發(fā)是地震發(fā)生后以EQIM系統(tǒng)產(chǎn)出的正式報震級為準,當達到設(shè)定的震級,系統(tǒng)自動觸發(fā)快速計算預(yù)先設(shè)定好的相關(guān)地震活動性參數(shù)和信息,繪制震中附近一定范圍內(nèi)的長期、中期、短期前兆觀測數(shù)據(jù)曲線以及對應(yīng)的異常變化數(shù)據(jù)曲線情況說明,然后將計算結(jié)果和繪制好的曲線及異常曲線情況說明以PowerPoint、Word等報告形式產(chǎn)出。判斷自動觸發(fā)的要素為地點和震級,二者同時滿足條件即可觸發(fā)。例如當福建省行政區(qū)內(nèi)發(fā)生M3.0（這里不區(qū)分ML和MS）級及以上地震,自動觸發(fā)條件滿足。借用MATLAB函數(shù)庫timer函數(shù)定時刷新功能定時監(jiān)測EQIM系統(tǒng)正式報信息,假設(shè)LastTime為EQIM最新正式報地震的發(fā)震時間,BeforeTimew為本地存儲前一次地震發(fā)震時間,則監(jiān)測過程如圖2。

t=timer('TimerFcn','m文件','Period',循環(huán)間隔,'ExecutionMode','fixedSpacing','TasksTo Execute',inf);%創(chuàng)建自動觸發(fā)定時器;

start(t);%開始執(zhí)行

由于前兆儀器觀測測項較多,自動繪制長中短期需要一定的時間,因此需要考慮整個程序執(zhí)行完一次需要耗時多長,根據(jù)測試結(jié)果需要4.5 min左右,循環(huán)間隔要大于運行完成時間,同時兼顧時效性,因此設(shè)定為5 min。循環(huán)次數(shù)需設(shè)置為無限次循環(huán)inf。

圖2 會商資料自動觸發(fā)監(jiān)控過程Fig.2 Monitoring process automatically triggered by the consultation information

2 前兆觀測信息產(chǎn)出

當前繪制前兆觀測曲線并形成 PPT會商報告主要是通過EIS2000（地震前兆信息處理與軟件系統(tǒng)）和中國地震前兆臺網(wǎng)數(shù)據(jù)處理系統(tǒng)[9-10]篩選臺站和測項繪制曲線,并查找相應(yīng)時間段的觀測曲線數(shù)據(jù)異常干擾情況說明,該方法需耗費較多的人工精力,在產(chǎn)出圖件方面可能存在曲線坐標格式、美觀度及標準化等方面不統(tǒng)一。MATLAB在繪圖方面的強大功能可以較好的解決該問題。

2.1 數(shù)據(jù)庫基礎(chǔ)信息

福建前兆臺網(wǎng)基礎(chǔ)信息數(shù)據(jù)庫采用的是分布式Oracle數(shù)據(jù)庫,其遵循《地震前兆數(shù)據(jù)庫結(jié)構(gòu)》（DB/T 51-2012）設(shè)計規(guī)范[11]。MATLAB自動繪制觀測數(shù)據(jù)曲線及獲取異常干擾情況說明需要獲取Oracle數(shù)據(jù)庫相關(guān)觀測數(shù)據(jù)表、數(shù)據(jù)字典信息及日志表,周克昌等[11]對地震前兆數(shù)據(jù)庫系統(tǒng)及其相關(guān)表信息已經(jīng)做了較為詳細的介紹,這里不在贅述。

2.2 數(shù)據(jù)庫訪問機制

MATLAB訪問Oracle數(shù)據(jù)庫有三種方式：

（1）OLE-DB (Ob ject Lin king and Embedding DataBase)。

（2）ODBC(Open DataBase Connectivity);

（3）JDBC (Java Data Base Connectivity)。

本文采用MATLAB作為后臺程序供Java調(diào)用,因此采用JDBC訪問機制。Java的Ojdbc14.jar文件提供MATLAB橋接Oracle的作用,將該文件拷貝到$HOMER2018(因版本可能會有所差別)javajar oolbox中,然后打開classpath.txt文件添加語句$matlabroot/java/jar/toolbox/ojdbc14.jar。連接語句如下,連接成功狀態(tài)見圖3:conn=database(數(shù)據(jù)庫實例名,用戶名,密碼,'oracle.jdbc.driver.OracleDriver',['jdbc:oracle:thin:@',IP,':1521:'])。

圖3 Matlab成功連接Oracle數(shù)據(jù)庫狀態(tài)Fig.3 Status of MATLAB successfully connecting to Oracle Database

2.3 曲線自動成圖

參數(shù)信息：繪制曲線需要臺站名稱、臺站代碼、測項代碼、采樣率、測項分量名稱以及測項單位等六個參數(shù)信息,形成單獨M文件。曲線自動繪制是通過對這些參數(shù)信息的循環(huán)讀取依次獲取觀測數(shù)據(jù)自動實現(xiàn)繪制。MATLAB從Oracle獲取觀測數(shù)據(jù)需要解決兩個問題：分鐘值和秒值觀測數(shù)據(jù)字段是CLOB類,無法直接獲取;原始觀測數(shù)據(jù)和預(yù)處理數(shù)據(jù)個別時段有出現(xiàn)缺記,數(shù)據(jù)庫中數(shù)值缺記是用“NULL”或“NULLALL”字符表示,MATLAB是無法識別除NAN字符外其他表示數(shù)值的字符。對應(yīng)的解決辦法是利用MATLAB提供的getSubString子類將oracle.sql.CLOB轉(zhuǎn)換為String類型,用regexprep函數(shù)將“NULL”和“NULLALL”全部替換成“NAN”。

曲線成圖：自動繪制后的曲線是系統(tǒng)默認模式,需要按我們需求對橫坐標和縱坐標刻度值、縱坐標單位以及橫坐標標題的標注進行規(guī)范,規(guī)范后曲線如圖4。

（1）橫坐標刻度：根據(jù)數(shù)據(jù)類型標注時間格式,原始數(shù)據(jù)和預(yù)處理數(shù)據(jù)時間精確到日,時間格式為“YYYY-MM-DD”,整點值和日值數(shù)據(jù)時間精確到月,時間格式統(tǒng)一為“YYYY/MM”,根據(jù)figure大小選擇自適應(yīng)模式等間隔標注刻度。

（2）橫縱表標題：標題統(tǒng)一為“臺站名稱-觀測測項分量名稱（起始日期-結(jié)束日期）”,其中日期格式為“YYYYMMDD”。

（3）縱坐標刻度：刻度標注根據(jù)figure大小自適應(yīng)模式自動標注,需要值得注意的是靜水位數(shù)值表示的是井口距離探頭的數(shù)據(jù),靜水位數(shù)值上升說明水位下降,反之水位上升,體應(yīng)變觀測同水位觀測原理一致,觀測值上升表示向壓縮狀態(tài)轉(zhuǎn)變,觀測值越小表示向拉張狀態(tài)轉(zhuǎn)變。因此,除了體應(yīng)變觀測和靜水位測項縱坐標刻度值為反方向,其他測項的縱坐標刻度值離原點越遠值越大。

（4）縱坐標單位：MATLAB縱坐標單位默認顯示在左側(cè)中間,規(guī)范統(tǒng)一為縱坐標正上方,語句為：text(min(xdata),axis(4),unit,'HorizontalAlign ment','center','fontsize',15)。

圖4 規(guī)范化觀測曲線圖件Fig.4 Curve diagram of standardized observation

3 結(jié)果產(chǎn)出

根據(jù)對地震活動性參數(shù)的計算結(jié)果以及曲線自動繪制成圖結(jié)果分別產(chǎn)出初步的震情專報和前兆各學科觀測曲線及其異常干擾情況說明分析報告。Office應(yīng)用程序套件具有豐富的功能,調(diào)用 Office VBA方法可以快速實現(xiàn)對文檔、電子表格和演示文稿進行創(chuàng)作、格式設(shè)置和操縱。通過 MATLAB調(diào)用 office VBA的對象模型、屬性和方向?qū)崿F(xiàn)對 PPT和 WORD模塊的設(shè)計,實現(xiàn)報告的自動化產(chǎn)出。

3.1 震情專報產(chǎn)出

初步的震情專報自動產(chǎn)出是通過MATLAB調(diào)用本地WORD服務(wù)器,以預(yù)先設(shè)定好的相對固定格式和內(nèi)容形成的模板產(chǎn)出得到專報內(nèi)容,假設(shè)2019年10月25日福建莆田仙游縣發(fā)生3.7級地震,模擬自動快速產(chǎn)出初步震情專報,產(chǎn)出結(jié)果見圖5。

圖5 模擬震情專報自動產(chǎn)出結(jié)果—以2019年10月25日福建莆田仙游縣3.7級地震為例Fig.5 Automatic output results of simulated earthquake situation special report --taking the Xianyou M3.7 earthquake in Putian,Fujian Province on October 25th,2019 as an example

3.2 觀測曲線及異常干擾情況說明產(chǎn)出

MATLAB調(diào)用PowerPoint服務(wù)器,將繪制完成觀測曲線及曲線異常干擾質(zhì)量情況說明自動導入PPT文件,一張幻燈片版式由標題和繪制好的觀測曲線組成,標題以“臺站名稱-儀器名稱-短（中或長）期”命名,其中短期為曲線時長小于3個月,數(shù)據(jù)為原始采樣值（在沒有特殊要求下默認采用觀測儀器原始產(chǎn)出值：秒值、分鐘值、整點值、日值,下同。）,中期為時長大于三個月小于6個月,數(shù)據(jù)一般采用整點值,長期為時長大于6個月以上,數(shù)據(jù)一般采用日值,備注頁說明曲線異常干擾情況,導出結(jié)果見圖6：

圖6 觀測曲線及異常干擾情況說明自動導入PPT結(jié)果Fig.6 Observation curve and abnormal interference description automatically imported into PPT

4 結(jié)語

本文主要利用MATLAB軟件在數(shù)據(jù)處理和計算能力方面的強大功能,能夠快速計算相關(guān)參數(shù)結(jié)果及占用較少的內(nèi)存獲取較大字節(jié)的地震前兆觀測數(shù)據(jù),保證部署在服務(wù)器端實時監(jiān)測地震發(fā)生后快速判斷是否滿足條件并自動產(chǎn)出相應(yīng)可靠的資料,爭取應(yīng)急時間,采用JAVA對程序內(nèi)容進行可視化,滿足日常會商對圖件和資料個性化產(chǎn)出的需求,極大的提高工作效率。