國立臺灣大學理學院大氣科學研究所 碩士論文
Department of Atmospheric Sciences College of Science
National Taiwan University Master Thesis
颱風路徑與結構同化研究─系集卡爾曼濾波器 Assimilation of Tropical Cyclone Track and Structure
Based on the Ensemble Kalman Filter
連國淵 Guo-Yuan Lien
指導教授： 吳俊傑 博士 Advisor: Chun-Chieh Wu, Ph.D.
中華民國 98 年 6 月
June, 2009
致 謝致 謝 致 謝致 謝
感謝我的指導老師：吳俊傑教授，提供我良好的研究資源，與充分自由的研究環 境，耐心指導我們的研究，並關心我們的生活。
感謝台大物理系與大氣系的教授們在課堂上授予我各方面的知識。感謝 Prof.
Fuqing Zhang 與 Zhiyong Meng 無私提供我 EnKF 資料同化程式，使我能順利進行此研 究。感謝我的口試委員：曾忠一、林沛練、林博雄、楊舒芝教授與鄭明典主任，謝謝 你們費心閱覽我的論文，並提出許多寶貴的建議與指教。
在研究室的三年，感謝研究室的夥伴們。感謝昆炫、椿喜、忠權學長、占慧學 姐，熱心指導初踏入大氣研究領域的我，不時關心我的研究，並且在你們身上看到對 學術的熱忱，是我的好榜樣。感謝新淦學長、怡瑄、苡珊學姐，謝謝你們的陪伴，我 們共同完成一些事、經歷一段過程，並且你們領先我一步，讓我可以跟著你們的腳步 走。感謝紹良、呂易學弟、詩潔、啟雲學妹，與我一同學習，分享生活經驗。感謝助 理如茵和婉禎默默的付出。
感謝女友宛真，陪我一路走來。感謝登山時同行的夥伴們，與我一起感受這片土 地的美好。感謝佔盡我人生大部分時光的國中好友們，讓我永遠有可以依賴的地方。
感謝愛我的父母與家人，撫養我長大，並讓我可以無憂無慮的學習。
感謝其他未提及的人們，以及這塊土地。
摘 要摘 要 摘 要摘 要
近 30 年來，颱風的路徑預報有穩定的進展，但要在模式中建構出具有正確中心 位置、移速、以及合理渦旋結構的颱風初始場，一直是颱風數值模擬的一大挑戰，渦 旋植入、虛擬渦旋資料同化、或是渦旋重新移位等方法皆設計用以改善颱風的初始 化。隨著 EnKF 的發展，以 EnKF 直接同化渦旋位置的方法亦被提出，透過一計算渦 旋中心位置的觀測算符，可使模式中的渦旋保持在觀測路徑上。
在本研究中，我們進一步定義更多與颱風路徑和結構相關的觀測算符，包括中心 位置、渦旋移速與海表面軸對稱風速等。前二者的觀測量可使用作業單位依據衛星資 料所做的颱風定位，軸對稱風速剖面則可由經驗公式擬合衛星與如 DOTSTAR 或 T- PARC 的飛機觀測資料而來。
我們在配置 EnKF 資料同化系統的 WRF 模式中進行了颱風渦旋初始化和快速更 新週期同化分析等兩類實驗。在初始化實驗中僅同化颱風特殊觀測量，未採用任何現 有的渦旋植入方案。颱風的路徑與軸對稱平均風速剖面可在初始化時段末期與給定的 觀測資料相符，並且在僅同化最低層風速剖面的狀況下，整個垂直結構可被完善建 立。本方法最大的優點在於 EnKF 得到之分析場平衡性佳且相容於使用的模式，因此 後續預報的颱風強度得以穩定維持，不會在預報初期有劇烈調整的狀況。而在快速更 新週期同化分析實驗中，除了同化颱風特殊觀測量外，常態性探空儀與投落送資料也 被一併同化，以進行數天的長時段模擬分析。海表面軸對稱風速結構主要由 2008 年 T-PARC 實驗中 C-130 飛機連續 4 次穿越颱風中心偵察任務的資料而來，此資料對颱 風軸對稱風速同化有相當大的助益。實驗結果顯示，包括眼牆置換的過程在內，颱風 的路徑與結構的演變可被完善掌握。
由此實驗結果，我們認為這個新概念與技術提供了一個颱風初始化的有效方法，
並有潛力用以設計探討颱風結構的高解析度實驗，以及改善作業颱風模式的路徑與強 度預報。本研究已成功運用此方法，有效同化 2008 年 T-PARC 四架飛機所聯合觀測 之寶貴颱風資料，可為颱風動力探討帶來新的突破契機。
關鍵詞：颱風初始化、資料同化、系集卡爾曼濾波器、T-PARC
Abstract
Over the past 30 years, track forecasts of tropical cyclones (TC) have been in steady progress, but initializing a realistic vortex in the correct location and with the correct storm motion and structure remains a challenging task. Some techniques, including vortex bogusing, bogus data assimilation, and relocation, have been designed to improve the TC initialization. With the progress of ensemble Kalman filter (EnKF), assimilating vortex position based on the EnKF given an “observation operator” that computes the vortex position has also demonstrated capability in keeping a vortex along the observed track.
In this study, some new and effective observation operators related to the TC track and structure are proposed, including center position, velocity of storm motion, and sea surface axisymmetric wind structure. The observational quantities of first two parameters can be available from operational centers mainly based on the satellite analysis, and the radial wind profile can be evaluated through curve fitting using empirical formula, along with the aircraft surveillance data such as from DOTSTAR and T-PARC.
By assimilating these parameters based on the EnKF, we carry out two types of experiments in high-resolution WRF model: one for TC initialization, and the other for update cycle analysis. In initialization experiments, only special parameters for TCs are assimilated in a 24-hour period, without any extant bogus scheme. The TC track and axisymmetric wind profile well follow the specified observation data at final time of the initialization period. The overall vertical structure can be suitably constructed by assimilating only one-level wind profile. Moreover, one important benefit of this method is that almost no subsequent adjustment follows, indicating that the initial condition for forecast simulation after initialization period is dynamically balanced, as well as model- compatible. In update cycle experiments, both special parameters for TCs and conventional radiosonde and dropwindsonde data available from GTS are continuously assimilated, in order to perform an update cycle analysis of several days. The sea surface axisymmetric wind structure is determined from 4 continuous reconnaissance flights by C-130 aircraft during T-PARC in 2008. The result shows that the track and structure evolution of TCs, including the eyewall replacement cycle, can be captured in this simulation, indicating the usefulness of observations from reconnaissance missions in this method.
Our results suggest that this new technique provides an effective means of improving TC initialization and has good potential to help conducting some high-resolution numerical experiments to better understand the dynamics of TC structure, and to improve the operational TC model forecast.
Key words: tropical cyclone initialization, data assimilation, ensemble Kalman filter, T-PARC
目 錄目 錄 目 錄目 錄
致謝...i
摘要...ii
Abstract...iii
目錄...iv
圖表目錄...vii
第一章 前言...1
1.1 颱風初始化方法回顧...1
1.1.1 渦旋植入與重新移位...2
1.1.2 虛擬渦旋資料同化...2
1.1.3 同化渦旋中心位置...3
1.2 研究動機與目的...4
第二章 研究工具與方法...6
2.1 動力模式簡介...6
2.2 系集卡爾曼濾波器...6
2.2.1 卡爾曼濾波器與擴張卡爾曼濾波器...7
2.2.2 系集與樣本協方差矩陣...9
2.2.3 分析與預報方程...10
2.2.4 系集平方根濾波器...10
2.2.5 協方差擴張...11
2.2.6 協方差局地化...12
2.2.7 應用於中尺度區域模式...13
2.2.7.1 狀態變數...13
2.2.7.2 初始系集產生方式...13
2.2.7.3 擾動側邊界條件...14
2.2.7.4 巢狀網格...14
2.3 描述颱風渦旋的特殊觀測量...14
2.3.1 中心位置...15
2.3.2 移動速度...16
2.3.3 海表面軸對稱風速結構...17
2.3.3.1 颱風軸對稱風速結構經驗公式...17
2.3.3.2 觀測資料使用方式...20
第三章 颱風渦旋初始化實驗...22
3.1 模式設定...22
3.2 觀測資料與實驗設計...23
3.3 初始化時段...24
3.3.1 同化路徑與軸對稱風速結構的結果（TK-MS）...24
3.3.2 不同化任何資料的結果（NONE）...26
3.3.3 僅同化路徑的結果（TK）...26
3.3.4 颱風垂直結構的建立...27
3.4 預報表現...27
3.5 敏感度測試...28
3.5.1 模式解析度...28
3.5.2 系集規模...29
3.5.3 協方差擴張...30
3.6 討論...31
第四章 快速更新週期同化分析實驗...34
4.1 模式設定...34
4.2 觀測資料與實驗設計...35
4.2.1 例行性無線電探空儀...35
4.2.2 投落送...35
4.2.3 颱風路徑與軸對稱風速結構...36
4.2.4 實驗設計...37
4.3 同化分析時段...38
4.4 各組實驗的預報表現...39
4.4.1 TK-MS-TP-ALL 實驗的颱風路徑與結構預報...39
4.4.2 額外飛機觀測資料對同化路徑與軸對稱風速結構實驗的影響...40
4.4.3 額外飛機觀測資料對僅同化路徑實驗的影響...41
4.5 討論...41
第五章 總結...43
5.1 結論...43
5.2 未來展望...44
附錄 A 本研究的 EnKF 執行步驟...46
附錄 B 模式的垂直層設定與海表面風求取方式...47
參考文獻...48
附表...52
附圖...55
圖表目錄圖表目錄 圖表目錄圖表目錄
表 1 應用 WRF 模式上的 EnKF 資料同化系統之狀態變數列表。...52
表 2 本研究定義之三種描述颱風渦旋的特殊觀測量。...52
表 3 鳳凰颱風初始化實驗的各組實驗設計。...53
表 4 以同化颱風特殊觀測量做颱風渦旋初始化與傳統渦旋植入方案的比較。...53
表 5 辛樂克颱風同化分析時段內可用的觀測資料數目。C-130 的內核資料指距離 颱風中心 80 公里內的投落送資料，其他資料指 80 公里外的投落送資料。...53
表 6 辛樂克颱風同化分析時段內的 9 架次飛行觀測任務概況。...54
表 7 辛樂克颱風快速更新週期同化分析實驗的各組實驗設計。...54
圖 1 資料同化更新週期流程示意圖。...55
圖 2 Willoughby 片段連續剖面公式的示意圖以及觀測資料擬合範例。 (a) 上圖粗 線 Vi 為內部風速曲線（最大風速半徑以內）， Vo 為外部風速曲線（最大 風速半徑以外），兩區域交接處以下圖之權重函數 w 平滑化，最後所得風速 剖面為上圖著色區。 (b) 著色區為 Diana 颶風在 1984 年 9 月 11 日的觀測風速 剖面，粗線為公式擬合結果，其中外部風速曲線採雙指數遞減疊合。（取自 Willoughby et al. 2006，Fig. 1、Fig. 2）...55
圖 3 鳳凰颱風初始化實驗的模式範圍設定。總共 2 層網格，其中最內層網格會追 隨颱風渦旋中心移動，圖中顯示了起始時間（7 月 25 日 21 時）和結束時間 （7 月 26 日 12 時）的位置。...56
圖 4 鳳凰颱風初始化實驗的時間流程。...56
圖 5 鳳凰颱風在 2008 年 7 月 26 日 12 時由觀測資料決定出來的海表面軸對稱風速 曲線。圖中綠色倒三角型標記為 DOTSTAR 任務之投落送測得的海表面風速， 橘色正方形標記為 JTWC 的 8、10 級風暴風半徑估計，淺紫色正方形標記為 JTWC 的近中心最大風速估計，紅色粗線為 Willoughby 片段連續剖面公式， 即為實際同化至模式中的風速曲線。...57
圖 6 在 2008 年 7 月 26 日 12 時執行的 DOTSTAR 偵察觀測任務飛行路線與，疊加 上當時的低軌道衛星影像。黑色實心圓圈為投落送投擲地點，風標為投落送 測得的 925 hPa 高度風（單位：knots）。...57 圖 7 鳳凰颱風 TK-MS 實驗的系集路徑。灰色粗線為觀測颱風路徑，路徑上時間標
記的前 2 位數字為日期，後 2 位數字為小時。彩色粗線為系集平均路徑，細 線為每個系集成員的路徑，依颱風中心最低海平面氣壓以不同顏色顯示（如 圖例）。...58 圖 8 鳳凰颱風 TK-MS 實驗中海表面風場與軸對稱平均風速結構的變化。上、中、
下圖各為 7 月 25 日 15 時、26 日 0 時、26 日 12 時。 (a) 皆為系集平均。色階 為海表面風速量值（單位：m/s），等值線為海平面氣壓（單位：hPa），黑 色粗線為初始化時段起始時間至當時的觀測路徑（同化至模式中的資料），
藍色粗線為同一時段的模式颱風路徑。 (b) 黑色粗線為觀測的颱風海表面軸 對稱風速剖面（同化至模式的資料），紅色粗線為模式系集平均的海表面軸 對稱風速剖面，灰色細線為個別系集成員的風速剖面（單位：m/s）。...59 圖 9 鳳凰颱風 TK-MS 實驗中海表面軸對稱風速和中心最低海平面氣壓隨時間的演 變。 (a) 色階為系集平均的海表面環狀平均切向風速（單位：m/s）。 (b) 紅 色粗線為系集平均颱風中心海平面氣壓，灰色細線為個別系集成員的中心氣 壓（單位：hPa）。...60 圖 10 鳳凰颱風 TK-MS 實驗 (a) 颱風中心位置（單位：km）與 (b) 軸對稱風速徑
向剖面（單位：m/s）的平均誤差與系集散布隨時間的演變。細線為系集平均 和同化觀測資料的差值，粗線為系集成員間的標準差，其中實線表示預報系 集的結果，虛線表示分析系集的結果。粗點線為觀測誤差。（詳細定義請參 考內文中 (a) (39)、(40) 式與 (b) (41)、(42) 式。）...60 圖 11 Chen et al.（2007）在二維正壓模式中同化渦旋中心位置的表現。細線為颱風
中心位置的平均誤差，粗線為系集散布（單位：km）。其中虛線表示預報系 集的結果，實線表示分析系集的結果（實線、虛線意義與圖 10 相反）。（取 自 Chen et al. 2007，Fig. 10(b)）...61 圖 12 同圖 7，但為鳳凰颱風 NONE 實驗。...61 圖 13 同圖 8，但為鳳凰颱風 NONE 實驗。...62
圖 14 同圖 9，但為鳳凰颱風 NONE 實驗。...63
圖 15 同圖 10，但為鳳凰颱風 NONE 實驗，並少了預報系集（原圖虛線）的結果。 ...63
圖 16 同圖 8，但為鳳凰颱風 TK 實驗。...64
圖 17 同圖 9，但為鳳凰颱風 TK 實驗。...65
圖 18 同圖 10，但為鳳凰颱風 TK 實驗。...65
圖 19 鳳凰颱風 (a) TK-MS、(b) TK、(c) NONE 實驗初始化時段結束時（2008 年 7 月 26 日 12 時）通過颱風中心的東西方向垂直剖面結構。色階為水平風速 （單位：m/s），等值線綠色實線為位渦（單位：PVU），藍色虛線為位溫 （單位：K），黑色箭頭為垂直環流向量。...66
圖 20 鳳凰颱風 TK-MS 實驗的系集路徑預報。灰色粗線為觀測颱風路徑，路徑上時 間標記的前 2 位數字為日期，後 2 位數字為小時。彩色粗線為系集平均路徑， 細線為每個系集成員的路徑，依颱風中心最低海平面氣壓以不同顏色顯示 （如圖例）。紅色叉號為預報起始時間的系集平均颱風中心位置。...67
圖 21 鳳凰颱風 TK-MS 實驗的系集強度預報。 (a) 颱風中心最低海平面氣壓。綠線 為模式預報結果，紅線為中央氣象局給定的觀測值。 (b) 颱風近中心最大風 速。綠色實線為模式中出現的最大風速，綠色虛線為最大軸對稱平均切向風 速，皆取模式最低層（海表面）資料，紅線為中央氣象局給定的觀測值。..67
圖 22 同圖 8，但為鳳凰颱風 LOW 實驗。...68
圖 23 同圖 10，但為鳳凰颱風 LOW 實驗。...68
圖 24 同圖 8，但為鳳凰颱風 HIGH 實驗。...69
圖 25 同圖 10，但為鳳凰颱風 HIGH 實驗。...69
圖 26 同圖 8，但為鳳凰颱風 SMALL 實驗。...70
圖 27 同圖 10，但為鳳凰颱風 SMALL 實驗。...70
圖 28 同圖 8，但為鳳凰颱風 BIG 實驗。...71
圖 29 同圖 10，但為鳳凰颱風 BIG 實驗。...71
圖 30 同圖 8，但為鳳凰颱風 INFLA-0.5 實驗。...72
圖 31 同圖 10，但為鳳凰颱風 INFLA-0.5 實驗。...72
圖 32 同圖 8，但為鳳凰颱風 INFLA-0.95 實驗。...73
圖 33 同圖 10，但為鳳凰颱風 INFLA-0.95 實驗。...73
圖 34 辛樂克颱風同化分析實驗的模式範圍設定。總共 3 層網格，其中最內層網格 會追隨颱風渦旋中心移動，圖中顯示了起始時間（9 月 9 日 3 時）和結束時間 （9 月 13 日 3 時）的位置。...74
圖 35 辛樂克颱風同化分析實驗的時間流程。...74
圖 36 辛樂克颱風同化分析時段內常態性探空與投落送資料的空間分布。紅色圓圈 為常態性探空儀觀測位置，帶有底色的標記是投落送觀測位置。各飛機的代 表色見圖例，其中 C130 飛機投擲的投落送分為兩種：「C130 outer」代表位 於颱風中心 80 公里外的投落送，「C130 inner」代表位於颱風中心 80 公里內 的投落送。同一種底色不同標記代表不同的飛行架次。...75
圖 37 由美軍 C-130 飛機穿越偵察與 DOTSTAR 之投落送決定出來的辛樂克颱風海 表面軸對稱風速結構。左欄 (a)、(c)、(e)、(g) 為 C-130 連續 4 天偵察任務測 得的風速剖面，黃色細線為飛行高度（700 hPa）風速，橘色細線為機載 SFMR 測得的海表面風速，紅色粗線為同化至模式中的風速曲線，綠色倒三 角型標記為 DOTSTAR 任務之投落送測得的海表面風速，橘色正方形標記為 JTWC 的 8、10、12 級風暴風半徑估計，淺紫色正方形標記為 JTWC 的近中 心最大風速估計。右欄 (b)、(d)、(f)、(h) 為最接近 C-130 觀測時間的低軌道 衛星影像，疊加上 C-130 與 DOTSTAR 偵察任務的飛行路線，黑色實心圓圈 為 DOTSTAR 的投落送投擲地點，風標為投落送測得的 925 hPa 高度風（單 位：knots）。...76
圖 38 辛樂克颱風 TK-MS-TP-ALL 實驗的系集路徑。灰色粗線為觀測颱風路徑，路 徑上時間標記的前 2 位數字為日期，後 2 位數字為小時。彩色粗線為系集平 均路徑，細線為每個系集成員的路徑，依颱風中心最低海平面氣壓以不同顏 色顯示（如圖例）。...77
圖 39 辛樂克颱風 TK-MS-TP-ALL 實驗中海表面風速與海平面氣壓的演變，各小圖 的時間如左上角所示。色階為海表面風速量值（單位：m/s），等值線為海平 面氣壓（單位：hPa）。...78 圖 40 同圖 39，但為前 30 分鐘累積降雨量（單位：mm）。...79 圖 41 辛樂克颱風 TK-MS-TP-ALL 實驗中 9 月 11 日 5 時的 (a) 模式最大回波值（單
位：dBZ）與 (b) 同一時間低軌道衛星影像的比較。...80 圖 42 辛樂克颱風 TK-MS-TP-ALL 實驗中海表面軸對稱風速、軸對稱 30 分鐘降雨
量和中心最低海平面氣壓隨時間的演變。 (a) 色階為系集平均的海表面環狀 平均切向風速（單位：m/s）。 (b) 色階為系集平均的環狀平均 30 分鐘降雨 量（單位：mm）。 (c) 紅色粗線為系集平均颱風中心海平面氣壓，灰色細線 為個別系集成員的中心氣壓（單位：hPa）。...80 圖 43 辛樂克颱風 TK-MS、TK-MS-TP-OUT 與 TK-MS-TP-ALL 實驗 (a) 颱風中心
位置（單位：km）與 (b) 軸對稱風速徑向剖面（單位：m/s）的平均誤差與 系集散布隨時間的演變。紅色線為 TK-MS 實驗，藍色線為 TK-MS-TP-OUT 實驗，綠色線為 TK-MS-TP-ALL 實驗。細線為系集平均和同化觀測資料的差 值，粗線為系集成員間的標準差，皆取預報系集的結果。粗點線為同化入模 式中的觀測資料。（詳細定義請參考內文中 (a) (39)、(40) 式與 (b) (41)、
(42) 式。）...81 圖 44 辛樂克颱風 TK-MS-TP-ALL 實驗以 9 月 10 日 3 時為預報起始時間的預報
（TK-MS-TP-ALL_1003）。紅色橫向分隔線代表預報起始時間，在此以前皆 顯示分析實驗結果，以後皆顯示預報實驗結果（但亦繪出海平面氣壓的分析 結果以供比較）。 (a) 色階為系集平均的海表面環狀平均切向風速（單位：
m/s）。 (b) 色階為系集平均的環狀平均 30 分鐘降雨量（單位：mm）。 (c) 黑色與紅色粗線各為分析實驗與預報實驗的系集平均颱風中心海平面氣壓，
灰色細線為個別系集成員的中心氣壓（單位：hPa）。...82 圖 45 同圖 44，但為 9 月 11 日 3 時的預報（TK-MS-TP-ALL_1103）。...83 圖 46 同圖 44，但為 9 月 12 日 3 時的預報（TK-MS-TP-ALL_1203）。...83
圖 47 辛樂克颱風 TK-MS、TK-MS-DOT、TK-MS-TP-OUT、TK-MS-TP-ALL 四組 實驗以 9 月 10 日 3 時為預報起始時間的系集路徑預報比較。灰色粗線為觀測 颱風路徑，路徑上時間標記的前 2 位數字為日期，後 2 位數字為小時。彩色 粗線為系集平均路徑，依颱風中心最低海平面氣壓以不同顏色顯示（如圖 例），淺色的細線為個別系集成員的路徑。紅色叉號為預報起始時間的颱風 中心位置，在此以前為分析實驗結果，以後為預報實驗的結果，各組預報實 驗以不同顏色顯示（系集平均路徑以底色區分，系集成員路徑則直接以線條 顏色區分），見圖上標示。...84 圖 48 同圖 47，但為 9 月 11 日 3 時的預報比較。...85 圖 49 同圖 47，但為 9 月 12 日 3 時的預報比較。...85 圖 50 同圖 47，但為辛樂克颱風 TK、TK-DOT 兩組實驗 9 月 10 日 3 時的預報比較。
...86 圖 51 同圖 47，但為辛樂克颱風 TK、TK-DOT 兩組實驗 9 月 11 日 3 時的預報比較。
...86 圖 52 本研究採用的模式垂直層設定與分佈。圖中的數字為 WRF 模式的垂直座標
。...87
第一章 前言 第一章 前言 第一章 前言 第一章 前言
1.1 颱風初始化颱風初始化颱風初始化方颱風初始化方方法回顧方法回顧法回顧法回顧
近年來由於數值模式與資料同化方法的持續精進，數值模式產品已在天氣預報中 扮演相當重要的角色，在熱帶氣旋（tropical cyclones；TCs；後文皆統稱為颱風）的 預報方面也不例外，模式預報是重要的參考項目。但和一般天氣系統不同的是，颱風 渦旋的尺度較一般的綜觀天氣系統小且結構較為複雜（尤其是在眼牆附近），並且多 生存在缺乏觀測資料的海面上，使得現今作業上以測站、探空氣球、飛機與衛星觀測 資料為主的資料同化方式，在颱風渦旋的初始化上效果有限，常常無法僅用這些傳統 資料捕捉到颱風渦旋的真正位置與強度。這不但導致路徑預報的準確度降低，更使得 以動力模式進行颱風強度預報困難重重，一直無法獲得顯著的改善。
要解決這個問題，已有很多專為颱風渦旋設計的初始化方案被提出，如渦旋植入
（ vortex bogusing ） 與 重 新 移 位 （ relocation ） 、 虛 擬 渦 旋 資 料 同 化 （ bogus data assimilation；BDA）等。這些方案的主要精神是利用颱風高度旋轉與軸對稱的特性，
設法將我們認知上「合理」的渦旋結構放入模式初始場中的正確位置。但隨意改變初 始場經常會導致模式本身的不平衡，使得模擬初期會遭遇一段調整過程，整個模擬的 品質因而降低。又當我們植入軸對稱的渦旋結構時，非軸對稱結構與駛流的資訊不易 保留，亦會導致路徑預報變差。現今有相當多複雜的設計用來減少上述問題，但尚未 完全克服颱風渦旋初始化的問題。
也由於颱風高度旋轉與軸對稱的特性，傳統的颱風觀測通常以中心位置、移動速 度、中心氣壓、近中心最大風速、暴風半徑等資訊來表達，這些資料主要是藉由衛星 觀測推演而得，可大致描述一個颱風渦旋的狀態，在主觀預報上佔有相當重要的角 色。可惜模式通常無法直接使用這些資訊，僅能在渦旋初始化的過程中間接採用，對 模式預報的幫助有限。
以下對目前常用的幾種資料同化方式 進行簡要的回顧，並且也介紹 Chen and Snyder（2007）所提出的以系集卡爾曼濾波器（ensemble Kalman filter；EnKF）同化 渦旋中心位置的新方法。本研究即以此方法為基礎，做進一步深入探討與更詳細的檢 驗。
1.1.1 渦旋植入與渦旋植入與渦旋植入與渦旋植入與重新重新重新重新移位移位移位移位
Kurihara 等人（Kurihara et al. 1993、1995、1998；Bender et al. 1993）研發「地球 物理流體動力實驗室」（Geophysical Fluid Dynamics Laboratory；GFDL）颱風預報模 式中的颱風渦旋初始化方案，對渦旋植入有深入的研究與探討。他們的主要概念是透 過颱風模式積分來產生軸對稱風場，並設法與非軸對稱流場結合，可得到一平衡且具 有非軸對稱流場特徵的颱風結構，實際用於預報時，濾除原先分析場中較差的颱風結 構，以前述高解析度虛擬資料取代之。他們的研究結果顯示，透過這一系列方法可以 有效重建初始場中的颱風結構，進而減少 GFDL 颱風模式的路徑與強度預報誤差。
另一方面，若不做複雜的颱風初始化，我們可只將初始場中位於錯誤位置的颱風 結構移位至正確位置上，也可以使模式的颱風預報有部份改善。例如美國國家環境預 報中心（National Centers for Environmental Prediction；NCEP）的全球模式就有採用 此渦旋重新移位方案（Liu et al. 2002）。
現今透過這類颱風初始化方案已能獲得不錯的颱風模擬與預報，但即使經過許多 複雜的處理，模擬初期模式不平衡的問題仍然可能存在。並且，這些人為主觀決定的 程序也留存著很大的不確定性，包括如何適當地重建颱風的風場結構，如何給定溫 度、壓力、水汽等其他相關氣象變數，以及如何有效地同化使用額外觀測資料，如飛 機觀測等，不同的選擇可能會導致截然不同的路徑模擬結果，值得加以探討（ Wu et al. 2000；Chou and Wu 2008）。
1.1.2 虛擬虛擬虛擬虛擬渦旋資料同化渦旋資料同化渦旋資料同化渦旋資料同化
為解決渦旋植入方案可能造成模式不平衡的問題，也有人不做直接資料取代的動 作，改以變分資料同化（variational data assimilation）的技術將若干筆虛擬渦旋資料 同化至模式中，稱為虛擬渦旋資料同化。其中 Zou and Xiao（2000）首先提出以四維 變分資料同化將數小時時間窗口（time window）中的虛擬海平面氣壓資料同化至模 式中，藉著模式動力方程對其他相應氣象變數的限制，可以得到一較符合模式動力的 初始場。隨後 Xiao et al.（2000）亦提出同化虛擬三維風場資訊的方案，而 Wu 等人
（周昆炫 2003；Wu et al. 2006）則對多種虛擬渦旋資料同化的方式有詳盡的比較，並
以渦旋尺度的地轉調整（geostrophic adjustment）的概念強調了風場資訊對於颱風初 始化的重要性。
雖然以虛擬渦旋資料同化得到的初始場較無平衡性的問題，但四維變分資料同化 的計算相當耗費電腦時間，虛擬渦旋資料的建立也包含著諸多不確定性，缺乏簡單一 致的方法。
1.1.3 同化渦旋中心位置同化渦旋中心位置同化渦旋中心位置同化渦旋中心位置
Chen and Snyder（2007）提出衛星觀測得到的渦旋位置可視為一種特殊的觀測量
（observable，意為「可觀測量」，指在資料同化中可以被直接觀測的變數，後文皆 簡稱為觀測量），直接以資料同化的技術將其同化至模式中，就像同化由都卜勒雷達 觀測得到的徑向風速或是衛星觀測得到的亮度溫度一般。以資料同化的語言來說，就 是要定義一用以求取渦旋中心位置的觀測算符（observation operator），再藉此將觀 測得到的實際渦旋位置同化入模式中。這類特殊觀測量的觀測算符相當複雜，並且是 非線性的，例如要計算渦旋中心位置，通常需要先尋得氣壓極小值或渦度極大值的網 格點，再計算此點附近特定範圍內的某種權重平均位置。傳統的變分資料同化方式對 此處理不易，尤其是在四維變分資料同化中，更難以加入如此複雜的觀測算符；而在 EnKF 的資料同化系統上，則易於同化帶有複雜觀測算符的觀測量，並容許非線性的 觀測算符。
Chen 等人在二維模式上檢驗此方法的可行性。他們發現當初始猜測場的渦旋位 置誤差較渦旋本身空間尺度來得小時，EnKF 可使渦旋平順地移位至更接近觀測位置 的地方，而得到一個相當不錯的分析場，效果與傳統渦旋重新移位的方法類似。但若 初始猜測場的渦旋位置誤差與渦旋本身空間尺度相當甚至是更大時，由於物理空間上 誤差的特性會偏離高斯分布（Gaussian distribution），不符合 EnKF 的假設，因此分 析結果會不盡理想。但他們也指出，只要給予密集且精確的渦旋位置觀測資訊，進行 週期性（cycling）同化，模式中渦旋位置就不會產生太大的偏差，而這些密集、精確 的渦旋位置正是我們由衛星觀測資料中所能獲得的。因此，以 EnKF 同化渦旋中心位 置來使模式中渦旋保持在觀測路徑上依然是一個相當可行的方法。
他們並且發現連續同化渦旋中心位置的過程中，除了系集平均的渦旋中心位置保 持在觀測路徑上外，系集成員間的位置與結構差異亦可控制在一定的範圍內，即背景 場的不確定性可維持在較小的狀態。而且若進一步定義描述渦旋強度與形狀的觀測算 符，在同化渦旋位置的同時亦同化渦旋強度與形狀等資訊，可增加此方法的可用性。
他們並在連續同化這些渦旋參數的時段後再進行預報實驗，預報準確度比起未進行此 同化的對照組有所改善。
值得一提的是，此方法的效用雖然與過去渦旋植入與重新移位、虛擬渦旋資料同 化等方法類似，但由於 EnKF 與所搭配的模式動力相容的特性，得到的分析場的平衡 性相當好，可望改善過去隨意改變初始場導致模擬初期會遭遇一段調整過程的問題，
此方法的可行性與效用值得大家更進一步探討檢驗。
1.2 研究動機研究動機研究動機與目的研究動機與目的與目的與目的
如上述，Chen 等人已在二維模式上顯示 EnKF 應用在渦旋初始化之成效，而在 他們的相關工作（Chen and Snyder 2006）中，亦使用天氣研究與預報模式（Weather Research and Forecasting model ； WRF model ） 結 合 資 料 同 化 研 究 測 試 平 台 （ Data Assimilation Research Testbed；DART）對 2003 年至 2005 年間的 5 個颱風進行同化中 心位置與最低海平面氣壓的測試，證實在複雜的三維模式中此方法依然可行。但他們 使用的解析度最高只有 36 公里，無法解析眼牆附近的結構，並且未對各種參數的敏 感度進行較詳細的試驗。
另一方面，傳統的颱風報告上常有幾項描述颱風渦旋的特殊觀測量，除了中心位 置和最低海平面氣壓之外，尚有颱風移動速度、近中心最大風速與各風速值對應的暴 風半徑等。以往，這些我們熟知的參數幾乎僅能用於主觀預報上，對於模式模擬及預 報的用處有限，但若仿照上述同化渦旋中心位置的方法，將和這些參數相關的觀測算 符都建立在 EnKF 資料同化系統中，就可直接同化這些資料，或許對模式模擬及預報 可產生若干助益。
因此，本研究將使用配置有 EnKF 資料同化系統的 WRF 模式進行高解析度的模 擬，重新檢驗在高解析度、包含完整物理過程（full-physics）的區域模式上同化颱風
中心位置與其他特殊觀測量的效用，並檢驗各項敏感度。後文內容概述如下：第二章 介紹研究工具與方法，包括 WRF 模式簡介與 EnKF 的基本原理和執行方式，以及本 研究欲同化的描述颱風渦旋的特殊觀測量。第三章為本研究的第一部份實驗，將以此 方法進行短時段之颱風初始化，藉此在個案研究上可代替傳統渦旋植入與重新移位的 方案，成為颱風渦旋初始化的另一項選擇。第四章則為第二部份實驗，試著在學術環 境有限的資源中進行涵蓋一颱風生命週期的更長時間模擬，探索將此方法應用於作業 預報上的潛力，並探討 T-PARC（the Observation System Research and Predictability Experiment - Pacific-Asian Regional Campaign ； THORPEX-PARC ； Elsberry and Harr 2008）計畫期間額外飛機觀測資料與本方法結合的效益。第五章為結論，並提出值得 後續進行的延伸研究方向。
第第
第第二二二二章 研究工具與方法章 研究工具與方法章 研究工具與方法章 研究工具與方法
本 研 究 所 使 用 的 動 力 模 式 是 Advanced Research WRF （ ARW ） 模 式 2.2.1 版
（Skamarock et al. 2005），搭配上 EnKF 資料同化系統。此資料同化模組的最初是由 Zhang 等 人 （ Zhang et al. 2006 ； Meng et al. 2007 ） 所 發 展 ， 採 用 Whitaker and Hamill （ 2002 ） 的 方 法 ， 屬 於 系 集 平 方 根 濾 波 器 （ ensemble square root filter ； EnSRF）的一種。為了進行本研究的實驗，我們修改了部份程式，最主要的差異是多 加入了實驗 中將用到的描 述颱風渦旋的 特殊觀測算符 ，以及加入處 理巢狀網格
（nesting grid）架構的能力，並可使內層網格的位置追隨颱風渦旋中心移動，提昇數 值計算的效率。
2.1 動力模式簡介動力模式簡介動力模式簡介動力模式簡介
Advanced Research WRF（ARW）模式主要由美國國家大氣研究中心（National Center for Atmospheric Research ； NCAR ） 的 中 小 尺 度 氣 象 部 門 （ Mesoscale and Microscale Meteorology division；MMM）發展與維護，用以作為一有彈性的、先進的 大氣模擬系統，廣泛適用於數公尺到數千公里的尺度，並可運作於各種電腦平台以及 支援高速平行運算。其主要的應用包括教學、理想模擬、參數化研究、資料同化研 究、預報研究、即時天氣預報與耦合模式（coupled model）應用等。詳細資訊請見 ARW 的使用者網頁：http://www.mmm.ucar.edu/wrf/users/ 。
2.2 系集卡爾曼濾波器系集卡爾曼濾波器系集卡爾曼濾波器系集卡爾曼濾波器
系集卡爾曼濾波器（ensemble Kalman filter；EnKF；Evensen 1994、2003；Zhang et al. 2007；曾忠一 2006）最初由 Evensen 在 1994 年所提出，用以解決傳統卡爾曼濾 波器應用於非線性模式上的限制，以及過份龐大而不可及的儲存量與計算量。由於 EnKF 的核心概念簡單易懂，又不需要如四維變分資料同化中切線性算符（tangent linear operator）、伴隨方程（adjoint equation）與反向積分等的複雜處理，計算上相 對容易，因此很快就在大氣與海洋方面獲得廣泛的應用。EnKF 以完全相同的模式用 於分析與預報，因此分析場的平衡性相當好，而且 相容於該模式的動力性質。EnKF
也可說是傳統的系集預報（ensemble forecast）與資料同化兩者的巧妙結合，使系集 預報的成員更具代表性。後文將從傳統的卡爾曼濾波器談起，介紹本研究所使用的 EnKF 資料同化系統的理論基礎與實作方式，而逐步的執行流程則詳列於附錄 A。
2.2.1 卡爾曼濾波器卡爾曼濾波器卡爾曼濾波器卡爾曼濾波器與擴張卡爾曼濾波器與擴張卡爾曼濾波器與擴張卡爾曼濾波器與擴張卡爾曼濾波器
卡爾曼濾波器（Kalman filter；Kalman 1960）是一種序列資料同化（sequential data assimilation）的方法，在模式持續向前積分的過程中，只要有新的觀測資料加 入 ， 就 綜 合 當 時 的 預 報 場 和 這 些 觀 測 資 料 做 極 小 方 差 估 計 （ minimum variance estimate），得到新的分析場。為了極小方差估計的需要，我們除了要處理模式場在 時間上的推演外，也需要知道誤差的統計特性，以誤差協方差（error covariance）來 表示，並同時對誤差協方差也進行分析和預報。
若以向量與矩陣的形式來表達，模式場中的預報變數（prognostic variables）構成 的向量稱為狀態向量（state vector），以 xf 和 xa 表示，其中上標 f 和 a 分別代表 預報場和分析場。而狀態向量中任兩個元素的誤差協方差即構成誤差協方差矩陣
（error covariance matrix），以 Pf 和 Pa 表示，上標的意義與前面一致。由於 Pf 代 表 的 是 加 入 觀 測 資 料 更 新 前 的 背 景 狀 態 ， 又 或 可 稱 之 為 背 景 誤 差 協 方 差 矩 陣
（background-error covariance matrix）。
序列資料同化的執行流程如圖 1 所示。在一個更新週期中，包含有分析步驟和預 報步驟，分析步驟以上一次分析的短時段預報結果作為初始猜測值，加入觀測資料得 到新的分析場，此分析場也就是當次預報的初始場；而預報步驟則由模式將此分析場 向前積分直到下一個要同化觀測資料的時間，完成一次資料同化更新週期。一般大氣 模式的更新週期常為 6 小時或 12 小時（但本研究為 30 分鐘）。
分析步驟的方程式如下：
K=PfHTH PfHTR−1 (1)
xa= xfK yo−H xf (2)
Pa=I−K HPf (3)
在這裡 H 稱為觀測算符，用以將模式的狀態向量轉換成觀測量 ； R 為觀測誤差的 協方差矩陣，常為對角矩陣（diagonal matrix），即不同觀測間的誤差無相關性；由 (1) 式求得的 K 稱為卡爾曼增益矩陣（Kalman gain matrix），代表的是模式預報和觀 測資料之間的權重關係，用以在 (2) 式與 (3) 式中分別更新狀態向量與誤差協方差矩 陣。 yo 是觀測量構成的向量， I 是單位矩陣（identity matrix）。
而預報步驟包括模式場和誤差統計特性的預報，模式狀態向量的預報如下：
xn1f =M xn
a (4)
誤差協方差矩陣的預報則為：
Pn1f =M Pn
aMTQ (5)
其中下標 n 為時間指標（本文在與時間不相關的其他方程式中省略此指標）， M 為所使用的線性模式， Q 為模式誤差的協方差矩陣，亦即由於模式本身不完美所造 成的預報誤差。
值得注意的是，在卡爾曼濾波器中，觀測算符 H 和所使用的模式 M 都是線性 的，因為只有線性的模式在運算時能保持狀態向量誤差的高斯分布特性，而高斯分布 是卡爾曼濾波器的一個重要假設。對於輕微非線性的模式和觀測算符，我們仍能以類 似 卡 爾 曼 濾 波 器 的 方 程 式 來 處 理 ， 稱 之 為 擴 張 卡 爾 曼 濾 波 器 （ extended Kalman filter；EKF），其分析和預報方程分別為：
K=PfHTH PfHTR−1 (6)
xa= xfK[ yo−h xf] (7)
Pa=I−K HPf (8)
xn1f =mxn
a (9)
Pn1f =M Pn
aMTQ (10)
與卡爾曼濾波器的差異之處有： m 為所使用的非線性模式， h 為非線性觀測算符，
而 原 先 的 M 和 H 在 這 裡 分 別 為 m 和 h 的 切 線 性 算 符 ， 也 稱 為 Jacobi 矩 陣
（Jacobian）。但對於非線性程度較大的模式而言，擴張卡爾曼濾波器的適用性並不
一定良好，在某些非線性模式中會導致協方差演變的不穩定。
另外，在複雜的大氣模式上要執行擴張卡爾曼濾波器，除了面對非線性模式的問 題外，更大的困難是 Pf 和 Pa 的大小遠超出現今電腦可以處理的範圍。通常大氣模 式的維度可達 108，因此誤差協方差矩陣的維度將達到 1016，光是儲存就有很大的問 題，更不可能對其進行分析與預報。直到 EnKF 的提出，大大節省了卡爾曼濾波器的 計算量，才讓上述的資料同化流程有機會直接用在複雜的大氣模式上。
2.2.2 系集與樣本協方差矩陣系集與樣本協方差矩陣系集與樣本協方差矩陣系集與樣本協方差矩陣
EnKF 必須建構於系集（ensemble）模式上，也就是說要同時模擬很多份有著小 差異的模式場，稱之為系集成員（ensemble members）。系集成員的平均稱為系集平 均（ensemble mean），其意義相當於傳統確定性預報中的單一模式場；而系集成員 散布的情形則可代表模式場的機率密度（probability density）分布，亦即以有限的隨 機樣本來估計實際連續的誤差分布。
有了作為樣本的系集，EnKF 直接以預報樣本 xk
f 的協方差來求取 Pf ，亦即：
Pf= cov [xf, xf] = 1 K−1
∑
k=1 K
[x k f −xf][x kf −xf]T (11)
xf= 1 K
∑
k=1 K
x k f (12)
其中下標 k 為系集編號的指標， K 為系集成員的總個數，頂線表示系集平均。分 析誤差協方差矩陣 Pa 與分析樣本 xa f 的關係亦同。
通常觀測資料的數量都比模式狀態向量的維度還要小很多，我們亦可將 (1) 式中 的 Pf HT 與 H Pf HT 由狀態空間與觀測空間中預報誤差的互協方差直接計算：
Pf HT= cov [xf, h xf] = 1 K−1
∑
k=1 K
[x k f −xf][h x k f −hxf]T (13)
H Pf HT= cov [hxf, h xf] = 1 K−1
∑
k=1 K
[hx kf −h xf][h x k f −hxf]T (14)
而不用直接求取龐大的 Pf ，可大幅減少計算時間，並使用非線性的觀測算符 h 。
2.2.3 分析與預報方程分析與預報方程分析與預報方程分析與預報方程
傳統 EnKF 的分析步驟與一般卡爾曼濾波器類似，即以先前介紹過的分析方程分 別更新每一個系集成員：
K=PfHTH PfHTR−1 (15) xk a = x k
f K [y k
o −h x k
f ] (16)
但經 Burgers et al.（1998）證明，這裡用於不同系集成員的觀測量 y ko 須由原先的觀 測量加入協方差為 R 並呈高斯分布的隨機誤差而得：
y ko =yo k , cov [ yo, yo] = 1 K−1
∑
k=1 K
k
k
T = R (17)
才能使更新後的分析樣本具有正確的誤差統計特性。值得一提的是， (16) 式對每一 個系集成員的運算也將同時使系集平均如同一般卡爾曼濾波器般更新。
而 EnKF 的預報步驟，在完美模式（無模式誤差）的假設下，則完全等同於一般 的系集預報，也就是直接對各個系集成員做預報，以未經修改的非線性模式在有限樣 本上做誤差的推演。與擴張卡爾曼濾波器相比，不但計算量縮小至可以負荷的程度，
也免除了切線性模式未能如原始非線性模式正確傳遞誤差分布的疑慮。
為了使樣本具有足夠的代表性，系集成員的數量不可以太少（至少需要數十個以 上），否則 很可能會低估 實際上的誤差 ，甚至得到錯 誤的資訊，稱 為取樣誤差
（sampling error）。要減少取樣誤差的影響，除了盡量增加系集成員的數量（但這也 代表了計算負擔的增加）外，一些人為的處理方式也經常被使用，如協方差擴張
（covariance inflation）或協方差局地化（covariance localization）等，這些處理方式 將於後文中介紹。
2.2.4 系集平方根濾波器系集平方根濾波器系集平方根濾波器系集平方根濾波器
實際執行 EnKF 的分析步驟時，由於諸多好處，我們常不直接處理系集成員本 身，而是以各系集成員與其平均之間的偏差（ensemble perturbation）替代之：
x ' k f = x k
f − xf , x 'k a = x k
a − xa (18)
並可免除對觀測資料加入隨機誤差的必要，得到確定性的分析方程。由於這一類的分 析方案都涉及誤差協方差矩陣的平方根，因此通稱為系集平方根濾波器（ ensemble square root filter；EnSRF）。
本研究所使用之 EnKF 系統的分析方案如下（Whitaker and Hamill 2002）：
xa= xfK [yo−hxf] (19)
x ' k a = x ' k
f −K h' x k f , h ' xk f = h x k
f − hxf (20)
K=PfHTZdT−1 ZdZr−1 (21)
其中 Zd 和 Zr 分別是 H Pf HTR 與 R 的平方根：
H Pf HTR=ZdZTd (22)
R=ZrZrT (23)
它們都是 L×L 矩陣（ L 表示要一起同化的觀測量之個數），並且有非唯一解，可 透 過 Cholesky 分 解 （ Cholesky factorization ） 或 奇 異 值 分 解 （ singular value decomposition）等方法求得。此分析方案不需對觀測資料加入隨機誤差，在更新系集 平均時所使用的 K 與 (15) 式同，但以 (20) 式更新每一個系集成員的偏差時則需採 用由 (21) 式計算而得的 K 。
若採序列處理法，即一次僅同化一個觀測量，則 H Pf HT 與 R 皆變成純量，
K 的表示式可進一步化簡為：
K= K , =
1
H PfRHTR
−1 (24)2.2.5 協方差擴張協方差擴張協方差擴張協方差擴張
由前面的介紹可知，觀測量在 EnKF 的分析步驟中有兩方面的功能：其一是更新 模式預報的系集平均，使其接近於實際觀測結果，如 (19) 式；其二是更新每個系集
成員，使系集成員間的散布程度降低（縮減誤差協方差），亦即得到更加可靠（誤差 分布更小）的分析場，如 (20) 式。但以上的分析方程都是在無限系集大小與完美模 式假設下推導而來，使用較少的系集成員、較差的初始猜測，以及忽略模式誤差的結 果，都會導致 EnKF 的分析低估了應有的不確定性，進而使系集成員越來越相似，最 終可能喪失 EnKF 的可用度（Burgers et al. 1998）。強迫增加 EnKF 分析所得到的誤 差協方差可改善上述問題，稱為協方差擴張（covariance inflation）。
本研究做協方差擴張的方法為（Zhang et al. 2004）：
x ' k newa = 1− x ' k
a x ' k
f (25)
其中 為介於 0 ~ 1 之間的未定參數，若取 0 代表無協方差擴張，選擇的數值越大則 會使分析的系集偏差越接近原先預報的系集偏差（較大），亦即越強烈的協方差擴 張。此方法的優點為誤差協方差的人為調整只會在有觀測資料加入時才發生，並且可 以對不同的觀測資料類型設定各自的協方差擴張程度（選取不同的 值）。
2.2.6 協方差局地化協方差局地化協方差局地化協方差局地化
原始的 EnKF 分析方案會同時更新模式範圍內所有的狀態變數，更新的依據就是 每個狀態變數之間的誤差協方差。在一般的認知中，物理空間中相距過遠的兩個狀態 變數之間的相關性應該是微小，但在小系集的情況下，卻很容易得到錯誤的估計，將 使距離觀測位置遙遠的狀態變數有不合理的更新。為避免此問題，我們可人為定義一 觀測資料的影響範圍，在此範圍中使由有限樣本估計而得的誤差協方差隨距離遞減，
過遠的誤差協方差即設為零（不更新此狀態變數），稱為協方差局地化（ covariance localization）。
本研究使用 Gaspari et al.（1999）的五階函數做協方差局地化，誤差協方差依據 狀態變數在物理空間中與觀測位置的距離 r 乘上一係數 C ：
C r =
{
−12114
rrrr00
55−1212
rrrr00
445858
rrrr00
33−5353
rrrr00
22−51
rr0
4−23
rr0
−1 , 0, rr≤r≤r0 0 (26)由此式可知在 r=r0 時誤差協方差就減少為原估計的 0.208 倍，而在 r≃1.8 r0 以外就 幾乎等於 0。
在這邊要特別強調，協方差擴張和協方差局地化都是在有限計算資源下很人為化 的處理方式，雖然對 EnKF 的效用有立即的改善，但仍需仔細檢驗其是否帶來預期之 外的影響。本文將在 3.5 節中對此兩者做敏感度測試。
2.2.7 應用於中尺度區域模式應用於中尺度區域模式應用於中尺度區域模式應用於中尺度區域模式
2.2.7.1 狀態變數狀態變數狀態變數狀態變數
ARW 的預報變數包括網格點座標（經地圖投影後的正交座標，與經緯度座標在 方向上和比例上皆不完全一樣）上的水平風速分量 u 、 v 、垂直速度 w 、擾動位 溫（perturbation potential temperature） ' 、擾動重力位（perturbation geopotential）
' 、擾動乾空氣柱質量（perturbation dry air mass in column） 'd ，以及水象相關 成 份 的 混 和 比 ， 其 中 水 象 粒 子 種 類 的 多 寡 會 與 所 選 定 的 雲 微 物 理 過 程 （ cloud microphysics）有關（Skamarock et al. 2005）。本研究選取的雲微物理過程為 WSM 6- class graupel scheme，預報的水象相關成份有以下 6 種：水汽（water vapor； qv ）、
雲水（cloud water； qc ）、雲冰（cloud ice； qi ）、雨（rain； qr ）、雪（snow；
qs ）、霰（graupel； qg ）。以上這些預報變數就是我們在 EnKF 更新中會處理的 狀態變數，詳列於表 1。
2.2.7.2 初始系集產生方式初始系集產生方式初始系集產生方式初始系集產生方式
要進行 EnKF 資料同化，首先一定需要具備初始系集。初始系集可以很簡單地透 過在背景分析場加入隨機擾動而得到，擾動的大小要接近合理的背景誤差，以期能正 確代表誤差的統計分布。實作上加入隨機擾動的方式可以有非常多種，可擾動部份或 全部的變數，結果有可能產生平衡的場，也有可能產生不平衡的場。研究顯示通常不 理想的擾動所造成的不良影響不會持續太久，模式經過一段時間的積分後會漸趨平 衡，不同變數之間的相關性也會逐漸確立，配合上一定數量的觀測資料同化進入模式 中後，EnKF 即可正常表現（Evensen 2003）。
本研究透過 WRF 的三維變分資料同化系統（3D variational data assimilation tools for the WRF model；WRF-Var）來產生初始系集（僅採用部份程式，不需實做任何變 分資料同化），得到的初始系集將具有和 WRF-Var 中相同的背景誤差協方差（Zhang et al. 2006）。擾動是加入在 WRF-Var 的變換流函數場中（transformed streamfunction field；Barker et al. 2003），再反演回水平風、溫度與氣壓場，因此初始系集在大尺度 上是平衡的。其他的預報變數，包括垂直運動場和水汽場等，並未加以擾動。
2.2.7.3 擾動側邊界條件擾動側邊界條件擾動側邊界條件擾動側邊界條件
由 於 WRF 是 區 域 模 式 ， 運 作 時 需 要 事 先 給 定 側 邊 界 條 件 （ lateral boundary condition），若我們直接給予每一個系集成員完全相同的側邊界條件，會導致接近邊 界處的不確定性在模式積分中被不當移除，對後續的 EnKF 資料同化有不良影響。因 此，本系統在每一次的預報步驟中，會對每個系集成員的側邊界條件給予適當的隨機 擾動，有助於維持 EnKF 的穩定性。
2.2.7.4 巢狀網格巢狀網格巢狀網格巢狀網格
WRF 模 式 允 許 使 用 巢 狀 網 格 （ nested grids ） ， 其 中 雙 向 回 饋 （ two-way feedback）的設計可使內、外層網格保有一致性，亦即各層網格間對應網格點的資料 永遠是近似相同的。依此特性我們可將 EnKF 的分析步驟同時作用於每一層網格上，
對應的網格點將有一模一樣的更新，便可達成巢狀網格模式中的資料同化。
由於本研究關注於颱風的高解析度模擬，因此還需要採用追隨渦旋移動的巢狀網 格（vortex-following nested grids），將有限的計算資源做最有效的利用。此項特殊設 計包含一自動尋找颱風中心的演算法，會持續將內層網格移動至對齊颱風中心的位 置，維持眼牆附近有最高的解析度。我們將此項能力建構至這套 EnKF 系統中，使其 相容於資料同化的流程。
2.3 描述颱風渦旋的特殊觀測描述颱風渦旋的特殊觀測描述颱風渦旋的特殊觀測量描述颱風渦旋的特殊觀測量量量
本研究欲同化的颱風參數包括中心位置、移動速度、以及海表面軸對稱風速等，
皆為一般描述颱風渦旋常用的參數，但要從模式變數求取的運算較為複雜且具高度非 線性，因此鮮少看到有人將其直接應用於資料同化上。我們先將這三種颱風特殊觀測 量的幾項屬性列於表 2，後文將仔細說明相應觀測算符的詳細定義方式和觀測資料來 源及使用流程。
之所以選擇這三種參數，主要是想涵蓋過去各種颱風渦旋初始化方案的共同理 念，期能以同化這些參數來達成類似颱風渦旋初始化的效果，成為它們的替代方案。
中心位置和移動速度兩者共同代表的就是颱風的路徑，與渦旋植入或是重新移位至正 確位置的意義相仿；而在颱風結構方面，過去有許多方法嘗試決定出詳盡的三維渦旋 結構，但實際上由觀測所能獲得的資訊卻極為不足，使得結構的重建相當困難，不確 定性非常大，因此在本研究中我們僅選取近海表（即 10 米高度）的軸對稱風速結構 來同化，這算是實務上最容易取得，且具代表性的颱風結構資訊。在我們後面的實驗 中會證明（參見 3.3.4 颱風垂直結構的建立），由於 EnKF 的更新是架構在完整的模 式動力之上，藉著同化這樣的一維資料（單層軸對稱風速剖面）就有能力產生出合理 的三維颱風渦旋結構。
2.3.1 中心位置中心位置中心位置中心位置
颱風中心位置的決定是以低層氣壓中心為主，因為低層氣壓中心較容易定出穩定 且平順的颱風路徑，有利於我們以資料同化的方式來修正颱風路徑。實際上颱風中心 位置要以兩個觀測量來表示，即東西向與南北向的網格點座標值，計算方式是先找到 700~900 hPa 氣壓面上垂直平均重力位高度（geopotential height）極小值的網格點，
計算此網格點附近一定區域內重力位高度的平均值，再以「該地重力位高度低於此區 域平均值的量值」對鄰近網格點做權重平均（僅採計此量值為正值的網格點），得到 精確且路徑較為穩定的颱風中心位置：
wn=
{
〈 Φ〉 − 〈Φ0 n〉 , 〈Φ, 〈 Φnn〉 〈Φ〉〉 ≥ 〈Φ〉 (27)iTC=
∑
n=1 N
wn in
∑
n=1 N
wn
, jTC=
∑
n=1 N
wn jn
∑
n=1 N
wn
(28)
其中 Φ 為重力位高度，下標 n 作為此區域內各網格點的指標，角括號 〈 〉 表示低層 垂直平均，頂線表示區域平均， i 、 j 各為東西向與南北向的網格點座標，得到的 iTC 和 jTC 就是颱風中心位置的網格點座標，也就是我們要同化入模式中的觀測量。
以上的計算流程雖然有不少可任意選取的變數，但實際上當颱風有一定強度的時候，
各種選擇所得到的颱風中心位置不會有太大的差異。
颱風中心位置的觀測值可直接取自各作業單位提供的中心定位資料，這些定位主 要是依據多種衛星觀測而得。可取得的定位資料時間解析度通常為 3 小時或 6 小時，
我們視實際需要可把這些定位資料做 cubic spline 內插到 EnKF 同化的更新週期上。
在觀測誤差的設定方面，考量作業上合理的定位誤差（Edson et al. 2003），並觀 察實際進行多次 EnKF 同化測試的效用，本實驗中颱風中心位置的觀測誤差在單一方 向座標（東西向與南北向）上的值 Ri ,TC= Rj , TC 設定為 10 ~ 20 公里，因此合成的中 心位置觀測誤差值相當於：
Rpos=
Ri , TC2 R2j , TC=2 Ri ,TC (29)約在 14 ~ 28 公里之間。
2.3.2 移動速度移動速度移動速度移動速度
定義了颱風中心位置的計算方式後，移動速度的計算相對簡單，就只是由目前颱 風位置與一段時間前颱風位置的差異來求得，代表過去這段時間內的平均移速，與颱 風中心位置的同化相仿，移動速度的觀測量也要用網格點座標上東西向與南北向的分 量共同表示：
i 'TC= iTC−iTC, prev
∆ t , j 'TC= jTC− jTC, prev
∆ t (30)
這段時間間隔 ∆ t 不宜太短，因為短時間內微小的路徑擺動可能會算出劇烈變化 的颱風移速，會導致同化過程中劇烈且不合理的更新，對結果有不良影響；但也不宜 太長，以免忽略掉有意義的路徑轉折，合理的時間間隔可取在約 2 ~ 6 小時之間，本 研究皆取 3 小時，即計算過去 3 小時的平均移速。此運算雖然簡單，但實際上代表的 卻是牽涉到模式動力（介於這段時間內的模式積分）的相當複雜的觀測算符。
颱風移動速度的觀測值一樣可由各作業單位提供的中心定位時間序列計算而得，
計算方式同 (30) 式，只是將式中模式的颱風位置改成觀測的颱風位置。而颱風移速 的觀測誤差，經反覆測試，在單一方向分量上的值 Ri ' ,TC=Rj ' ,TC ，較合適的設定約 在 1.0 ~ 1.5 m/s（3.6 ~ 5.4 km/h）之間。
2.3.3 海表面海表面海表面海表面軸對稱風速結構軸對稱風速結構軸對稱風速結構軸對稱風速結構
海表面軸對稱風指的是海表面以颱風中心為圓心的環狀平均切向（tangential）風 速，是一個一維的剖面，每個半徑對應到一個風速值。由於環狀平均的動作，此數值 與颱風的移動和非軸對稱結構無關，因此同化此資料只是調整模式中颱風的軸對稱平 均結構，並不會使原先非軸對稱的資訊喪失，亦不會有影響駛流強度的顧慮。
實際同化這條風速剖面時，颱風中心位置依然由前述方法求得，但為了系集成員 間的一致性，對所有的系集成員皆採用系集平均的中心位置為原點。表面風速的定義 取 標 準 的 10 米 高 風 速 ， 在 海 洋 上 大 約 相 當 於 WRF 模 式 中 跟 隨 地 形 （ terrain- following）的垂直座標 = 0.9988 的地方，因此我們特別將垂直方向上最底層（含 水平風速的）網格直接設定在這個高度，也就是以最底層網格點的風速作為海表面風
（詳見附錄 B）。半徑 0 ~ 400 公里間的風速剖面會選取數十個點同時同化入模式 中，最外圍大約每隔 20 公里取一點，越靠近中心處取點的密度則逐漸加大。
此項資料的「觀測值」要從何而來是個有趣且複雜的問題。對大部分的颱風來 講，這種觀測是不存在的，但我們還是可以由若干間接的觀測得到部份資訊，再以某 些被檢驗過的參數化公式擬合之，建立整個剖面的軸對稱風速曲線。下面將分別說明 幾種颱風軸對稱風速結構的經驗公式，以及由觀測資料決定出經驗公式中各項參數的 方法。
2.3.3.1 颱風軸對稱風速結構經驗公式颱風軸對稱風速結構經驗公式颱風軸對稱風速結構經驗公式颱風軸對稱風速結構經驗公式
一般的理想實驗常使用 Rankine 渦旋（Rankin vortex）或是修正 Rankine 渦旋
（modified Rankine vortex）結構來建立理想颱風渦旋，部份渦旋植入方案也使用類似 的理想結構來做颱風渦旋的初始化：
