OpenACC 平行化語言用於網格式定位程式之開發及其在地震預警之應用

國立臺灣大學理學院地質科學研究所 碩士論文

Department of Geosciences College of Science

National Taiwan University Master Thesis

OpenACC 平行化語言用於網格式定位程式之開發及 其在地震預警之應用

A development of grid search earthquake location program using OpenACC parallel computing toolkit and

its application to earthquake early warning

陳璽安 Hsi-An Chen

指導教授：吳逸民 博士

Advisor: Yih-Min Wu, Ph.D.

中文摘要

對於地震預警系統來說，在地震發生後能快速且準確的估計震源位置和規模 是相當重要的項目。在地震定位方面，使用簡化的一維速度構造模型對於地質構 造較複雜之區域無法得到良好的震源位置。因此本研究之目的在於利用三維速度 構造模型建立一套準確且計算過程簡單的地震定位程式。首先將台灣地區切成密 集的三維空間格點，格點間距在經度、緯度和深度上皆為 1 公里，總格點數為 6.88 百萬，並假設各點皆為可能之震源位置，利用三維速度構造模型與三維空 間的波線追跡法來計算理論 P 波走時，並建立各格點至中央氣象局地震速報觀測 網 (Real-Time Data stream, RTD) 之地震測站的理論走時資料庫。當地震發生 時，將接收到地震 P 波訊號的測站兩兩配對建立等時差 (Equal differential time, EDT) 曲面，再以逐點計算殘差 (residual) 總和，找出最小殘差值的格 點即為震源位置。為了處理運算過程中產生的大量數據，本研究利用 OpenACC 平行化套件發揮圖形運算處理器 (GPU) 的運算優勢來減少計算過程需要花費的 時間。利用 TAIGER 炸測資料及 115 個歷史地震事件對於定位程式的可靠度進行 測試，誤差分析顯示本研究之定位方法能在不增加運算時間的情況下，有效提升 地震定位的精度，因此能進一步改善地震預警系統之執行效率。相較於中央氣象 局地震目錄之定位結果，震央平均差異為 4.77 公里，深度差異為 5.11 公里。

關鍵詞:地震、地震即時定位、地震預警

Abstract

The rapid and reliable hypocenter and magnitude estimation are two of the most issue for earthquake early warning (EEW) system after an earthquake occurs. In earthquake location aspect, for early warning purpose, the simplified one-dimensional (1-D) velocity model cannot reflect the complex structures in certain area and result in large location uncertainties. Therefore, an accurate location program based on three-dimensional (3-D) velocity model with a simple-calculating process is developed for earthquake location. First, we divided Taiwan region into 3-D grids, 1 km interval in longitude, latitude and depth, which are assumed as possible hypocenters. Thus, the study area consists of 6,880,000 grids. P wave travel times between each grid and seismic station are stored in a pre-calculated database using the 3-D velocity model and 3-D ray tracing method. When an earthquake occurs, an equal differential time (EDT) surface is then conducted in each different triggered station pairs. By comparing the residuals between EDT surfaces and each grid, the grid with the minimum residual will be considered as the hypocenter. In order to deal with massive computation load during the comparing process, we apply an OpenACC parallel computing toolkit to GPU (Graphic Processing Unit) to reduce the processing time. Records of 115 earthquakes and 10 TAIGER explosion events are used to test our location program performance. Results show our location program can provide a reliable hypocenter location in short times and thus it could improve the efficiency of the EEW application. The average location difference between result of this study and Central Weather Bureau (CWB) earthquake catalog for epicenter and depth are 4.77 km and 5.11 km respectively.

Keywords: earthquake, earthquake real time location, EEW

目錄

論文口試委員審定書 ... I 中文摘要 ... II

Abstract ... III

目錄 ... IV 圖目錄 ... VI 表目錄 ... VII

第一章 緒論... 1

1.1 前言 ... 1

1.2 動機和目的 ... 2

第二章 文獻回顧 ... 4

2.1 地震預警歷史與台灣地震預警系統發展 ... 4

2.2 地震預警之地震定位方法 ... 7

2.3 現行中央氣象局地震預警系統之地震定位流程 ... 12

第三章 研究方法 ... 13

3.1 本研究之地震定位流程 ... 13

3.2 測站理論走時資料庫建立 ... 14

3.3 OpenACC 程式平行化 ... 16

第四章 結果與討論 ... 18

4.1 TAIGER 人工炸測震源資料測試 ... 18

4.1.1 地震定位結果 ... 19

4.1.2 定位結果分析 ... 19

4.2.1 地震定位結果 ... 26

4.2.2 定位結果分析 ... 26

4.3 案例討論：2016 年 2 月 6 日高雄市美濃區芮氏規模 (ML) 6.6 地震 ... 32

4.4 地震定位時效討論 ... 34

第五章 結論 ... 35

參考文獻 ... 36

附錄 A 本研究所使用之 GPU 規格 ... 39

圖目錄

圖 1.1 台灣歷史災害地震震央分布 (修改自鄭世楠等，1999) ... 3

圖 2.1 美國科學家 Coopper 提出之地震預警概念圖 (Cooper, 1868) ... 5

圖 2.2 台灣第一個裝設在花蓮地區的地震預警系統雛形 (Wu et al., 1999) ... 6

圖 2.3 虛擬子網地震預警系統虛擬子網地震預警系統 (Wu and Teng, 2002) ... 6

圖 2.4 等時差 (EDT) 曲面示意圖 ... 11

圖 2.5 以平面圖呈現等時差 (EDT) 曲面定位概念 ... 11

圖 3.1 臺灣地區的三維網格點範圍與 RTD 測站 (106 站) 分布圖 ... 15

圖 3.2 程式平行化計算流程 ... 17

圖 4.1 人工炸測震源位置 ... 20

圖 4.2 各個炸測點的震央距、深度和震源距誤差的重定位結果 ... 21

圖 4.3 各個炸測點的震央距、深度和震源距誤差的重定位結果 ... 22

圖 4.4 測站最大空餘角 (GAP) 值與震央距、深度和震源距誤差的比較 ... 23

圖 4.5 N2 和 S1 炸測點位的重定位結果以及使用之 RTD 測站分布 ... 24

圖 4.6 N2 和 S1 炸測點位的重定位結果以及使用之 RTD 測站分布，將 S1 炸測點位的測 站數目減少至與 N2 炸測點位相同 ... 25

圖 4.7 歷史地震資料重定位結果之比較 ... 28

圖 4.8 震央位置和深度之誤差比較 ... 29

圖 4.9 震源位置誤差與測站空餘角 (GAP) 之比較 ... 30

圖 4.10 測站空餘角 (GAP) 介於 180^o至 230^o之間的島外地震事件 ... 31

圖 4.11 2016 年 2 月 6 日美濃地震報告 ... 32

圖 4.12 美濃地震事件的定位結果 (第 1 報至第 12 報) ... 33 圖 4.13 本研究之定位程式於 GPU 平行化前後在不同的測站數目下得出震源位置所需要

表目錄

表 4-1 TAIGER 人工炸測震源位置、施測時間、紀錄規模和測站 (RTD) 數目 ... 20

第一章 緒論

1.1 前言

台灣本島位處於歐亞大陸板塊與菲律賓海板塊交接處，板塊碰撞劇烈，造就 了複雜的地質環境以及頻繁的地震活動。自從台灣 20 世紀初期有地震儀器觀測 以來，已經歷過多次災害性地震，如 1906 年梅山地震、1935 年新竹-台中地震、

1964 年白河地震、1986 年花蓮地震、1999 年 921 集集地震和 2016 年的美濃地 震等(如圖 1.1) ，均在各地造成不同的人員傷亡和財產損失。迄今的地震學研究 對於孕震機制 (如震源破裂過程) 已有相當多的探討，但是地震發生的原因仍受 到許多複雜的因素影響。在未來，由於板塊聚合作用，應力持續累積，具有災害 性規模的地震仍有機會發生，加上現代化發展造成人口分布密集，無形中增加地 震災害的威脅性，因此地震預警 (Earthquake Early Warning, EEW)成為地震防災 中相當重要的項目之一。

地震預警的目的係為當地震發生後，在強烈地表振動到達之前，透過初達波 資訊快速的評估地震資訊來對目標地區發出警報。一般地震發生時，地震波依照 傳播路徑可分為兩類：可在地球內部傳遞的體波 (body wave)和只能在地球表面 傳遞的表面波 (surface wave)。而體波又可依其波動性質分為介質振動方向平行 波傳方向的 P 波 (primary wave) 和振動方向垂直波傳方向的 S 波 (shear wave or secondary wave)，其中以 P 波速度最快，其次為 S 波，最後是表面波。經斷層錯 動產生的地震波其大部分的能量是藉由 S 波和表面波傳遞，因此 S 波和表面波 為造成災害的主要原因。由於地震 S 波的震幅資訊可由初達 P 波數秒訊號估計 求得 (Kanamori, 2005；Wu and Kanamori, 2005)。透過此特性，地震預警系統可 在地震發生後，一旦確認此地震具有致災特性便可對外發出警報，讓民眾或是大 眾運輸系統和工廠等設施有數秒至數十秒的反應時間採取相關措施以減少災

E(𝐫, T₀) = ∑(w_jdt_j)²

J−1

j=0

, dt_j = [T_j− T₀] − τ_j(𝐫)

1.2 動機與目的

決定地震規模和發生位置是地震預警系統中當重要的兩項參數 (Allen and Kanamori, 2003)，因此，快速且有效的地震定位程序是系統中不可或缺的元素。

本研究將針對地震定位做討論。在地震定位方面，目前中央氣象局 (Central Weather Bureau, CWB) 的地震預警系統所採用的定位方法參考自 HYPO71 程式 (Lee and Lahr, 1972) 其 使 用 線 性 逆 推 方 法 (Geiger, 1912) 和 一 維 半 空 間 (half-space) 連續速度模型 (Chen et al, 2015)。然而此方法存有三項缺點會造成地 震定位誤差，首先，台灣島位於兩個板塊隱沒系統之間，地質條件複雜，進而反 映出台灣的地下速度構造有相當程度的側向不均勻性。因此，簡化的一維速度模 型無法準確的反映地震波從震源位置傳遞至測站所需要的走時時間，進而影響地 震定位精準度。

其次，為震源深度和發震時間這兩個參數的解析度，由於傳統線性定位方法 係為求解觀測值與理論值之間的殘差值目標函數 E(𝐫, T₀)的最小值 (式 1-1)，

其中dt_j為走時殘差，w_j 為權重函數，r 為假設之震源位置，τ_j(𝐫) 為地震波理論 走時，Tj和 To分別為地震波觀測到時和假設的初始發震時間。由於 E(𝐫, T₀) 為走 時殘差的分散度，因此震源位置和發震時間這兩個參數之間存有互補效應 (trade-off)。假使各個測站與震源之間的震央距約略相等時，一個較深的震源位 置所造成的走時殘差可藉由減少發震初始時間來彌補，進而造成定位誤差來源。

最後，走時殘差的變化與震源位置之間並不是線性關係，因此在觀測資料不 足夠的情況下，式 1-1 所得出的解可能會落在區域最小值 (local minimum) 而不 是全域最小值 (global minimum) 而產生誤差。

為了改善現存地震預警系統之定位成果，本研究以中央氣象局地震速報系統 觀測網之 RTD (Real-Time Data stream)測站為架構，利用等時差 (Equal

(1-1)

differential time, EDT) 曲面定位法，並結合臺灣地區三維速度構造模型和格點搜 尋法進行快速的震源位置估算。

圖 1.1 台灣歷史災害地震震央分布

(資料來源 http://scweb.cwb.gov.tw/GraphicContent.aspx?ItemId=21)

1 2

8

9 16

19 7

4

18

5

3 6

14

11 17

12 10

15 13

1 斗六地震 2 梅山地震 3-6 南投地震系列

7 新竹-台中地震 8 中埔地震 9 新化地震 10-14 縱谷地震系列

15 恆春地震 16 白河地震 17 花蓮地震 18 921 地震 19 美濃地震

第二章 文獻回顧

2.1 地震預警發展歷史與臺灣地震預警系統發展

地震預警概念的雛形起源於 1868 年，由美國的科學家 J.D. Cooper 於 1868 年提出，其概念為利用電磁波傳遞速度比地震波波速還快的特性，當地震儀監測 到地震後即可以將地震警報利用電磁波的形式傳送至較遠方的城市警告當地居 民以達到預警效果 (圖 2.1) ，不過此概念因當時環境的諸多因素限制下而無法 實現。

日本的國道鐵路公司(Japanese Railway, JR)在 1964 年的 Niigata 地震發生之 後，著手在鐵路沿線上設置地震警報系統，並以最大地動加速度(Peak Ground Acceleration, PGA)為判斷標準，當地震測站偵測到 PGA 大於 40 gal (1 cm/s² =1 gal) 便 會 對 其 沿 線 上 的 列 車 發 出 警 告 。 隨 後 該 公 司 發 展 出 新 的 地 震 預 警 系 統 UrEDAS (Urgent Earthquake Detection and Alarm System) (Nakamura, 1988)，為最 早的商業化地震預警系統，此系統利用單一測站所接收到的 P 波資訊來進一步推 估地震發生的位置和大小。

墨西哥自 1985 年的 Michoacan 地震後，在 1991 年開始建立地震預警系統 SAS (Seismic Alarm System)，此系統在 1995 年成功的偵測到 Copala 地震，並於 地震發生後，在具有破壞性的剪力波到達墨西哥市前的 72 秒前成功發佈警報 (Espinosa-Aranda et al., 1995)。

臺灣發展地震預警的動機來自於 1986 年 11 月 15 日發生於花蓮外海芮氏規 模 (ML) 6.8 的地震，雖然該地震震央位於花蓮外海，卻在 120 公里遠的台北盆 地造成災情，根據地震波走時資料，S 波由花蓮傳遞到台北地區約需 30 秒的時 間，因此如果能在 30 秒內將地震資料處理完畢便能在破壞性剪力波到達前爭取 到預警時間。

中央氣象局於 1994 在花蓮地區建立台灣第一個地震預警系統雛形 (Wu et

al., 1999) (圖 2.2)，其目的為當花蓮地區發生地震時，可以透過近震源區的地震

觀測站將地震資訊迅速的傳回中央氣象局。到了 1997 年，鑒於先前預警系統的 經驗，中央氣象局挑選數十個台灣強地動觀測系統 (Taiwan Strong Motion Instrumental Program, TSMIP) 之測站，利用數據專線即時傳送地震儀器訊號至中 央 氣 象 局 ， 組 成 臺 灣 即 時 強 地 動 觀 測 網 (real-time strong-motion monitoring network)，並且利用測站接收到 P 波 10 秒窗口內的地動訊號預估規模 (M_L10) (Wu

et al., 1998)，並透過軟體開發引入區域子網的概念有效縮短資料處理時間，建立

起虛擬子網 (Virtual Sub Network, VSN) 地震預警系統 (Wu and Teng, 2002) (圖 2.3)。

圖 2.1 美國科學家 Coopper 提出之地震預警概念圖 (Cooper, 1868)。

圖 2.2 台灣第一個裝設在花蓮地區的地震預警系統雛形 (Wu et al., 1999)

圖 2.3 虛擬子網地震預警系統 (Wu and Teng, 2002)

(2-2)

2.2 地震預警之地震定位方法

對於地震減災來說，良好的震央位置可以有效提供救災單位評估災害大小。

此外，部分地震規模的預估方法亦受震源位置影響。例如，估算各地區的地表加 速度值 (PGA) 時，其衰減公式內的震源距參數便與震源位置相關 (Lin and Wu, 2010)。以下將介紹幾個用於地震預警研究中的定位方法:

(1) 線性化逆推定位法

此方法為一經典的地震定位方法，由於地震波 (P 波及 S 波) 傳達到各個接收 測站的時間會與測站和震源之間的距離遠近而有所變化，藉此特性，便可利用各 測站的地震波到時紀錄來反推震源位置，原理如下。

地震波傳達到測站的到時可用數學式來表示，

其中 𝑇_𝑜𝑏𝑠^𝑖 為第 i 測站接收到的地震波到時，𝑡_𝑜為地震發生之時間，r 為地震 波的傳播路徑，u 為地震波所經過的路徑上的慢度 (slowness) (即速度的倒數)。

因式 2-1 中的積分項為一非線性方程式無法直接求解，因此 Geiger 於 1912 年提 出藉由泰勒展開式將其近似化成線性關係 (式 2-2)，並利用假設初始震源和發震 時間進行疊代運算得出最佳解，

其中 𝑡^𝑖為地震波之理論走時，(∆x, ∆y, ∆z, ∆τ) 為疊代過程中震源位置 (x, y, z) 及發震時間 τ 的調整量，走時殘差 𝑟^𝑖為地震波的觀測走時與理論走時之差值。

由於有四項未知數 (震源位置和發震時間)，故至少要有四組觀測到時資料 才能進行定位，每經過一次疊代運算便會得出一組震源位置及發震時間的調整量，

(2-1) T_obsⁱ = t_o+ ∫ u(𝐫)ds

𝑟_{𝑠𝑡𝑎𝑡𝑖𝑜𝑛}

r_{𝑠𝑜𝑢𝑟𝑐𝑒}

𝜕𝑡^𝑖

𝜕𝑥 ∆x +𝜕𝑡^𝑖

𝜕𝑦∆x +𝜕𝑡^𝑖

𝜕𝑧 ∆x + ∆τ = 𝑟^𝑖

(2-3)

(2-5) (2-4)

𝜃_𝑖 = tan⁻¹(R_ud,ew)_i (R_ud,ns)_i (2) 單站定位法

Nakamura (1988) 提出透過單一測站的地震波形紀錄中水平方向與垂直方向 分量振幅之間的交互相關(cross-correlation) 函數來推估地震波入射方位角，

(R_ud,ew)_i = α(R_ud,ew)_i−1+ (X_ud)_i(X_ew)_i (R_ud,ns)_i = α(R_ud,ns)_i−1+ (X_ud)_i(X_ns)_i

其中 (𝑅_{𝑢𝑑,𝑒𝑤})_𝑖 和 (R_ud,ns)𝑖 分別為地震波振幅垂直方向分量與水平東西和南北

方向分量的交互相關函數，α 為平滑係數 (介於 0 到 1 之間)，(X_ud)_i、(X_ew)_i 和 (X_ns)_i為第 i 點在垂直、東西和南北方向分量的振幅值。求得交互相關函數後 (式 2-3 及式 2-4) 便可計算地震波入射方位角，

此外，Nakamura (1988) 利用歷史紀錄中地震波波形振幅初動的垂直分量與水平 分量的比值分別與震央距、震源深度和震源距作出經驗衰減關係式，藉由此關係 式便可進一步推估震央距離和震源深度。

(3) 測站幾何分布定位法

除了透過 P 波到時和波形資訊來進行地震定位之外，也可以利用觸發測站的幾何 分布型態及配合接收到 P 波的時間先後順序進行定位。Sambridge (1999) 提出 Voronoi cells 方法，其概念為假設 P 波波速固定，依照測站幾何分布將所有測站 兩兩之間畫出中垂線，最後會在各個測站周圍畫出多邊型區域 (Voronoi cell)，如 果有地震發生在某一測站所屬的多邊形範圍內則該測站會最先收到地震波訊號，

因此可以推估地震位置。測站的分布密度越高則 Voronoi cells 的面積越小，對於 地震位置的約束性就越好，但是此方法僅能求出地震震央的粗略位置，無法得知 震源深度。

(4) 等時差 (EDT) 曲面定位法

此方法最早可追朔至 Milne (1886) 使用作圖法進行二維空間的地震定位，其 中包括利用各測站的 P 波和 S 波之間的時間差 (S-P differential times)，配合假設 的地震波波速，即可在地圖上畫出以測站為圓心的同心圓，圓形的交點便是震央 位置。另外，將平面上對於兩測站有相同到時差的點連結起來會形成一雙曲線，

多個測站之間所構成的雙曲線的交點即為震央位置。

延續平面空間雙曲線定位法的概念，假設一個三維均質環境，一個點震源產 生的地震波傳達到兩個不同測站的路徑可以想像等同於從兩測站分別產生震波 (如圖 2.4 中的橘色同心圓和綠色同心圓) 到達震源位置，而兩個震波波前的交會 點則代表兩測站接收到地震波之到時差，因此將具有相同到時差的波前相連起來 便 會 在 空 間 中 畫 出 一 個 以 兩 測 站 為 焦 點 的 雙 曲 面 S_ij， 稱 為 等 時 差 (Equal differential time, EDT) 曲面，所有落在此曲面上的點即為可能的震源位置，並滿 足下列關係式，

其中t_i、 t_j 為測站 i 和 j 之 P 波觀測走時，T₀ 為假設之地震發震時間，T_i、 T_j 為 第 i 和第 j 站之 P 波觀測到時，由於正確的震源位置應落在 EDT 曲面 S_ij 上，因 此在至少有三組 EDT 曲面以上的情況下其交集即為震源位置 (圖 2.5)。在現實 情況下，EDT 曲面的幾何形狀取決於速度模型和使用的地震波相，當使用相同 的波相 (同為 P 波或是 S 波) 時，EDT 曲面為雙曲面。如果T_i、T_j 分別為同一測 站的 P 波和 S 波到時，則 EDT 曲面為一半球面。由於地下構造為非均質環境，

因此 EDT 曲面並不會是簡單平整的幾何形狀，而會因速度構造的不同而變形。

實際應用到地震定位上時，假設震源位置 q ， τ_i(𝐪) 、τ_j(𝐪) 為 q 點至第 i 和第 j 測站的理論走時，則相對到時差dt_ij (即觀測到時差與理論走時差之間的差 S_ij = t_i− t_j = (T_i− T₀) − (T_j− T₀) = T_i− T_j (2-6)

當 q 與真實震源位置非常接近時時，則 dt_ij 趨近於零。由 (式 2-6) 可知 EDT 曲 面定位法不需要發震時間這一項參數，因此可避免震源位置與發震時間之間的互 補效應。

Zhou (1994) 提出的主站法 (Master Station method, MSM) 為首先將 EDT 曲 面結合三維速度模型進行地震定位之研究。此方法分為兩個部分，首先將研究區 域劃分出數個格點，接著使用三維速度模型預先計算各個格點至地震測站的理論 地震波走時並存成參考檔。當地震發生後，以第一個觸發的測站當作參考值，若 有 N 個測站，則共可建立 N-1 個 EDT 曲面，接著對每個格點搜尋參考檔內的理 論地震波走時資料判斷是否被 EDT 曲面通過，有最多 EDT 曲面通過次數的格點 即為震源位置可能所在的區域。

Satriano et al., (2008) 進一步將 EDT 曲面結合 Voronoi cells 的地震定位概念應 用於地震預警上，把未觸發的測站加入考量，利用尚未接收到地震波之測站的幾 何位置分布，幫助最先觸發的測站縮小震央位置的範圍，加快定位速度。

本研究即採用 EDT 曲面定位法並結合 GPU 平行運算於格點搜尋的方法，詳 細流程與步驟於下一章節介紹。

圖 2.4 等時差 (EDT) 曲面示意圖

圖 2.5 以平面圖呈現等時差(EDT) 曲面定位概念 (a) 地震發生後，在最先觸 發的兩測站之間可以建立一曲面 (b) 隨著測站增加便可建立多組曲面，由於 空間中需要至少三組曲面其共同交集才可能是一點，因此最少需要三個測站才

(a) (b)

2.3 現行中央氣象局地震預警系統之地震定位流程

現行的中央氣象局地震預警系統的地震定位流程分為兩步驟 (Chen, 2015)，

先求出震央位置後再用格點搜尋法找出最佳的震源深度。首先，震央位置的估算 採用線性化逆推方法 (Geiger, 1912)，所使用的速度構造模型是由三維速度構造 模型 (Wu et al., 2009) 經平均化後得出的一維半空間 (half-space) 連續速度模型，

此速度構造模型以深度 40 公里當作邊界分為上下兩層，深度小於 40 公里所對應 的理論地震波波速 V(D) 為，

其中 D 為地殼深度。而深度大於 40 公里的理論地震波速則為，

當震央位置決定好後，再利用格點搜尋法的方式從深度 10 公里開始至 100 公里 結束，以 10 公里為間隔計算不同深度的走時殘差值 (travel time residual)，最後，

有最小走時殘差值的深度和其對應的震央位置即為該事件的震源位置。

此地震定位過程會在地震預警系統內持續運作，而定位程序的觸發條件為至 少要有六個測站偵測到 P 波到時才會開始進行運算。

V(D) = 5.103 + 0.067 × D

V(D) = 7.805 + 0.005 × D

(2-8)

(2-9)

第三章 研究方法

3.1 本研究之地震定位流程

本研究所使用的地震定位法為 EDT 曲面定位法，基於地震預警需要快速決 定地震位置的概念，本研究僅使用 P 波資料進行計算，整體的定位流程與 Zhou (1994) 提出的主站法類似，但是對其中的細項進行修改。

首先，將臺灣地區劃分為數個三維網格點，每個網格點皆為假設之震源位置 q ， 利 用 三 維 速 度 模 型 (Hunag et al., 2014) 和 虛 擬 曲 折 波 線 追 跡 法

(Pseudo-bending ray tracing method) ( Um and Thurber, 1987；Koketsu and Sekine, 1997) 預先計算出各格點至 RTD 測站之理論地震波走時並儲存成資料庫。

接著當地震發生後，將各測站的觀測到時兩兩配對為一組形成 EDT 曲面，

由式 2-7 可知 dt_ij 代表著震源位置對於 EDT 曲面的擬和程度，因此我們可以利 用最小方均根 (Root Mean Square, RMS) 值來找出最符合所有 EDT 曲面交集的 震源位置 qhypo，

其中 M 為 EDT 曲面的數量，J 為使用的測站數目。因此對於空間中所有的格點 q 都會有自己的 RMS 值，當 RMS(𝒒) 為所有格點之最小值 (即 global minimum)

時對應的 q = q_hypo即為最佳解。與主站法 (Zhou, 1994) 不同的地方為不把第一 個觸發測站的觀測到時作為唯一的參考值，而是使用所有測站觀測到時的配對組 合，減少對於單一觀測到時誤差的相依性。相對於 Satriano et al., (2008) 利用 Voronoi cells 加快震央定位速度，而本研究則利用 GPU 平行化計算概念加快格點 搜尋 (式 3-1) 過程。

RMS(𝐪) = √1

M∑(dt_ij)²

i,j

；M = C₂^J (3-1)

3.2 測站走時資料庫建立

本研究將臺灣地區範圍：120^o E-122.5 ^o E；21.5^o N-25.8 ^o N；深度範圍：1-64 km，切成 1x1x1 km 之三維網格點，總計有 6880000 個網格點 (圖 3.3)。接著需 計算各個網格點至 RTD 測站 (共 106 站) 之理論地震波走時並依據測站名稱將對 應的走時資料分別存檔。

波線追跡法 (ray tracing) 為計算單一震源之震波傳達到一測站的路徑和所 需 要 之 走 時 的 方 法 。 本 研 究 所 使 用 之 波 線 追 跡 法 為 虛 擬 曲 折 波 線 追 跡 法 (Pseudo-bending ray tracing method) ( Um and Thurber, 1987)，此波線追跡法概念 為假設兩點(震源和測站) 之間的波線路徑為一多節點之折線。首先，從一直線 開始，接著利用速度構造求取該線通過區域 x, y, z 三方向的速度梯度，藉由此梯 度方向調整節點位置以符合走時最短之路徑，接著便不斷增加路徑上的節點與調 整節點位置，當走時的變化率小於門檻值後便停止計算，此時會得到一多節點折 線，此折線會近似於理論波線之路徑，而沿著此折線路徑上的走時即為地震波理 論走時。Koketsu and Sekine (1997) 進一步將此波線追跡法改良成三維球面座標 系統，並考慮速度構造中的不連續面，以便更符合真實環境，本研究之理論地震 波走時資料即使用此改良版本進行計算。

圖 3.1 臺灣地區的三維網格點範圍與 RTD 測站 (106 站) 分布圖。此圖僅呈現 二維空間格點分布範圍，每個格點間距為 1 公里，共有 107500 個假設震源位 置。紅色三角行為 RTD 測站位置，綠色格點為假設震源位置。

3.3 OpenACC 程式平行化

當一地震發生，根據測站觸發的先後順序便可以在兩兩測站之間建立 EDT 曲面，接著本研究開發之定位程式便會根據觸發之測站名稱，從其對應的走時資 料庫一一計算網格點與各個 EDT 曲面的殘差值，有 N 個測站便會有 N(N-1) 組 EDT 曲面，因此將會產生龐大的計算量。託現代電腦科技發展之福，平行化運 算已成為分析大量數據不可或缺的工具。

早期 GPU (Graphic Processing Unit) 為一種圖形處理專用的加速硬體，無法 處理在 CPU 上執行的運算任務，而隨著程式架構的進步，程式開發者能藉由控 制繪圖晶片來進行資料運算。

OpenACC (Open Accelerators) (https://www.openacc.org/) 是由 Cray、CAPS、

Nvidia 和 PGI 共同開發的加速器程式標準，於 2012 年推出，其開發目的是為了 降低 GPU 程式平行化編程的難度，讓各領域的研究人員能輕鬆地運用 GPU 強大 的平行化能力改善現有程式碼 (C/C++/Fortran) 的執行速度，其中有兩項特點：

(1) 只要將編譯器提示(compiler directives) 嵌入程式碼中需要平行運算的區域，

編譯器便會將該區域的程式碼交給 GPU 進行運算。

(2) 搭載 OpenACC 的 GPU 加速程式碼可以在各種不同的 GPU 平台上執行。

本研究利用 OpenACC 平行化套件將網格計算過程 (式 3-1) 移植到電腦顯示 卡上 (圖 3.2)，藉由圖形運算處理器 (GPU) 在平行計算上的優勢來增加格點計 算過程的速度。

圖 3.2 程式平行化計算流程

第四章 結果與討論

4.1 TAIGER 人工炸測震源資料測試

為了研究台灣的地下構造，臺灣和國外學者共同合作於 2004 - 2009 年間執 行 TAIGER 計畫 (Taiwan Integrated Geodynamic Research) (Wu et al., 2013)，利用 大規模人工爆炸波施測所獲得的震波走時資料來探討臺灣三維地體構造，其中於 2008 年執行的的陸上炸測主動式震源實驗便是利用炸藥產生人工震源進行震測 剖面和速度構造反演之研究。由於我們無法得知自然界地震發生的真正位置，因 此可藉由已知位置的人工震源資料用於評估地震定位程式的結果表現。

本研究利用 TAIGER 計畫中 10 個人工炸測震源資料 (表 1) (圖 4.1) 來比較 兩種定位程式對於人工震源位置的回復程度。其中一種為本研究開發之網格式定 位程式，另一種為 3DLOC (Wu et al., 2003)，3DLOC 為使用三維速度模型之線性 化逆推定位程式。在本次測試中，這兩種定位程式都使用相同的三維速度模型 (Huang et al., 2014) 和虛擬曲折波線追跡法 (Um and Thurber, 1987；Koketsu and Sekine, 1997)，因此主要的差異在於得出震源位置的方式不同。

由於 3DLOC 定位程式是使用線性化逆推方法，經由反覆疊代運算得出地震 位置最佳解，因此需要輸入一個初始值 (震源位置和發震時間) 當作起點，如果 初始位置不夠接近真實震源位置，則最後經程式逆推反演得出的結果可能無法回 復到正確的位置。相對於線性化逆推方法，本研究開發之定位程式採用網格式計 算方法，直接計算空間中各格點的理論值與觀測值的擬合程度，因此不需要給定 一初始值，減少造成誤差的原因。

以下為對 3DLOC 定位程序使用兩種初始值條件與本研究之定位程序進行定 位結果比較，使用的測站為 RTD 之測站 (圖 3.1)。第一種初始值條件為使用原始 的炸測位置和施測時間；第二種則參考中央氣象局地震預警 (CWBEEW) 系統 定位程序的初始值條件，即初始震央位置為接收到 P 波之測站幾何中心，初始深 度為 10 公里，初始發震時間為第一個接收到 P 波訊號之測站的 P 波標定時間。

4.1.1 地震定位結果

重定位結果如圖 4.2 和圖 4.3，在震央誤差方面，3DLOC 在兩種測試條件下 和本研究之定位結果有相當類似的誤差分布，但是在 N2 炸測點位則有明顯的差 異，3DLOC 在使用 CWBEEW 的初始值條件下所得出定位結果有較大的誤差 (>

20 km)。和震央誤差相比，可以看出此兩種定位方法之結果整體上有誤差增加的 趨勢，而且對於各個炸測點位的誤差值有明顯的不同，另外 3DLOC 在 CWBEEW 的初始值條件下對於 N2 點位仍有較大的誤差產生 (約 10 km)。

4.1.2 定位結果分析

在重定位的測試結果中可以看到這兩種定位方法對於震央位置的回復能力 比震源深度還來得好，此情形符合普遍在進行地震定位工作時對於震源深度較難 掌控的狀況。從定位結果來看，即使是使用原始的震源位置和發震時間當作初始 值的線性逆推定位程式 (3DLOC) 或是少掉發震時間這項變數的本研究之定位 程式，皆無法有效改善深度誤差較大的問題，此現象可能是因為速度模型在淺層 部分的解析能力較差的關係所導致。

另外我們將 3DLOC 在使用 CWBEEW 的初始值條件下和本研究之定位程式 進行測站最大空餘角 (GAP) 對於定位結果的影響比較 (圖 4.4)。整體來看，

3DLOC 和本研究之定位程式的定位誤差和 GAP 值大致有正相關性。

對於 3DLOC 在改變初始值後對於 N2 炸測點有較大定位誤差的現象，可能 是因為 N2 炸測事件所使用的站測站數目較少或是測站分布接近線性的關係所造 成。圖 4.5 為 3DLOC 和本研究之定位程式對於 N2、S1 炸測點位的定位結果以 及分別使用的 RTD 測站分布位置，可以看出在有相似的 GAP 值情況下，N2 炸 測點位的測站分布較為線性。

為了檢驗是否為測站數目較少而造成定位誤差，我們將 S1 炸測點位使用的

S1 S2

S3P S3

N1 N2

S4

N3 & N4 N3P

變，因此對於線性逆推定位方法來說，初始值的改變很可能會造成定位結果的不 穩定，而本研究之定位程式則可以避免此狀況。

表 4-1 TAIGER 人工炸測震源位置、施測時間、紀錄規模和測站 (RTD) 數目

Site Code Longitude Latitude Elev. (m) Shooting Time (UT) M_L (CWB) Stations GAP S1 120.2299 23.5106 6 2017/02/27 17:01:49.486 3.02 11 176 S2 120.4131 23.42513 20 2017/02/27 17:30:57.350 2.62 13 118 S3 120.7066 23.28189 650 2017/02/26 18:02:45.860 1.89 10 84 S3P 120.6187 23.13875 335 2017/02/26 17:32:43.296 1.78 10 95 S4 121.1336 23.13466 390 20170/2/28 17:02:08.560 2.30 9 88 N1 121.0377 24.79534 55 2017/03/06 17:03:23.068 2.56 7 114 N2 121.2276 24.67742 590 2017/03/06 17:34:56.364 1.55 5 187 N3 121.4797 24.57874 405 20170/3/04 18:01:17.988 2.10 10 82 N3P 121.4797 24.57874 405 2017/03/04 17:01:17.824 2.67 11 82 N4 121.8035 24.44872 6 2017/03/05 17:03:40.008 2.43 8 188

圖 4.1 人工炸測震源位置

圖 4.2 各個炸測點的震央距、深度和震源距誤差的重定位結果。藍色星號為 3DLOC 使用原始炸測資料的結果，紅色叉號為本研究之定位程式結果

3DLOC

This Study

圖 4.3 各個炸測點的震央距、深度和震源距的定位誤差結果。藍色星號為 3DLOC 使用 CWBEEW 定位初始值得出的結果，紅色叉號為本研究之定位程式結果

3DLOC

This Study

3DLOC This Study

圖 4.4 測站最大空餘角(GAP) 值與震央距、深度和震源距 誤差的比較。藍色星號為 3DLOC 使用 CWBEEW 定位初 始值得出的結果，紅色叉號為本研究之定位程式結果

圖 4.5 N2 和 S1 炸測點位的重定位結果以及使用之 RTD 測站分布。實心 三角形為測站位置，黑色星號為原始炸測點位，藍色星號為 3DLOC 使 用 CWBEEW 定位初始值所得出之震央定位結果，紅色叉為本研究之定 位程式所得出之結果

N2

S1

N2

S1

圖 4.6 N2 和 S1 炸測點位的重定位結果以及使用之 RTD 測站分布，將 S1 炸測點位的測站數目減少至與 N2 炸測點位相同。實心三角形為測站位 置，黑色星號為原始炸測點位，藍色星號為 3DLOC 使用 CWBEEW 定 位初始值所得出之震央定位結果，紅色叉為本研究之定位程式所得出之 結果

4.2 歷史地震資料離線測試

本研究挑選 2014 至 2017 年，震源深度小於 60 公里，芮氏規模 (ML) 大於 4.5，島內和島外共 105 個地震事件。進行本研究開發之定位程式測試，並將定 位結果與中央氣象局地震預警 (CWBEEW) 系統之結果分別和中央氣象局地震 目錄做比較。氣象局的地震目錄是經由人工重新標定 P 波到時並重新定位過後的 結果，所使用的地震觀測網包括 RTD、短週期和寬頻觀測網。在氣象局地震預 警系統的定位過程中，當有至少六個測站觸發便會進行地震定位並產生第一報結 果，接著依據觸發的測站數目增加而產生後續的報數，其中為了避免誤報 (false alarm) ，第一報和第二報會由系統內存，第三報才會對外發布，因此在本測試 中將針對第三報的資訊進行比較，而本研究中的定位程式也只使用第三報的測站 紀錄進行定位。

4.2.1 地震定位結果

本研究開發之定位程式和 CWBEEW 之定位結果與氣象局地震目錄的比較 如圖 4.7，震央位置和深度平均差值分別為 9.26 km、7.18 km (CWBEEW) (圖 4.7a)，

4.77 km、5.11km (本研究之結果) (圖 4.7b)，從震央和深度的比較結果可以看出本 研究之定位方法能有效改善地震定位的精準度。

4.2.2 定位結果分析

圖 4.7 為島內和島外定位結果比較，對於台灣島東部外海的地震事件，受限 於地震測站的分布影響，造成在定位結果上有較大的不確定性，因此我們進一步 探討島內 (36 個地震事件) 和島外 (79 個地震事件) 的定位結果。首先，圖 4.8a 和 4.8b 的上圖為相對震央差距 (CWBEEW 之震央誤差 - 本研究之震央誤差)，

正值代表本研究得出的震央位置比 CWBEEW 更接近氣象局地震目錄的結果，對 於約八成的島內及島外的地震事件來說，本研究之定位程式能有效改善其定位結

出來的震源深度大於重定位目錄之結果，在此深度差距比較圖中，平均值越 接近零和標準差越小表示該方法有較好的定位結果，結果顯示出本研究之定位程 式對於島內和島外的地震定位結果皆有良好的改善。

在先前 TAIGER 炸測資料測試結果中，顯示出測站空餘角 (GAP) 會影響定 位結果。在此我們也測試 GAP 值對於本研究開發之程式和 CWBEEW 對於歷史 地震資料的定位結果之影響程度，圖 4.9 為 GAP 值對震源距誤差的比較，從結 果來看，本研究之定位程式和 CWBEEW 所得出的震源位置與中央氣象局地震目 錄的誤差都有隨著 GAP 值增加而變大的趨勢。但整體上而言，本研究之定位程 式的定位誤差小於 CWBEEW 的定位誤差，震源誤差分別為 7.79 km 和 13.13 km。

另外可以發現本研究之定位程式對於 GAP 值約在 180^o到 230^o 之間的島外地震 誤差有特別明顯的改善，圖 4.10 為坐落在此 GAP 值範圍的地震震央分布，對於 發生在花蓮外海的地震事件，如先前文獻回顧所提到，此區域的地震事件對於台 北盆地有較大的影響，因此定位結果的改善能增加地震預警系統評估結果的可靠 性。

圖 4.7 歷史地震資料重定位結果之比較 (a) 氣象局預警系統之結果 (藍色圓

(a)

(b)

Δddepth : 5.11 ± 6.02 km Δddepth : 7.18 ± 7.1 km

Δdepicenter : 9.26 ± 12.9 km

Δdepicenter : 4.77 ± 6.15km

CWBEEW

This Study

CWBEEW This Study

22% 78%

Δd_depth= -1.65± 8.43 km

圖 4.8 震央位置和深度之誤差比較。(a) 島內地震事件 (b) 島外地震事件

23% 77%

Δd_depth= -2.20± 8.18 km

Δ

Δ

Δd_depth= -3.82± 10.09 km

(a)

(b)

Δ

Δ

CWBEEW

This Study

(a)

(b)

圖 4.9 震源位置誤差與測站空餘角 (GAP) 之比較 (a) 本研究之定位結果 (b) CWBEEW 之定位結果

CWBEEW

This Study

圖 4.10 測站空餘角(GAP) 介於 180^o至 230^o之間的島外地震事件 (a) CWBEEW 之定位結果(b) 本研究之定位結果

(a) (b)

4.3 案例討論：2016 年 2 月 6 日高雄區美濃市芮氏規模 (M_L) 6.6 地震

在本節中，我們針對 2016 年 2 月 6 日發生的美濃地震事件進行定位結果的 比較。

美濃地震發生於 2016 年 2 月 6 日凌晨 3 時 57 分，根據中央氣象局的定位結 果，其主震位置位於北緯 22.92 度，東經 120.54 度，地震深度為 14.6 公里，屬 於淺層地震。主震芮氏規模為 6.6，最大震度 (7 級) 地區位於台南市新化區 (圖 4.11)。圖 4.12 為本研究之定位程式根據 CWBEEW 定位系統產生的第 1 至第 12 報定位結果中的測站資訊進行定位所得出的定位結果，並與 CWBEEW 的定位結 果比較，從定位結果可以看出本研究之定位程式的定位誤差皆小於 CWBEEW 的 定位誤差，其中 CWBEEW 的定位誤差變化從第 1 至第 12 報都相當穩定，而本 研究之定位程式地定位誤差則有擺盪的現象。

圖 4.11 2016 年 2 月 6 日美濃地震報告

(引用至 http://scweb.cwb.gov.tw/GraphicContent.aspx?ItemId=20&loc=tw)

(a)

(b)

CWBEEW This Study

圖 4.12 美濃地震事件的定位結果 (第 1 報至第 12 報) (a) 震央誤差 (b) 深 度誤差

4.4 地震定位時效性討論

圖 4.11 為本研究之定位程式在進行網格式搜尋的程式碼部分於 GPU 平行化 前後在執行上所花費的時間比較，所使用之 GPU 規格列於附錄 A。

由圖 4.11 中可以看出相對於僅使用 CPU 單核心進行運算所花費的時間，本 研究之定位程式經由 GPU 平行化後只需約十分之一的時間便能得出地震定位結 果，顯示出藉由程式的平行化技術可以有效的改善行大量數據資料運算的執行效 率。

在目前現行的 CWBEEW 系統中，對於每一報的地震定位結果皆需要至少 10 秒鐘以上的時間才能完成計算，因此利用本研究之定位程式能有效降低在計 算地震震源位置的過程中所需要花費的時間。

CPU (single core) -- i7 3770 GPU -- NVIDIA GTX 1060

圖 4.13 本研究之定位程式於 GPU 平行化前後在不同的測站數目下得出震源 位置所需要的計算時間。

第五章 結論

1. 本研究利用 EDT 曲面定位法，透過格點搜尋方式，只需要簡單的四則運算 即可定位出震源位置，此外還結合了台灣地區三維速度構造模型及虛擬曲折 波線追跡法 (Pseudo-bending ray tracing method) ，能計算出接近真實環境的 地震波傳播路徑，進而提升走時精確度，讓定位結果更加精準。

2. 由於本研究將需要複雜計算過程的理論地震波走時資料預先計算並建立成 資料庫，因此當地震發生後，藉由 GPU 高速平型化計算的幫助下大幅減少 格點搜尋過程中需要花費的時間，提升地震定位速度。

3. 在本研究之定位程式所得出的結果中，對於震央誤差的改善程度還不錯，原 因可能在於利用了各個測站對 (station pair) 的關係，將測站的空間分布幾 何考慮進去，進而對於震央位置有較好的控制能力，但是對於深度部分則沒 有太明顯的改善情況。

4. 本研究之定位程式的定位精度受限於在建立理論地震波走時資料庫的過程 中所設置的格點大小，如果要劃分成更精細的格點會需要龐大的記憶體容量 來儲存資料庫中的理論地震波走時資料，進而增加讀取過程中的時間，較好 的解決方式為結合本研究之定位程式與 3DLOC 定位程式，先由本研究之定 位程式提供較為可靠的地震震源位置，再由 3DLOC 進行更精細的地震定 位。

參考文獻

1. Allen, R. M. and Kanamori, H. (2003), the potential for earthquake early warning in southern California. Science, 300(5620), 786-789.

2. Cooper, J. D. (1868), Earthquake Indicator. San Francisco Bull., November 3, 1868.

3. Chen, D. Y., Wu, Y. M., and Chin, T. L. (2015), Incorporating low-cost seismometers into the Central Weather Bureau seismic network for Earthquake Early Warning in Taiwan. Terr. Atmos. Ocean. Sci., 26, 503-513, doi:

10.3319/TAO.2015.04.17.01(T).

4. Chen, D. Y. (2015), Development and Study of Earthworm Platform for Earthquake Early Warning in Taiwan, Doctoral Dissertation, Development of Geosciences, College of Science, National Taiwan University, Taipei, Taiwan.

124 p.

5. Espinosa-Aranda, J. E., Jimenez, A., Ibarrola, G., Alcantar, F., Aguilar, A., Inostroza, M., and Maldonado, S. (1995), Mexico City seismic alert system. Seis.

Res. Lett., 66(6), 42-53.

6. Geiger, L. (1992), Probability method for the determination of earthquake epicenters from the arrival time only, Bull. St. Louis Univ., 8, 60-71.

7. Huang, H. H., Wu, Y. M., Song, X., Chang, C. H., Lee, S. J., Chang, T. M., and Hsieh, H. H. (2014), Joint Vp and Vs tomography of Taiwan: Implications for subduction-collision orogeny, Earth Planet. Sci. Lett., 392, 177-191.

8. Koketsu, K., and Sekine, S. (1998), Pseudo-bending method for three-dimensional seismic ray tracing in a spherical earth with discontinuities, Geophys. J. Int., 132, 339-346.

9. Kanamori, H. (2005), Real-time seismology and earthquake damage mitigation.

Annu. Rev. Earth Planet. Sci., 33, pp. 195-214.

10. Lee W.H.K., and Lahr J.C. (1975), HYPO71 (revised): a computer program for determining hypocenter, magnitude, and first pattern of local earthquakes, US

Geological Survey Open-File Report 75–311, pp. 1–116.

11. Lin T. L. and Wu, Y. M. (2010), Magnitude determination using strong ground-motion attenuation in earthquake early warning, Geophys. Res. Lett., 37, L07304, doi:10.1029/2010GL042502.

12. Milne, J. (1886), Earthquake and other earth movements, Appelton, New York.

13. Nakamura, Y. (1988), On the urgent earthquake detection and alarm system (UREDAS). Proc. Of the 9^th

World Conference on Earthquake Engineering,

14. Sambridge, M. (1999), Geophysical inversion with a neighbourhood algorithm-I.

Searching a parameter space, Geophys. J. Int., 122, 108-124.

15. Satriano, S., Lomax, A., and Zollo, A. (2008), Real-Time Evolutionary Earthquake Location for Seismic Early Warning, Bull. Seism. Soc. Am.,98(3), 1482-1494.

16. Um, J., and Thunber, C. (1987), A fast algorithm for two-point seismic ray tracing,

Bull. Seismol. Soc. Amer., 77, 3, 972-986.

17. Wu, Y. M., Shin, T. C., and Tsai, Y. B. (1998), Quick and reliable determination of magnitude for seismic early warning, Bull. Seism. Soc. Am., 88(5), 1254-1259.

18. Wu, Y. M., Chung, J. K., Shin, T. C., Hsiao, N. C., Tsai, Y. B., Lee, W. H., and Teng, T. L. (1999), Development of an integrated earthquake early warning system in Taiwan-Case for the Hualien area earthquakes. Terr. Atmos. Ocean. Sci., 10(4), 719-736.

19. Wu, Y. M., and Teng, T. L. (2002), A virtual subnetwork approach to earthquake early warnig. Bull. Seism. Soc. Am., 92(5), 2008-2018.

20. Wu, Y. M., Chang, C. H., Hsiao, N. C., and Wu, F. T. (2003), Relationship of the 1998 Rueyli, Taiwan, earthquake sequence using three-dimensions velocity structure with stations corrections, Terr. Atmos. Ocean. Sci., 14, 4, 421-430.

21. Wu, Y. M., and Kanamori, H. (2005), Rapid Assessment of Damaging Potential of Earthquakes in Taiwan from the Beginning of P Waves, Bull. Seism. Soc. Am., 95, 1181-1185.

22. Wu, Y. M., Shyu, J. B. H., Chang, C. H., Zhao, L., Nakamura, M. and S. K. Hsu (2009), Improved seismic tomography offshore northeastern Taiwan: Implications for subduction and collision processes between Taiwan and the southernmost Ryukyu, Geophys. J. Int. 178, 1042-1054.

23. Wu, Y. M., Chang, C. H., Kuo-Chen, H., Huang, H. H. and Wang, C. Y. (2013), On the use of explosion records for examining earthquake location uncertainty in Taiwan, accepted by Terr.

Atmos. Ocean. Sci.

10.3319/TAO.2013.01.31.01(T).

24. Zhou, H. W. (1994), Rapid three-dimensional hypocentral determination using a master station method, J. geophys. Res., 99(15), 439–15, 455.

附錄 A 本研究所使用之 GPU 規格

產品名稱 GeForce GTX 1060 產品型號 06G-P4-6161-KR CUDA 核心數目 1280

預設 GPU 運作時脈 1506 MHz

記憶體規格 6144 MB GDDR5 記憶體時脈 8008 MHz

記憶體頻寬 192 GB/s

資料來源：

https://www.evga.com/Products/Specs/GPU.aspx?pn=f1aefe7f-3b97-444d-9859-ce b772861f3f