朱偉剛 馬 東 鄭 宇 楊興海

1 長春工程學(xué)院勘查與測繪工程學(xué)院，長春市同志街3066號，130021 2 西南交通大學(xué)地球科學(xué)與環(huán)境工程學(xué)院，成都市犀安路999號，611756

干旱指數(shù)作為監(jiān)測和評價旱情的重要工具，可以直觀、準(zhǔn)確、定量地表征干旱強度、時空分布等特征[1]。傳統(tǒng)干旱指數(shù)的計算依賴于地面單基站觀測與遙感反演等手段，如以降水與地表徑流為基礎(chǔ)的代表性指標(biāo)Z指數(shù)、標(biāo)準(zhǔn)化降雨指數(shù)SPI(standardized precipitation index)和標(biāo)準(zhǔn)化徑流指數(shù)SRI(standardized runoff index)等[2-3]。但基于單基站觀測數(shù)據(jù)計算的干旱指數(shù)受數(shù)據(jù)質(zhì)量和站點覆蓋面限制，難以監(jiān)測大尺度的干旱事件；基于傳統(tǒng)光學(xué)或微波遙感數(shù)據(jù)反演的干旱指數(shù)僅包含地表水、樹冠水與淺層地下水變化導(dǎo)致的氣象干旱事件，難以表征陸地水總儲量變化引起的水文干旱事件。

GRACE/GRACE-FO系列衛(wèi)星可觀測地球重力場變化，監(jiān)測地表水資源儲量變化情況，估計陸地水儲量盈虧。因此，可利用GRACE衛(wèi)星數(shù)據(jù)反演等效水高并進一步計算干旱指數(shù)，識別并定量表征水文干旱特征[4-5]。然而，由于GRACE實際空間分辨率較低(300 km×300 km)，對小區(qū)域的陸地水儲量變化不夠敏感，不適合小區(qū)域水文干旱監(jiān)測；同時，GRACE任務(wù)結(jié)束和GRACE-FO任務(wù)開始之間存在數(shù)據(jù)空白，會影響對水文干旱的評估效果。因此，需要尋找一種新的干旱監(jiān)測手段作為補充。

GNSS垂向位移時間序列包含環(huán)境負(fù)載(非潮汐大氣負(fù)載、非潮汐海洋負(fù)載、水文負(fù)載等)引起的地球彈性形變信息，并且對小區(qū)域的環(huán)境負(fù)載變化敏感[6]?？鄢浅毕髿庳?fù)載和非潮汐海洋負(fù)載等影響得到僅包含陸地水文負(fù)載信號的GNSS垂向位移時間序列，可以用來反演陸地水儲量變化及監(jiān)測因陸地水儲量變化引起的干旱事件[7]。國內(nèi)學(xué)者在GNSS垂向位移監(jiān)測干旱方面開展大量研究[8-9]。但這些研究僅將GNSS垂向位移與GRACE-DSI、CI(comprehensive meteorological drought index)、PDSI(Palmer drought severity index)等干旱指數(shù)進行對比分析，未將其轉(zhuǎn)化為干旱指數(shù)并提供干旱嚴(yán)重程度的定量表征。

隨著中國大陸構(gòu)造環(huán)境監(jiān)測網(wǎng)絡(luò)的建立，10多年來連續(xù)運行的GNSS基準(zhǔn)站為水文大地測量提供了寶貴數(shù)據(jù)。本文從GNSS原始時間序列中提取僅包含水文負(fù)載信號的垂向位移時間序列，然后計算GNSS垂向位移的標(biāo)準(zhǔn)化異常數(shù)據(jù)GNSS-DSI，并將GNSS-DSI與GRACE-DSI進行對比，最后利用GNSS-DSI分析云南地區(qū)的干旱特征，為干旱監(jiān)測與評估提供新的工具。

1 研究區(qū)域與數(shù)據(jù)來源

1.1 研究區(qū)域

云南省位于我國西南地區(qū)，地勢北高南低，海拔落差懸殊，氣候兼具低緯氣候、季風(fēng)氣候、山原、高原氣候等特點。全省包含長江、珠江、元江、瀾滄江、怒江、大盈江6大水系，年降水量充沛。雖然云南水資源量很大，但在時空分布上極不均衡，并且受大氣環(huán)流異常、水汽輸送減少、青藏高原和社會環(huán)境等因素影響，云南省干旱頻發(fā)[10]。因此，開展干旱監(jiān)測對主動抗旱減災(zāi)、旱情風(fēng)險管理與生態(tài)環(huán)境保護具有重要作用。

1.2 數(shù)據(jù)來源

1.2.1 GNSS坐標(biāo)時間序列

從國家地震科學(xué)數(shù)據(jù)中心獲取云南地區(qū)26個GNSS臺站2011～2019年的時間序列進行分析，站點分布如圖1所示，其結(jié)果由GAMIT/GLOBK軟件解算得到，該產(chǎn)品已在數(shù)據(jù)解算過程中扣除固體潮、海潮、極潮等因素的影響[11]。GNSS時間序列中包含由構(gòu)造運動引起的趨勢性變化、由環(huán)境負(fù)載引起的季節(jié)性變化和地殼運動與儀器更換引起的階躍。其中，季節(jié)性變化主要體現(xiàn)在垂直方向，主要影響因素為大氣、非潮汐海洋和陸地水儲量負(fù)載[12]。因此，本文僅采用GNSS時間序列中的垂直分量進行分析。

圖1 云南地區(qū)GNSS站垂向位移周年振幅及臺站分布Fig.1 Annual vertical displacement amplitude and distribution of GNSS stations in Yunnan

1.2.2 GRACE及降雨數(shù)據(jù)

從美國得克薩斯大學(xué)空間研究中心CSR(http:∥www2.csr.utexas.edu/grace/)下載GRACE/GRACE-FO RL06 V02 Mascon數(shù)據(jù)集。該數(shù)據(jù)產(chǎn)品已考慮橢球效應(yīng)，并應(yīng)用所有適當(dāng)?shù)男拚?，包含C20、C30替換、地心一階項改正、冰川均衡改正(GIA)、加入大氣海洋模型(GAD)、扣除冰后回彈、非潮汐海洋負(fù)載和大氣負(fù)載影響等[13]。本文直接應(yīng)用該數(shù)據(jù)集，以研究區(qū)每個站點為圓心、300 km為半徑的格網(wǎng)值計算站點處等效水高(EWH)時間序列，并將其轉(zhuǎn)換為月尺度GRACE-DSI[14]。除2017-06～2018-06外，缺失的GRACE數(shù)據(jù)采用樣條插值方法進行補充，該方法通過對鄰近點數(shù)值進行三次插值得到。降雨數(shù)據(jù)從全球降水氣候?qū)W中心(http:∥gpcc.dwd.de)獲取，時間分辨率為1個月，空間分辨率為0.25°×0.25°。

2 研究方法

2.1 GNSS水文負(fù)載信號提取

為提取僅包含陸地水文負(fù)載變化引起的地表彈性形變信息，使用德國地學(xué)研究中心GFZ(https:∥www.gfz-potsdam.de/en/home/)提供的環(huán)境荷載產(chǎn)品消除非潮汐海洋負(fù)載(non-tidal ocean loading,NTOL)和非潮汐大氣負(fù)載(non-tidal atmospheric loading,NTAL)效應(yīng)的影響。從圖2可以看出，GNSS垂向位移時間序列受非潮汐大氣負(fù)載效應(yīng)影響較大，周年振幅在2～3 mm之間；受非潮汐海洋負(fù)載效應(yīng)影響較小，周年振幅最大為0.5 mm。最后對垂向位移時間序列進行建模，基于最小二乘法估計所有參數(shù)并去除時間序列中的長期線性趨勢及跳變，以減小構(gòu)造運動(如地震、火山噴發(fā)等)等其他因素的影響，得到以水文負(fù)載信號為主導(dǎo)的垂向位移時間序列[15]。模型公式為：

圖2 GNSS站點NTAL與NTOL垂向位移周年振幅Fig.2 Annual vertical displacement amplitude of NTAL and NTOL at GNSS stations

y(ti)=a+bti+csin(2πti)+dcos(2πti)+

esin(4πti)+fcos(4πti)+

(1)

為保證計算干旱指數(shù)時數(shù)據(jù)的完整性，使用GMIS軟件補齊缺失的GNSS數(shù)據(jù)。該軟件可以很好地恢復(fù)GNSS時間序列中缺失的數(shù)據(jù)，并且減弱時間序列中白噪聲的影響[16]。以YNCX、YNHZ、YNMH站為例，使用GMIS軟件得到濾波與插值補齊后的時間序列(圖3)。圖中黑色為扣除長期趨勢、階躍、非潮汐海洋負(fù)載和非潮汐大氣負(fù)載等影響的GNSS垂向位移時間序列，紅色為濾波與插值補齊后的垂向位移時間序列。

圖3 部分GNSS站濾波與插值補齊后的垂向位移時間序列Fig.3 Vertical displacement time series of some GNSS stations after filtering and interpolation

2.2 水文干旱指數(shù)計算

GRACE對大尺度的陸地水儲量變化十分敏感，該優(yōu)勢促進了其在干旱監(jiān)測領(lǐng)域的發(fā)展。Zhao等[17]通過GRACE數(shù)據(jù)反演等效水高(EWH)，計算EWH的標(biāo)準(zhǔn)化異常數(shù)據(jù)GRACE-DSI，結(jié)果表明，GRACE-DSI可以獨立描述干旱事件。但GRACE存在長期數(shù)據(jù)空白，且受空間分辨率限制，會影響評估精度。受上述干旱評估方法啟發(fā)，本文僅基于GNSS垂向位移數(shù)據(jù)計算干旱強度指數(shù)GNSS-DSI，計算方法與GRACE-DSI相同，且不存在數(shù)據(jù)缺失。

水文負(fù)載與其引起的地球彈性形變呈負(fù)相關(guān)關(guān)系，即水文負(fù)載增加，GNSS臺站垂向位置向下運動；水文負(fù)載減小，GNSS臺站垂向位置向上運動。為方便與GRACE-EWH、GRACE-DSI進行對比分析，本文對所有垂向位移數(shù)據(jù)進行反符號處理，利用處理后的垂向位移數(shù)據(jù)(vertical crustal displacement, VCD)進一步計算GNSS-DSI。GNSS-DSI定義為GNSS-VCD的標(biāo)準(zhǔn)化異常，計算方法是通過月內(nèi)平均將每日VCD時間序列重新采樣得到每月VCD時間序列，減去每月(1～12月)VCD平均值(氣候平均值)得到VCD異常，再除以VCD的標(biāo)準(zhǔn)差。計算公式如下：

(2)

為評估GNSS-DSI的可靠性，結(jié)合GRACE數(shù)據(jù)進行對比，按照上述方法計算GRACE-DSI，并將干旱指數(shù)等級劃分成11個類別(表1)。

表1 干旱指數(shù)等級劃分Tab.1 Classification of drought index

3 結(jié)果與分析

3.1 VCD、EWH與DSI及其相關(guān)性分析

為評價GNSS量化水文干旱程度的能力，本文對GNSS和GRACE觀測的均勻分布在云南不同地理位置(中部、東北部、東部、東南部、南部、西南部、西部、西北部)的8個臺站的VCD、EWH以及進一步計算得到的DSI進行相關(guān)性分析。 圖4為研究區(qū)26個GNSS站點處GNSS-VCD、GRACE-EWH和2種DSI的相關(guān)性系數(shù)以及位于不同相關(guān)性系數(shù)區(qū)間的臺站數(shù)量統(tǒng)計(無相關(guān)[0～0.2]、弱相關(guān)(0.2～0.4]、中等相關(guān)(0.4～0.6]、強相關(guān)(0.6～0.8]、極強相關(guān)(0.8～1])。從圖4(a)可以看出，整個云南地區(qū)的GNSS-VCD與GRACE-EWH表現(xiàn)出高度的一致性，兩者的相關(guān)性系數(shù)在0.58～0.94之間，位于中-極強相關(guān)區(qū)間，除東部地區(qū)存在2個臺站呈現(xiàn)中等相關(guān)外，其他地區(qū)的24個臺站均呈現(xiàn)強-極強相關(guān)。從圖4(b)、4(c)看出，云南地區(qū)整體的GNSS-DSI與GRACE-DSI具有較好的相關(guān)性，除呈無相關(guān)的3個臺站和呈弱相關(guān)的2個臺站外，其余21個臺站均呈中-強相關(guān)。結(jié)合圖4和圖1中水文負(fù)載位移周年振幅可以看出，云南東部、東南部、東北部地區(qū)陸地水儲量變化較弱，位于這3個區(qū)域的GNSS臺站季節(jié)性變化還包含其他地球物理因素[18]，導(dǎo)致這3個區(qū)域的GNSS-DSI與GRACE-DSI相關(guān)性較弱，并存在不穩(wěn)定性。

圖4 GNSS-VCD與GRACE-EWH、 GNSS-DSI與GRACE-DSI相關(guān)系數(shù)和不同相關(guān)系數(shù)區(qū)間的臺站數(shù)量統(tǒng)計Fig.4 Correlation coefficient between GNSS-VCD and GRACE-EWH, GNSS-DSI and GRACE-DSI and the number of stations in different correlation coefficient intervals

圖5為8個區(qū)域的GNSS-VCD、GRACE-EWH和2種干旱指數(shù)的時間序列以及相關(guān)性系數(shù)(藍色和黃色分別表示GRACE/GRACE-FO衛(wèi)星在2016-11～2017-06期間的數(shù)據(jù)質(zhì)量下降以及2017-07～2018-05期間的數(shù)據(jù)缺失)。從圖中可以看出，云南中部、西南部地區(qū)GNSS-VCD與GRACE-EWH、GNSS-DSI與GRACE-DSI具有很強的相關(guān)性，GNSS-VCD與GRACE-EWH的相關(guān)系數(shù)R值分別為0.89、0.92(圖5(a)、5(f))，GNSS-DSI與GRACE-DSI的R值分別為0.78、0.74(圖5(i)、5(n))；陸地水儲量年振幅較小的云南東北部、東部、東南部地區(qū)，GNSS-VCD與GRACE-EWH的R值分別為0.76、0.74、0.68(圖5(b)、5(c)、5(d))，GNSS-DSI與GRACE-DSI的R值分別為0.58、0.58、0.45(圖5(j)、5(k)、5(l))；其他3個地區(qū)的GNSS-VCD與GRACE-EWH均表現(xiàn)出良好的一致性，R值分別為0.83、0.93、0.91(圖5(e)、5(g)、5(h))，GNSS-DSI與GRACE-DSI也具有一定的相關(guān)性，R值分別為0.64、0.65、0.64(圖5(m)、5(o)、5(p))。從DSI變化趨勢來看，8個區(qū)域的2種DSI分別在2012、2013和2019年出現(xiàn)較低值，最低達到-2以下，顯示出3次較為嚴(yán)重的水文干旱事件。

圖5 8個區(qū)域的垂向位移-等效水高和干旱指數(shù)相關(guān)性分析Fig.5 Correlation analysis of vertical displacement-equivalent water height and drought index in 8 regions

由于GNSS臺站在接收信號過程中會產(chǎn)生部分隨機噪聲和虛假的季節(jié)性信號，而GRACE存在長期的數(shù)據(jù)質(zhì)量惡化及數(shù)據(jù)空白，并且GRACE空間分辨率較低(300 km)，對小區(qū)域水文負(fù)載變化不及GNSS-VCD與GNSS-DSI敏感，從而影響其對水文干旱事件的評估效果。因此，上述因素導(dǎo)致2種干旱指數(shù)存在一定差異，結(jié)合2種數(shù)據(jù)更有助于分析水文干旱事件。

3.2 GNSS-DSI表征水文干旱事件

為證明GNSS-DSI表征水文干旱事件的能力，本文調(diào)查云南地區(qū)的水文干旱特征。圖6為8個區(qū)域的GNSS-VCD、GNSS-VCD氣候平均值以及GNSS-VCD異常和GNSS-DSI與降水(PRECI)異常時間序列的比較。表2為云南地區(qū)2011～2019年通過GNSS-DSI探測到的水文干旱事件及其特征。結(jié)合圖6和表2可以看出，GNSS-DSI探測到了云南地區(qū)發(fā)生的持續(xù)時間和嚴(yán)重程度不同的水文干旱事件，其中，2011-08～2012-03、2012-11～2013-04、2019-04～2019-12這3個時期的干旱事件最為顯著，尤其是從2019-04開始持續(xù)9個月的特大干旱(D4)，各區(qū)域的VCD異常與PRECI異常均出現(xiàn)大幅下降，VCD平均虧損值和虧損峰值分別為-4.5 mm、-16.1 mm，DSI值最小為-2.43，期間東北部、東部、東南部區(qū)域從2019-08左右開始，西部、西北部從2019-10左右開始，干旱嚴(yán)重程度逐漸減弱。2011-08～2012-03干旱事件嚴(yán)重程度相對較弱，除2012-02左右DSI出現(xiàn)波谷，VCD虧損峰值最低為-14.0 mm，DSI值最小為-1.75，等級為極端干旱(D3)，其他時期均為中度-嚴(yán)重干旱(D1-D2)。2011-08～2011-11干旱主要集中在東北部、東部、東南部區(qū)域，2011-12左右開始擴散到整個云南地區(qū)，于2012-03逐漸恢復(fù)。曹言等[19]基于TRMM降雨數(shù)據(jù)也探測到本次干旱事件。2012-11～2013-04干旱事件與2011年比較類似，但持續(xù)時間較短，DSI于2013-03左右出現(xiàn)波谷，VCD虧損峰值最低為-15.1 mm，DSI值最小為-1.96，等級為極端干旱(D3)，同時期南部、東南部的干旱強度相對其他區(qū)域較弱。金燕等[20]也識別到本次干旱事件。

表2 云南地區(qū)干旱事件總結(jié)Tab.2 Summary of drought events in Yunnan

除3次比較嚴(yán)重的干旱外，GNSS-DSI還探測到4次持續(xù)時間短、強度較弱的干旱。第1次為2015-06～09左右，主要分布在云南東部、西部、西北部(圖6(c)、6(g)、6(h)、6(k)、6(o)、6(p))，VCD異常均為負(fù)值，DSI在2015-08左右出現(xiàn)波谷，3個區(qū)域的DSI最低值分別為-1.29、-1.4、-1.29，等級為中度-嚴(yán)重干旱(D1-D2)。第2次為2016-07～10左右，主要分布在云南東部和東南部(圖6(c)、6(d)、6(k)、6(l))，本次干旱與2015年相似，DSI在2016-09左右出現(xiàn)波谷，DSI最低值分別為-1.42、-1.11，等級為中度-嚴(yán)重干旱(D1-D2)，通過查閱中國氣象局國家氣候中心干旱事件監(jiān)測信息，證實了本次干旱事件。第3次為2018-03～07左右，干旱范圍較大，除云南西南部以外其他7個區(qū)域均有干旱，DSI在2018-04出現(xiàn)波谷，其中東南部強度最大，VCD平均虧損值和虧損峰值分別為-6.1 mm、-12.7 mm，DSI最低值為-1.81，強度為極端干旱(D3)，其他地區(qū)顯示為中度-嚴(yán)重干旱(D1-D2)，于維等[21]也發(fā)現(xiàn)了本次干旱事件。第4次為2018-10～12左右，主要分布在云南東北部，2018-11 DSI達到低值-1.41，等級為嚴(yán)重干旱(D2)，VCD平均虧損值與虧損峰值分別為-1.1 mm、7.4 mm，同時期東部、東南、西北也出現(xiàn)中度干旱(D1)。Jiang等[22]利用GNSS垂向位移監(jiān)測云南地區(qū)陸地水儲量變化，同樣發(fā)現(xiàn)了本文除2016年以外的6次干旱事件，并指出期間陸地水儲量均呈虧損狀態(tài)，說明該時期個別月份的PRECI異常雖為正值，但少量降水并未起到緩解干旱的作用，導(dǎo)致GNSS-DSI仍呈下降趨勢。

上述研究結(jié)果表明，基于GNSS-VCD的標(biāo)準(zhǔn)化異常GNSS-DSI可成功識別和量化因陸地水儲量虧損引起的水文干旱事件，說明GNSS-DSI可作為水文干旱程度的評價指標(biāo)應(yīng)用于實際水文干旱監(jiān)測。

4 結(jié) 語

本文將GNSS-VCD與GRACE反演的GRACE-EWH進行對比分析，發(fā)現(xiàn)GNSS-VCD可以很好地體現(xiàn)陸地水的季節(jié)性變化?；诖?，本文利用GNSS-VCD進一步計算GNSS-DSI，然后對比分析GNSS-DSI與GRACE-DSI，并利用GNSS-DSI研究云南地區(qū)水文干旱事件，得到以下結(jié)論：

1)云南省整體的GNSS-DSI與GRACE-DSI具有較好的相關(guān)性，共21個臺站呈中強相關(guān)，中部、南部、西南部、西部、西北部5個區(qū)域的GNSS臺站GNSS-DSI與GRACE-DSI相關(guān)性最好。東部、東南部、東北部個別GNSS臺站季節(jié)性變化包含其他地球物理因素，導(dǎo)致這些區(qū)域的GNSS-DSI與GRACE-DSI相關(guān)性較弱。

2)GNSS-DSI探測到云南地區(qū)7次水文干旱事件并可提供定量表征，期間陸地水儲量均呈虧損狀態(tài)，最為嚴(yán)重的是2019年干旱事件，強度大且持續(xù)時間長。7次干旱期間降雨異常多為負(fù)值或呈下降趨勢，個別月份降水異常雖為正值，但由于降水量不足、氣候變化和人類活動等因素影響，陸地水儲量仍處于虧損狀態(tài)，導(dǎo)致GNSS-DSI為負(fù)值。

3)GNSS-DSI可以作為GRACE-DSI的補充，將GNSS-DSI與GRACE-DSI相結(jié)合，可以提高水文干旱監(jiān)測的時空分辨率，更全面地了解干旱特征，改善對水文干旱的預(yù)測、評估與管理。