東海海域海底地形研究：以衛星測高所反演的重力場推估

國 立 交 通 大 學

土 木 工 程 學 系

碩 士 論 文

東海海域海底地形研究：以衛星測高所反演的重力場推估

An Investigation of Bathymetric Model over the East China Sea:

Prediction by Altimetry-Derived Gravity Field

研 究 生：鄭 詠 升

指導教授：黃 金 維

東海海域海底地形研究：以衛星測高所反演的重力場推估

研究生：鄭 詠 升 指導教授：黃 金 維 國立交通大學土木工程學系 摘 要 本研究結合數種衛星測高資料反演東海重力場，進而推估研究區域之海深模 型。在波形重定方面，採用次波形門檻值演算法並搭配高斯濾波去除雜訊；在反 演重力場方面，分別比較了逆 Vening Meinesz 公式 與最小二乘配置法，並分析 其精度。在計算海底地形方面，採用重力地質法計算出最終的海深模型，並分析 參考場與地質密度差對海深計算的影響。重力成果方面，與船測重力異常檢核後 ，發現波形重定技術能有效提升約 3~6 mgal 的精度；於海深成果方面，用大陸 礁層船測海深值評估精度，可發現本研究計算之模型較 ETOPO1 模型獲得更高 之精度，於淺海區域與深海區域精度約比百分之四。 關鍵字： 測高；東海；波形重定；逆 Vening Meinesz 公式；最小二乘配置法； 重力地質法

An Investigation of Bathymetric Model over the East China Sea:

Prediction by Altimetry-Derived Gravity Field

Student：Yung Sheng Cheng Advisor：Dr. Cheinway Hwang

Department of Civil Engineering National Chiao Tung University

Abstract

We investigated the bathymetric model over the East China Sea by altimetry-derived gravity. A sub-waveform threshold is used and a Gaussian filter is

applied in the retracking processing. Two methods, the Inverse Vening Meinesz

(IVM) formula and least square collocation (LSC), are used to derive gravity

anomalies from retracked altimeter data. The bathymetric models are constructed by

gravity-geological method (GGM) using the altimetry-derived gravity. Different

reference models and density contrasts for GGM are taken into account. The use of

retracking technique improves the gravity estimation by 3-6 mgal, based on the

evaluation with shipborne gravity. The accuracy of the predicted bathymetric model

in this study outperforms ETOPO1 DEM model. In shallow sea and abyssal region,

the accuracy of the predicted bathymetric model in this study is improved by 4 %.

Keywords: Altimeter, East China Sea, Retracking, Inverse Vening-Meneisz formula, least-squares collocation, Gravity-Geological Method.

誌謝

轉眼間，已在交大度過了兩年研究生涯，從原本剛大學畢業，從家鄉台南到 新竹來求學，一眨眼，已經通過碩士班口試，最後完成了一本碩士論文。 回想碩士生涯這段時間，除了專業領域的精進，也學習到很多的人生經驗與 生活的經驗，相信這對於我的未來是有很大的助益，而在這段求學過程，要感謝 的人真的太多了。首先，要感謝黃金維教授這兩年來悉心的指導，讓我學習到很 多知識，也使本論文能順利完成。感謝我的口試委員提供諸多寶貴意見和指正， 使本論文更臻完備。感謝博士班學長姐蕭董、蟳蟳哥、亨利、小悶、乃特、宜珊、 偉尊、彥杕、雅琪等在我有問題不懂時候，願意費心指導我與幫助我解決難題。 感謝碩士班學長姐元旎、登登、信瑜姐的指導。感謝研究所的同學貓兄、宅少爺、 蓉哥與其它碩士班同學和學弟妹在課業上與生活上的幫助。 最後，感謝我家人，爸媽和姐，使我在求學過程中讓我無後顧之憂。

目 錄

中文摘要……… i 英文摘要……… ii 致謝……… iii 目錄……… iv 圖目錄……… vi 表目錄……… viii 第一章 緒論... 1 1-1 研究動機 ... 1 1-2 文獻回顧... 4 1-3 論文架構 ... 6 第二章 測高衛星資料介紹與處理 ... 7 2-1 測高衛星介紹 ... 7 2-2 測高衛星原理 ... 11 2-3 波形重定 ... 14 2-4 次波形門檻值演算法 ... 15 第三章 船測資料介紹與處理 ... 19

3-1 船測資料介紹 ... 19

3-2 系統誤差改正與粗差剔除 ... 22

第四章 重力場與海深 ... 29

4-1 計算海深方法 ... 29

4-2 Inverse Vening Meinesz 公式(IVM) ... 34

4-3 最小二乘配置法(LSC) ... 36 4-4 高斯濾波 ... 37 4-5 重力地質法(GGM)... 38 4-6 向下延續法(DWC) ... 40 第五章 重力場計算結果與討論 ... 42 第六章 海深測量計算結果與討論 ... 50 第七章 結論與建議 ... 61 7-1 結論 ... 61 7-2 建議 ... 62 參考文獻... 63

圖 目 錄

圖 1-1 東海海域海底地形圖 ... 3 圖 2-1 各測高衛星在研究區域之重複軌道軌跡分佈圖 ... 10 圖 2-2 衛星測高幾何圖 ... 13 圖 2-3 波形概要圖 ... 14 圖 3-1 本研究所蒐集之我國大陸礁層船測海深資料軌跡圖 ... 20 圖 3-2-a 本研究所蒐集之 NGDC 船測海深資料軌跡圖 ... 21 圖 3-2-b 本研究所蒐集之 NGDC 船測重力資料軌跡圖... 21 圖 3-3 釣魚台列嶼周圍之船測資料 ... 22 圖 3-4 船測編號 hm951026 之軌跡圖 ... 27 圖 3-5 船測編號 hm951026 之重力比較圖 ... 28 圖 3-6 編號 hm951026 航次(a)平差前與(b)平差後重力與 EGM2008 重力 差之直方圖... 28 圖 4-1 研究流程圖 ... 31 圖 4-2 IVM 計算流程圖 ... 32 圖 4-3 LSC 計算流程圖 ... 33

圖 4-4 Inverse Vening Meinesz 公式的幾何變化圖 ... 35

圖 4-5 重力地質法的幾何架構圖 ... 39 圖 4-6 利用向下延續法，所計算出的海深各層的地質密度差圖 ... 41 圖 5-1 檢核重力網格精度的 NGDC 船測重力軌跡圖 ... 47 圖 5-2 Case 4 減去 Case 2 之差異圖 ... 47 圖 5-3 Case 3 減去 Case 4 之差異圖 ... 48 圖 5-4 Case 3 減去 Case 5 之差異圖 ... 48 圖 5-5 Case 3 所計算的重力網格圖 ... 49

圖 6-2-a Case A 模型與大陸礁層海深值之差異圖 ... 56

圖 6-2-b Case B 模型與大陸礁層海深值之差異圖 ... 57

圖 6-2-c Case C 模型與大陸礁層海深值之差異圖 ... 58

圖 6-2-d Case D 模型與大陸礁層海深值之差異圖 ... 59

表 目 錄

表 2-1 各測高衛星主要相關資訊 ... 11 表 3-1 利用二次多項式平差船測重力之結果 ... 24 表 5-1 本研究之六種不同重力模型 ... 45 表 5-2 各種重力網格與 NGDC 船測重力的比較差值 ... 45 表 6-1 各種海深網格計算方式 ... 53 表 6-2 各種海深網格與大陸礁層船測海深值的比較差值 ... 54 表 6-3 各種海深網格與 NGDC 船測海深值的比較差值 ... 55 表 6-4 重力模型影響海深比較分類 ... 55 表 6-5 不同重力模型影響海深比較統計 ... 55

第一章 緒論

1-1 研究動機

東海(East China Sea)海域(圖 1-1)遍佈著許多小島，除了琉球群島之外，極具

爭議的釣魚台群島也位於此海域(圖 1-1 之星號)，而東海也位於全世界地震發生 頻率最頻繁的東亞地震帶上，除此之外，此海域東側亦為黑潮流經之地，且黑潮 對台灣的氣候、漁業等影響極大。因此東海海域為一具有豐富實驗性之區域。 海深測量對於很多領域來說，具有很重要的意義。對大地測量領域來說，海 深測量可以精確求定地球的形狀與建立更精密的重力場；對地球科學領域來說， 海深測量可幫助了解地球內部物質與板塊移動的情形；對海洋學領域來說，海深 測量有助於研究洋流、潮流與海底地形的構造；此外，海深測量也可應用於幫助 國家邊界的界定。目前全世界海深測量發展至今，大致可分為三種方式來測量海 深，(一)船測聲納測量：也是目前求定海深精度最高的測量方法，其技術隨著科 技發展，目前主要有單音束測深系統與多音束測深系統，其中以多音束測深系統 之精度較高，但船測測量對於大範圍的區域來說，此方法有其限制。例如：測量 費用較高，測量時間過長，甚至還有受到主權問題或是地形之限制等。(二)利用 衛星影像來計算海深：此法雖然可以計算海深，但受限於光學的反射特性，只能 計算水深約 20~50 公尺的區域。(三)為利用重力測量推算海深：目前有利用船測 重力、衛星重力或是衛星測高重力來求算海深，於船測重力方面，其單點精度不 錯，但是若測大範圍區域，其經費與時間都很可觀，且其船測軌跡分佈不均；而

衛星重力方面，多數應用於求定長波長重力場，故較不能應用於海深研究；在衛 星測高重力方面，其觀測量為海水面高度，實際應用時必須先轉化為重力值，其 技術發展今已有 40 多年歷史，針對觀測量進行改正方法也日漸成熟，且測高衛 星能提供大量的高精度資料，因此本文將以衛星測高資料求算重力並結合船測重 力來反演東海海域之海深。 由於研究區域東海海域範圍廣(經度122°E−132°E、緯度22°N−32°N)，且 散布許多小島，而衛星測高在接近陸地，其足跡(footprint)會因為碰到陸地進而 破壞反射之脈衝波形，導致衛星測距之精度不佳，因此本文將利用波形重定 (waveform retracking)求得近岸衛星測高改正之海水面高，進而用於計算研究區域 之海深，也利用測高衛星重力所計算之重力異常來評估波形重定對於精度上的貢 獻。除此之外，因為本研究區域之釣魚台列嶼，受其主權問題，所以其周遭沒有 船測資料，因此，本研究也將以衛星測高重力反演之海深來研究釣魚台列嶼周遭 海底地形。最後，將本文計算之各種海深模型與 ETOPO1 海深模型進行比較， 期許成果能優於 ETOPO1 海深模型。

圖 1-1 東海海域海底地形圖， 黑色框線區域為研究區域範圍，

紅色星星為釣魚台列嶼所在。

1-2 文獻回顧

衛星測高技術最早是為了監測海水面變化而設計，自從 Brown (1977)提出了

平均脈衝回傳理論，也帶動了測高衛星反射之波形的研究與發展，如 Martin et al.

(1983)發展之β 參數演算法、Wingham et al. (1986)發展之 OCOG(off center of − gravity)演算法、Davis (1997)所發展之門檻值(threshold)演算法、Hwang (2006)提

出了改進門檻值(improved threshold)演算法與 Yang et al. (2012)提出的次波形

(Subwaveform)門檻值演算法。除了波形重定的理論發展之外，國內外學者也將

此技術應用於各方面，如陳松安 (2001)使用β −5參數演算法於 ERS-1 測高衛星

上，Deng (2003)對於 ERS-1 及 ERS-2 測高衛星資料進行波形重定，並且於台灣

沿海和澳洲沿岸船測重力資料進行比較，文獻中提到進行波形重定較未進行波形

精度提升了約 4 mgal。

在 70 年代衛星搭載測高儀執行大地任務以後，海洋重力場的計算變成為研

究的主題，最早期是於 Rummel and Rapp (1977)提出相關研究文章。Moritz (1980)

提到最小二乘配置法(Least Squares Collocation, LSC)。Rapp (1985)結合了 Geoos-3

及 Seasat 測高衛星資料進行全球海洋重力異常與大地起伏之計算。Hwang (1989)

亦使用 Geoos-3 及 Seasat 測高衛星資料做高精度的重力異常預估。Hwang and

Parsons (1995)結合 Seasat、Geosat/ERM、ERS1/ERM、TOPEX/POSEIDON 四種

測高衛星資料和船測資料求出雷克雅内斯海脊 (Reykjanes Ridge)的詳細重力。

量轉換為重力異常，於文獻中提到 IVM 計算殘餘重力異常於當時為一新技術，

與船測重力進行比較，精度約可達 5~14 mgal。

自從 1978 年 Seasat 衛星發射以來，利用衛星測高資料逆推海深研究就以展

開，再這數十年間發展中，根據文獻均有不錯的成果，例如 Dixon et al. (1983)、

Smith and Sandwell (1994)皆提出利用測高資料來求得海深的方法；Hwang (1998)

將測高衛星資料求得之重力異常以快速傅立葉變換(Fast Fourier Transfer, FFT)及

其逆變換來求得海深模型，最後將所計算的模型與其他已知的全球海深模型做比

較，並利用船測海深的資料進行精度評估；Hsiao et al. (2011)將測高資料以重力

1-3 論文架構 本文共分七章節進行討論，各章節內容介紹如下： 第一章：主要分為研究動機、文獻回顧與本論文章節架構之說明。 第二章：主要分為兩大主題，主題一為介紹本研究所使用之測高衛星資料，並說 明其原理。主題二為簡單介紹波形重定技術並詳述本研究所採用之波形重定技術 之理論。 第三章：介紹本研究所使用之船測資料，並說明資料處理之過程。 第四章：介紹本研究所採用之理論方法，將詳述如何從衛星測高資料求得重力異 常，並結合船測資料或是其他網格模型資料來求解本研究區海深模型之理論與方 法。 第五章：分析本研究區求得之重力場成果，並探討波形重定有無對測區之影響。 第六章：分析本研究區計算之海深成果，並比較全球知名之海深模型，並將本研 究之海深成果進行釣魚台群島之討論。 第七章：對本論文做結論與建議。

第二章 測高衛星資料介紹與處理

2-1 測高衛星介紹

本研究所收集到的衛星測高資料包含 Geosat/GM，ERS-1/GM 衛星之大地任

務衛星( Geodetic Mission, GM)與重複軌道任務衛星(Exact Repeat Mission, ERM)

Geosat/ERM、ERS-1/ERM、Jason1、ERS-2/ERM、TOPEX/POSEIDON 及 ENVISAT 等測高衛星資料。大地任務衛星資料均無重複軌跡，只有個別週期資料，因此大 地任務衛星資料全球覆蓋率高，故可精確求得海洋重力場與大地水準面。重複軌 道任務衛星資料由於在重複軌道上，有週期性的重複觀測量，故可求得較高之海 水面高精度，但覆蓋率並不高。以下將分別介紹對本研究所收集各衛星。圖 2-1 為各測高衛星於本研究的地面軌跡圖；表 2-1 為本研究收集知各衛星主要相關資 訊分別。由於大地任務資料分佈於研究區覆蓋率較高，因此本研究主要以大地任 務資料為主要觀測量，而對此 GM 資料觀測量進行波形重定改正。 (1) ERS-1/ERS-2 衛星 歐洲太空總署於 1991 年 7 月 17 日發射了歐洲第一顆遙感衛星 ERS-1，是一個太 陽同步衛星，高度約 785 公里。ERS-1 的飛行任務亦分為 GM 與 ERM 兩種，且 任務結束於 2000 年 3 月。其中 ERM 任務可以分為 3 天與 35 天兩種重複軌跡任 務，於本研究中採用的 ERS-1 和-2 衛星資料皆為 35 天的重複軌跡資料。ERS-2 衛星則於 1995 年 4 月 21 日發射，與 ERS-1 的軌道相同，也是太陽同步衛星。 ERS-1 和 ERS-2 雖然為遙感探測衛星，但均載有 Ku 波段（13.8GHz）的雷達測

高儀，而 ERS-1 和 ERS-2 衛星對於海洋相關活動有顯著貢獻，例如:冰河監測、

海洋表面地形和海洋環流等等。

(2) Geosat 衛星

在 1985 年 3 月 12 日美國海軍發射了 Geosat 衛星，衛星高度為 800 公里。其

任務分為 GM 和 ERM 兩種。Seeber (1993)提到 Geosat 測距儀觀測精度為±3.5公

分，主要的大地任務在於獲取高密度的海洋數據，以改進現有地球重力場與海洋 大地水準面的精度。飛行其間共獲得全球海洋上 2.7 億個觀測數據，飛行軌跡總 長度為 2 億公里長，地面軌跡之平均間距為 4 公里。大地任務於 1986 年 9 月 30 日結束後，衛星調整至固定軌道，開始執行重複軌道任務，一直到 1990 年 1 月 因能源衰退而終止其功能，Geosat/ERM 的重複軌跡週期為 17 天。Geosat/GM 加 上 Geosat/ERM 的數據成果，是當時最長的連續全球海面高度觀測資料，有十分 重要的科學價值。 (3) TOPEX/POSEIDON 衛星 TOPEX/POSEIDON（簡稱 T/P）衛星為美國太空總署（NASA）與法國太空 中心（CENS）共同合作發射。於 1992 年 8 月 11 日發射，T/P 衛星及其測高系 統是當今測高技術發展上精度水準最高的產品，觀測量（海水面高）精度約為 3.3 公分。早期測高衛星軌道高度均為 800 公里左右，軌道傾角均在 100 度左右， 而 T/P 衛星軌道高度達 1336 公里，軌道傾角為 66 度，重複週期為 10 天。而低 傾角的目的主要是為了適用全球主要大洋的在地球上分佈而設計的。T/P 衛星於

2006 年 1 月結束其測高任務，而當 Jason-1 衛星於 2001 年底發射之後進入原本 T/P 衛星軌道，T/P 衛星則平移到兩 Jason-1 衛星軌道之間，以增加測高衛星之空 間解析度。 (4) Jason-1 衛星 Jason-1 衛星為美國國家太空總署為延續 TOPEX/POSEIDON 衛星之地球 觀測任務，於 2001 年 12 月所發射之衛星，目前仍持續蒐集資料中。T/P 衛星海 水面精度已達 3.3 cm，而 Jason-1 衛星海水面精度更達到 2.5 cm，因此其精度相

較於 T/P 衛星改善許多，但 Frappart et al.（2006）提到 Jason-1 衛星在陸地上的

資料精度較差，其衛星軌道高度 1336 公里，軌道傾角為 66 度，重複週期為 10

天，衛星任務接續著 T/P 衛星的任務，持續量測地球之海平面高度。於 2012 年

4 月 12 日將其軌道調整至大地軌道，開始執行大地任務。

(5) Envisat 衛星

ENVISAT 衛星為歐洲太空總署 (ESA) 為延續 ERS 之地球觀測任務，於

2002 年 3 月所發射之衛星。ENVISAT 為一太陽同步衛星，飛行高度約 800 公里，

重覆週期為 35 天，搭載了十種感測器，其中包含主動式微波感測器─先進合成

孔徑雷達(SAR)系統，簡稱為 ASAR，飛行任務方面均與 ERS 任務相同，應用於

海洋相關活動。由於 ENVISAT 衛星於 2012 年 4 月 8 日發生訊號中斷，歷經一

個月的聯繫未果，因此 ESA 於 2012 年 5 月 9 日正式宣布終止 ENVISAT 衛星的

(a)

(b) (c)

圖 2-1 各測高衛星在研究區域之重複軌道軌跡分佈圖。(a)ERM 任務，其中紅色 為 Geosat/ERM 軌跡；黑色為 ERS-1/ERM、ERS-2 與 ENVISAT 軌跡；藍色為 TOPEX/POSEIDON 和 Jason-1 軌跡。 (b)Geosat/GM 任務。 (c)ERS-1/GM 任務。

表 2-1 各測高衛星主要相關資訊 衛星 研制單位 衛星 傾角(度) 重複 週期(天) 跨赤道 間距(Km) 波形重定 有無 Geosat/GM U.S.Navy 108 不重複 約 2 有 Geosat/ERM U.S.Navy 108 17 160 無 ERS-1/GM ESA 98.5 不重複 約 4 有 ERS-1/ERM ESA 98.5 35 80 無 ERS-2 ESA 98.5 35 80 無 T/P NASA/CNES 66 10 320 無 Jason-1 NASA 66 10 320 無 Envisat ESA 98.5 35 80 無 2-2 測高衛星原理 衛星測高的基本原理為經由衛星搭載的雷達測高儀(radar altimeter)並向地球 表面發射一脈衝(pulse)，隨著地球表面反射脈衝，然後由衛星天線接收此反射信 號，並計算衛星發射脈衝至接收訊號所經過的時間進而求得衛星至海水面的距 離。如公式 2-1: H₌c×∆t/2 (2-1) 其中 H 是衛星至地球表面的距離，c 是光速， t∆ 為雷達波所經過的時間。

測高衛星的基本觀測量為海水面高(sea surface height, SSH)，此乃由 2-1 式求 出之 H 值與衛星高度值 h 在海水面上之高度差，通常 h 為衛星運動參考體計算 同時刻衛星與參考橢球面間之高度。如公式 2-2: SSH =h−H (2-2) 其中衛星測高 h 之詳細幾何關係如圖 2-2 所示，可以用下列公式來表達: h =N +R+∆H +H+d (2-3) N 為大地起伏(geoidal height)

R 為海水面地形(sea surface topography, SST)

H ∆ 為瞬間效應，例如:潮汐效應、對流層效應、電離層效應等等 H為衛星測高儀觀測量 D 為徑向軌道誤差 海水面地形為大地水準面與平均海水面之差異量，其大小可以達 1 至 2 公 尺，通常平均海水面是指海水面沒有因時間而產生變化，因此也稱為固定海水 面。而大地起伏則是由橢球面至大地水準面之距離。

圖 2-2 衛星測高幾何圖 H ∆ h d N R H 參考橢求面 大地水準面 平均 海水面 瞬時 海水面 真軌道 計算軌道

2-3 波形重定

測高波形是由測高儀天線發射微波脈衝，經海水面反射脈衝能 量後，由測高

儀上裝置測高儀的自動增益系統(automatic gain control, AGC)，以不連續方式接

收並以值域 0-255 儲存，將接收到之值 連線則形成波形。波形(圖 2-3)所示，分為 三部分：熱雜訊(Thermal Noise)部份、上升區部份與下降區部份，上升區通常稱 之為前緣(Leading Edge)。圖 2-3 中上升區之一半，通常稱之為前緣中點，通常為 測高衛星距離觀測為發射時間至回傳脈衝最大振幅之一半，並定該時間點為預設 閥門(Tracking Gate)，因此時間間隔乘上閥門個數與光速即為距離：

2

Ga

c

R

=

×

∆

∆

(2-4) 其中， Ga∆ 為一個閥門時間間隔，ΔR 為一個閥門代表之距離。 圖 2-3 波形概要圖 紅色線為熱雜訊區，藍色線為上升區，黑色線為下降區 熱雜訊 上升區 下降區

於深海中，海水面可以近似為一平面，所以經海水面反射而回傳之脈衝波形 非常良好，故測距精度較高，但是當衛星接近陸地約 20 公里時，其足跡會由於 接觸到陸地或是較為複雜之海水面而導致反射波形被污染，使得測高衛星的預設 閥門不在前緣中點，進而影響衛星測距之 精度，所以需藉由波形重定求得近岸衛 星測高改正後之海水面高。常用的波形重定方法：OCOG 演算法、β 參數演算法、

門檻值演算法(Threshold Algorithm)、改良型門檻值演算法(Improved Threshold

Algorithm)與次波形門檻值演算法(Subwaveform)等。其中次波形門檻值演算法與 改良型門檻值演算法一樣是一種以門檻值演算法為基礎而發展的演算法，是參考 改良型門檻值演算法使用次波形的概念，使用更具物理意義的方法來搜尋理想次 波形，Yang et al. ( 2012)已證實次波形門檻值演算法能有效改善南極附近海域的 測高觀測量精度，進而提升測高重力的品質，故本論文採用此演算法來改善衛星 測高的資料。 2-4 次波形門檻值演算法 次波形門檻值演算法主要分成兩個部分。第一部分為從測高衛星的反射波形 中選取有效的次波形。第二部分為使用門檻值演算法對有效次波形進行波形重 定。針對第一部分，對於反射波形中有效的次波形之選取，將以典型的散射波形 Brown 模型為出發點，表示如下:

   ≥ < − + = _{− −} τ τ σ τ α τ t t e t erf A t P )] _(t1 ₎ 2 ( 1 [ 2 ) ( (2-5) 其中 A 表示波形振幅，t 表示採樣時間，τ 表示波形前緣中點，erf 為誤差方程， σ 則與有效波高與衰減參數有關。而在測高衛星穩定之下，α 可視為常數。因此， (2-5)式表明反射波形主要分成波形前緣和波形後緣兩部分，而波形前緣受σ 和 A 影響，波形後緣則受σ，A 和 α 影響。 不同有效波高之波形需作分析，因此在這演算法中選取相關分析作定量分 析。相關公式如下: r Sr r Sr r ' ' = (2-6) 2 1(Pr'( ) Pr') 1 1 '

∑

₌ − − = k i i k Sr (2-7) ) Pr ) )(Pr( Pr' ) ( (Pr' 1 1 ' ₁ − − − =

∑

₌ i i k r Sr _ik (2-8) 2 1(Pr() Pr) 1 1

∑

= − − = k i i k Sr (2-9) 其中，P_r'(i)和P_r(i)，i=1，．．．，k 分別為參考波形前緣的反射波與任意波形

之反射波，Pr'與Pr分別為其平均波形，Sr'，Sr分別為任意波形與參考波形之 標準差，_{Sr' 為任意波形與參考波形之協變方。相關系數 r 表示兩波形之間相關}r 性，r = 0 表示兩者不相關，0 < r < 1 兩者正相關，–1 < r < 0 負相關，|r| = 1 兩者 完全線性相關。 在選取參考次波形時，考慮到反射波波形之採樣值存在誤差，因此在進行波 形相關性分析時，不適合把 Brown 模型的理想散射波形的全部採樣點當作參考 波形，其可能使得波形相關係數變小。因而選取其中誤差較小的波形前緣之採樣 點當作參考波形前緣，共 22 個採樣點。之後再利用此參考波形前緣與實際觀測 之反射波形進行移動相關性分析，以 ERS-1 每一實際觀測波形有 64 個閥門值為 例，可得到 43 個相關係數，從中選取相關係數極大值對應之次波形為有效次波 形。而第二部分，當決定了有效次波形後，利用門檻值演算法進行波形重定，公 式如下所列。

∑

∑

= = = is i is i i i t P t P A 1 1 2 4 ) ( / ) ( (2-10)

∑

= = 5 1 5 1 i i N P P (2-11) N N l A P Th P T =( − ) + (2-12)

first k k k l k k k r i P P P T G G G G + − − − + = − − − − 1 1 1 1 ( ) (2-13) is 為反射波波形前緣個數，A 為波形前緣振幅， P_i(t)為常態化後之第 i 點波形， N P 為前五個常態化後波形之平均值，Th 為門檻值，G 為第 k 點門檻且需優於常_k 態化後波形 T_l ，G 為波形重定門檻值。 _r 於本研究中，實驗了不同的門檻值，最後以門檻值百分之二十的成果最佳。

第三章 船測資料介紹與處理

3-1 船測資料介紹

用於本研究之船測資料來源分為我國調查船和國際上所有實測之船測重力

與船測海深資料，我國調查船資料為內政部民國 96 至 98 年大陸礁層計畫所測得

的海深值(船測軌跡圖 3-1 )(參考於 97 年度我國大陸礁層調查計畫)。而國際資料

主要來自於美國海洋及大氣管理局(NOAA)所屬的 National Geophysical Data

Center(NGDC)資料庫，資料內容包括了船測重力與海深值(船測軌跡圖 3-2)。 (http://www.ngdc.noaa.gov/mgg/gdas/ims/trk_cri.html ) NGDC 船測重力與海深值在本研究中除了作為海深計算的參考場之外，船測 重力值亦作為評估測高衛星重力模型精度的檢核資料。在大陸礁層海深值方面， 因為測深儀設備較精密，且有利用檢核測線交錯的方式，做水深資料交錯比對， 以確認定位與水深量測之正確性，故此海深值的精度較高，因此於本研究中主要 是當作檢核最終的計算海深模型。但是在釣魚台列嶼周遭，由於 NGDC 船測資 料比我國大陸礁層船測資料量還要多，因此檢核釣魚台周遭之海深模型將以 NGDC 船測資料為主(如圖 3-3 所示)。

(a) 圖 3-2-a 本研究所蒐集之 NGDC 船測海深資料軌跡圖

(b) 圖 3-2-b 本研究所蒐集之 NGDC 船測重力資料軌跡圖

圖 3-3 釣魚台列嶼周圍之船測資料 (來源:100 年度東海與南海島礁資料收集及建置工作期末報告書) 3-2 系統誤差改正與粗差剔除 NGDC 資料庫蒐集 1930 年以後全球各機構之船測資料，內容包括重力、磁 力、海深等資訊。但由於早期量測船的定位精度較差，進而導致船載重力儀器量 測精度相對較差，直到了 1980 年以後，隨著 GPS 定位系統的發展，船載定位的 精度因此逐漸改善。此外，NGDC 的誤差來源還包括了重力儀的漂移、重力參考 場的不一致、基站偏移等。因此為了得到較可靠的船測重力數據資料，於本研究

中使用二次多項式來擬合 EGM2008 重力(Pavlis et al., 2012)與船測所得之重力異

常間的差值，最後再將擬合後的差值與原差值進行比較，並將超過三倍中誤差數

可以由正常重力公式(Heiskanen and Moritz, 1967)為其基準，表示成下式: φ δ 2 1 0 a sin a g_e = + (3-1) 其中δ 為正常重力，φ 為大地緯度，g_e a 與₀ a 為此公式中的係數。 ₁ 在船測重力中，假設量測船的船速和航向都能維持固定，這樣基站偏移值和 基準偏移皆可視為 DC shift，如果假設量測船以固定速率沿著子午線航行，然後 我們以 Maclaurin 級數將 cosine 函數展開至二階，那我們可將(3-1)式改寫成時間 的多項式(3-2)式，且儀器飄移所造成之誤差也可利用(3-2)式來進行擬合。因此本 研究利用下式擬合衛星重力與船測重力間的差值。 2 2 1 0 ct c t c g_t = + + δ (3-2) 其中c 為常差，₀ c 與₁ c 為此公式中的系數，t 為自某一參考時刻起算之時間。最₂ 後利用最小二乘法求解後，再以平差後的差值來改正原有的船測重力值。 本研究利用(3-2)式來擬合船測重力和 EGM2008 重力(展開至 2160 階)的差 值。表 3-1 為經平差後獲得各航線c 、₀ c 和₁ c 之係數值。本研究選取 NGDC 船₂ 測重力中編號 hm951026 航線重力為基準(軌跡如圖 3-4 之紅色線)，利用 EGM2008 重力異常網格內插出船測重力位置的重力異常，圖 3-5 為平差前後船

測重力與 EGM2008 重力之關係，可發現兩者之間有一偏移量存在。所以若要一 起使用船測重力和 EGM2008 資料，則需改正此偏移量。圖 3-6 為編號 hm951026 航線平差前和平差後船測重力與 EGM2008 重力間差值的直方圖，從圖中可以發 現平差前其差值偏向正值（差值平均值為 8.6 mgal），也表示著 EGM2008 重力與 NGDC 船測重力之間存在著系統誤差，而平差後顯示其差值平均值為零且呈現常 態分佈。 表 3-1 利用二次多項式平差船測重力之結果 航次 c 0 c 1 c 2 84001311 -343.0403 71 22.370872 -0.735421 84003111 -197.8517 02 1.506723 0.871886 84003211 2.163079 -0.013879 -0.026283 84004211 4.977842 -0.880892 -1.215139 cook08mv 1.085584 0.408563 0.032856 cook09mv 3.824027 -1.621727 -2.089390 cook11mv 0.611860 -0.827811 2.842034 dme06 1.707851 -6.370574 -11.11592 8 dme07 217.27881 1 1.892987 -0.082336 dme10 3.460532 6.074065 -3.571870 dme24 53.242791 -3.823391 -0.048114 gh7407 -0.509721 -0.415944 0.019855 gh7501 1.244237 -0.273187 -0.001758 gh7504 -6.380317 12.285740 191.30020 4 gh7505 -0.018273 0.019659 0.010580 gh7702 3.577856 -0.661268 -0.202348 hm910613 5.932938 -69.223473 2163.073401

97 hm940217 -0.312659 -0.507566 0.141894 hm951026 8.512429 0.727957 0.097260 hs00t551 0.462885 -0.118551 -0.001312 hs7303 -8.422756 0.066541 0.005032 hs7304 -0.553957 -0.157137 -0.013754 hs7403 3.666843 -0.068526 -0.002767 hs7405 -4.038178 -0.204497 -0.000097 hs7501 -0.378923 0.043585 0.025080 hs7602 -1.899107 -0.133901 -0.003010 hs7603 -2.937129 -0.150387 0.007098 hs8102 1.520271 -0.005633 0.004660 hs8202 -8.543534 0.158420 0.069275 hs8203 -4.582956 0.123156 -0.009188 hs98tis1 6.462436 -1.535300 0.678081 hs99tis2 4.458984 -0.176126 0.006001 ht8301 0.849748 0.150576 0.019997 ht8401 0.635402 -2.095575 1.292226 ht8402 6.354468 0.160042 -0.077067 ht8403 5.077652 0.129470 0.007088 ht8502 2.825485 0.000279 0.002055 ht8510 1.830816 -0.133561 -0.274086 ht8605 2.116995 0.103719 -0.033687 ht861016 3.610312 0.104716 -0.048008 ht861102 3.965043 0.139434 -0.041579 ht8612 4.282286 -0.091953 -0.031223 ht8701 5.162246 0.122208 -0.022369 ht871110 2.341797 0.049004 -0.036354 ht880211 1.747715 0.171625 0.134496 ht880521 2.858505 0.009731 -0.050582 ht89t201 -0.326654 -30.866158 -207.842883 ht89t202 0.986975 0.611148 -0.008191 ht96t361 3.013343 0.758438 -0.037894 ht96t362 5.838448 0.720726 -0.025581 ht96t363 5.190033 0.028211 0.006744 ht97t373 6.006801 0.095080 -0.013653 ht97t374 6.375858 -3.011126 -0.424033

ht97t375 6.093875 1.183224 0.031028 ht97t376 5.382970 0.299847 0.021269 ht97t382 2.325528 1.071583 26.761016 ht97t383 3.936528 1.080320 3.563850 ht97t384 6.083109 -0.012442 -0.004198 indp05wt 18.073434 16.027554 -9.208000 jare27l5 1.185912 -0.536720 -4.161367 jare28l3 0.492962 -3.369432 4.327528 kh7003 2.424790 0.485427 -0.153602 kh7004 6.719925 5.944769 32.035439 kh7104 11.551908 2.467410 9.365644 kh7201 2.286420 15.241379 -208.429120 kh7202 -0.870644 0.010841 0.005145 kh7304 -9.016261 -3.377166 6.239237 kh7403 16.490431 -0.333843 -0.155044 kh7501 12.328660 0.432187 -0.060074 kh7502 4.988963 -0.509275 -0.058392 kh7505 16.762371 0.162449 -0.156198 kh7602 1.614401 0.050684 1.647971 kh7604 -5.340889 -6.181572 9.564644 kh7605 28.982474 -4.325373 -0.068972 kh7801 -4.533915 -0.949488 0.072376 kh7802 -5.189978 -0.784352 1.077192 kh79 -54.35125 2 -0.007054 0.110076 pol6202 49.467748 -13.070004 -334.571608 pol7201 23.790957 -0.596896 -0.034289 rc1206 2.904479 -0.049611 0.105664 rc1217 1.151362 -0.231665 -0.371459 rc1404 -0.487206 0.117638 0.118080 rc1710 0.212192 16.242900 69.625197 rc1711 -1.807808 -0.272465 0.670098 rc2006 -2.539064 0.570000 2.120642 rc2007 -1.509043 -0.159559 0.045232 mw9006 0.813513 -13.559318 81.141410 rc2008 -2.496085 1.990835 15.103490

rc2009 -1.601850 -3.297858 2.901898 v2008 4.085503 11.563009 -6.504647 v2107 -11.04880 5 -0.221654 2.530694 v2108 -6.816892 2.370434 1.379245 v2110 -11.85829 2 -46.369863 834.87359 5 v2814 7.936683 -0.533399 5.419354 v3309 2.069485 -149.204092 -704.179446 v3310 2.511896 0.102777 -0.013036 vit51 171.86464 1 -1.296134 -0.076946

圖 3-4 船測編號 hm951026 之軌跡圖

(紅色線為船測編號 hm951026 之軌跡)

圖 3-5 船測編號 hm951026 之重力比較圖

(b)

圖 3-6 編號 hm951026 航次(a)平差前與(b)平差後重力與 EGM2008 重力差之直方 圖

第四章 重力場與海深

4-1

計算海深方法 本文研究流程如參考圖 4-1，將波形重定後之測高資料利用兩種方法求得重 力異常，第一種為 IVM 搭配去除-回複法計算得重力異常(流程如參考圖 4-2)，第 二種為最小二乘配置法搭配去除-回複法求得重力異常(流程如參考圖 4-3)，再來 將兩種重力異常成果與 NGDC 船測重力比較，較佳之重力成果進行反演海深模 型。在計算海深模型方面，採用重力地質密度法計算本研究區之海深，本研究總 共有四種模型。第一種為將 DTU10 作為參考場，考慮不同地質密度差求得之海 深；第二種為將 DTU10 作為參考場，採用單一地質密度差(1.64 g cm−3 )求得之 海深；第三種為將 NGDC 船測資料作為參考場，考慮不同地質密度差求得之海 深；第四種為將 NGDC 船測資料，採用單一地質密度差(1.64 g cm−3 )求得之海 深。最後將四種海深成果與我國大陸礁層計畫之船測海深值進行比較與分析。本 章將會逐一介紹本文計算重力異常與海深之方法。

圖 4-1 研究流程圖

進行波形重定 計算殘餘大地水準梯度 利用 IVM 計算 重力異常 利用 LSC 計算 重力異常 與船測資料比較 建立最佳重力模型 測高衛星 觀測量(SSH) DTU10 model 濾波 & 重新取樣 資料編輯 船測資料 利用 GGM 建立海深模型 利用 GGM 建立海深模型

圖 4-2 IVM 計算流程圖 參考 重力異常 Add 計算最內 圈效應 重力異常 波形重定後 測高衛星觀測量 計算殘餘大地水準梯度 計算南北向、東西向 垂線偏差 IVM 計算殘餘重力異常

圖 4-3 LSC 計算流程圖 波形重定後 測高衛星觀測量 計算殘餘大地水準梯度 LSC 計算殘餘重力異 常 Add 參考 重力異常 重力異常

4-2 Inverse Vening Meinesz 公式(IVM)

IVM 求得重力異常的流程是將垂線偏差以最小二乘配置法求得南北向與東

西向之垂線偏差分量，然後利用一維快速傅立葉轉換計算出殘餘重力異常，最後

考慮重力場和最內圈效應影響後進而求得重力異常。詳細方法如下︰

根據 Vening Meinesz 公式(Hwang, 1998)，在 p 點的重力異常值可以經由垂

線偏差(DOV)的全球積分求得(圖 4-4)， q qpd H p g ε σ π γ σ

∫∫

′ = ∆ 4 ) ( 0

(4-1) 其中，(φ, λ)為緯度和經度，ε 為從 q 點到 p 點的垂線偏差組成，_pq γ 為參考橢球₀ 上的正常重力，ψ 為 p 到 q 的球面距離。_pq ψ 可藉由下面公式求得， _pq

(

q p

)

q p q p pq φ φ φ φ λ λ

ψ =sin sin +cos cos cos −

cos

(4-2)

在(4-1)式中的 kernel function 如下定義：       +       + + − = ′ 2 sin 1 2 sin 2 2 sin 2 3 2 cos 2 sin 2 2 cos pq pq pq pq pq pq H ψ ψ ψ ψ ψ ψ

(4-3)

圖 4-4 Inverse Vening Meinesz 公式的幾何變化 而(4-1)式中離散形式的數值積分為

( )

( )(

)

( )

(

)

[

]

( )

( )

(

)

[

sin

]

( )

} cos { 4 sin cos 4 cos 1 1 cos 1 1 1 1 0 cos cos 0 1 1 1 η α λ ξ α λ π λ φ γ α η α ξ λ π λ φ γ λ φ φ φ φ φ φ λ λ λ φ F H F F H F F H g qp qp qp qp qp qp qp p n q n q n q p ∆ ′ + ∆ ′ ∆ ∆ = + ∆ ′ ∆ ∆ = ∆

∑

∑ ∑

= − = =

(4-4)

其中ξcos =ξcosφ、ηcos =ηcosφ、∆λqp =λq−λp、∆φ和∆λ為緯度與經度的網格

間格，F 為一維的傅立葉轉換。使用 (4-4)式需要南北向與東西向垂線偏差分₁ 量，該些分量已提供在一規則網格上，這些格點可以由以下描述獲得。一個延軌 跡的垂線偏差是藉由兩個連續海面高度(SSH)的斜率一次近似。 d H H g =−( 2 − 1)

(4-5)

d 為兩點間的距離。而這需要南北向與東西向垂線偏差分量，然而我們可以利用 最小二乘配置法(LSC)(Moritz, 1980)來求得這些向量。如下公式：

l ) C (C C s sl ll n 1 − + =       = η ξ

(4-6)

向量 l 為延軌跡的垂線偏差，向量 s 為南北向 (ξ) 和 東西向 (η)分量，C 為 s_sl 和 l 的協變方矩陣， C_ll為 l 和 l 的協變方矩陣，C_n為雜訊矩陣，C_n為對角矩陣， 理論上為延軌跡的垂線偏差的變方。當使用參考重力模型時，重力模型的誤差須 考慮協方差函數。 4-3 最小二乘配置法(LSC) 此方法為根據先驗精度或是觀測值與估計值之間關係來計算協變方矩陣以 估算未知值，也是一種可以結合不同種類之大地觀測量來決定地球形狀和引力場 的方法(Moriz,1980)。相關公式如下︰

( )

C

(

C D

)

e g_res = _{ge e}+ −1 ∆ (4-7) 其中∆g_res為殘餘重力異常值。C_gε,C 分別為殘餘重力異常與殘餘大地梯度協變_ε 方矩陣、殘餘大地梯度之協變方矩陣。D 為雜訊值。 e 為海水面梯度觀測量。 求得殘餘重力異常後，可以利用下列公式(4-8)計算重力異常。

res short long g g g g =∆ +∆ +∆ ∆

(4-8)

其中∆g為重力異常，∆glong為長波長重力異常，通常由全球擾動位模式計算求 得，於本研究中採用 EGM2008 展到 2160 階，∆g_short為短波長重力異常，由於本 研究資料位於海水面而趨近於零，所以忽略不予計算，∆g_res為殘餘重力異常。 4-4 高斯濾波 由於測高衛星觀測量經波形重定後易含有高頻雜訊，所以於本研究中以高斯 濾波拿去高頻雜訊，下列為高斯濾波公式： 2 ) ( 2 1 σs d e W − = |ds|≤ Rs (4-9) 其中 2 A R_s = ， 6 A = σ 而 A 為搜尋罩窗大小 高斯濾波為低通濾波器（low-pass filter），在進行濾波時，搜尋罩窗大小 A 會影響最後的成果，若搜尋罩窗太小則無法偵測到粗差也無法除去原始觀測量的 高頻雜訊，至於搜尋罩窗太大的話，將會使資料太過平滑，如此一來，可能會使 觀測量失真。因此本研究在實驗階段嘗試不同的搜尋罩窗大小，最後以搜尋罩窗 大小 15 公里的效果最佳。

4-5 重力地質法(GGM)

本文算海深的方法為重力地質法，此法在 1972 年由 Ibrahim 和 Hinze 所提

出。在這理論方法中，所觀測之重力值可分為短波長和長波長的重力值

)

(

)

(

)

(

i

g

i

g

i

g

_obs

=

∆

_short

+

∆

_long

∆

(4-10)

其中∆g_obs(i)為觀測重力值，於本研究中觀測重力值為測高重力，∆g_short(i)與 ) (i glong ∆ 分別為短波長和長波長之重力值。短波長重力值由地形效應所引起，因 此反映出淺層海底地形的變化，而長波長重力則反映出深層的海底地形質量變 化。圖 4-5 為重力地質法的幾何架構，從圖中指出船測海深或是模型的海深可以 做為控制點，並利用布格平板公式推導短波長重力值 ) ) ( ( 2 ) ( ref short j G E j E g = ∆ − ∆ π ρ

(4-11) ) ( j g_short ∆ 表 示 第 j 點 的 短 波 長 重 力 ， G 為 萬 有 引 力 常 數 (=6.672×10-8 3

cm /g·sec2)，∆ρ為地質密度差(density contrast)(海水與岩石地質密度差，約 1.64 g −3

cm )，E(j)為第 j 點的海深值，E_ref為參考基高程值，通常以該區域最深的海深

值為參考高程值。將測高重力值減去短波長重力值後，可得到長波長重力值。將

波長重力網格，並利用下列公式，將其轉換為海深值。 ref short E G i g i E + ∆ ∆ = ρ π 2 ) ( ) (

(4-12)

E(i)為利用 GGM 計算之海深模型網格。(4-12)式的 E(j)，本文分別以 NGDC 船測 海深與 DTU10

2

′

×

2

′

網格海深值，搭配測高重力網格 ，進而得到最終的海深 模型。 圖 4-5 重力地質法的幾何架構 67 . 2 = ρ _g −3 cm 海平面 海底地形 海水 1 j j 2 P j ) (j₁ E ) (j₂ E ref E ) (jn E

4-6 向下延續法(DWC)

影響 GGM 海深模型的因素除了

∆glong( j),∆gobs( j)與∆gshort( j)的精度之

外，地質密度∆ρ也是決定海深結果的重要因素，本文利用向下延續法(downward continuation)，求定海底各層的重力異常值，然後比較各層重力後，計算出各層 的地質密度差，以取代傳統利用 1.64 g cm−3_。 Parker(1977)提到利用向下延續法求定不同層重力異常值，表達公式如下： ) , ( ) , ( 2 2 2 12 1 2 y x h f f h y x h f f e G f f G = π∆ x+ y (4-13) 其中，Gh2(fx,fy)與Gh1(fx,fy)分別為在利用二維傅立葉轉換之重力場，f 與x fy 分別為元素 x 與 y 之頻率，∆h₁₂ =h₂ −h₁為不同層之高程差。 為了降低雜訊的影響，所以在向下延續法中使用了高斯濾波，高斯方程式如下： 2 2 2 2_{( )} ) , ( k fx fy y x Gau f f e F = + (4-14) 其中，k 為高斯方程式頻率域中的半頻寬，通常 k 的大小決定了向下延續法成果 的解析度。 當向下延續法使用了高斯方程式，可以表達如下式： 2 2 2 2 2 2 12 2 1 ) ( 2 ) , ( ) , ( h fx fy k fx fy y x h y x h f f G f f e G = π∆ + − + (4-15)

通常在使用高斯濾波時後，雖然可以去除高頻的雜訊，但也可能因此忽略了真實 的訊號，因此在決定 k 之大小時，需要謹慎考慮。 本文研究區域最深約 7300 公尺，故利用 DWC 實際計算地質密度時，以每 730 公尺為一層，則總共分為 10 層，再分別將海平面的重力值，約化至各層， 計算出每一層所對應的地質密度差，結果如圖 4-6 所示。 圖 4-6 利用向下延續法，所計算出的海深各層的地質密度差。

第五章 重力場計算結果與討論

本章將分析利用衛星資料計算之六種不同的重力模型(如表 5-1)，分別為：

Case 1：利用 LSC 法計算重力網格(無濾波、無波形重定)。 Case 2：利用 IVM 法計算重力網格(無濾波、無波形重定) 。 Case 3：利用 IVM 法計算重力網格(使用高斯濾波、使用波形重定) 。 Case 4：利用 IVM 法計算重力網格(無濾波、使用波形重定)。 Case 5：利用 LSC 法計算重力網格(使用高斯濾波、使用波形重定)。 Case 6：利用 LSC 法計算重力網格(無濾波、使用波形重定)等。 本研究計算六種模型之重力異常成果，如表 5-2 所示。 本文用於檢核重力成果的資料為 NGDC 船測重力，由於本研究區域東海海 域的大陸棚約在海深 500m 以內，因此為了使檢核成果更加清楚與利於分析重力 成果，本文將檢核的船測重力資料分為淺海區(海深值小於 500m)與深海區(海深 值大於 500m)(圖 5-1)。除了本研究利用測高資料計算的重異常成果之外，亦同

時檢核目前著名的全球重力網格 DTU10 與 EGM2008(展開至 2160 階)，DTU

(Danmarks Tekniske Universitet)為 Anderson et al (2010)製作之 1 分全球重力網

格，同時 DTU10 模型也包含全球的海深資料(亦為 1 分解析度)，而 EGM2008(Earth

Graviational Model 2008)為 Pavlis et al. (2012)製作之全球重力網格，而此網格模

型的空間解析度為 5 分。所有比較結果一樣列於表 5-2。一般而言，波形重定技

訊，所以有無濾波也會有影響，因此，本章節將從各種不同的重力成果來探討淺 海與深海區域波形重定有無的比較、濾波有無的差別、IVM 與 LSC 的差別等主 題來作介紹。 (1) 淺海與深海波形重定比較 從表 5-2 中之 Case1 至 Case 6 中，深海地區精度明顯比淺海地區精度好，這 也顯示了測高衛星技術在近海地區精度較差，其反射波形會受到陸地或複雜的海 水面而受到汙染。從 Case 1 與 Case 2 的結果顯示，無論在深海地區或是淺海地 區，若測高資料無使用波形重定技術處理，則所計算出的重力網格精度遠不如其 他有進行波形重定的網格。整體來說，在本研究區域波形重定的技術可以有效提

升約 3~6 mgal 的重力精度。圖 5-2 為 Case 4 減去 Case 2 之差異圖，可以發現波

形重定的有無對於東海海域重力的成果有著全面性的差異，幾乎全研究區域存在 約±10mgal 的差值(對比無波形重定)

，且越靠近陸地的地方其差異越大。因此，

波形重定技術對重力值的變化影響最大，除了有效改善淺海地區的重力精度之 外，也能改善深海地區之重力精度。 (2) 濾波有無的差別 圖 5-3 為 Case 3 減去 Case 4 之差異圖，由此圖可發現，高斯濾波的使用造 成局部地區的重力值發生明顯的 3~5 mgal 的變化。由圖可發現琉球海溝周遭其

差異較大，這是因為此區域重力異常變化起伏較大。所以經由高斯濾波的處理

後，可以平滑波形重定後的觀測量，使得最後重力成果能夠移除雜訊。整體而言，

高斯濾波的使用，則可提升重力精度約 0.5 mgal。

(3) IVM 與 LSC 的差別

從 Case 1 與 Case 2 來看，在無使用波形重定技術下，IVM 法比 LSC 法計算

的重力網格精度改進約 3 mgal；而若是使用波形重定後，則僅改善約 0.5 mgal。

圖 5-4 為 Case 3 減去 Case 5 之差異圖，結果顯示使用使用 IVM 及 LSC 的差別，

僅僅在岸邊有較大的差值，其餘地區幾乎沒有變化。

由表 5-2 可知，於所有模型中，Case 3 的重力網格於本研究中，無論是在

淺海區或是深海區，都擁有最佳的精度，其精度也比 EGM2008 和 DTU10 等全

球網格更佳(約提升 0.5~1 mgal)。因此在計算海深模型所需要的重力網格部分，

表 5-1 本研究之六種不同重力模型

Case 計算方法 採用波形重定 使用濾波

Case 1 LSC NO NO

Case 2 IVM NO NO

Case 3 IVM YES YES

Case 4 IVM YES NO

Case 5 LSC YES YES

Case 6 LSC YES NO

表 5-2 各種重力網格與 NGDC 船測重力的比較差值 (mgal) Model depths Mean Std ev Max Min Case 1 LSC (no filter) (no retracking) All -0.3 13.4 94.7 -106.2 > -500m (shallow) -2.4 14.6 75.6 -92.4 < -500 (deep) 0.2 13.1 94.7 -106.2 Case 2

IVM (no filter) (No retracking)

All -0.4 10.0 94.6 -79.0 > -500m (shallow) -2.4 10.9 74.1 -67.7 < -500 (deep) 0.2 9.7 94.6 -79.0

Case 3 IVM (filter used) (retracking used)

All -0.4 7.0 86.3 -76.5 > -500m (shallow) -3.1 7.8 77.4 -76.5 < -500 (deep) 0.3 6.6 86.3 -57.3 Case 4

IVM (no filter) (retracking used) All -0.5 7.2 87.5 -77.6 > -500m (shallow) -2.7 8.0 79.8 -77.6 < -500 (deep) 0.1 6.9 87.5 -58.9 Case 5 LSC (filter used) (retracking used) All -0.5 7.1 86.5 -77.1 > -500m (shallow) -3.2 7.8 76.4 -77.1 < -500 (deep) 0.1 6.6 86.4 -56.6 Case 6 LSC (no filter) (retracking used) All -0.6 7.2 88.5 -78.4 > -500m (shallow) -2.8 8.2 88.5 -78.4 < -500 (deep) -0.1 6.9 87.5 -58.5 DTU10 All -0.4 7.1 87.6 -77.7 > -500m (shallow) -3.3 7.9 78.5 -77.7 < -500 (deep) 0.3 6.7 87.6 -59.5 EGM2008 All -0.4 7.3 92.3 -76.9 > -500m (shallow) -2.9 8.1 86.2 -76.9 < -500 (deep) 0.2 6.9 92.3 -58.4

圖 5-1 用於檢核重力網格精度的 NGDC 船測重力軌跡(紅色點)，淺藍色為淺海區 (海深值小於 500m)，深藍色為深海區(海深值大於 500m)。

圖 5-3 Case 3 減去 Case 4 之差異圖

第六章 海深測量計算結果與討論

本研究利用上一章節 Case3 重力網格來計算海深網格，因考量不同的參考

場、單一地質密度差與不同層地質密度差，所以本計算過程分為 4 個海深模型(表

6-1)，其計算海深之方法皆採用重力地質法(GGM)。在參考場方面，Case A 與

Case B 採用 DTU10(取樣至兩分網格)，而 Case C 與 Case D 採用 NGDC 船測重

力與海深值。在地質密度差方面，Case A 與 Case C 考慮不同層的地質密度差， 而 Case B 與 Case D 採用單一地質密度差(1.64 g −3 cm )。本文評估海深網格精度 方面，則利用內政部提供之大陸礁層計畫之船測海深值(軌跡參考圖 6-1)作為本 研究所計算之海深模型。由於研究區域位於東海海域，考量其地形特性，海深值 小於 2000m 多半位於平坦的大陸棚與沖繩海潮，介於 2000 至 4000m 之間則是琉 球海溝之海深變化梯度最大區域，大於 4000m 則為琉球海溝末端和深海海床。 因此本研究將檢核船測海深值分為三個部分，分別為海深值小於 2000m、2000m 到 4000m 和大於 4000m。除了分為以上三部分之外，釣魚台也位於本研究區域 之內，所以本研究也會將所計算之四個模型檢核釣魚台周遭之海深值，但因為我 國大陸礁層計畫無釣魚台周遭之海深，而 NGDC 有幾條航線軌跡較於接近釣魚 台列嶼，故在檢核釣魚台周遭海深值方面，本文採用 NGDC 船測值作為檢核。 本研究除了檢核自行計算之四個海深模型之外，也與目前全球較多人採用之一分 網格 ETOPO1 與一分 DTU10 進行比較。除此之外，本研究也探討不同的重力參 考模型對於海深成果的影響。

檢核結果如表 6-2、表 6-3 、圖 6-2-a~圖 6-2-d 所示。表 6-2 為海深網格與大 陸礁層船測海深值在研究區域之差值統計。表 6-3 為海深網格與 NGDC 海深值 在釣魚台群島區的差值統計。圖 6-2-a 至圖 6-2-d 則為各海深模型與大陸礁層海 深值之差異圖。本章節為了清楚的分析海深成果將分為三大部分，第一部分為整 個研究區域，將分別探討不同海深對於各海深模型之精度分析、不同參考場對海 深模型精度影響與地質密度差對海深模型精度之影響；第二部分則分析各海深模 型於釣魚台列嶼之精度；第三部分則探討討不同的重力參考模型對於海深成果的 影響，表 6-4 為不同重力模型的分類，其中重力模型選擇了表 5-2 中成果最好的 Case 3 與未做波形重定的重力模型 Case 2 來計算海深模型，而表 6-5 則為不同重 力模型求得之海深模型與大陸礁層海深值之差值統計。 (1) 東海海域

由表 6-2 可發現，在小於 2000m 的區域，Case A 與 Case B 的結果均較 ETOPO1

網格佳，於精度(std ev)提升了約 5~7 公尺，且平均值也比 ETOPO1 更接近 0。在 海深介於 2000 至 4000m 的區域，Case A 的平均值較 ETOPO1 佳，精度方面亦提 升了約 8 公尺。在海深大於 4000 公尺的區域，精度評估結果顯示 Case A 與 ETOPO1 模型成果差不多。而 ETOPO1 模型在不同海深成果會優於其他海深模 型的部分，吾人認為原因可能為 ETOPO1 模型並非單單是利用測高衛星資料所 計算之模型，其結合了所有全球船測資料的混合模型，故在計算海深方面的資料 來源比本研究的資料來的多。而海深越深的地方精度會比淺海地方來的差，也顯

示了純重力計算之海深模型，在離海平面越遠的地方，因為重力訊號比較不能反

映出真實的海底地形變化，所以其推算出的海深精度會較差。從 Case A 與 Case B

可以發現 Case A 整體成果優於 Case B，由此可見，採用不同層地質密度差在反

應深層海底地形方面會比採用單一地質密度差(1.64 g −3

cm )來的理想。而 Case C 與 Case D 之計算成果均不如 Case A 與 Case B，因此於本研究區域中，以 NGDC

船測資料作為參考場並不適合東海海域，另外從圖 6-2-c 與 6-2-d 中 Case C 與 Case D 的差異圖中，可以發現局部地區出現很大的差異，這也顯示以 NGDC 船 測資料作為參考場時，由於船測軌跡分佈的不均，若部分地區重力異常變化較 大，也會影響了最後海深的結果。所以在 4 個模型中，由 Case A 所計算的海深 成果最佳，其海深成果圖如圖 6-3 所示。 (2)釣魚台群島 在釣魚台區域的結果方面(表 6-3)，可以發現原先在整個研究區域中成果最

佳的 Case A 卻不如預期，而 ETOPO1 與 Case D 的成果較佳，其中 Case D 的成

果無論標準差或是平均值均優於 ETOPO1，原因為檢核釣魚台群島的海深資料為 NGDC 船測資料，而 ETOPO1 模型本身就包含了 NGDC 船測資料，而 Case D 則 是以 NGDC 船測資料作為參考場並搭配單一地質密度差所計算之海深模型。但 從上述可以得知，重力地質密度法在局部地區，以 NGDC 船測資料作為參考場 時，其成果會較佳，且在淺海區域，採用採用單一地質密度差(1.64 g −3 cm )會比 考慮不同層地質密度差來的適合。

(3)不同重力模型對海深成果的影響

從表 6-5 成果所示，Case 2 與 Case 3 於東海海域在重力異常方面差了 3 mgal，

在海深成果的反應上，於淺海區域差了約 60 多公尺，於深海區域差了約 30~40 公尺。 圖 6-1 用於檢核重力網格精度的大陸礁層海深值軌跡(紅色點)，淺藍色海深值為 小於 2000 公尺區域，中藍色海深值為介於 2000 至 4000 公尺區域，深藍色海深 值為大於 4000 公尺區域。 表 6-1 各種海深網格計算方式 Model 參考場模型 地質密度差 計算方法 case A DTU10 兩分重 力與海深模型 考慮不同層地 質密度差 GGM case B DTU10 兩分重 力與海深模型 單一地質密度 差(1.64 g cm-3) GGM case C NGDC 船測與 海深值 考慮不同層地 質密度差 GGM case D NGDC 船測與 海深值 單一地質密度 差(1.64 g cm-3) GGM

表 2-2 各種海深網格與大陸礁層船測海深值的比較差值 (m) (東海海域) Model depths Mean Std Max Min Case A > -2000m 9.2 107.3 1092.3 -820.8 -2000m~-4000m 9.7 201.7 1086.4 -755.6 < -4000m -92.2 196.6 1119.6 -1137.4 Case B > -2000m 6.6 105.8 1101.5 -870.4 -2000m~-4000m 3.9 218.3 1113.5 -864.4 < -4000m -112.6 199.5 1214.9 -1117.4 Case C > -2000m 23.1 149.3 1163.9 -995.2 -2000m~-4000m 17.5 296.4 1260.9 -1159.8 < -4000m -38.3 232.7 1211.7 -1377.8 Case D > -2000m 9.7 120.0 935.7 -804.5 -2000m~-4000m 11.2 301.7 1493.1 -1209.7 < -4000m -79.0 262.8 1373.4 -1497.7 ETOPO1 > -2000m 11.8 112.0 1006.1 -949.0 -2000m~-4000m 48.8 209.1 1128.9 -890.0 < -4000m -12.8 196.8 1101.7 -1343.0 DTU10 > -2000m -6.0 105.4 882.0 1201.9 -2000m~-4000m -9.9 224.7 857.6 1170.7 < -4000m 110.6 200.9 1174.8 -1223.3

表 6-3 各種海深網格與 NGDC 船測海深值的比較差值 (m) (釣魚台區域)

model Mean Std max min

Case A 10.3 108.1 243.8 -240.5 Case B 20.5 124.7 358.5 -314.9 Case C 20.1 101.9 629.7 -267.9 Case D 18.3 67.9 411.5 -172.2 ETOPO1 -46.3 69.9 166.7 -249.9 表 6-4 重力模型影響海深比較分類 Model 重力場模型 參考場模型 地質密度差 計算方法

Case I Case 2 DTU10

兩分重力與 海深模型

考慮不同層地 質密度差

GGM

Case II Case 3 DTU10

兩分重力與 海深模型 考慮不同層地 質密度差 GGM 表 6-5 不同重力模型影響海深比較統計(m)

Model depths Mean Std Max Min Case I > -2000m -26.4 171.2 956.3 -1635.3 -2000m~-4000m -15.5 246.9 1186.9 -1233.1 < -4000m 89.4 221.6 1121.1 -1187.7 Case II > -2000m 23.1 149.3 1163.9 -995.2 -2000m~-4000m 17.5 296.4 1260.9 -1159.8 < -4000m -38.3 232.7 1211.7 -1377.8

第七章 結論與建議

7-1 結論 本文利用數種測高衛星資料，經由波形重定後，以兩種計算重力異常的方 法，計算東海海域(經度122E−132E、緯度22N−32N)的重力網格，並進而 利用此網格搭配重力地質密度法推求出此區域的海深模型，綜合本研究的之方法 與計算成果，可以歸納出下列幾點結論： （1）對於測高衛星觀測數據不佳部分，可以經由波形重定進行修正，而在接近 陸地的測高衛星觀測量可以經由波形重定的技術有效的改正，如此一來，也可以 讓可用的資料變多，避免粗差產生，使其觀測量精度更好，進而有效提升最後的 重力模型精度。因此於本研究中，波形重定技術能在東海海域有效提升 3~6 mgal 的重力模型精度。 （2）本文所使用的兩種重力異常方法，由於觀測資料量龐大，在計算時間方面 上皆需要一段時間，在 IVM 與 LSC 計算時間上，IVM 明顯較 LSC 短，在重力 成果方面，IVM 法（使用高斯濾波）於東海海域能獲取最佳的重力網格。 （3）本文利用 IVM 法（使用高斯濾波） 於東海海域內，也能獲得比目前著名 的全球重力網格 DTU10 與 EGM2008（展開至 2160 階）更高的精度。 （4）重力地質密度法於東海海域內，能獲得比 ETOPO1 更高精度的海深模型， 即使是海深梯度變化極大的琉球海溝也是如此。而從不同模型的成果中，亦發現

採用重力地質密度法，參考場的不同與地質密度的考量皆會影響最後的海深成 果。 （5）在釣魚台區域，可以發現 Case D 的成果較 ETOPO1 佳，這也顯示了重力 地質密度法，若在局部地區將可以採用船測資料作為參考場來計算海深模型。 7-2 建議 綜合了上述的結論，於本章節可以歸納出以下幾點建議： （1）若未來要利用衛星測高資料來計算當地區域之海深，採用波形重定的技術 將可使計算之海深成果更佳，由於波形重定後，其觀測量易含有高頻的雜訊，故 在使用高斯濾波時，搜尋罩窗的大小將視不同區域進行調整，以東海區域中為例 子，本研究實驗之搜尋罩窗以 15 公里之結果最佳，但在不同研究區域中，可以 試驗最適合於該區域之搜尋罩窗大小，但在海底地形變化大的地區，如海溝或中 洋脊等部分，也可能存在高頻的雜訊，因此使用高斯濾波時需要謹慎考慮這一部 分。 （2）在採用重力地質密度法反演海深時，不同的研究區域可以實驗不同的參考 場，或是實驗不同的地質密度差，皆可得到不同的海深成果。 由於未來有更多更高精度的測高衛星正在發展（例如 Jason1 的 GM 任務）， 若能結合這些測高衛星觀測量，將可提高海洋重力異常之精度，進而改善利用重 力推算海深之成果。

參考文獻

陳松安，2001。衛星測高資料分析南海海潮、海水位、環流及渦漩，博士論文， 國立交通大學土木工程學系，新竹。

Andersen, O.B., 2010. The DTU10 Gravity field and Mean sea surface Second international symposium of the gravity field of the Earth (IGFS2), Fairbanks, Alaska.

Brown, G.S., 1977. The average impulse response of a rough surface and its applications, IEEE Transactions on Antennas and Propagation, Vol. AP-25, No.1, pp. 67-74.

Davis, C.H., 1997. A robust threshold retracking algorithm for measuring ice-sheet surface elevation change from satellite radar altimeter, IEEE Transactions on Geosciences and Remote Sensing, Vol. 35, No.4, pp. 974-979.

Deng, X., 2003. Improvement of geodetic parameter estimation in coastal regions from satellite radar altimetry, Ph.D. Dissertation, Curtin University of Technology, Australia.

Dixon, T.H., M. Naraghi, M.K. McNutt, and S.M. Smith, 1983. Bathymetric prediction from Seasat altimeter data,Journal of Geophysical Research., Vol. 88, pp. 1563-1571.

Frappart, F., S. Calmant, M. Cauhope, F. Seyler, and A. Cazenave, 2006. Preliminary results of Envisat RA-2 derived water levels validation over the Amazon basin, Remote Sensing of Environment, Vol. 100, pp. 252-264.

Heiskanen, W.A. and H. Moritz. 1967. Physical Geodesy, W.H. Freedman and Company, San Francisco and London.

Hsiao Y.S., J.W. Kim, B.Y. Lee, and C. Hwang, 2011. Bathymetry estimation using the gravity-geologic method: An investigation of density contrast predicated by the downward continuation method, Terrestrial, Atmospheric and Oceanic Sciences, Vol. 22, pp. 347-358.

Hwang, C., 1989. Height precision gravity anomaly and sea surface height estimate from Geos-3/Seasat satellite altimeter data, Rep.399, Department of Geodetic Science and Surveying, Ohio State University, Columbus, OH.

Hwang, C. and B. Parsons, 1995. Gravity anomalies derived from Seasat, Geosat, ERS-1 and TOPEX/POSEIDON altimetry and ship gravity: a case study over the Reykjanes Ridge, Geophysical Journal International, Vol.122, pp. 551-568

Hwang, C., 1998. Inverse Vening Meinesz formula and deflection-geoid formula: application to prediction of gravity and geoid over the South China Sea, Journal of Geodesy, Vol. 72, pp. 304-312.

Hwang, C., J. Guo, X.L. Deng, 2006. Coastal gravity anomalies from retracked Geosat/GM altimetry: improvement, limitation and the role of airborne gravity data. Journal of Geodesy, Vol.80, pp. 204–216.

Ibrahim, A., and W.J. Hinze, 1972. Mapping buried bedrock topography with gravity, Ground Water, Vol.10, pp. 18–23.

Martin, T.V., H.J. Zwally, A.C. Brenner, and R.A. Bindschadler, 1983. Analysis and Retracking of Continental Ice Sheet Radar Altimeter Waveforms, Journal of Geophysical Research., Vol. 88, pp. 1608-1616.

Moritz, H., 1980. Advanced Physical Geodesy, Herbert Wichmann, Karlsruhe.

Parker, R. L., 1977. Understanding inverse theory. Annual Review of Earth and Planetary Sciences, Vol. 5, pp. 35-64.

Pavlis, N. K., S. A. Holmes, S. C. Kenyon, and J. K. Factor, 2012.The development and evaluation of the Earth Gravitational Model 2008 (EGM2008), Journal of Geophysical Research., Vol. 117, pp. B04406.

Rapp,R.H.,1985.Detailed gravity anomalies and sea surface heights derived from Geos-3/Seasat altimeter data,Rep.365,Department of Geodetic Science and Surveying, Ohio State University, Columbus, OH.

Rummel, R., and R. H. Rapp, 1977. Undulation and anomaly estimation using Geos-3 altimeter data without precise satellite orbit, Bull. Geodesique, Vol. 51, pp. 73-88.

Seeber, G., 1983.Satellite Geodesy, De Gruyter, Berlin.

Smith, W.H.F., and D.T. Sandwell, 1994. Bathymetric prediction from dense satellite altimetry and sparse shipboard bathymetry, Journal of Geophysical Research., Vol. 99, pp. 21803-21824.

Wingham, D.J., C.G. Raplex, and H. Griffiths, 1986. New techniques in satellite tracking system, Proceedings of IGARSS’88 Symposium, Zurich, Switzerland, September, pp. 1339-1344.

Yang, Y., C. Hwang, H.J. Hsu, E.D. Chen, and H. Wang, 2012. A sub-waveform threshold retracker for ERS-1 altimetry :a case study in the Antarctic Ocean, Computers & Geosciences, Vol. 41, pp. 88–98.