尹金方 梁旭東 陳鋒 劉英 何會中 梁釗明 鄒海波 徐建軍 郝世峰 謝衍新

（1 中國氣象科學(xué)研究院災(zāi)害天氣國家重點實驗室，北京 100081；2 浙江省氣象局，杭州 310017；3 江西省氣象災(zāi)害應(yīng)急預(yù)警中心, 南昌 330046；4 廣東海洋大學(xué)海洋與氣象學(xué)院，湛江 524088）

0　引言

大氣再分析資料在天氣、氣候研究以及各種應(yīng)用中已經(jīng)得到廣泛的應(yīng)用。針對區(qū)域天氣、氣候研究以及能源、環(huán)境等應(yīng)用需求，高時空分辨率的區(qū)域再分析得到快速發(fā)展。美國NCEP于2005年開始建立北美區(qū)域再分析系統(tǒng)（North American Regional Reanalysis，NARR）[1]，ECMWF開發(fā)了覆蓋整個歐洲的高分辨率再分析系統(tǒng)（High Resolution Regional Reanalysis for Europe，HRRRE）[2]，美國NCAR等機(jī)構(gòu)共同建立了北極區(qū)域的再分析系統(tǒng)（Arctic System Reanalysis，ASR）[3]，美國阿拉斯加州立大學(xué)等機(jī)構(gòu)建立了白令海峽大氣再分析系統(tǒng)（Chukchi–Beaufort High-Resolution Atmospheric Reanalysis，CBHAR）[4]，芬蘭和瑞典在ECMWF的組織下完成了波羅的海地區(qū)的再分析系統(tǒng)[5]，德國氣象局組織建立了德國及歐洲區(qū)域的再分析系統(tǒng)（CORDEX：High-resolution Regional Reanalysis for the European CORDEX Domain）[6]，日本在全球再分析的基礎(chǔ)上建立了東亞區(qū)域高分辨率（25 km）再分析系統(tǒng)。

近30年來，中國氣象局積累了豐富的觀測資料，為東亞區(qū)域大氣再分析提供了很好的基礎(chǔ)。針對再分析資料的廣泛需求，中國氣象局設(shè)立了“氣象資料質(zhì)量控制多源數(shù)據(jù)融合與再分析”攻關(guān)任務(wù)（2015—2020 年），東亞區(qū)域再分析是該攻關(guān)任務(wù)的一個攻關(guān)方向。通過攻關(guān)任務(wù)的實施，將用 6年時間，初步建立我國的東亞區(qū)域再分析資料集并建立再分析業(yè)務(wù)技術(shù)體系。

東亞區(qū)域再分析的研究任務(wù)包括：同化系統(tǒng)優(yōu)化及觀測資料同化效果研究；適合東亞區(qū)域的模式物理過程試驗研究；東亞區(qū)域再分析系統(tǒng)集成與再分析數(shù)據(jù)集的研制；東亞區(qū)域再分析數(shù)據(jù)集檢驗與評估。研究目標(biāo)是發(fā)展有效同化多源觀測資料的東亞區(qū)域再分析關(guān)鍵技術(shù),建成我國第一代東亞區(qū)域大氣再分析系統(tǒng)，制作10年（2008—2018年）的較高分辨率（東亞區(qū)區(qū)域10～15 km，局地3～4 km）的東亞區(qū)域再分析數(shù)據(jù)產(chǎn)品。

1　東亞區(qū)域再分析系統(tǒng)建設(shè)

建立集成的東亞區(qū)域再分析系統(tǒng)，包括資料預(yù)處理、資料同化、數(shù)值模式、檢驗與評估四個部分。

1.1　資料收集及預(yù)處理

目前可利用的觀測資料包括常規(guī)地面和探空、飛機(jī)觀測、GPS、雷達(dá)、風(fēng)廓線、衛(wèi)星導(dǎo)風(fēng)、野外試驗等觀測資料。區(qū)域再分析的特色之一是采用更多的全球再分析中沒有使用的資料。因此，東亞區(qū)域再分析中針對探空、雷達(dá)、地面觀測、地基GPS/Met、飛機(jī)報、風(fēng)廓線等資料進(jìn)行進(jìn)一步的收集，以便得到當(dāng)前東亞地區(qū)最為完善的觀測資料集。以探空資料為例（圖1），黑色圓點表示探空觀測的國際交換站空間分布，紅色圓點表示未參加交換的探空觀測空間分布。紅色圓點這部分觀測資料即為新增部分，在其它再分析資料中是沒有使用到的。此外，區(qū)域再分析中探空資料采還包括時間加密觀測，從而一定程度上增加了探空資料的空間分辨率和時間分辨率。需要指出的，雷達(dá)觀測資料是全部未參與國際交換。

圖1　探空觀測的空間分布（紅色點為中國未參加國際交換的探站點，填色代表地形）Fig. 1 Spatial distribution of sounding stations (The red dots represent the observational data which are not shared with the Global Communications System (GTS), and the shaded with colors is terrain)

1.2　同化系統(tǒng)

基于Gridpoint Statistical Interpolation（GSI）資料同化系統(tǒng)[7]，開展了GSI同化系統(tǒng)的本地化移植和優(yōu)化，實現(xiàn)了地面、探空、雷達(dá)、風(fēng)廓線、GPS/MET等多源資料的同化，從而建立了基于GSI的循環(huán)同化模塊。通過長時期模擬試驗測試了GSI同化系統(tǒng)的穩(wěn)定性。針對不同的觀測資料，開展了觀測資料的同化效果分析。例如，同化風(fēng)廓線資料后，對低層風(fēng)場改進(jìn)明顯，從而激發(fā)上升運(yùn)動。

GSI中原始雷達(dá)徑向風(fēng)資料同化是以徑向風(fēng)作為觀測算子進(jìn)行同化，缺乏有效的切向風(fēng)信息，存在一定不足。為了改善GSI同化系統(tǒng)對雷達(dá)資料的同化效果，重點改進(jìn)了雷達(dá)徑向風(fēng)同化算子。GSI中的原始徑向風(fēng)同化算子為問題

當(dāng)雷達(dá)仰角較小時（小于5°），垂直分量對雷達(dá)徑向風(fēng)影響較小，將徑向風(fēng)中垂直分量扣除后式（1）變?yōu)椋?/p>

GSI中原始徑向風(fēng)同化方法是將式（2）的右端作為觀測算子。通過式（2）可以發(fā)現(xiàn)，原始徑向風(fēng)觀測算子僅能引入單一的徑向風(fēng)增量，通過同化系統(tǒng)后會產(chǎn)生風(fēng)場兩維增量，即同化系統(tǒng)在吸收徑向風(fēng)信息后根據(jù)背景誤差協(xié)方差矩陣調(diào)節(jié)得到風(fēng)場的切向分量。也就是說，分析風(fēng)場在徑向上吸收了觀測資料而體現(xiàn)出觀測場的特征，而切向分量則無法反映觀測場的切向特征。

為了改進(jìn)雷達(dá)徑向風(fēng)的同化效果，引入了以IVAP反演方程組為基礎(chǔ)的同化算子[8-9]。該方法采用提取單位分析單元內(nèi)的平均風(fēng)信息在同化系統(tǒng)背景場和物理條件的約束下，求取式（4）左右兩端極小化得到分析風(fēng)場。

根據(jù)上述同化理論，改進(jìn)了的GSI雷達(dá)徑向風(fēng)同化算子，并開展了個例模擬試驗分析。針對2013年菲特臺風(fēng)，選取2013年10月6日06—12時刻雷達(dá)徑向風(fēng)資料進(jìn)行。采用原始方法（RW）和新方法（IVAP）分別同化雷達(dá)資料后可以看到，同化雷達(dá)徑向風(fēng)在一定程度上增強(qiáng)了近臺風(fēng)中心的風(fēng)場結(jié)構(gòu)，尤其是在臺風(fēng)第一、四象限的風(fēng)速增加較為明顯。從圖2中可以看出，改進(jìn)算子試驗?zāi)M出的臺風(fēng)風(fēng)場增量更為顯著，尤其是在臺風(fēng)第四象限模擬出了約10～15 m/s的西南風(fēng)增量，使整個風(fēng)場增量場閉合形成了較為完整的逆時針渦旋結(jié)構(gòu)，也增強(qiáng)了整個臺風(fēng)的結(jié)構(gòu)和強(qiáng)度。這主要得益于IVAP方法將部分切向風(fēng)信息的引入，增強(qiáng)了對整個風(fēng)場的調(diào)整程度，得到了更為準(zhǔn)確的風(fēng)場。

圖2　同化后3 km高度層風(fēng)場增量和海平面氣壓：（a）原方法和（b）改進(jìn)方法Fig. 2 Horizontal wind increments at 3 km level above ground after radar data assimilation with (a) the original method and (b) the improved method. The shaded in colors is wind speed and the contours are sea level pressure

圖3 給出了三組試驗對臺風(fēng)未來12 h的路徑和強(qiáng)度預(yù)報。由圖可見，得益于雷達(dá)徑向風(fēng)同化對臺風(fēng)初始結(jié)構(gòu)和風(fēng)場特征的改進(jìn)，RW試驗和IVAP試驗預(yù)報的臺風(fēng)路徑均比對照（CTL）試驗有明顯改進(jìn)，三組試驗的平均路徑誤差分別為77.2、61.3和44.7 km。RW試驗在1～5 h內(nèi)路徑誤差波動明顯，而后均優(yōu)于CTL試驗。IVAP試驗得益于同化對1～5 km有雷達(dá)資料區(qū)域的風(fēng)場的改進(jìn)，其預(yù)報的前8 h一直維持在40 km范圍內(nèi)，而后逐漸增加到80 km，比RW試驗有明顯優(yōu)勢。從強(qiáng)度模擬上看，RW試驗和IVAP試驗?zāi)M的中心氣壓強(qiáng)度誤差均小于CTL試驗，雖然其強(qiáng)度數(shù)值上離觀測還有較大差距，但I(xiàn)VAP試驗給出了最接近的結(jié)果。從近中心最大風(fēng)力上來看，CTL試驗在登陸前模擬的風(fēng)場強(qiáng)度明顯偏弱，RW試驗減小了模擬偏差，而IVAP試驗?zāi)M出與實況相當(dāng)?shù)娘L(fēng)速。但是對于臺風(fēng)登陸后IVAP模擬的風(fēng)場卻明顯偏強(qiáng)。整體而言，對臺風(fēng)12 h的路徑和強(qiáng)度預(yù)報表明，雷達(dá)徑向風(fēng)同化通過對臺風(fēng)初始結(jié)構(gòu)和風(fēng)場的改進(jìn)，使預(yù)報的臺風(fēng)路徑誤差比不同化時有明顯減小，而改進(jìn)同化算子后，預(yù)報精度有進(jìn)一步的提高。

圖3　臺風(fēng)路徑和強(qiáng)度預(yù)報對比：（a）三種預(yù)報路徑與觀測路徑；（b）預(yù)報中心最低海平面氣壓誤差；（c）預(yù)報最大風(fēng)速。其中，OBS為觀測，CTL為不同化雷達(dá)資料的比較試驗，RW為GSI原方法同化雷達(dá)資料，IVAP為基于IVAP方法的改進(jìn)同化算子同化試驗Fig. 3 Comparisons of predicted (a) tracks, (b) MSLP, and(c) MWS of Typhoon Fitow to observations with no radar data assimilation (CTL), radar data assimilation with the original method (RW) and the improved method (IVAP)

1.3　數(shù)值模式系統(tǒng)

東亞區(qū)域大氣再分析系統(tǒng)采用Weather Research and Forecasting Model（WRF-ARWv3.8.1）中尺度模式[10]。在模式配置過程中，綜合考慮了東亞區(qū)域的天氣和氣候特征，確定了大氣再分析系統(tǒng)的流程和模式參數(shù)配置。模式區(qū)域的設(shè)置主要考慮以下五個方面：1）模式區(qū)域空間范圍，模式空間區(qū)域需包括全中國的范圍，尤其是南?！熬哦尉€”范圍（圖4）；2）主要天氣系統(tǒng)，利用近10年的ERA-Interim再分析數(shù)據(jù)分析了高中低3個主要層次的天氣形勢[11]。850 hPa考慮了索馬里急流到中南半島的西南季風(fēng)水汽輸送帶、太平洋近赤道東北信風(fēng)水汽輸送高值區(qū)的影響；500 hPa主要捕捉東北亞大槽、北歐高壓脊和烏拉爾山高壓脊的位置演變的影響（圖5）；200 hPa層次包括中緯度高空急流區(qū)的主體位置，側(cè)重西邊界效應(yīng)；3）模式邊界地形，模式邊界應(yīng)處在地形變化緩慢的區(qū)域，同時避免在海陸交界的地區(qū)；4）地圖投影及投影因子，由于區(qū)域范圍較大，在設(shè)置模式中心點時，應(yīng)控制地圖投影因子在1.0左右；5）模式區(qū)域地表特征，下墊面屬性對模式的影響較明顯，范圍的設(shè)置需同時兼顧地表屬性。

圖4　模式區(qū)域設(shè)置（填色表示地形高度）Fig. 4 Spatial coverage of the model domains (shade represents terrain)

根據(jù)本攻關(guān)任務(wù)的需求，結(jié)合東亞區(qū)域的地理特征，模式分辨率設(shè)置為12 km，垂直層次為74層，模式頂高10 hPa。為了使模式能夠更好地分辨復(fù)雜地形，在模式低層進(jìn)行層次加密，即模式的厚度層從低層到高層，厚度逐漸增加?；趯|亞地區(qū)天氣系統(tǒng)的模擬研究的調(diào)研結(jié)果，開展了不同物理過程組合的系列敏感性試驗。根據(jù)試驗結(jié)果，云微物理參數(shù)化方案選取Thompson方案[12]，邊界層參數(shù)化方案采用YSU方案[13]，積云對流參數(shù)化方案采用新Kain-Fritsch方案[14]，陸面參數(shù)化方案為Noah-MP方案[15]，長短波輻射參數(shù)化方案采用rrtmg方案[16]。

圖5　東亞區(qū)域500 hPa層主要天氣系統(tǒng)Fig. 5 Main synoptic systems at 500 hPa level over East Asia

東亞區(qū)域再分析系統(tǒng)設(shè)計中采用了循環(huán)同化，為了確定循環(huán)同化的時間間隔，采用地表氣壓隨的時間變化來確定模式起轉(zhuǎn)（spin-up）時間。圖6顯示的是模式區(qū)域平均地表氣壓隨時間的變化。從圖中可以看出，模式積分300 min后，模式地表氣壓傾向變化趨于穩(wěn)定。因此，本系統(tǒng)中spin-up時間設(shè)為6 h。

東亞區(qū)域大氣再分析系統(tǒng)整體框架如圖7所示，系統(tǒng)用hybrid-3Dvar同化方法。系統(tǒng)采用ERA-Interim再分析資料作為背景場，每天18 UTC時冷啟動一次，spin-up時間為6 h。隨后模式每3 h一次循環(huán)，分別在同化時刻生成再分析數(shù)據(jù)，并輸出每小時一次預(yù)報場。對于模式的集合成員的產(chǎn)生采用兩種方式，其中第一組10個成員疊加隨機(jī)擾動，第二組10個成員采用不同的物理參數(shù)化方案，第三組10個成員采用第二組的參數(shù)化方案并疊加隨機(jī)擾動。

1.4　初步試驗評估

圖6　模式區(qū)域平均地表氣壓傾向隨時間的變化Fig. 6 Evolution of the domain-averaged surface pressure tendency

圖7　東亞區(qū)域再分析系統(tǒng)流程Fig. 7 Flowchart of the reanalysis system

為了測試本再分析系統(tǒng)的性能，首先開展了為期一個月（2014年7月）的控制試驗（不同化任何資料）。利用模式區(qū)域內(nèi)的常規(guī)站點觀測資料對模式結(jié)果進(jìn)行了檢驗。從高空要素（表1）和地面要素（表2）評估結(jié)果來看，控制試驗所得的高空要素的均方根誤差均小于原有ERA-Interim，地面要素的均方根誤差略大于ERA-Interim資料。進(jìn)一步分析地面要素場發(fā)現(xiàn)，地面要素場能夠反映出時間、空間變化的中小尺度特征。此外，對1，4和10月做了同樣的試驗和評估，結(jié)果一致。整體而言，本系統(tǒng)能夠保留原有ERA-Interim大尺度信息的同時，又能夠反映出局地的細(xì)微特征。因此，可將本再分析系統(tǒng)用于今后再分析數(shù)據(jù)的研制。

表1　降尺度試驗高空要素檢驗結(jié)果Table 1 Vertical levels verification of the downscaling experiments

表2　降尺度試驗地面要素檢驗結(jié)果Table 2 Surface verification of the downscaling experiments

2　結(jié)論及展望

本文介紹了中國氣象局“氣象資料質(zhì)量控制多源數(shù)據(jù)融合與再分析”攻關(guān)任務(wù)的第五攻關(guān)方向——東亞區(qū)域大氣再分析項目取得的進(jìn)展。圍繞著東亞區(qū)域大氣再分析需求，收集整理各類觀測資料，基于WRF模式及GSI同化系統(tǒng)建立了東亞區(qū)域大氣再分析系統(tǒng)。該系統(tǒng)由資料預(yù)處理、資料同化、數(shù)值模式和檢驗與評估四部分組成。模式設(shè)置綜合考慮了東亞區(qū)域的地理地貌特點、天氣氣候特征等因素。改進(jìn)了GSI同化系統(tǒng)對雷達(dá)資料徑向風(fēng)的同化算子，使得同化徑向風(fēng)能夠同時引入更加合理的風(fēng)場信息。基于已搭建的東亞區(qū)域大氣再分析系統(tǒng)，開展了個例模擬和批量降尺度試驗。個例試驗結(jié)果表明，該系統(tǒng)能夠較好的再現(xiàn)典型天氣過程特征。降尺度試驗結(jié)果檢驗顯示，高空要素的均方根誤差均小于原有ERA-Interim，地面要素的均方根誤差略大于ERA-Interim資料。整體而言，該系統(tǒng)初步具有在全球再分析的基礎(chǔ)上提高區(qū)域再分析資料性能的能力，可將本系統(tǒng)用于今后再分析數(shù)據(jù)的研制。今后將針對常規(guī)觀測、雷達(dá)、地基GPS、飛機(jī)報等資料的同化及進(jìn)一步優(yōu)化完善，在此基礎(chǔ)上構(gòu)建東亞區(qū)域大氣再分析數(shù)據(jù)集。

[1]Mesinger F, Dimego G, Kalnay E, et al. North American Regional Reanalysis. Bull Amer Meteor Soc, 2006, 3(3): 343-360.

[2]Dahlgren P, Landelius T, K?llberg P et al. A high-resolution regional reanalysis for Europe. Part 1: Three-dimensional reanalysis with the regional HIgh-Resolution Limited-Area Model (HIRLAM). Q J R Meteorol Soc, 2016, 142: 2119–2131.

[3]Bromwich D H, Wilson A B, Bai L, et al. A comparison of the regional arctic system reanalysis and the global ERA-Interim reanalysis for the Arctic. Q J R Meteorol Soc, 2016, 142, 644-658.

[4]Krieger J R, Zhang J, Liu F, et al. Regional reanalysis over the Beaufort and Chukchi Seas. 2012, AGU Fall Meeting.

[5]Fortelius C, Andr? U, Forsblom M. The BALTEX regional reanalysis project. Boreal Environment Research, 2002, 7: 193–201.

[6]Bollmeyer C, Keller J D, Ohlwein C, et al. 2015. Towards a highresolution regional reanalysis for the European cordex domain. Q J R Meteorol Soc, 2015, 141: 1–15.

[7]Shao H, Derber J, Huang X Y, et al. Bridging research to operations transitions: Status and plans of community GSI. Bull Amer Meteor Soc, 2016, 97: 1427-1440.

[8]Liang X D. An integrating velocity–azimuth process single-doppler radar wind retrieval method. J Atmos Oceanic Technol, 2007, 24: 658–665.

[9]Chen F, Liang X D, Ma H. Application of IVAP-based observation operator in radar radial velocity assimilation: The case of typhoon Fitow. Mon Wea Rev, 2017, 145: 4187-4203.

[10]Skamarock W C, Klemp J B, Dudhia J, et al. A description of the advanced research WRF version 3. NCAR Tech. Note NCAR/TN-475+STR, 2008.

[11]Dee D P, Uppala S M, Simmons A J, et al. The ERA-Interim reanalysis: Con figuration and performance of the data assimilation system. Q J R Meteorol Soc, 2011, 137(656): 553-597.

[12]Thompson G, Field P R, Rasmussen R M, et al. Explicit forecasts of winter precipitation using an improved bulk microphysics scheme. Part II: Implementation of a new snow parameterization.Mon Wea Re., 2008, 136(12): 5095-5115.

[13]Hong S-Y, Noh Y, Dudhia J. A new vertical diffusion package with an explicit treatment of entrainment processes. Mon Weather Rev,2006, 134(9): 2318-2341.

[14]Kain J S. The kain-fritsch convective parameterization: An update.J Appl Meteorol, 2004, 43(1): 170-181.

[15]Niu G-Y, Yang Z-L, Mitchell K E, et al. The community Noah land surface model with multiparameterization options (Noah-MP)I: Model description and evaluation with local-scale measurements.J Geophys Res: Atmos, 2011, 116(D12): D12109.

[16]Mlawer E J, Taubman S J, Brown P D, et al. Radiative transfer for inhomogeneous atmospheres: RRTM, a validated correlated-k model for the longwave. J Geophys Res: Atmos, 1997, 102(D14):16663-16682.