基于OpenCV的模擬地震記錄矢量化軟件

史鵬飛+張翔+豐繼林+王茂發(fā)+韓定良

摘 要：隨著地震觀測技術(shù)與手段的不斷更新，數(shù)字化地震觀測逐漸取代傳統(tǒng)的模擬地震觀測。為了保存歷史數(shù)據(jù)，如何有效提取模擬地震記錄中的珍貴信息是地震工作者的研究重點。該項目提出了基于OpenCV的模擬地震記錄矢量化算法，使用C#開發(fā)相應(yīng)軟件。該軟件利用有震波形區(qū)域提取方式可以快速、準(zhǔn)確地將模擬地震記錄矢量化。

關(guān)鍵詞：模擬地震記錄；數(shù)字化；OpenCV；有震波形區(qū)域提??；波形追蹤

0 引言

世界上第一張地震圖要追溯到1889年，英國人米爾恩和尤因利用安置在德國波茨坦的現(xiàn)代地震儀記錄下了發(fā)生在日本的一次地震。此后70年的時間里，人類觀測的地震信息都是以模擬地震記錄的方式存在。模擬地震記錄蘊(yùn)含著大量寶貴的地震信息，具有不可估量的價值。直到19世紀(jì)70年代，美國率先建成了由十個數(shù)字化臺站組成的數(shù)字化臺網(wǎng)。隨后數(shù)字地震記錄技術(shù)逐漸取代模擬地震記錄技術(shù)，一直發(fā)展沿用至今。這造成了地震數(shù)據(jù)應(yīng)用研究中，模擬地震信息與數(shù)字地震信息之間有一條明顯的數(shù)據(jù)斷裂。如今大多數(shù)模擬地震記錄存放在倉庫里經(jīng)受著自然的侵蝕，蘊(yùn)藏著寶貴信息的模擬地震記錄無法被計算機(jī)處理。若能有效提取這些模擬地震信息，將會推動地震研究工作快速發(fā)展。

1 研究現(xiàn)狀

目前為止，對于模擬地震記錄的研究主要分為兩個階段：存儲與矢量化。

存儲是將模擬地震記錄中的信息由易損壞的紙質(zhì)信息轉(zhuǎn)化為方便存儲與備份的數(shù)字圖片信息。2012年開始，哈佛大學(xué)投入大量精力進(jìn)行相關(guān)模擬地震記錄電子存檔工作，該項計劃被命名為哈佛大學(xué)地震波形歸檔計劃。2015年，河北地震局蔣宏毅等人研發(fā)了一套模擬地震記錄數(shù)字化存儲管理系統(tǒng)，實現(xiàn)了模擬地震記錄的數(shù)字化存儲與查詢。但數(shù)字存儲只是將紙質(zhì)模擬地震記錄掃描為數(shù)字圖片，沒有進(jìn)行矢量化，計算機(jī)仍無法識別和處理模擬地震記錄，也無法挖掘記錄中的寶貴信息‘5]。

許多科研工作者進(jìn)行模擬地震記錄矢量化的研究：Teves-Costa（ 1999）使用CADcore和AutoCAD開展模擬地震記錄的矢量化工作；Pintore（2005）基于人工神經(jīng)網(wǎng)絡(luò)提出了Teseo方法，可以手動和自動矢量化模擬地震記錄；徐濤等（2014）利用Matlab GUI開發(fā)了一個交互的模擬地震矢量化程序；王茂發(fā)等（2016）提出了CSF算法，并開發(fā)了相應(yīng)的程序。雖然矢量化研究工作取得了顯著成績，但許多矢量化工作只是將一小部分含有復(fù)雜波形的圖片截取出來單獨(dú)進(jìn)行矢量化研究，脫離了整張圖片使得后期波形拼接與時間拼接變得難以實現(xiàn)。 鑒于跨平臺計算機(jī)視覺庫OpenCV在人臉識別、車牌識別、動態(tài)手勢檢測等方面的準(zhǔn)確性與高效性，因此，本文基于OpenCV設(shè)計了新的模擬地震記錄矢量化算法，并使用程序設(shè)計語言C#開發(fā)軟件。其中包含分塊提取波動區(qū)域、自動提取復(fù)雜波等特有的功能模塊，能夠快速、準(zhǔn)確地處理整張模擬地震記錄。

2 軟件矢量化過程

2.1 流程

本軟件實現(xiàn)模擬地震記錄矢量化的過程共分為8個步驟（圖1所示）：加載地震圖紙、圖紙二值化、搜索起點、有震波形區(qū)域提取、掃描波動區(qū)域、判斷波形類型、拼接和反演。下面給出每一個步驟的具體描述。

2.2加載地震圖紙

將整張圖紙加載至矢量化軟件中是模擬地震記錄矢量化最基本的工作。早期的模擬監(jiān)測是將一段時間（一般為1天）內(nèi)3個方向（東西、南北、上下）的地震波動記錄在一張圖紙上。觀察整張圖（如圖2）紙可知，波形密集，包含的信息量巨大，該項目使用C#常用控件PictureBox實現(xiàn)了軟件的圖片加載功能，可以加載、處理整張模擬地震記錄。

2.3 圖紙二值化

許多因素如光線、圖紙的新舊程度等，會對模擬地震記錄掃描而成的圖片產(chǎn)生影響，所以必須對圖紙進(jìn)行二值化處理，它是在圖片預(yù)處理階段減少圖片噪聲非常有效的方式。一張灰度的光柵圖片類似一個二維的數(shù)組，每一個元素在數(shù)組中的值在0-255之間。利用OpenCV的Threshold函數(shù)對圖片進(jìn)行二值化處理：設(shè)置一個閾值（目前實驗獲得經(jīng)驗數(shù)值為155），將大于閾值的像素點置為255，小于閾值的置為0，如圖3所示：

2.4搜索起點

為了確保數(shù)據(jù)的一致性，軟件在提取模擬地震記錄中的信息時，始終在同一個坐標(biāo)系下處理。以圖片的左上角為原點，圖片的上邊緣為x軸正值方向，圖片的左邊緣為y軸正值方向。取原點附近的一個點c（x，0）為起點，沿著垂直方向向下搜索，每次步長為一個像素點。搜素算法具體描述如下：

（I）水平方向

水平方向x值保持不變（公式1），為了避免圖片左側(cè)空白邊緣對搜索結(jié)果的影響，一般x取值稍大于0（實驗取值為10）。 X_i=x（i=l，2…m）

（1）

（2）垂直方向

垂直方向y值由0不斷變大，每次以一個步長即一個像素向下搜索，并判定像素值的變化。當(dāng)遇到相鄰兩個像素點值跳變時記錄這個像素點坐標(biāo)到數(shù)組中（x，Y_n1），然后繼續(xù)向下搜索，再次遇到跳變時繼續(xù)記錄（x，Y_n2）…。

則第一條波的起點坐標(biāo)為公式（2）：

直到搜索y值到圖紙的下邊緣，數(shù)組中記錄的點（x，Y_n1），（x，Y_n2），……，（x，_nm）。（點的編號為1，2…m）這樣每一條波形起點的信息都可以確定。第h（l≤h≤m/2）條波的起點坐標(biāo)信息為（公式3）：

2.5有震波形區(qū)域提取

模擬地震記錄中最重要的數(shù)據(jù)是記錄地震發(fā)生時的復(fù)雜波，矢量化工作的難點是矢量化這震波形區(qū)域的波形。通過觀察發(fā)現(xiàn)這些復(fù)雜波只是整張圖紙很小的一個局部，所以采取有震波形區(qū)域提取的方式：先找到復(fù)雜波所在的第一個區(qū)域，在這個波動區(qū)域內(nèi)進(jìn)行矢量化操作，接著尋找復(fù)雜波所在的第二個區(qū)域進(jìn)行矢量化操作，直至完成所有波動區(qū)域的波形矢量化（圖4）。所有的矢量化工作基于圖紙的局部卻沒有脫離整體圖紙，確保了波形信息與時間信息拼接操作的可行性。如果在整張模擬地震圖紙上進(jìn)行矢量化復(fù)雜波的操作，這樣既增加難度，也降低效率和準(zhǔn)確率。endprint

2.6掃描波動區(qū)域

2.6.1 判斷波形類別

首先，使用搜索起點算法將波動區(qū)域的所有波形起點找到。然后，按照Y值由小變大的順序以每條波的起點水平掃描這條波。掃描算法具體描述如下：

（1）選取在此區(qū)域內(nèi)第一條波的起點（X_c1，y_c2）。

（2）水平方向，X_C1不變。垂直方向，以一個像素點為步長，先沿著v軸正方向搜索至（ X_C1，Y_C1+△y） （Ay最大為波形寬度的平均值，若在此區(qū)域內(nèi)出現(xiàn)像素跳變點則記錄其坐標(biāo)（ X_c1.y_D1），若不出現(xiàn)跳變點則記錄點（ X_C1，y_D1）=（X_C1，Y_C1+△y）；后沿著y軸負(fù)方向搜索至（X_C1，Y_C1-△y），若在此區(qū)域內(nèi)出現(xiàn)像素跳變點則記錄其坐標(biāo)（X_C1，y_D2），若不出現(xiàn)跳變點則（X_C1，y_D2）=（X_C1，Y_C1-△y）。然后記錄在（X_C1，Y_C1±△y）范圍內(nèi)中點坐標(biāo)同時記錄該點的像素值。

（3）水平方向：X_C2=X_C1+1。

（4）重復(fù)（2-3）操作直至波動區(qū)域邊界。

使用掃描算法掃描平滑波時，由于Ay的限制，可以基本消除與復(fù)雜波交叉的影響，追蹤到平滑波的中線，所以記錄下來平滑波中線的像素基本都為黑點。但由于復(fù)雜波是有振幅的，而△_y值最大為線寬，所以掃描復(fù)雜波時許多中點的像素點為白色（圖5）。因此，當(dāng)掃描一條波形記錄下的白色像素點多于其他波形時，則判定該波為復(fù)雜波同時進(jìn)行標(biāo)記。

2.6.2 追蹤平滑波

判斷波形類型后，根據(jù)標(biāo)記信息進(jìn)行平滑波的跟蹤。利用改進(jìn)的掃描算法追蹤平滑波：掃描時，若在區(qū)域（y+△y，y-△y）}H現(xiàn)像素跳變點，記錄此范圍內(nèi)的中點并將此范圍的所有黑色像素點置為白色像素點，反之則不變。這樣既保證記錄平滑坡中線坐標(biāo)的正確性，也不會破壞復(fù)雜波的波形，追蹤效果如圖5所示：

2.6.3 去除噪聲

平滑波追蹤完成后，復(fù)雜波追蹤變得更加簡單。但儀器記錄和保存時各種因素可能會導(dǎo)致圖片存在干擾信息，為了避免這些噪聲對追蹤復(fù)雜波的影響，必須再次進(jìn)行去噪處理。首先，利用OpenCV對圖片進(jìn)行SmoothMedian處理，消除椒煙噪聲與斑點噪聲。其次，利用FindContours函數(shù)提取所有圖形的輪廓，記錄每一個輪廓的id，計算每一個輪廓的周長。最后，通過比較保留周長最長的輪廓，其余全部去除——這樣就保存了復(fù)雜波信息而去除了噪聲，圖6是保留的復(fù)雜波輪廓。

2.6.4追蹤復(fù)雜波

如同上面一樣，可以計算出復(fù)雜波起始處的中點坐標(biāo)，然后以此點為起點追蹤復(fù)雜波。在追蹤過程中，值得注意的是，在波峰和波谷處需要記錄此處的v值上邊緣值與下邊緣值，這樣可以正確反映地震的能量釋放情況。追蹤復(fù)雜波的算法與追蹤平滑波的算法相似，不同之處在于：沒有△y的限制，當(dāng)中間點的y值大于相鄰兩點的v值時（波谷，v軸正方向向下），記錄其下邊緣的v值；當(dāng)中間點的v值小于相鄰兩點的v值日寸（波峰），記錄其上邊緣的v值；搜素到區(qū)域邊界，圖7中紅線為波形的中線，藍(lán)點為提取的關(guān)鍵點。

2.7 拼接和反演

在前期工作中已經(jīng)實現(xiàn)了拼接與反演，這里只做簡單描述。每條波形上的點都以坐標(biāo)的形式存儲在數(shù)組中并且標(biāo)有唯一的id，首先根據(jù)點的坐標(biāo)信息進(jìn)行波形的拼接，然后根據(jù)圖紙上標(biāo)注的時間信息進(jìn)行時間的拼接，最后將波形信息反演，如圖8所示：

3 結(jié)束語

該項目設(shè)計了基于OpenCV的模擬地震記錄矢量化軟件，分8個步驟矢量化模擬地震記錄的柵格圖片，通過給出的測試結(jié)果可以看出，反演出來的波形較好的符合原始圖紙中的波形。下一步，將進(jìn)一步提升算法效率和改善軟件，期望開發(fā)出自動校驗?zāi)K以確保矢量化的質(zhì)量。endprint