國立臺灣大學理學院海洋研究所 碩士論文
Institute of Oceanography College of Sciences National Taiwan University
Master Thesis
台灣北部短週期噪訊研究
1. 周遭噪訊層析成像 2. 噪訊來源研究
On Short Period Ambient Noise of Northern Taiwan
1. Ambient Noise Tomography
2.
Probing the Source of Ambient Noise
陳映年 Ying-Nien Chen
指導教授：龔源成 博士 喬凌雲 博士
Advisor : YuanCheng Gung, Ph.D.
Ling-Yun Chiao, Ph.D.
中華民國 98 年 7 月
July, 2008
摘要 摘要摘要 摘要
由於操作簡單以及相對於傳統地震學在表面波層析成像法上佔有許多優勢，
因此藉由計算測站之間連續紀錄的交互相關函數而得到測站之間格林函數的技 術，現在已被廣泛應用在地震學的研究上。本研究應用此技術於分析 2006 年北台 灣三個地震網三方向的連續紀錄，其中包含中央氣象局地震觀測網北部測站以及 中央研究院在新竹和陽明山所架設的微震觀測網。針對每條波徑，我們分別計算 垂直、徑向與橫向三個方向的交互相關函數，垂直與徑向的交互相關函數代表雷 利波的格林函數，而橫向的則代表洛夫波的格林函數。接著我們測量交互相關函 數 1 秒至 5 秒之間的群速度與相速度值。經過資料篩選之後，我們使用品質穩定 的頻散曲線，並配合多重尺度參數法，反演北台灣雷利波與洛夫波的二維相速度 與群速度速度構造，再利用二維速度構造的結果，建構北台灣淺層地殼三維速度 構造。結果顯示，高解析度的速度模型與地質構造單元十分吻合。
我們也藉著以下方法來研究周遭噪訊的特性：（1）分析交互相關函數正負方向 訊號的相對強度；（2）測量交互相關函數相對於一整年之平均值的振幅變化，建立 背景能量起伏與時間和方位角之間的關係；（3）計算沿海測站連續紀錄的頻譜密度 以及方均根值隨時間的變化。結果顯示：（1）近岸的海浪可能是主要的噪訊來源，
而當能量從海洋傳到陸地時，海底地形可能扮演著重要的角色；（2）不同方向之交 互相關函數在時空變化的特性上十分相似，意味著不同方向的噪訊在淺層地殼中藉 由散射而被充分混和，達到近似散射場的環境；(3) 大氣的擾動可能是影響交互相 關函數時空變化特性的原因。
關鍵字 ： 周遭噪訊、多重尺度參數法、背景能量變化
Abstract
Retrieving Green functions between stations by cross-correlating continuous seismic records has quickly become a popular technique in seismology for its operational simplicity and various advantages over traditional surface wave tomography.
We apply this technique to three component continuous seismic data recorded from three networks in northern Taiwan, including Tatun Volcanic Area array, Hsinchu array and northern part of Central Weather Bureau Seismic Network, for the time period from Jan, 2006 to Dec, 2006. For each station pairs, we derive Love waves from T–T (transverse) component cross-correlation functions (CCF), and Rayleigh waves from Z-Z (vertical) and R-R (radial) component CCF respectively. We measure group and phase velocities for the period range from 1 to 5 seconds. With careful data selection, the qualified dispersion curves are used to derive two dimensional (2-D) phase and group velocity maps for both Rayleigh and Love waves with multi-scale inversion technique. The 2D maps are then used to develop a 3-D shallow velocity structure of the northern Taiwan.
We also attempt to probe the sources of ambient noise by several approaches: (1) analyzing the relative strength between the causal and acausal CCF; (2) measuring the relative strength of CCF amplitudes with respect to their own annual average as a function of time and azimuth to determine the background energy flow; and (3) computing power spectra density and rms amplitudes as a function of time for
representative costal stations. The results show that (1) offshore ocean waves are likely
observed temporal variations of CCF.
Keywords: seismic ambient noise、 multi-scale parameterization、background energy flow
目錄 目錄 目錄 目錄
口試委員會審定書 i
中文摘要 ii
英文摘要 iv
目錄 v
圖目錄 vii
表目錄 xi
第一章 緒論 1
第二章 理論背景與資料處理流程 4
第三章 頻散分析與資料篩選 10
3-1 CCF Whitening 10
3-2 群速度頻散分析 12
3-3 相速度頻散分析-影像轉換技術 14
3-4 資料篩選 18
第四章 表面波層析成像 25
4-2 三維速度構造 33
第五章 北部周遭噪訊分析 39
5-1 CCF 與北部噪訊特性 39
5-1-1 CCF 非對稱性與北部噪訊來源 39
5-1-2 CCF 非稱性之際節變化 41
5-1-3 周遭噪訊背景能量變化 47
5-1-4 CCF 強度與測站間距離的關係 49
5-2 單一測站連續紀錄分析 51
5-2-1 單一測站譜密度之機率函數 51
5-2-2 單站連續紀錄 rms 分析 54
第六章 結論 62
參考文獻 64
附錄 測站儀器修正 68
圖目錄 圖目錄 圖目錄 圖目錄
圖 1-1 測站分佈圖與有效資料時間。 3
圖 2-1 原始資料與 3D-Normalization 之後的結果。 7 圖 2-2 測站 HC04 垂直方向分別使用三種 Normalization 之後，相
對於原始資料的扭曲程度。 8
圖 2-3 資料經過 1-bit 正規化與 3D 正規化之後，計算 CCF 的結果。 9 圖 3-1 三方向(Z-Z、R-R 與 T-T) CCF 頻譜圖。 11 圖 3-2 測站 ENA 與 HC23 之 CCF，whitening 前後的波形以及頻譜比
較。 11
圖 3-3 以 HC24-HC11 為例，比較訊號 Whitening 前後頻散分析的差
異。 12
圖 3-4 頻散分析示意圖。 13
圖 3-5 走時-週期影像圖。 15
圖 3-6 速度-週期影像圖。 16
圖 3-7 挑選合理相速度頻散曲線示意圖。 17
圖 3-8 波徑 TAP-NSK 之 CCF 篩選過程。 19 圖 3-9 資料篩選前，三方向 CCF 的訊噪比跟 STD 與一整年群速度比
20
圖 3-13 測站相對位置以及經篩選之後的 CCF 之 STD。 24 圖 4-1 利用多重尺度參數法所產生的不同等級的模型參數與模型。 28 圖 4-2 運用多重尺度參數法所得到的模型解，其擬合程度與模型變
異數之間的消長曲線。 29
圖 4-3 不同 damping（λ）值之 0.35Hz 雷利波群速度速度構造圖。 30 圖 4-4 利用多重尺度參數法逆推之北台灣二維速度構造。 32 圖 4-5 台灣北部地質圖。（經濟部地質調查所） 33
圖 4-6 雷利波的敏感算核。 35
圖 4-7 初始模型與逆推之後的模型符合頻散曲線的能力圖。 36
圖 4-8 逆推前後符合資料的能力的差異圖。 37
圖 4-9 以 Kim 為初始模型以及三維速度構造逆推結果。 38 圖 4-10 以 Wu ㄧ為初始模型以及三維速度構造逆推結果。 40 圖 5-1 測站間一整年之相同方向的 CCF 非對稱性，以及研究區域
周圍 12 個測站相對於測站 NSK 之 CCF 的非對稱性。 40 圖 5-2 各測站不同方向（R-Z、Z-R）CCF 的非對稱性以及測站 NSK
對於周圍 12 個測站之非對稱圖。 41
圖 5-3 各測站以及周圍測站對測站 NSK 季節性 CCF(Z-Z)之非對稱
圖。 43
圖 5-4 各測站以及周圍測站對測站 NSK 季節性 CCF(R-R)之非對稱
圖。 44
圖 5-5 各測站以及周圍測站對測站 NSK 季節性 CCF（T-T）之非對稱
圖。 45
圖 5-6 2005～2006 年 HC16- HC05 之 Z-Z 與 R-R 的 CCF 相對強度
變化。 46
圖 5-7 不同方角之 BEF 隨時間變化圖 48
圖 5-8 最大振幅與測站距離的關係 50
圖 5-9 測站 TAP 垂直方向連續紀錄之 PDF 分析與 Z-Z 方向 CCF 的
結果。 52
圖 5-10 沿海 6 測站全年之垂直方向連續紀錄之 PDF 以及其相對地
理位置。 53
圖 5-11 測站 TWB1 rms 未經處理前之連續紀錄以及分析結果。 56
圖 5-12 測站 TWB1 去除儀器突波後，三方向連續紀錄 rms 分析結果。 57
圖 5-13 測站 rms 紀錄與其頻譜能量分析圖。 58
圖 5-14 北部 7 測站與 TWA 之間的 CCF 波包振幅變化分析。 59
圖 5-15 三方向 CCF 振幅變化分析。 60
圖 5-16 測站 TWB1 七月份 rms 變化圖。 60 圖 5-17 2006 年 7 月份測站 NST 與 HC23 之間 Z-Z 與 R-R 方向之 CCF
分析。 61
圖 3 每天的 CCF 對一整年疊加的 CCF 之相關係數圖。 70
表目錄 表目錄 表目錄 表目錄
表 4-1 逆推二維速度構造所使用之波徑數 25
表 5-1 分潮週期 57
表 A 儀器極軸反轉時間表 70
第一章 緒論
系統中的噪訊(noise)並非無用的，早期統計物理學的研究顯示，藉由計算線性 系統中任兩點噪訊的連續紀錄之交互相關函數(Cross-correlation Function, CCF)可 以得到兩點之間的格林函數(Green’s function)(R.Kubo, 1966)，相當於介質的脈衝反 應函數（impulse response function）；換言之，研究噪訊可以得到系統的物理特性。
此概念建立在 20 世紀早期布朗運動的研究(Einstein, 1905 )以及日後通訊儀器雜訊 的研究上(Lee et al., 1950)，近年來則廣泛地應用在日震學(Duvall et al., 1993)、聲 學(e.g. Weaver and Lobkis, 2001;Derode et al., 2003 )、流體力學(Godin, 2006)及熱力 學(e.g. Weaver and Lobkis, 2002; Godin, 2007)的研究，而 Campillo 和 Paul (2003) 首先將此方法應用在地震學的研究上。在地震學的應用上，這是一門相當年輕的 嶄新領域，充滿著各種可能性及挑戰。
在地表上兩測站之間的格林函數主要反映在表面波的能量上(Shapiro et al., 2005)，因此藉由計算測站間連續紀錄之 CCF 可以得到兩測站間表面波之格林函 數，再藉由頻散分析(Dispersion analysis)，可以得到兩測站間表面波傳遞之群速度 與相速度，藉以了解地殼及上部地函之速度構造，這種方法稱為周遭噪訊層析成 像法(Ambient seismic Noise Tomography, ANT )。相對於傳統的分析方法，由於此 方法並不需要天然地震，故而沒有震源不確定性（source location, mechanism and fault finiteness）的問題，其應用範圍亦不受限於天然地震的時空分布。ANT 的解 析度與精確性乃取決於測站分布的密度以及訊號的長度，即使在地震發生數目極 少的區域，依然能夠藉此方法大幅改善傳統表面波層析的解析度。
除了在表面波層析成像法的應用之外，分析周遭噪訊也可以提供大地構造的 即時訊息。例如火山噴發前常有岩漿入侵或是岩漿庫膨脹等現象，這都可能造成 地表震波格林函數的變動。因此，藉由分析測站間的 CCF 所提供的即時訊息，即 可作為評估火山活動的依據（Brenguier et al., 2008 ; Courtland, 2008）。
周遭噪訊除了可以提供地震學上的應用之外，造成地表噪訊的自然機制也是 目前相當熱門的科學問題。目前的若干研究顯示，周遭噪訊的產生主要源自於海 洋、大氣與陸地間的交互作用(Rhie and Romanowicz, 2004；Kedar and Webb, 2005；
Rhie and Romanowicz, 2008)。根據研究顯示，波浪碎波時拍擊水體而引發的長波 振盪，即所謂的亞重力波(infragravity waves)，可能是地球連續震盪(hum)的主要來 源，其能量藉由表面波的型式傳遞到陸地上(Henderson et al., 2006；Uchiyama and McWilliams, 2008)。因此研究周遭噪訊的時空變化與特性也提供了沿岸地形與海水 運動的訊息。當噪訊來源是均勻分布在空間時，所得到的交互相關函數在正負兩 邊會有相當好的對稱性；反之，不對稱性則暗示了噪訊來源之方向性。若干研究 指出，較強的噪訊源通常指向與測站－海岸間距較短的方向(e.g. Stehly et al., 2007 ;
Gu et al., 2007； Yang et al., 2008； Kraeva et al., 2009；Keith et al., 2009; .尤, 2008)。
在本研究中，我們使用北臺灣 34 個測站 2006 年之三方向（three component）
短週期連續地震資料，其中包括了中央氣象局地震觀測網北部的 17 個測站，以及 中央研究院在新竹(2004.5~2006.12)與陽明山(2003.5~)所架設的微震觀測網的 17 個 測站(圖 1-1)。藉由計算其 CCF 所得到短週期(2s~5s)的表面波格林函數，經過頻散 分析與適當的資料篩選之後，配合多重尺度有限參數法，分別反演了不同週期之 北台灣淺層地殼的二維剪力波(shear wave)速度構造，再配合不同的初始模型，反 演三維速度構造。
此外，我們也藉由分析 CCF 之非對稱性與測站之間方位角的關係、其隨時間 不同所產生的變化情形以及在頻率域的能量分布狀況，探討北台灣短週期噪訊的 特性。
(b)
圖 1-1，(a)測站分佈圖，圖中黃色三角形為新竹微震觀測網之測站，紅色三角形 為陽明山微震觀測網之測站，藍色三角形則為中央氣象局地震監測網之測站，紫 色代表測站間之大圓路徑。（b）有效資料時間。
月份 月份 月份 月份 (a)
第二章 理論背景與資料處理流程
在線性系統中，藉由計算均勻散射波場中測站間的連續紀錄的交互相關函數
（cross-correlation function ,CCF），可以得到測站間的格林函數。在不考慮訊號能 量會隨傳遞距離增加而衰減的情況下，此關係可以藉由簡單的推導證明(Weaver and Lobkis, 2001):
有限空間中的彈性波場可以用常態模式展開(Normal mode expansion)表示:
( )
( , ) m( ) msin m mcos m
m m
u t u a ω t b ω t
=
∑
ωr +r (式 2-1) 將式 2-1 對時間微分可以得到速度場:
( )
( , ) m( ) mcos m msin m
m
v r t =
∑
u r a ω t−b ω t (式 2-2) 而空間中兩點之間的理論格林函數:( ) ( ')
( , ', ) m m sin m ( )
m m
u u
G t ω t H t
=
∑
rω r ×r r (式 2-3)
0B0B0B0B0B
其中，H t 為 Heaviside step function ( )
對時間微分之後，可以得到速度場的格林函數
( )v ( , ', ) m( ) m( ') cos m ( )
m
G r r t =
∑
u r u r ω t×H t (式 2-4)在式 2-2 中，a 與m b 是由波場來源所決定的模式係數，在均勻散射場中，m a 與m bm 為無相關的隨機係數，因此
2
n m nm
a a =F δ
2
n m nm
b b =F δ
配合式 2-1 以及式 2-5 計算散射波場中兩點之間的 CCF:
( )v ( ) ( , ) ( , )
AB A B
C τ = v r t v r t+τ
( )
{
( )
( )
( ) }
( )
,
( ) ( ) ( )
cos cos cos sin sin cos sin sin
v
AB n A m B
n m
n m n m
n m n m
n m n m
n m n m
C u u
a a t t
a b t t
b a t t
b b t t
τ
ω ω τ
ω ω τ
ω ω τ
ω ω τ
=
× +
− +
− +
+ +
∑
r r
2 m( )A m( ) cosB m
m
F u u ω τ
=
∑
r r . （式 2-6）比較式 2-6 與式 2-4，可以得到
{ }
( ) 2 ( ) ( )
( ) ( , , ) ( , , )
v v v
AB A B A B
C τ =F G r r τ +G r r −τ （式 2 -7）
因此，藉由計算兩測站間連續記錄的 CCF，我們可以得到兩測站間的格林函 數。
然而，地表並非理想的均勻散射場；地震訊號、當地的特殊訊號、儀器的突 波…等，都會影響我們計算 CCF 的結果，因此在進行交互相關計算前，必須先將 資料進行處理，讓它比較符合均勻散射場的假設。
在計算各個測站間地震連續紀錄的 CCF 之前，我們分別參考 Bensen et al.(2007) 所提出來的資料處理流程，以及尤(2008)研究北台灣周遭噪訊雷利波，發現其能量 主要分佈在 2 秒至 5 秒之間的特性，先對原始訊號進行移除儀器響應(instrument response)，再降低取樣率(down sampling)至每秒兩點以增加運算速度，接著藉由帶 通濾波(bandpass filtering)，將訊號頻段限制在 1 秒至 10 秒之間。最後，為了避免 進行 CCF 運算時受到儀器突波(glitches)、地震訊號或是近測站雜訊的影響，我們 將資料進行正規化(Normalizatoin)處理。
ㄧ般而言，常用的正規化方式有兩種，分別是 running-absolute-mean 正規化 (Bensen et al., 2007)與 1-bit 正規化 (Shapiro et al., 2005)。根據 Bensen et al.(2007) 的研究指出，經由 running-absolute-mean 正規化之後的資料所計算出的 CCF 其訊
噪比(signal to noise ratio, SNR）比 1-bit 正規化來得高。1-bit 正規化的優點在於其 運 算 較 為快 速 ，而 running-absolute-mean 正 規化則較能保持原始訊號的相位 (phase)；然而，使用此兩種方法都抹除了原始訊號在空間上的能量分佈情形。因 此 ， 本 研究 使 用另 一種 正 規化 的 方式 （ 3C-Normalization），將三個方向(three components)的振幅大小除上每一時間點空間上的長度（式 2-8，圖 2-1）。這種正規 法的方式除了可以保有訊號空間上的特性之外，在保存波相的能力方面則介於 1-bit 與 running-absolute-mean 這兩種正規化之間（圖 2-2）；在運算上，由於此正 規化方法不受座標旋轉的影響，因此僅須對於原始訊號進行一次正規化處理，大 幅提昇水平方向訊號處理的速度。
( , )
( , ) 1, 2,3
| ( ) |
i i
S i t
S i t i
= t =
S （式 2-8）
另外，我們也比較了 1-bit 與 3C 正規化方法之後所計算出的 CCF。圖（2-3）
顯示單日 CCF 運算的結果：就波形而言兩個方法並沒有太大的差異，然而使用 3C 正規化之後的資料依然保存原始訊號在三方向的能量比例，因此 CCF 的振幅變化 更能表現出能量在不同方向傳遞的相對大小。
在 進 行 正 規 化 之 後 ， 我 們 分 別 計 算 測 站 間 三 方 向 (Vertical 、 Radial 與 Transverse，之後分別以 Z-Z、R-R 與 T-T 表示)的 CCF。為了得到測站間 R-R 與 T-T 方向的 CCF，需先將水平方向的訊號依測站之間的方位角進行座標旋轉(式 2-9)，並投影成與波徑方向垂直(橫向,Transverse)與平行(徑向，Radial)的振動，再 分別計算兩測站間 T-T 與 R-R 的 CCF。最後，每條波徑分別得到 Z-Z、R-R 與 T-T 三組的 CCF，而 Z-Z 與 R-R 的 CCF 代表雷利波的格林函數，而 T-T 的 CCF 則代 表洛夫波的格林函數。
圖 2-1 原始資料（a~c）與 3C Normalization (d~f)之後的結果。在 3C Normalization 之後，並不會改變各個方向的能量比例；而 1-bit Normalization 之後，各個方向的 能量相同。
（a）
（b）
（c）
（d）
（e）
（f）
圖 2-2 測站 HC04 垂直方向分別使用三種 Normalization 之後，相對於原始資料的 扭曲程度。每 3 小時計算一次資料之頻譜圖，將其最大值定為 1；並在計算正規 化之後的波相與原始波相的差值後，統計一整年之資料在頻率域之分布機率圖以 及波相改變的機率圖。(a)原始資料的頻譜分布機率(b)1-bit normalization 的頻譜分 布 機 率 (c)3C-normalization 的 頻 譜 分 布 機 率 (d) Running-absolute mean Normalization 的頻譜分布機率 (e) 1-bit normalization 之後波相與原始資料的差 異 (f) 3C-normalizatio 之後波相與原始資料的差異(g) Running-absolute mean Normalization 之後波相與原始資料的差異。在圖 a~d 中，藍色細線表示各頻率密 度最大值之連線，而在圖 e~g 中，藍色細線表示各頻率平均值之連線，黑色的線 表示與平均值相差一個標準差之連線，可以清楚看到 3C-normalization 保存波相的 能力介於 1-bit normalization 與 Running-absolute mean Normalization 之間。
(a) (b) (c) (d)
(f) (g) (e)
圖 2-3 資料分別經過 1-bit Normalization 與 3C Normalization 之後，計算 CCF 的結 果。紅線與黑線分別表示 1-bit Normalization 與 3C Normalization 的 CCF，並將振 幅限制在正負 1 之間；兩者的波形差異極小。
第三章 頻散分析與資料篩選
為了增強有效訊號的品質，我們先將每天的 CCF 疊加成以月為單位和一整年 為單位的資料，並將 CCF 正負方向的訊號進行平均。再根據群速度頻散分析的結 果，設定篩選資料的標準。最後將篩選過後的資料，進行相速度的頻散分析。以 下將分別介紹資料頻散分析的方法以及資料篩選的流程。
3-1 CCF Whitening
在均勻的散射場中，兩測站之間連續記錄的 CCF 就是兩測站間的格林函數，
亦即介質的脈衝反應函數。然而，地表並非均勻的散射場。噪訊產生的機制、其 強度在空間中的分佈差異以及儀器響應等因素，使得本研究中 CCF 的能量主要集 中在三秒附近，且能量頻譜也隨著測站組的方位而改變 (圖 3-1)。
由於 CCF 在各頻率的能量並非一致，因此在頻散分析時可能會造成所量測的 速度值無法代表帶通濾波頻段的中央頻率(central frequency)。為了確保所量測的速 度值足以代表所選取的中央頻率，我們將 CCF 在頻率域進行 whitening(圖 3-2)，亦 即將訊號頻率域的能量拉平；然而，whitening 也有可能放大了 CCF 中不可靠的訊 號，造成測量的誤差，所以，適當的資料篩選流程是必須的。
圖 3-1 三方向(Z-Z、R-R 與 T-T)CCF 頻譜圖。圖中縱軸為方位角，正北為零，順 時針為正。橫軸為週期。圖中顯示 CCF 的能量集中在 2~5 秒之間，而隨著方位角 的不同，頻譜也隨之改變。
圖 3-2 測站 HC24 與 HC11 之間的 CCF，whitening 前後的波形以及頻譜。左圖為 whitening 前的結果，右圖則為 whitening 之後的結果。
3-2 群速度頻散分析
我們使用時頻譜分析法（Frequency-time analysis， Levshin et al. 1989）測量 訊號不同週期的群速度值。首先，利用帶通濾波(頻寬 0.2Hz)針對不同週期進行濾 波，再將濾波後的訊號藉由希爾伯特轉換(Hilbert transform)計算該訊號的波包，並 將波包的最大值視為訊號到達的時間。兩測站之間的距離除以波包到時，即可得 到該週期訊號的群速度值，示意圖如圖(3-3) 。
我們比較 Whitening 對於分析結果的影響，從圖 3-4 中發現其分析結果並不相 同。頻譜分析顯示原始 CCF 在 0.1Hz 和 0.3Hz 之間有兩能量峰值(圖 3-2)，是以帶 通濾波的訊號並無法確實代表所選取的中央頻率。在頻散分析前將 CCF 進行 Whitening，將可改善此問題。
圖 3-3 群速度頻散分析示意圖，以 NSY-NSK 為例，(a)將原始訊號依不同的中央
圖 3-4 以 HC24-HC11 為例，比較訊號 Whitening 前後頻散分析的差異。左圖為 Whitening 前頻散分析結果，右圖為 Whitening 後頻散分析結果；由於原始訊號 在 0.1Hz 和 0.3Hz 之間有兩能量峰值（圖 3-2），因此無法保證濾波之後的訊號 能 夠 代 表 所 選 取 的 中 央 濾 波 頻 率 ， 故 本 研 究 在 頻 散 分 析 前 先 將 訊 號 進 行 Whitening。
3-3 相速度頻散分析-影像轉換技術(Image transformation technique)
在遠場（far field）的前提下，也就是測站與震源之間的距離必須大於表面波 波長的三倍，fundamental mode 表面波的格林函數可以寫成（Dahlen and Tromp, 1998）
Re{ ( ) exp( )} (8 )1/ 2cos( )
AB AB 4
G i t kS k t π
ω −ω ≈ π −ω + ， （式 3-1）
其中 AB
AB
k C
= ω ，為 AB 之間的平均波數（wave number），t 為表面波走時，CAB為
AB 之間的平均相速度， 為 AB 之間的距離，S sin( )
= R ，代表表面波的幾何衰 減（geometrical spreading），R 則為地球半徑。
當
AB 4
k t π
ω
− +
=0 時，相當於格林函數的波峰，那麼我們可以求得在頻率ω之 下的平均相速度
( ) / 8
CAB T
t T
= ∆
− ， （式 3-2）
T 為表面波之週期。
運用上述的原理，我們求出格林函數的相速度主要分成下面幾個步驟：
1. 我們先利用帶通濾波(頻寬 0.1Hz)針對不同週期進行濾波，接著將濾波後的 CCF 的振幅正規化，亦即將振幅限制在正負 1 之間。最後，我們得以建立走 時-週期影像（time-period image），振幅大小則用顏色表示(圖 3-5)。
2. 利用式 3-2 將走時-週期影像中的走時 t 轉換成相速度（圖 3-6）；圖中每一條 由不同週期的峰值所組成的最平滑曲線，都有可能是真正的相速度的頻散曲
3. 我們把利用時頻譜分析法所量得的群速度頻散曲線當作參考值，試圖將所有 可能的相速度頻散曲線中，其相對應的群速度頻散曲線與參考值相差最少 的，當作最合理的相速度頻散曲線。然而，理論上 cycle skip 所造成的每一條 相速度頻散曲線的群度度值都是一樣的。為了降低 cycle skip 的效應，增加此 方法的可行性，我們只將振幅大於 0.6 的峰值連線；換言之，每一條頻散曲線 不僅在各頻率的能量皆必須大於 0.6，而且必須是最平滑的曲線。
4. 最後，我們算出每一條頻散曲線相對應的群速度值，並計算週期 2 秒到 5 秒 之間該群速度值與用時頻譜分析法量得的群速度值的差值，並將差值最小 的，當成最合理的相速度頻散曲線（圖 3-7）。
圖 3-5 走時-週期影像圖。橫軸為週期，縱軸為訊號走時，而顏色則代表訊號經過 正規化之後的振幅大小，藍色的線表示中央濾波頻率為 3.5Hz 的 CCF。
圖 3-6 速度-週期影像圖。將圖 3-2 的走時經由公式 3-2，轉換成速度。圖中六條 紅線代表各個週期中峰值（local maxima）的最平滑連線，它們都是可能的相速度 頻散曲線。
圖 3-7，圖中棕色曲線為時頻譜分析計算群速度的結果，也就是『參考頻散曲線』。 在選定的週期（2～5s）之內，分別計算各個頻散曲線（綠色與藍色，圖 3-4 中各 個頻散曲線相對應的群速度頻散曲線）與參考頻散曲線的差。計算結果顯示，藍 色曲線相差最少，而紅線則是相對應的相速度頻散曲線。
3-4 資料篩選
首先，我們先剔除不符合遠場假設（∆ ≤3 λ ）的資料。接著，為了有效增強 訊號的品質，我們將每天的 CCF 分別疊加成以月以及一整年為單位的訊號，並計 算其群速度值。通常，疊加的時間越長，所得到的訊號品質越好；然而，我們發 現某些月份的群速度值與一整年疊加的結果所計算出來的相差甚多，造成該組訊 號的群速度的標準差(standard deviation, STD)變大。
為了避免不好的訊號在疊加的過程中影響最後的訊號品質，我們先剔除該月 份群速度值與一整年疊加的結果的群速度值相差大於 0.5km/s 的資料（圖 3-8），若 剩餘的月份超過 3 個月，便將剩餘的資料重新疊加，重新測量其群速度值並計算 其標準差。若該組資料的標準差大於群速度值的 10％，就將該組資料剔除(圖 3-9)。
同時，為了衡量篩選前後整體訊號的強度變化，我們計算了訊號的訊噪比：
將訊號振幅最大值除上整條訊號(0~300s)之 rms。我們發現篩選之後的資料，不但 標準差明顯下降，而且訊噪比提高，顯示整體而言訊號變得比較穩定(圖 3-9、圖 3-10)。以標準差作為資料篩選的門檻，我們發現部份保留的 CCF 訊噪比偏低，但 其月頻散結果仍非常穩定，甚至 CCF 中主要訊號外的較低震幅訊號的一致性也相 當高，顯示 CCF 中除了 fundamental mode 的表面波外，應有更多值得進一步探討 的有用訊號（圖 3-11）。
最後，我們將篩選之後的訊號還原成原始訊號(左右訊號未平均)，並把測站依 其地理位置分成三組，計算訊號正負方向的訊噪比(圖 3-12)以及標準差(圖 3-13)。
結果發現：1.同組的訊號因其站距較短，所以訊號強度較強(平均 SNR >8)；2.不同 組間的訊號強度較弱，這可能是因為能量隨站距增加而衰減以及噪訊方向性所造
圖 3-8 波徑 TAP-NSK 之 CCF 篩選過程。(a)最上方之 CCF 為經過篩選之後疊加的 結果，下方為每個月份分別疊加的結果。黑色之 CCF 表示符合資料篩選門檻而被 保留，紅色的 CCF 則表示不符合門檻而被剔除的月份，藍色的 CCF 則是一整年疊 加的結果。（b）測站間相對位置 （c）一整年（藍色）以及每個月（紅色、黑色）
份群速度頻散分析之結果。在分析單站連續紀錄時發現測站 TAP 在 2～6 月以及 10～12 月出現了問題，藉由我們的資料篩選流程，可以有效地把有問題的資料剔 除。
圖 3-9 資料篩選前，三方向 CCF 的訊噪比跟 STD 與一整年群速度比值的關係 圖。圖中紅線表示挑選門檻，而 STD 與一整年群速度的比值高於此門檻的資料將 會被剔除。圖中 AVG<STD,SNR>則表示 CCF 之 STD 與 SNR 的平均值。
圖 3-10 經過資料篩選之後，CCF 的 SNR 與 STD 作圖。與資料篩選前相比，各方 向 CCF 平均之 STD 都明顯下降，且平均 SNR 值上升，顯示資料經過篩選之後變 得更為穩定。平均而言，Z-Z 方向的 CCF 之 STD 最小，SNR 也最大，表示資料 的品質最好。從圖中顯示 STD 與 SNR 並沒有明顯的相關性，也就是說，即使 SNR 較低，仍然是很穩定的訊號。
22
圖 3-12 測站相對位置以及經篩選之後的 CCF，還原成原始訊號後（左右訊號未 平均）正負訊號之 SNR。(a)本研究區域之測站依地理位置分成三組，以不同的 顏色加以區別。圖 b~d 表示不同方向之 CCF 的 SNR，縱軸為噪訊之來源測站，
橫軸為接受訊號之測站，字的顏色表示組別。
（a）
（b）
（c） （d）
圖 3-13 測站相對位置以及經篩選之後的 CCF，還原成原始訊號後（左右訊號未 平均）正負訊號之 STD。(a)本研究區域之測站依地理位置分成三組，以不同的顏 色加以區別。圖 b~d 表示不同方向之 CCF 的 STD，縱軸為噪訊之來源測站，橫 軸為接受訊號之測站，字的顏色表示組別。跟圖 3-12 比較發現，大致上訊號的 SNR 越高，其 STD 越小，但是即使訊號的 SNR 偏低，其 STD 也可能相當小。
(a)
(b)
(c) (d)
第四章 表面波層析成像
本研究針對由CCF表面波頻散結果，使用多重尺度有限參數法（Chiao and Kuo, 2001;Chiao and Liang, 2003）進行二維速度構造的反演，並將不同週期的雷利波與 洛夫波的群速度與相速度的二維速度構造的結果，配合不同的初始模型，反演符 合各點頻散結果的一維速度構造（Hermann, 2004），從而建立北臺灣淺部的三維速 度構造。以下將分別介紹2維與二階段反演的3維速度構造的方法與結果。
4-1 運用多重尺度參數法反演二維速度構造
我們分別測量了 0.25、0.3、0.35、0.4、0.45 以及 0.5 Hz 之表面波之群速度與 相速度。在資料篩選之後，不同頻率可以使用的波徑數如表 4-1。其中，我們使用 T-T 方向的資料進行反演北台灣洛夫波的速度構造；在雷利波的資料選取上，我們 使用 Z 方向的資料進行相速度速度構造的反演，並結合 Z-Z 與 R-R 資料中標準差 較小的資料，反演群速度速度構造。
表 4-1 逆推二維速度構造所使用之波徑數
Frequency （Hz） Rayleigh (U) Rayleigh (C) Love (U) Love (C)
0.25 485 402 315 315
0.3 493 425 342 342
0.35 497 437 355 355
0.4 496 434 359 359
0.45 485 430 358 358
0.5 454 397 319 319
4-1-1 多重尺度參數法
在逆推問題中，我們常用有限的觀測資料來反演無限自由度的連續物理模 型，因此無法得到滿足資料的唯一模型解。在本研究中，我們使用阻尼最小平方 法（DLS），在滿足資料的基礎上同時也要求最小化模型的長度（model norm），也 就是求得式 4-1 的最小值的模型解。
F =||Gm d− ||2 +λ||m||2 式 4-1 =(Gm−d) (T Gm−d)+λm mT
=m G Gm m GdT T − T −d GmT +d dT +λm mT 最小化 F :
∂F/∂mT =G Gm G dT − T +λm= 0 ⇒(G GT +λI m) =G dT
m=(G GT +λI)−1G dT =G dg 式 4-2
接著我們使用多重尺度有限參數法，將模型空間進行離散小波轉換(Discrete Wavelet Transform, DWT，G' =GWT，W 為小波轉換運算子）T ，拆解成不同尺度的 基底，再根據逆小波轉換(Inverse Wavelet Transform, IWT)，先反演長波長訊號所提 供的大尺度的構造訊息，再依序往小尺度的構造反演，累積模型中各區域的細微 構造；而模型所能解析的最小尺度則取決於由各區域的波線覆蓋密度。
因此，在多重尺度有限元素參數法的運算之下，式子 4-2 可以改寫成：
( T T ) 1 T T
m=W W G GW+λI −W G d 式 4-3 m: 模型參數 （model）
W: 運算子（operator），這裡指的是逆小波轉換
我們先設定一個涵蓋北台灣的矩形空間當作我們的模型，而整個矩形設為等 級 1 的模型參數。接著將矩形經度方向的中點與緯度方向的中點連線，將模型分 成 2X2 的矩形區塊，其模型參數為 2，依此類推，一直到等級 7，總共有 4225 個
格點（(2(7 1)− +1)2＝4225）。不同等級的模型參數所反演出來的構造以及格點大小
如圖 4-1，隨著等級的增加，模型所能提供的細微構造也就越清楚；而相對於波徑 密度較低的地區，依然保有大尺度的構造特性。
接著，為了避免過份滿足資料中的誤差而造成模型解的錯誤，我們可以藉由 選取阻尼因子λ，來決定相信資料的程度；換言之，我們寧可保守選擇資料中可 靠的部份也不要過份滿足資料。在本研究中，我們分別使用擬合程度(variance reduction，v.r.)與模型變異數(model variance)來估計模型符合觀測資料的能力以及 模型的可信度
2 2
|| ||
. . (1 )
|| ||
d Gm
v r d
= − − 式 4-4
，而λ值可以用兩者之間的消長曲線（trade-off curve，圖 4-2）來決定：模型變 異數越小，代表模型越平滑，也越可靠；而擬合程度越高，則表示模型越能反應 資料（圖 4-3）；通常λ值會採用消長曲線的轉折點處。
圖 4-1 以雷利波 0.35Hz 群速度層析成像為例，damping 值為 0.1，利用多重尺度
圖 4-2 0.35Hz 雷利波群速度運用多重尺度參數法所得到的模型解，其擬合程度與 模型變異數之間的消長曲線，藍點上的數字為λ的大小。在此例中我們選取λ值 為 0.1。
4-1-2 二維速度反演結果
在資料篩選之後，我們分別反演了表面波 2～4s 之間六個頻率的群速度與相速 度之速度模型（圖 4-4）；雖然各模型間所提供週期相差不大（≤2s），但是由於量 測前所有的 CCF 已經經過 whitening，因此可以確定不同週期間的模型差異反應了 不同深度的速度構造。
整體而言，這些短週期表面波的二維速度分布與淺層地質構造的特性相當吻 合。在臺北盆地、蘭陽平原以及西部麓山帶屬於相對低速區，而在雪山山脈以及 中央山脈北部地區則屬於相對高速區（圖 4-5）：台北盆地以及蘭陽平原以沉積層 為主，西部麓山帶則之岩層以第三紀之碎屑岩為主，地表則覆蓋著由膠結鬆散的 砂頁岩互層所組成的頭嵙山層、紅土台地以及沖積層；而雪山山脈以及中央山脈 北部主要以硬頁岩-板岩相所構成的輕度變質岩組成。
圖 4-5 台灣北部地質圖，三角形為本研究使用之測站。（經濟部地質調查所）
4-2 三維速度構造
在解出不同週期的二維速度模型之後，可以得到每一個格點上的頻散曲線，
我們可進一步逆推各格點不同深度的速度分佈，並藉由這種兩階段的反演建構三 維速度構造。
我們分別採用 Kim（Kim et al, 2005）以及 Wu（Wu et al, 2007）使用體波資料 所建立的台灣三維速度構造當成初始模型。在深度的分層則為：前 0.5 公里分別為 0.2 與 0.3km，其餘的格點大小則固定為 0.5km，一直到 10 公里為止。藉由計算模 型格點上雷利波的敏感算核（sensitivity kernel），證實不同的週期的表面波確實能 夠反應不同深度的速度構造（圖 4-6），而本研究中所使用的訊號頻段，大致上反 應了地表淺層四公里以上的速度構造；從dc dVs/ 以及dc dVp 的曲線中發現，表面/ 波對於剪力波的速度變化較為敏感，因此，我們選擇在固定Vp Vs 比值（1.73）的/ 情形之下，分別反演每一格點上不同深度的剪力波波速。
在圖 4-7 中，我們分別計算了不同的初始模型以及進行逆推之後，滿足頻散曲 線的能力。結果顯示，兩初始模型在速度構造相對低速帶，滿足表面波頻散曲線 的能力明顯較差，而經過逆推之後，兩模型滿足頻散資料的能力皆大幅提升。我
們計算了逆推前後滿足資料的能力差異(圖 4-8)，發現在速度相對低速區改善最 多，這也顯現使用表面波資料進行反演淺層速度構造的優勢：初始模型所使用的 體波其行進方向幾乎垂直地表，相較於在地表傳遞的表面波，對於淺層區域的側 向速度變化無法提供很好的解析能力。整體而言，逆推前 Wu 的模型比 Kim 的模 型還要符合二維的模型結果(Wu : 74﹪; Kim:67﹪)，經過逆推之後，兩者皆符合 95
％以上，且 Kim 的模型又高了 Wu 的模型 0.3％，顯示各點的一維逆推方法在不同 的初始條件下，其對資料的滿足度均可達到高度的收斂。然而，進一步比較發現，
兩結果在不同深度的平均速度有著明顯的差異(圖 4-9、4-10)，顯示此一維反演對 初始模型有著相當的依賴度。
在圖 4-9 以及圖 4-10 中，我們分別展示了所使用的初始模型以及相對應的逆 推結果，接下來我們以逆推前後模型的變異數(model variance)以及平均速度變化進 行討論。兩初始模型的模型變異數都小於 0.25，而逆推之後模型變異數都大於 0.7，
顯示我們得到的模型較初始模型更具有異質性，亦即可以提供更多的細節。就速 度變化而言，逆推之後在相對低速區速度比初始模型來的低，而在相對高速區則 比初始模型來的高；整體而言，逆推之後的平均速度皆小於初始模型。而不同的 初始模型所得到的結果十分接近：在沖積層、紅土台地堆積以及蘭陽平原屬於相 對低速帶，而從淺層一公里的剖面中發現，在雪山山脈西部，與其走向平行處也 出現局部的低速帶，這與地質構造中在陸臺沈積環境下造成的含煤地層-木山層以 及淺海相沈積環境的南港層相符。
總言之，短週期的 ANT 確實可以提供淺層構造的更多細節，然而不同的初始 模型所得到的結果確有相當的差異，顯示初始模型的選擇會影響我們最後的結
圖 4-6 雷利波的敏感算核。圖 a~c 為不同週期的雷利波敏感算核，藍線與黑線分 別表示dc dVs/ 以及dc dVp 隨深度變化之曲線，而圖/ d 則為速度構造模型，圖中 藍線與黑線分別表示 Vs 與 Vp 隨深度變化的曲線。
（a） （b） （c） （d）
圖 4-7 初始模型與逆推之後的模型符合頻散曲線的能力圖。逆推之後兩模型滿足 頻散曲線的能力皆達到 95％以上，已看不到初始模型的影響，顯示此逆推方法可 得到收斂的結果。
圖 4-9 以 Kim 為初始模型以及三維速度構造逆推結果。左圖為深度 1 公里、2 公 里以及 3 公里之初始模型，右圖則為相對應之逆推結果。
第五章 北部周遭噪訊分析
在本研究中，周遭噪訊連續紀錄之 CCF 除了提供淺層地殼速度反演的素材之 外，藉由分析異常的 CCF 訊號亦可提供儀器修正的依據，而 CCF 訊號強度的非對 稱性以及其時空變化，更隱含著噪訊來源的時空分佈以及強度變化的訊息。本章 節將藉由分析 CCF 的時空變化以及單站的連續紀錄，探討北部周遭噪訊的時空特 性以及測站訊號特性對於 CCF 的影響，至於儀器修正的過程與結果，請參見附錄。
5-1 CCF 與北部噪訊特性
5-1-1 CCF 非對稱性與北部噪訊來源
當噪訊並非均勻分佈在空間中時，會造成 CCF 正負兩邊訊號強度的差異，因 此，藉由分析正負訊號振幅的比值，可以找出各方向噪訊的相對強度。根據尤(2008) 研究北台灣垂直方向 CCF 的結果顯示，北部噪訊主要來自北部與西北部海岸。在 此，我們分別計算了 R-R、T-T 、 Z-Z、R-Z 與 Z-R 方向 CCF 之正負方向最大振 幅的比值（圖 5-1、圖 5-2），發現雖然不同方向的 CCF 所計算出來得結果略有不 同，但是都呈現兩個特性：
1. 雖然產生噪訊(1~10s)的機制依然不明，但是根據研究顯示，主要與海浪跟沿岸 地形之間的交互作用有關，因此能量會隨距離海岸增加而衰減。本研究結果顯 示，北方來得訊號較南方來得強，這可能是由於北部的測站離海岸較近，而南方 的測站距離海岸較遠所致，這部份的結果符合主要噪訊來源地為海岸區域的推 論。
2. 西北方的噪訊比東南方的噪訊來的強，這可能是由於這兩區地沿海地形不同所 導致，但進一步的確認仍有待結合更多的東部沿海測站資料（北緯 24.5 度以南）。
圖 5-1 左圖是每個測站相同方向之 CCF 非對稱性（由上而下分別是 R-R、Z-Z 與 T-T），線段長度表示正負訊號最大振幅之比值，線段越長表示訊號相對較強的方 向。右圖是把位於本研究區域中心位置的測站 NSK 當作接收訊號的測站，分析從 周圍測站來的 CCF 訊號強度，以更清楚顯示不同地理位置相對強度的差異。非對 稱性強弱以長度表示，箭頭表示訊號來向。
圖 5-2 左圖為各測站不同方向 CCF 的非對稱性（R-Z & Z-R），右圖為測站 NSK 對於周圍測站之非對稱圖，其強弱以長度表示，箭頭表示訊號來向。。
5-1-2 CCF 非對稱性之季節性變化
我們發現 CCF 的振幅比值會隨著時間有所不同，顯示噪訊的相對強度除了隨 著方位改變而有所不同之外，也會隨著時間不同而有所變化（圖 5-3 ~ 圖 5-5）。
從 3 個方向不同月份的 CCF 相對強度變化圖中，我們發現 Z-Z 方向的季節性 相對強度變化最為明顯：較強的訊號在 1~3 月時，以北偏西 30 度至北偏東 50 度 之間為主，在 7~12 月時，西方的訊號強度則有增強的趨勢；而 R-R 方向的 CCF 季節性相對強度變化雖然沒有 Z-Z 方向來的明顯，較強的訊號來自北偏西 45 度至 北偏東 50 度之間，但是在 7~9 月之間北偏西 30 度方向的訊號有明顯增強的現象；
至於 T-T 方向的 CCF 相對強度變化則最不明顯，較強的訊號始終來自於北偏東 45 度左右，但是在 7~9 月之間，北偏西 90 度至 45 度之間的訊號明顯增強，隨即又 消失不見。
為了進一步研究 CCF 相對強度季節性變化的特性，我們增加新竹微震觀測網 2005 年的資料進行分析，發現在 9 月至 11 月之間，南北方向的相對強度發生變化，
(圖 5-6)。然而，為更全面的探討此現象，仍有待結合更多的南方測站以及更長時 間的資料。
圖 5-3 各測站以及周圍測站對測站 NSK 季節性 CCF(Z-Z)之非對稱圖。
圖 5-5 各測站以及周圍測站對測站 NSK 季節性 CCF（T-T）之非對稱圖。
圖 5-6 2005～2006 年 HC16- HC05 之 Z-Z 與 R-R 的 CCF 相對強度變化。每 10 天疊加一次，每次移動 5 天，將訊號最大振幅定為 1，振幅大小以顏色表示。正 方向的訊號只出現在每年的 9 月至 11 月，這可能是相對強度季節性的變化所造成。
5-1-3 周遭噪訊背景能量變化
除了相對強度的研究，我們也分析單邊 CCF 的能量變化，以探討噪訊背景能 量流變化（Background Energy Flow, BEF. Stehly et al., 2006）。由於測站 EGS(龜山 島)的資料不到半年，而測站 TAP 的訊號也有問題，因此先將此兩測站剔除。接著，
剩餘的 CCF 每十天疊加一次，每次移動五天，並計算其波包之最大振幅。最後，
算出各方位角各個時期波包振幅與整年平均的差值，並將最大值定為 1（圖 5-7）；
換言之，比年平均為大時，其值大於 0，反之，比年平均值小時，其值小於 0。
結果顯示，不同方向的 CCF（R-R、Z-Z & T-T）其 BEF 的變化十分類似，這 很可能是噪訊在淺層地殼可充分散射，導致各分量的能量交換度高，因此各分量 的 CCF 在時間上之強度變化具有高度的相關性。
若將不同方位角間能量變化特性，配合台灣季風交替的時節，則可以將 BEF 依方位角分成兩區：1.在方位角-30 度至 180 度之間，一年之中以春、冬二季能量 最強，對應到東北季風；2.在方位角-30 度至-180 度之間，能量則以夏、秋二季最 強，相對應的則是西南季風。
結合以上觀察，雖然產生噪訊的機制不明，但是噪訊能量變化的時空特性卻 與台灣季節性氣候變化相互吻合，顯示噪訊強度變化及可能與大氣擾動有關。
5-1-4 CCF 強度與測站間距離之間係
除了噪訊的時空分佈，測站間的距離遠近也是影響訊號強度的原因之一。由 於 CCF 的能量主要反應在表面波上，而表面波的振幅又會隨著傳遞的距離增加而 減少（geometrical spreading），在距離極小時（sin(△)～△），理論上表面波的振幅 大小與距離開根號成反比。我們將最大振幅對測站間距離作圖（圖 5-8），發現雖 然振幅隨著距離的增加而減少，但並非理論所預期，與距離開根號成反比。這是 因為 CCF 的強度還同時會受到噪訊在空間中分佈的影響(Yang et al., 2008)，此外，
沿岸地形也可能是影響 CCF 強度的原因之ㄧ，再加上各方向的背景能量(BEF)並非 保持定值，因此我們量測的只是一個一整年平均的結果；換言之，除非我們能夠 扣除噪訊來源強度的影響，否則難以將所有 CCF 的強度放在相同的基準點上做比 較。
圖 5-8 由上而下分別為 Z-Z、R-R 與 T-T 之 CCF 最大振幅與測站距離的關係，藍 色表示 SNR 較大的一方，紅色則為 SNR 較小的一方。
5-2 單一測站連續紀錄分析
在 5-1 節中我們藉由分析 CCF 的振幅變化來了解噪訊源的時空特性。接 下來，我們藉由計算單一測站連續紀錄的譜密度機率函數以及 rms (Root Mean Square)隨時間的變化，探討單一測站連續紀錄之特性與 CCF 之間的關係。
5-2-1 單一測站譜密度之機率函數
為了探討噪訊能量在各頻率分佈的情形，我們先將垂直方向的速度型連續紀 錄轉換成加速度型的連續紀錄，並每一小時計算其頻譜圖(power spectra density, PSD, webb, 1998 ; McNamara et al., 2004; Tanimoto, 2005)，每次移動半小時，最後 統計各週期其強度（db,10 log (10 m2/s4)）之分佈機率（Probability density function, PDF, McNamara et al., 2004）。
我們比較六個靠岸的測站的 PDF，發現在北部與西部的測站其 PDF 在 1~5s 有明顯的突起，顯示此頻段的強度相對較強；反之，東部的兩個測站在此頻段內 並沒有明顯的起伏。結合沿岸地形的變化，發現這與 5-1 節分析噪訊相對強度的結 果相符：水深淺，將有利於碎浪形成，進而成為噪訊來源（圖 5-9）。
另外，我們也發現 PDF 異常曲線與其 CCF 之間的關係。比較 TAP 單日 PDF 的結果，發現其 PDF 分佈曲線的形狀可概分成兩類（圖 5-10(b)）；當 PDF 呈現
『Type-A』分佈時(圖 5-10(b) 十月份)，其 CCF 是雜亂而無清楚訊號（圖 5-10(c)）； 反之，當 PDF 呈現『Type-B』分佈時(圖 5-10(b) 九月份)，可以在 CCF 中看到清 楚的訊號。雖然這些雜亂的 CCF 經過疊加之後會互相抵銷，並不會影響分析訊號 群速度與相速度的結果，但是會壓抑訊號的強度，因此在分析 BEF 時並不考慮與 此測站相關的波徑，以免影響分析結果。
圖 5-9 沿海 6 測站全年之垂直方向連續紀錄之 PDF 以及其相對地理位置。在台灣 西部與北部的測站，其 PDF 在 1~5s 之間(綠色方框)訊號有明顯突起，表示此頻段 訊號相對其它頻對訊號能量較強較強；而東部測站則沒有此現象。此結果與 CCF 之非對稱性十分吻合。
圖 5-10 測站 TAP 垂直方向連續紀錄之 PDF 分析與 Z-Z 方向 CCF 的結果，綠色 方框為 1~5 秒的頻率範圍。(a)單日 PDF 的結果，日期標示在圖的上方。9 月單日 分析結果明顯跟 10 月的結果不同。(b)9 月與 10 月 PDF 之結果。10 月的 PDF 呈現 U 字型，在此訂為 Type-A，而九月的 PDF 則訂為 Type-B。(c) 當 PDF 為 Type-A 時， CCF 無訊號；當 PDF 為 Type-B 時， CCF 訊號良好。
（a）
（b）
（c）
5-2-2 單站連續紀錄 RMS 分析與 CCF
藉由分析單站連續紀錄之 rms 隨時間的變化，可以了解當地周遭噪訊的能量 變化，並探討影響噪訊大小的可能原因。為了降低儀器的突波以及地震影響分析 的結果，我們每 2 分鐘計算一次連續紀錄之 rms，每次移動一分鐘，因此每天有 1440 筆資料；也就是說最高解析頻率是 720 cpd (cycles per day)。我們發現儀器的 突波以及地震訊號依然會對分析結果造成影響（圖 5-11），因此，在分析前有必要 對 rms 的紀錄進一步處理。
為了移除異常的大訊號影響分析結果，我們計算每一點資料前 50 分鐘 rms 記 錄的平均，以該平均值的正負百分之五十當作門檻，若超過該範圍則將該點資料 的值定為平均值。
移除儀器突波的訊號之後，從測站 TWB1 之 rms 的變化中(圖 5-12)，可以清 楚的看到每天有 1～2 次的起伏，而且能量以南北向最大。為了進一步了解噪訊能 量之變化週期特性，我們挑選測站中連續紀錄裡突波訊號較少以及沒有突發的天 氣變化（颱風、鋒面，http://www.cwb.gov.tw/V6/index.htm）的時間，進行頻譜分 析。結果發現，在週期 1.8~2.2 cpd 之間有極大值，與分潮週期（表 5-1）比較發 現，周遭噪訊的能量變化與半日潮（Semi-diurnal tide）的週期接近，顯示海水面 的垂直運動會造成噪訊的強度變化(圖 5-13)。
為了檢驗潮汐的變化使否會反應在 CCF 的結果上，我們分析北部 7 個測站
（TWS1、TWY、NWF、TWB1、ILA、EGS 與 NCU，其 rms 結果均有明顯之半 日潮週期）與測站 TWA 每 4.5 小時計算一次的 CCF 的訊號波包振幅大小（圖 5-14）。結果顯示，雖然大多數的 CCF 訊號大小變化符合半日潮的週期，然而其變
性並不高，這與之前 BEF 的分析結果並不相符，確切的原因尚待未來更進一步的 研究。
而當颱風侵台期間，rms 的值明顯增大且半日潮的訊號明顯被蓋過（圖 5-16）， 一直到颱風過後，半日潮的現象才又漸明顯。而此現象也反應在 CCF 上，我們分 析七月份的 CCF 變化(圖 5-17)，發現當颱風來時，所有的 CCF 訊號除了強度明顯 增強之外，其頻譜的能量分佈也比平常（2-5s）來得低頻（4-8s），反應了天氣現 象確實可改變噪訊源的頻譜。
綜合以上結果，雖然造成噪訊的機制依然不明，但是結果顯示在穩定氣候下，
噪訊的強度與海水面的變化有關；而當颱風時期，颱風所造成的劇烈天氣氣象(強 風、降雨、大氣擾動…等等)反而成為噪訊能量的主要來源。
圖 5-11 測站 TWB1 rms 未經處理前之連續紀錄以及分析結果。(a) rms 連續紀 錄，圖中的峰值為儀器的突波造成。（b）將 rms 之振幅限制在 70 um，可以發 現其變化有明顯的週期運動。（c）rms 頻譜分析結果，並未看到明顯的潮汐週期，
是因為突波的訊號影響了分析的結果，因此有必要將該訊號去除。
（a）
（b） （c）
圖 5-12 測站 TWB1 去除儀器突波後，三方向連續紀錄 rms 分析結果。從 rms 的分析中發現，三方向之能量並不相同，南北方向之能量最大；而頻譜分析的 結果顯示，半日潮週期為 1.9cpd，而全日潮週期為 0.98cpd，分別與分潮中的N2 與P 接近。 1
表 5-1 分潮週期 ，Melchior (1983)
潮別 代號 分潮名稱 週期（cpd）
K2 日月合成半日潮 2.00546
S2 主太陽半日潮 2
M2 主太陰半日潮 1.93227
N2 主太陰橢圓潮 1.89598
半日週期
2N 2 二次主太陰橢圓潮 1.85969 K1 日月合成全日潮 1.00274 P1 主太陽全日潮 0.997262 O1 主太陰全日潮 0.929536
全日週期
Q 主太陰橢圓潮 0.893244
圖 5-13 測站 rms 紀錄與其頻譜能量分析圖。中間的地圖為本研究區域，而紅色 三角形則表示分析的測站。在 rms 圖中，橫軸為時間，縱軸為 rms 值，黑線表示 rms 之結果，紅線表示 rms 經過帶通濾波（0.5～3 cpd）之曲線；而在頻譜能量 圖中，擁有最大能量的頻率則標示在圖的右上方。所有的測站其 rms 震盪週期大 約 2 cpd 左右，十分吻合半日潮之週期；而測站 NSK 可能因為遠離海岸造成頻 譜能量的峰值不明顯。
圖 5-14 北部 7 測站與 TWA 之間的 CCF 波包振幅變化分析。中間的地圖為本研 究區域，而藍色三角形為測站 TWA 之位置，周圍的紅色三角形表示七個對 TWA 進行 CCF 分析的測站。在 CCF 振幅的變化圖中，細黑線表示 CCF 之振幅，紅 線則為 CCF 振幅變化曲線經過代通濾波(0.5~3 cpd)之後的曲線，藍線則為各測站 rms 變化曲線經過帶通濾波（0.5～3 cpd）之後的曲線。在頻譜能量圖中，1.8～
2.2 cpd 之間的最大值標示在圖的右上方。除了測站 EGS(龜山島)之外，其餘的 測站與 TWA 之間的 CCF 訊號強度變化皆與半日潮的週期相近。
圖 5-15 三方向 CCF 振幅變化分析。圖中紅線表示 R-R 的 CCF，黑線表示 Z-Z 之 CCF，而藍線則表示 T-T 之 CCF。從測站 NWF 與 NSK 方向來得噪訊造成三 方向 CCF 的變化相當一致（in phase），以及從測站 TWS1 來得噪訊則造成 T-T 與 Z-Z 方向的 CCF 的變化相似之外，其餘的 CCF 變化曲線並不一致。
圖 5-17 2006 年 7 月份測站 NST 與 HC23 之間 Z-Z 與 R-R 方向之 CCF 分析。
(a)將整個月份 CCF 最大振幅訂為 1，振幅大小以顏色表示。當颱風來時，CCF 的振幅明顯增強。(b)單日 CCF 之譜能量分析，橫軸為週期，正負值分別表示正 負方向之 CCF 分析結果。當颱風發生時，頻譜發生變動。(c) 2006 年七月份三個 颱風的路徑示意圖，藍色圓圈為颱風中心位置。颱風名稱、編號以及發生時間標 示在圖右方。(d)測站位置圖。
(a) (b)
(c)
(d)
第六章 結論
根據以上的分析結果，我們可以歸納出以下幾點結論以及未來研究方向：
1. 藉由分析交互相關函數正負方向訊號的異常現象，可對儀器的感應器極性與 記錄器的時間偏移進行修正，除確保 ANT 結果不受儀器問題影響外，也具 體改善該儀器的資料品質。
2. 經過頻散分析與合理的資料篩選之後，配合多重尺度有限參數法，分別得到 北台灣雷利波與洛夫波的群速度與相速度的二維速度構造。反演所得速度模 型與地表地質構造十分吻合。
3. 結合不同週期的二維速度構造的結果，配合不同的初始模型，進行三維速度 構造的反演，發現初始模型的選擇會影響我們的結果，因此未來將會結合表 面波敏感算核以及 CCF 的資訊進行逆推，反演台灣三維的淺層速度構造，並 進一步探討台灣淺層地殼彈性夠造的非均向性。
4. 藉由分析 CCF 正負訊號的非對稱性，我們發現訊號強度可能與測站與海岸間 的距離以及海岸地形有關；而從非對稱性隨時間的變化中，發現疑似噪訊相 對強度季節性的變化。此外，我們也藉著分析單邊 CCF 訊號的強弱變化，探 討不同方位的噪訊的能量變化；而不同方向的 CCF 皆呈現了極為相似的 BEF 結果，並且其特性與台灣季風交替的時節相符，顯示大氣的擾動可能是控制 整體噪訊強弱的主要原因。
5. 沿海測站的地理位置配合連續紀錄的 PDF 結果，發現在淺水區的測站，該 PDF 在 1～5s 之間相對強度較強，可能是因為該地理環境較有利於碎浪形
7. 當颱風接近台灣時，所有連續紀錄的能量皆大幅增加，而且 CCF 的結果除了 振幅增強之外，訊號也比平常（2～5s）來的低頻（4～8s），顯示颱風所造成 的天氣現象改變了噪訊源的頻譜。
8. 未來將運用 Beamforming （Brzak et al. 2009）的技術，搜尋噪訊的來源分佈，
並分析交互相關函數的非對稱性與測站之間方位角的關係、其隨時間不同所 產生的變化情形以及在頻率域的能量分布狀況；並結合海洋物理相關的研究 以釐清海水波動能量傳遞至陸地的效率與近岸海底地形之間的關係，進一步 了解噪訊的主要來源與產生的機制。
參考資料參考資料 參考資料參考資料:
Brzak, K., Y.J. Gu, A. Okeler, M. stickler , and A. lerner-Lam(2009), Migration imaging and forward modeling of microseismic noise sources near southern Italy, Geochem.
Geophy. Geosyst., 10,Q01012,doi:10.1029/2008GC002234
Brenguier, F. , ShapiroHHHHH, N. M.HHHHH, Campillo, M. , Ferrazzini, V. , Duputel, Z., Coutant , O.
and Nercessian, A. (2008) Towards forecasting volcanic eruptions using seismic noise , Nature Geoscience 1, 126 – 130 , doi:10.1038/ngeo104
Chiao, L.-Y., and W.T. Liang, (2003) Multiresolution parameterization for geophysical inverse problems, Geophysics, 68(1), 1-11.
Campillo, M., and A. Paul (2003), Long-Range Correlations in the Diffuse Seismic Coda,Science,299,547
Courtland, R. (2008) Earth science: Harnessing the hum, Nature, 453, 146-148 (2008) doi:10.1038/453146a
Herrman, R.B. (2004) Computer Programs in Seismology,
http://www.eas.slu.edu/People/RBHerrmann/ComputerPrograms.html Derode, A., E. Larose, M. Tanter, J. de Rosny, A. Tourin, M. Campillo, B. A. van
Tiggelen, and M. Fink (2003) Recovering the Greens function from field-field correlation in an open scattering medium (L), J. Acoust. Soc. Am.,113
Dahlen, F.A & Tromp, J., 1998. Theoretical Global Seismology, Princeton Univ. Press, Princeton, New Jersey.
McNamara, D.E. and R.P. Buland (2004) HHHHHAmbient Noise Levels in the Continental
