GPS 追蹤站藉由長期連續接收衛星資料，與參考框架中控制點進行解算，進 而獲得地面點位之空間資訊，解算得到之 GPS 序列資料可提供地表點位的動態連 續行為，藉此估計速度場，獲得點位運動資訊。然而經長期連續觀測地表點位所 得到之 GPS 序列資料，可能會因為自然環境 (如：斷層活動導致地殼變動) 或是 人為因素 (如：改變儀器天線位置、抽取地下水導致地層下陷)，造成序列資料行 為趨勢具有不連續性。一般的 GPS 資料處理軟體，在處理資料時，並無考慮資料 的不連續性以及觀測量可能含有粗差等問題，直接估計速度場，如此會造成估計 的速度場會與實際地表點位的動態行為有所差距，無法正確描述點位運動的趨 勢。本章針對長期 GPS 序列資料可能存在之資料不連續性、觀測量含有粗差等問 題，提出改善之方法，藉此提升速度場的品質。
4-1 GPS 資料解算
本研究以 IGS ( International GNSS Service )追蹤站作為坐標參考框架，依據 IGS 亞洲框架為基礎進行每日成果解算，接著再將這些每日解算成果與 IGS 全球 框架控制網的成果進行合併，獲得測站於 IGS 全球框架下之點位空間資訊。選用 IGS 控制網作為參考框架的主要原因在於 IGS 提供分佈全球的 GPS 測站資訊，包 含四百多個測站，並且每週進行成果解算，更新其測站資訊( Dow et al., 2005 )，相 較於 ITRF 需經過較長的時間提供新的解算成果，以 IGS 測站作為參考框架的控制 網，其測站資訊較能夠符合目前地表之動態行為。
以本研究所處理之台灣地區資料為例，實際資料處理中所使用到的 IGS 亞洲 區域測站包括：GUAM、TSKB、DAEJ、SUWN、YMSM、FLNM、TWTF、PIMO、
NTUS、KUNM、IISC、KIT3、BAHR、NSSP、ZECK、ZWEN、TELA、NICO、
ANKR、TUBI、SOFI、MATE，各站之地理位置如圖 4-1-1。
圖 4-1-1 IGS 亞洲測站位置分佈圖
實驗中台灣追蹤站的資料來源由內政部地政司衛星測量中心所提供，衛星追 蹤站包含陽明山(YMSM)、鳳林(FLNM)、太麻里(TMAM)、墾丁(KDNM)、北港 (PKGM) 、 金 門 (KMNM) 和 馬 祖 (MZUM) 七 站 。 資 料 類 型 為 RINEX (Receiver INdependent EXchange) 格式，資料時間從 2001 年第 1 日到 2007 年第 116 日，資 料接收的時間為每 30 秒接收一筆衛星資料。台灣追蹤站的地理位置見圖 4-1-2，其 中 YMSM 站、FLNM 站、TMAM 站、PKGM 站、KDNM 站位於台灣本島，KMNM 站與 MZUM 站屬於離島部分，分別位於台灣的西方與西北方。
圖 4-1-2 台灣追蹤站地理位置分佈圖
GPS 資 料 的 處 理 軟 體 採 用 美 國 麻 省 理 工 學 院 ( Massachusetts Institute of Technology, MIT )所開發的 GPS 資料分析處理軟體 GAMIT ( GPS Analysis at MIT ) 和 GLOBK ( Global Kalman filter VLBI and GPS analysis program ) ( Herring et al., 2006 )進行 GPS 觀測資料處理與解算的工作。首先使用 GAMIT 軟體解算衛星追蹤 站於 IGS ( International GNSS Service )亞洲區域框架下的每日成果(包括台灣地區 七個追蹤站和 IGS 亞洲區域追蹤站)，接著將解算的成果與在 IGS 全球框架下測站 進行合併，之後藉由 GLOBK 軟體估計出代表該期間之最適坐標與速度場，進而 獲得在 IGS 全球框架的動態行為模式下，台灣追蹤站的坐標和速度場。本研究所 使用之資料介紹整理如表 4-1-1。
表 4-1-1 GPS 資料介紹
資料來源 IGS 資料中心 & 內政部地政司衛星測量中心 資料時間 2001 年第 1 日到 2007 年第 116 日
資料接收時間間隔 30 秒
資料格式 RINEX
GPS 資料處理軟體
GAMIT & GLOBK ( MIT )參考框架 IGS 全球參考框架
IGS 追蹤站
GUAM、TSKB、DAEJ、SUWN、TWTF、PIMO、TCMS、TNML、BJFS、WUHN、NTUS、KUNM、
IISC、KIT3、BAHR、NSSP、ZECK、ZWEN、
TELA、NICO、ANKR、TUBI、SOFI、MATE，
台灣追蹤站 YMSM、FLNM、TMAM、KDNM、PKGM、
KMNM、MZUM 4-2 GPS 序列資料品質提升
一般從 GPS 資料處理軟體得到之 GPS 序列資料，會發現資料本身可能存有粗 差，或是在經過長期連續觀測的地表點位得到之點位，其點位資料的行為趨勢並 不連續。為了確保所獲得之 GPS 序列資料是具有可靠性，本研究將透過不連續運 動偵測、粗差偵測與週期性訊號處理等方式，提升序列資料的品質。序列資料的 處理方式包含不連續運動偵測、直線擬合、粗差偵測、週期性訊號處理等，處理 流程如圖 4-2-1。
GPS 序列資料
不連續運動偵測
不連續 連續資料 資料
直線擬合
直線擬合 斷點偵測
斷點前 資料 斷點後
資料
直線擬合 直線擬合
直線擬合
直線擬合 速度場估計
週期性 訊號處理
判斷斜率 是否相等
不同段 速度場估計
合併二段 資料
粗差偵測
粗差偵測 粗差偵測
不同連續 資料
相同連續 資料
圖 4-2-1 長期連續 GPS 序列資料處理流程圖
4-2-1 不連續運動偵測
得到斜率資料的平均值與標準差後，即可根據所設定之門檻，決定斷點候選 點，找出斜率變化劇烈的資料所在位置。在此設定斷點的門檻為當所計算的斜率 數值超過三倍的標準差範圍以上，即視為可能之斷點，將這些斷點作為候選點，
並以人工輔助，判定斷點位置，將資料分段處理。
4-2-2 速度場估計
本研究對序列資料估計速度場，採用的直線擬合公式( Wolf & Ghilani, 1997, pp.
187-190 )為
b mx
y
= + (4.2.4)其中 y 為坐標值，
x
為觀測資料時間，m
為直線斜率，b
為截距。平差方式採用間接平差法( Adjustment of Indirect Observations ) ( Mikhail &
Ackermann, 1976, pp. 159-161 )。間接平差法之數學模式如下：
( ) l , x = 0
F
(4.2.5)其中
l
代表觀測量，x
代表參數。上式經過線性化後，可得
d B v
l + + Δ =
(4.2.6)其中
v
代表殘差，B 代表參數x
的偏微分係數矩陣，Δ 代表參數x
的改正數，d
代表常數向量。上式經過整理後可寫成
( l d ) f
B
v + Δ = − + =
(4.2.7)根據最小二乘平差原理，可得
( B t WB )
Δ=B t Wf
(4.2.8)其中
W
為觀測量的權矩陣。上式可簡化為
t
N Δ =
(4.2.9)其中
N
=B t WB
，t
=B t Wf
。 解算參數改正數Δt N − 1
=
Δ (4.2.10)
殘差計算為
Δ
−
=
f B
v
(4.2.11)平差後的觀測量 lˆ
v l
l ˆ = +
(4.2.12)後驗單位權中誤差為
r Wv v t
=
2
ˆ
0
σ
(4.2.13)其中 r 代表自由度。
4-2-3 粗差偵測
由於一般 GPS 資料處理軟體，並沒有考慮資料本身是否含有粗差( blunders ) 的問題，可能會造成所估計之速度場品質不可靠。因此在本研究中，透過粗差偵 測( Blunders Detection )，剔除序列資料中含有粗差的觀測量，藉此提升資料品質。
在此所設定之粗差門檻有二項，只要符合二項門檻其中之一，即視為粗差。
第一種是計算序列資料的殘差( residuals )後驗方差，經過直線擬合步驟，根據 平差後得到之殘差後驗方差作為判定粗差的門檻，當資料的殘差大於平差後之殘 差，視為粗差。殘差的餘因子矩陣計算如下：
l
vv Q Q l
Q
= −ˆ ˆ
(4.2.14)其中 Q 代表觀測量的餘因子矩陣，
Q l ˆ l ˆ
的計算式如(4.2.15)式t l
l BN B
Q ˆ ˆ = − 1
(4.2.15)其中 B 代表參數的偏徵分係數矩陣，
N 代表法方程式矩陣。 − 1
第二種是計算殘差資料中的平均值和標準差，在此設定當殘差數值位於三倍 標準差以外，則視為粗差。將含有粗差的觀測量剔除後，再重新擬合直線，估計 速度場。平均值與標準差的計算方式，可參考(4.2.2)式與(4.2.3)式。
4-2-4 週期訊號處理
地表變動的原因除了由於地球內部能量所引發的力學過程造成板塊運動之 外，地球亦會受到宇宙環境中各種因素的影響，例如：天文潮汐，即太陽、月球 對地球的引潮力。根據地球物理與地質的相關研究( Tanaka et al., 2006; Mentes, 2008 )表示，地震與潮汐之間的關係確實存在。因此藉由長期連續觀測地表點位，
其點位之動態資訊應含有天文潮汐因素對地表所造成之週期性影響。
在考慮觀測量中可能含有週期性訊號的情況下，對剔除粗差後的觀測量殘 差，視情況進行週期性訊號處理。對於 GPS 解算得到之時間序列資料，Nikolaidis (2002)的研究中以多個正弦波組成之多項式函數進行時間序列的資料擬合。在此採 用正弦波( sine wave )週期函數進行觀測量殘差的週期性訊號擬合，並將擬合得到 之週期訊號扣除，重新估計速度場。
所使用的正弦波函式為
( )
( x o p )
A
y = × sin 2
π+ /
(4.2.16)其中 A 為正弦波的振幅， p 為一個正弦波的週期，
o
代訊號的位移，x
為資料 時間。4-2-5 統計檢定
當
F
>F α / 2 , v 1 , v 2
，即拒絕虛無假設H 。 0
df
；若二個樣本的標準偏差顯著不同，則使用(4.2.22)式四捨五入取整 數計算自由度：若經檢定後，發現二條直線的斜率不相同，則判定為不連續資料；若經檢定
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.2
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.1
KDNM Y Residuals # 1790
圖 4-3-2 原始資料殘差
0 200 400 600 800 1000 1200 1400 1600 1800
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.2
KDNM Y Offset 5084820.6647 m Jumper Candidate # 38
Series data Jumper Candidate
圖 4-3-4 斷點偵測
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.1
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.1
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.1
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.2
KDNM Y Offset 5084820.6647 m # 1790
圖 4-3-8 分段估計速度場
2000 2001 2002 2003 2004 2005 2006 2007 2008
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.1
KDNM Residuals 1 Removing Blunders Y # 361
圖 4-3-10 剔除第一段資料粗差
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.1
KDNM Residuals 2 Removing Blunders Y # 489
圖 4-3-11 剔除第二段資料粗差
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.1
KDNM Residuals 3 Removing Blunders Y # 392
圖 4-3-12 剔除第三段資料粗差
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.2
KDNM Y Offset 5084820.6647 m Removing blunders # 1242
圖 4-3-13 分析處理後三段資料速度場 估計
2000 2001 2002 2003 2004 2005 2006 2007 2008
-0.1
KDNM Y Residuals Removing blunders # 1242
圖 4-3-14 分析處理後三段資料殘差
當序列資料經判定為不連續資料，以分段處理的方式，對各段進行粗差偵測 及速度場估計後，須決定分段後之資料所估計之速度場是否相同，在此根據估計
得到之斜率與斜率精度，依照 4-2-5 節所述之統計檢定方式，作為處理後之不同段
台灣地區七個追蹤站在 XYZ 三個方向的 GPS 序列資料分析處理成果，詳見附 錄 A.1；七個追蹤站在 XYZ 三個方向之原始資料時間序列，請見附錄 A.2；七個 追蹤站經分析處理後的時間序列，請見附錄 A.3。
表 4-3-2 台灣追蹤站 2001 至 2007 年 GPS 序列資料分析處理前後精度比較表 原始資料之速度場精度(mm/yr) 資料分析處理後速度場精度(mm/yr) Stations
X 方向 Y 方向 Z 方向 平均精度 X 方向 Y 方向 Z 方向 平均精度 FLNM ± 0.19 ± 0.31 ± 0.24 ± 0.25 ± 0.12 ± 0.35 ± 0.29 ± 0.25 YMSM ± 0.24 ± 0.49 ± 0.29 ± 0.34 ± 0.22 ± 0.43 ± 0.28 ± 0.31 KDNM ± 0.32 ± 0.54 ± 0.32 ± 0.39 ± 0.22 ± 0.53 ± 0.30 ± 0.35 KMNM ± 0.31 ± 0.72 ± 0.47 ± 0.50 ± 0.20 ± 0.32 ± 0.41 ± 0.31 MZUM ± 0.28 ± 0.62 ± 0.48 ± 0.46 ± 0.14 ± 0.32 ± 0.26 ± 0.24 PKGM ± 0.16 ± 0.33 ± 0.25 ± 0.25 ± 0.12 ± 0.27 ± 0.30 ± 0.23 TMAM ± 0.23 ± 0.27 ± 0.26 ± 0.25 ± 0.15 ± 0.27 ± 0.28 ± 0.23 表 4-3-2 為台灣追蹤站 2001 年至 2007 年 GPS 序列資料以原始資料所估計得 到之速度場精度，以及資料經過分析處理後之速度場精度成果。由數值成果顯示，
經過處理後之資料，在速度場估計方面，品質大部分都有提升，僅有少數幾個方 向的速度場精度未提升。如：FLNM 站的 Y 與 Z 方向，PKGM 站的 Z 方向及 TMAM 站的 Z 方向，其中精度未提升之原因可能在於原始資料數量過多，而資料分析後 因為採取分段處理的關係，造成該期間之資料量與原始資料差異太大，導致二者 的自由度有所差距。因為在進行直線擬合時，當觀測量數量愈多，所擬合出來的 直線斜率精度就愈高。由擬合直線的斜率與斜率精度的曲線圖(圖 4-3-15)即可得知 這樣的趨勢。圖中藍色線代表斜率，紅色線代表斜率精度，從圖上顯示隨著觀測 量個數增加，直線斜率愈趨穩定，斜率精度的曲線亦逐漸降低(意即精度愈高)。因 此當觀測量個數有所差異時，會以擁有較大自由度者獲得較高之精度。
從七個追蹤站之平均速度場精度顯示，速度場的精度範圍從原本的± 0.3 mm/yr~± 0.5 mm/yr 間，經過處理之後，精度範圍提升在 ± 0.2 mm/yr ~ ± 0.4 mm/yr 間，顯示 GPS 序列資料經分析處理後，能夠有效提升速度場的品質。
0 50 100 150 -12
-10 -8 -6 -4 -2 0 2 4 6
Observations
(m /y r)
FLNM Up Slope
slope std slope
圖 4-3-15 直線擬合斜率及其精度曲線圖
一般 GPS 資料處理軟體，直接根據解算得到之序列資料估計速度場，並沒有
一般 GPS 資料處理軟體，直接根據解算得到之序列資料估計速度場，並沒有