1972—2011年東昆侖山木孜塔格峰冰川面積變化與物質平衡遙感監測

蔣宗立，張俊麗，張 震，劉時銀，魏俊鋒，郭萬欽，祝傳廣，黃丹妮

(1.湖南科技大學資源環境與安全工程學院，湘潭 411201;2.安徽理工大學測繪學院，淮南 232001;3.云南大學國際河流與生態安全研究院，昆明 650500;4.中國科學院西北生態環境資源研究院冰凍圈科學國家重點實驗室，蘭州 730000)

0 引言

冰川作為冰凍圈的重要組成部分，其動態變化與海平面變化、氣候變化、水循環等有緊密聯系[1-4]。目前評估冰川變化主要通過監測冰川面積、長度、高程的變化或者估算物質平衡等。冰川面積變化的遙感監測研究已經日漸成熟，冰川高程變化或者物質平衡的監測主要使用大地測量法[5-8]。光學立體像對、歷史地形圖、合成孔徑雷達干涉測量(interferometry synthetic aperture Radar，InSAR)是獲取數字高程模型(digital elevation model ,DEM)的主要手段，能比較精確地用于冰川多年平均物質平衡的估算[9-12]。

木孜塔格峰位于青藏高原東昆侖山西段、阿爾金山自然保護區的西南角，是昆侖山第二大冰川分布區。木孜塔格峰西側冰川位于車臣河流域(流域代碼為5Y62)，東側冰川位于向陽湖和雪景湖等流域(5Z13)，少數冰川位于阿其格庫勒湖和鯨魚湖流域(5Z14)。Wei等[13]根據我國第一、二次冰川編目評估了1970s—2009年間的5Z13流域冰川面積，指出其縮減了0.14±0.07%/a，5Z14流域冰川面積縮減了0.18±0.08%/a；Ye等[14]根據Landsat影像認為1976—2013年間阿雅格庫木庫里湖和可可西里湖流域區(5Z1)冰川面積縮減了0.14%/a，塔里木河流域(5Y6)冰川面積縮減了0.12%/a。不過針對木孜塔格峰地區冰川面積變化缺乏詳細報道。Gardner等[2]利用ICESat GLAS數據評估了東昆侖山和青藏高原內陸地區冰川，得到其在2003—2009年間高程變化為-0.01±0.12 m/a；Brun等[15]評估了2000—2016年間的高亞洲冰川物質平衡，不過2000年以前的冰川物質平衡狀態未知；郭萬欽等[16]通過冰川表面流速與邊界變化分析了2008—2009年間該地區魚鱗川冰川的躍動過程，但沒有對躍動前后物質遷移進行分析。

本文采用Landsat遙感影像、歷史地形圖、航天飛機雷達地形測繪使命(shuttle Radar topography mission，SRTM)DEM和TerrSAR-X/TanDEM-X等資料分析了木孜塔格峰地區冰川面積變化、冰川物質平衡、躍動冰川和前進冰川，為認識該地區冰川變化規律及水資源與水循環提供數據支撐，為氣候變暖情景下冰川變化提供事實依據。

1 研究區概況及數據源

1.1 研究區概況

木孜塔格峰(圖1)位于東昆侖山西段，N36°16′～36°42′，E87°5′～87°39′之間，新疆與西藏的交界處，主峰海拔6 973 m，是東昆侖山最高峰，主脊線呈NE—SW走向。該區域現代冰川極為發育，根據我國第二次冰川編目[17]，區內共發育現代冰川214條，面積共662.2 km2。該地區冰川屬于夏季補給型，夏季受青藏高原熱低壓影響，暖空氣在西風環流和山地地形控制下抬升，產生降水，為冰川提供補給[18]。

圖1 研究區示意圖

1.2 數據源及預處理

本研究使用的主要數據有SRTM-X DEM、地形圖、TerraSAR-X/TanDEM-X(表1)、中國第一次(1972年)、第二次(2010年)冰川編目數據[19]、Landsat遙感影像和氣象數據。研究區冰川平均坡度為12.2°，地形圖冰川區高程精度優于8 m[20]。SRTM DEM數據的高程精度為16 m(95%置信區間)，并隨研究區域的不同而呈現不同精度[21]。考慮到SRTM-X對冰雪的輕微穿透和季節波動，SRTM-X DEM可代表1999年消融期末的表面高程[7-8]。TerrSAR-X/TanDEM-X與SRTM-X均為X波段，在1999—2011年間冰川高程變化和物質平衡中穿透差異影響可以忽略。冰川面積和基于冰川主流線的冰川長度通過冰川編目數據獲取[19]。

表1 TerrSAR-X/TanDEM-X數據Tab.1 TerrSAR-X/TanDEM-X data

此外，本研究使用了中國氣象局(http://data.cma.gov.cn)提供的中國地面氣溫、降水日值/月值 0.5°×0.5°格點數據集。

2 研究方法

2.1 基于差分合成孔徑雷達干涉測量的冰川高程提取

利用SRTM-X作為參考DEM，使用TerraSAR-X/TanDEM-X雙基站合成孔徑雷達數據進行差分干涉，不斷地通過迭代方法去除基線不準確導致的趨勢誤差，最后將優化的差值加在原先已配準的SRTM-X之上，得到TanDEM-X DEM，空間分辨率采樣成30 m，詳細流程如圖2[22]。

圖2 TerrSAR-X/TanDEM-X差分干涉獲取DEM流程Fig.2 Flowchart for DEM generation from TerrSAR-X/TanDEM-X differential Interferometry

2.2 DEM配準

基于地形圖DEM分別對SRTM DEM和TerrSAR-X/TanDEM-X DEM進行配準。其配準原理是根據2期DEM高程偏差dh與坡度α、坡向φ存在三角函數關系[23]，即

(1)

式中a，b，c為系數，可通過非冰川區坡度小于5°的dh/tanα與φ組成的散點圖經過回歸分析得到。為了消除高程偏差異常值的影響，利用高程偏差的5%和95%分位數剔除異常值。

最后利用非冰川區高程差殘差與最大曲率的線性關系，對冰川區的高程差殘差進行校正[24]。

2.3 誤差評估

冰川面積誤差以空間分辨率的一半為距離的緩沖區面積，1∶100 000地形圖相當于空間分辨率為27 m數據[7-8]。冰川面積變化的誤差e的公式為

(2)

式中ε1和ε2分別為2期冰川面積誤差。

冰川長度誤差取空間分辨率的1/2，冰川長度變化誤差也采用式(2)計算，但此時ε1和ε2分別為2期冰川長度誤差。

冰川高程變化的不確定性E用無冰區高程殘差的平均值Em和標準差σ計算[25]，即

(3)

式中N為空間去相關處理后的像元個數。本研究將所有DEM空間分辨率統一為30 m，取空間去相關距離為600 m[25]。校正后的DEM數據誤差特征分布如表2。

表2 校正后DEM數據誤差特征分布Tab.2 Statistics of errors of the adjusted DEMs

2.4 物質平衡評估

本研究基于Huss的研究[26]，采用850±60 kg/m3作為冰川體積—物質轉換的密度，其中60 kg/m3作為冰川物質平衡估算結果誤差進行計算。

3 結果與分析

3.1 冰川面積變化

1972—2011年間木孜塔格峰地區冰川總面積縮減了1.03±2.42%，年均縮減率為0.03±0.06%。其中1972—1999年間年均面積縮減率為0.03±0.09%，1999—2011年間年均面積縮減率為0.02±0.21%。大多數冰川趨于穩定，僅有47條冰川呈現不同程度的退縮(圖3和表3)。冰川面積萎縮率最高的是5Z141E32冰川(冰川面積為0.8±0.1 km2)，年均縮減率為0.21±0.41%；其次是蛇頭川冰川，年均縮減率為0.15±0.13%。2條冰川呈現前進現象，為木孜塔格冰川和淙流冰川，分別前進了224 m和80 m。

(a)1972—1999年 (b)1999—2011年 (c)1972—2011年

圖3 1972—1999年、1999—2011年、1972—2011年木孜塔格峰地區冰川高程變化
Fig.3GlacierelevationchangesofUlughMuztaghduring1972-1999,1999-2011,1972—2011

表3 木孜塔格峰地區冰川面積變化Tab.3 Glacier area changes of Ulugh Muztagh

3.2 冰川物質平衡

研究結果顯示(圖3和表4)，木孜塔格峰地區冰川1972—2011年間減薄2.63±0.22 m，冰川體積減少了1.17±0.15 km3，呈現微弱的負平衡(-0.06±0.01 m w.e./a)趨勢，其中1999年以前為-0.11±0.02 m w.e./a，1999年后幾乎接近于平衡狀態(0.02±0.04 m w.e./a)。該區域冰川物質平衡空間差異性大，不同冰川呈現不同的冰川物質平衡狀態。1972—1999年間大多數冰川呈現不同程度的負平衡狀態，也有少數冰川接近于平衡狀態或者呈現正平衡趨勢，比如冰鱗川冰川和龍頭溝冰川。1999—2011年間近2/3的冰川呈現正平衡趨勢。本研究按照1972—2011年間冰川面積變化將冰川分為3類：穩定型、退縮型和前進型。由于魚鱗川冰川期間發生過躍動，但魚鱗川冰川在1972—2011年間內總體上屬于退縮型，因此在排除了魚鱗川冰川后，退縮型冰川物質平衡為-0.03±0.01 m w.e./a，穩定型冰川物質平衡為-0.03±0.01 m w.e./a，前進型冰川物質平衡為-0.25±0.01 m w.e./a，魚鱗川冰川物質平衡也為-0.25±0.01 m w.e./a。可知前進型冰川的物質虧損大于穩定型和退縮型。從圖3可以看出，木孜塔格冰川和淙流冰川末端均出現明顯的增厚現象，而上游是不斷虧損的，冰川上游的物質通過冰崩、物質重分配轉移到下游，導致下游冰川末端增厚和末端前進。

表4 木孜塔格峰地區冰川物質平衡分布特征Tab.4 Glacier mass changes during the investigated periods

3.3 前進冰川與躍動冰川

通過對比1972—2011年間的Landsat影像可知，木孜塔格冰川和淙流冰川均發生過多次前進。木孜塔格冰川在1972—1986年間前進了100 m，1993年前進了25 m，1994—1999年間前進了33 m，1999—2004年間前進了46 m，2004—2008年間前進了20 m；淙流冰川1992—1993年間前進了50 m，1994—1999年間20 m，1999—2011年間前進了10 m。根據上文分析，冰川不斷前進的過程中，冰川物質呈現虧損狀態，不過1999年以后冰川物質是正積累的，這是降水量增加的結果。

魚鱗川冰川在2007—2011年間發生躍動，北側末端前進距離達548±34 m，從圖3(b)和圖4—5可以看出，冰川躍動是從中部開始的。由于1999—2011年間冰川表面高程變化較小，根據郭萬欽等[16]的分析，冰川躍動主要發生在2007—2011年間，而1999年后冰川表面高程變化不大，因此可將SRTM DEM代表躍動前的冰川表面高程信息，Tan-DEM-X代表躍動后的冰川表面高程信息。結果表明，積蓄區表面高程至少平均下降27 m，最大凹陷為59 m，體積減少了2.7×108m3，接收區表面高程至少平均升高43 m，最大抬升113 m，體積增加了2.1×108m3。大量冰體向下游轉移，冰川末端向前推移。

圖4 躍動冰川表面高程變化Fig.4 Elevation change of surged glacier

圖5 躍動冰川沿主流線高程變化Fig.5 Elevation change of surged glacier along the main flow line

4 討論

4.1 不確定性分析

基于大地測量法的冰川物質平衡估算受DEM精度的影響，不過本研究經配準后無冰區高程差殘差小于1 m，在長時間段內的冰川物質平衡估算時對結果影響較小。本研究中冰川物質平衡估算誤差最高不超過0.04 m w.e./a。SRTM和TanDEM-X波段的數據獲取時間都是2月份，本研究用來代表前一年消融期末的冰川表面高程，期間季節波動引起的高程差與X波段對冰雪的穿透可能存在抵消，降低了穿透引起的不確定性。本研究基于保守假設對物質平衡季節波動影響進行評估，10月至次年2月物質積累假設為2 a的全年降水量轉換(根據中國氣象局格網數據1999年平均值為129 mm；2000年為212 mm；2011年為217 mm；2012年為236 mm)，1972—1999年、1999—2011年、1972—2011年季節波動引起的誤差分別為0.01 m w.e./a，0.03 m w.e./a和0.01 m w.e./a，該誤差包含在不確定性評估結果中，一般可以忽略。

4.2 冰川變化影響因素

木孜塔格峰地區降水較少，影響冰川發育的主要因素是低溫。圖6為基于氣象格網數據的1971—2011年間氣溫和降水變化，其中P1—P4為氣象格網中心點。冰川積累和消融均發生在夏季(6—8月)，由圖6推測，1972—2011年間夏季氣溫升高是冰川退縮的主要原因，其中1971—1999年間氣溫明顯上升(0.17 ℃/a)，降水量雖有波動但總體上呈現微弱上升趨勢(0.2～0.7 mm/a)，由于降水量微弱增加不足以彌補氣溫上升帶來的冰川物質損失，因此期間冰川物質呈現負平衡趨勢。2000—2011年間氣溫呈微弱上升趨勢(0.02～0.04 ℃/a)，降水量上升趨勢比較顯著(3.4～6.3 mm/a)，盡管氣溫仍然處于較高的水平，但降水顯著增加，2005年后波動較小，因此期間冰川物質平衡負平衡趨勢比較微弱，接近于平衡狀態，說明降水的增加在一定程度彌補了冰川物質損失。

(a)氣溫變化

(b)降水量變化

4.3 冰川前進與躍動原因

1999—2011年間2條前進冰川均呈現正積累狀態，但1972—1999年間卻呈現強烈的物質虧損。從圖3可以看出木孜塔格冰川和淙流冰川1972—1999年間上游物質虧損嚴重，因此冰川上游強烈消融，冰川內部液態水的潤滑導致冰川前進，與正常的前進不同。此外，2條前進冰川分別位列同一山峰的南北坡，坡度比較平緩，平均坡度分別為6°和8°，冰體相對寬大，末端相對穩定，為前進創造了條件。

關于躍動冰川的原因解釋目前尚無定論，一般認為是液態水的潤滑作用，液態水來源主要有積累產生的壓融水、冰體內部儲存的融水或者液態降水等[8]。該地區降水較少，現有氣象資料表明2007年前后并未出現強降水現象。王寧練等[27]在昆侖山玉珠峰冰川內部發現富含水冰層，其富含水冰層的形成與氣候變暖相關。此外，郭萬欽等[16]分析表明魚鱗川冰川躍動前冰裂隙發育，冰川躍動主要開始于2008年10月，夏季液態降水可能通過冰裂隙進入冰體從而誘發冰川躍動。

4.4 與其他研究的對比

本研究區冰川面積相對穩定，結果小于其他研究結果[13-14]，主要是因為范圍比較小、時間尺度與其他研究不一致。不過，從鄰近區域的研究(表5)也可以看出，該地區冰川面積相對于其他地區比較穩定。阿爾金山[28]冰川面積萎縮嚴重，其次是各拉丹冬[29]和普若崗日[30]，西昆侖[31]、阿汝錯[8]、崗扎日[32]木孜塔格峰地區冰川面積變化較小，西風環流影響可能從西往東減弱，導致冰川面積減小率也存在從西往東增加的態勢。

表5 木孜塔格峰及毗鄰地區冰川面積變化Tab.5 Glacier area changes in Ulugh Muztagh and surroundings

Brun等[15]的結果顯示木孜塔格峰地區2000—2016年間冰川物質平衡估算為0.06 m w.e./a，與本研究1999—2011年間冰川物質平衡結果(0.02±0.04 m w.e./a)一致。Zhou等[33]利用Hexagon KH-9 和SRTM DEM(C波段)評估1973—2000年間的該地區物質平衡為-0.06±0.12 m w.e./a，與本文1972—1999年間冰川物質平衡結果(-0.11±0.02 m w.e./a)比較接近。兩者差異可能是由于Zhou等[33]使用的C波段SRTM DEM對冰雪穿透的存在一定不確定性及存在一定的數據空值(32%)，而本研究數據相對完整。

對比1999年以后臨近區域的冰川物質平衡結果(表6)，可以看出木孜塔格峰周圍1999年以來大多數冰川區域呈現正物質平衡趨勢。木孜塔格峰西南方向約150 km的聳峙嶺[34]地區冰川正平衡趨勢較為明顯，聳峙嶺往西的阿汝錯[8]、土則崗日[12]、西昆侖[12]及西昆侖西部邊緣地區的瓊木孜塔格[12]等區域冰川物質都呈現比較明顯的正平衡趨勢。受西風環流的影響顯著，降水增多，冰川物質呈現正平衡現象[35]。木孜塔格峰東南方向約200 km的崗扎日[32]地區冰川物質也呈現弱正平衡，但也比木孜塔格峰較為強烈，可能是因為冰川局部地形和氣候存在一定差異。崗扎日往南和往東方向的普若崗日[10]、西各拉丹冬[36]等地區冰川則呈現負的物質平衡趨勢，這些區域受西風環流的影響極小，在氣候變暖的影響下，物質虧損相對嚴重。

表6 木孜塔格峰及毗鄰地區1999年以來的物質平衡Tab.6 Glacier mass changes from 1999 in Ulugh Muztagh and surroundings

5 結論

采用歷史地形圖數據，Landsat，SRTM DEM數據和TerraSAR-X/TanDEM-X數據，對木孜塔格峰地區1972—1999年間與1999—2011年間的冰川面積變化和物質平衡進行了對比研究。結果表明：

1)1972—2011年間木孜塔格峰地區冰川面積年均縮減率為0.02±0.06%。其中1972—1999年間年均縮減率為0.03±0.09%，1999—2011年間年均縮減率為0.02±0.21%。47條冰川表現為退縮，2條冰川表現為前進。

2)1972—2011年間木孜塔格峰地區冰川物質呈現微弱的負平衡(-0.06±0.01 m w.e./a)趨勢，其中1999年以前為-0.11±0.02 m w.e./a，主要受氣溫升高所影響；1999年后幾乎接近于平衡狀態(0.02±0.04 m w.e./a)，因為降水的增加彌補了氣溫升高造成的物質虧損。該區域冰川物質平衡空間差異性大，不同冰川呈現不同的冰川物質平衡狀態。退縮型冰川和穩定型冰川物質平衡均為-0.03±0.01 m w.e./a，前進型冰川物質平衡呈現較為強烈的物質虧損(-0.25±0.01 m w.e./a)。

3)木孜塔格峰地區冰川前進與正常前進冰川不同，其前進可能是冰川內部消融引起的，液態消融水的潤滑導致冰川上游物質向下轉移。冰川躍動的原因也可能與內部消融有關，與液態降水進入冰川也可能相關。