2014—2020年科其喀爾巴西冰川運動速度年際與季節變化特征分析

趙晉彪， 張 震，， 許楊楊， 王榮軍， 蔣宗立

（1.安徽理工大學 空間信息與測繪工程學院，安徽 淮南 232001；2.安徽理工大學 礦山采動災害空天地協同監測與預警安徽普通高校重點實驗室，安徽 淮南 232001；3.安徽理工大學 礦區環境與災害協同監測煤炭行業工程研究中心，安徽 淮南 232001；4.中國科學院 西北生態環境資源研究院 冰凍圈科學國家重點實驗室，甘肅 蘭州 730000；5.湖南科技大學 地球科學與空間信息工程學院，湖南 湘潭 411201）

0 引言

運動是冰川的固有屬性，是冰川區別于其他冰體的主要標志［1］。冰川運動在其動力學系統的控制下，通過改變物質的分布，影響冰川幾何形狀、物質平衡、水熱及邊界的變化［1］。冰川運動速度作為冰川動力學模型的重要參數，不僅能響應物質平衡變化，同時也能反應冰川的動態信息［2-3］。在全球變暖的背景下，研究冰川運動及其變化規律，對區域經濟發展、資源合理利用和自然災害預警有著重要意義。

冰川運動速度的提取方法主要包括實地測量和遙感監測［4-5］。實地測量是在冰面上布設花桿，通過花桿的位置變化測量冰川運動速度。盡管實地測量能獲得精度高、連續性好的數據，但由于山地環境惡劣和成本的限制，無法進行長周期、大范圍的觀測［6-7］。遙感監測具有獲取速度快、監測范圍廣的優勢，可彌補實地測量的不足［8-9］?；谶b感影像的冰川運動速度提取方法主要有雷達差分干涉測量、SAR 偏移量跟蹤和光學影像互相關［10］。由于山地冰川地形復雜且消融較為強烈，獲取高相干影像較為困難，因此很少采用雷達差分干涉測量［11-12］?；赟AR 數據的偏移量跟蹤方法受天氣影響小，但長時間冰川消融導致較大表面差異，較適合短時間尺度的冰川運動速度提?。?，13］。而光學影像互相關方法雖受天氣影響，但歷史資料豐富，已被廣泛應用于年尺度的冰川運動速度研究中［14-15］。

隨著遙感數據時空分辨率的不斷提高，基于遙感數據提取冰川運動速度已成為當前研究的熱點［16-18］。部分學者利用遙感數據對大區域尺度的冰川運動特征進行分析，比如Sam 等［19］基于GoLIVE數據分析了喜馬拉雅山冰川運動特征；Dehecq 等［20］基于Landsat 等數據分析了亞洲高山區冰川運動變化特征。也有部分研究者針對單條冰川進行時空變化特征分析，比如Das 等［21］對米亞冰川1992—2019 年的運動速度時空分布特征及影響因素進行分析；周中正等［22］基于Sentinel-1A 影像分析了2016年崗納樓冰川表面運動速度變化。

科其喀爾巴西冰川表面分布有大量表磧，對其運動速度進行研究，有助于了解表磧覆蓋型冰川對氣候變化的復雜響應。以往的研究更側重于冰川運動速度的空間分布以及區域性變化［23-24］，對于冰川運動速度年內特征及年際變化特征的認識仍然不足。因此，本研究利用Landsat 8影像提取2014—2020 年科其喀爾巴西冰川的運動速度，分析其年內波動特征和年際變化特征，為冰川動力學模型及物質平衡變化提供基礎支撐。

1 研究區概況

科其喀爾巴西冰川（冰川編號：5Y674A0005），位于新疆溫宿縣，屬亞大陸型冰川。冰川表面高程范圍為3 020～6 342 m，冰川總長約25.1 km，面積83.56 km2，冰儲量15.79 km3，物質平衡線海拔高度為4 350 m［25］。該冰川消融期較長［26］，從4 月中下旬開始消融，9 月中下旬消融結束。冰川區受潮濕氣流影響［27］，降水豐富，5—9 月降水量約占全年降水量的75%。據冰川編目資料顯示［25］，該冰川消融區分布有大量表磧，表磧面積為17.5 km2，占整個冰川面積的21%，其表磧厚度隨海拔降低逐漸增加，在冰川末端厚度達到200 cm［28］。

2 數據與方法

2.1 數據

2.1.1 Landsat 8

本研究使用Landsat 8 L1T 影像數據，該級別影像是經過輻射校正和幾何精校正處理的數據。根據影像質量及數據可獲取性，結合表面積雪狀況及研究時段間隔，主要采用冰川表面積雪較少，冰川紋理清晰的數據用于獲取冰川運動速度。共選取2014—2021 年間的16景Landsat 影像（行列號：147，31），用于冰川運動速度提取。并按照影像成像時間間隔組成影像對（表1），時段1～7為年運動速度影像對，時段8～21為季節運動速度影像對。數據獲取網站（https：//earthexplorer.usgs.gov/）。

圖1 科其喀爾巴西冰川位置示意圖Fig.1 Location of Koxkar Baxi Glacier

表1 用于提取科其喀爾巴西冰川運動速度的Landsat 8影像對Table 1 The Landsat 8 image pairs for extracting Koxkar Baxi Glacier flow velocity

2.1.2 ASTER GDEM

ASTER GDEM（ASTER Globe Digital Elevational Model）由美國NASA 和日本METI 共同研制并面向全球免費發布。該產品基于ASTER 近紅外影像數據生成，其水平空間分辨率為30 m［29］，數據覆蓋了全球99%的陸地表面（83° N～83° S）。本研究使用的是第三版數據（ASTER GDEM V3），主要用于獲取研究區高程信息。該產品的發布網站為NASA EARTH DATA（https：//search.earthdata.nasa.gov/search）。

2.1.3 冰川編目數據

冰川邊界數據來源于中國第二次冰川編目數據集（http：//www.ncdc.ac.cn），該數據是以Landsat TM/ETM+和SRTM V4 為數據源，采用高效波段比值閾值分割方法，并結合人工修訂所獲得。

2.1.4 其他數據

氣候數據是來自歐洲中期天氣預報中心（ECMWF）的第五代大氣再分析產品（ERA5），是一種綜合性的再分析數據，本文中利用大氣再分析產品獲取氣溫和降水數據（數據獲取網站為https：//cds.climate.copernicus.eu/#！/home）。

ITS_LIVE（The inter-mission Time Series of Land Ice Velocity and Elevation）數據集［30］提取自Landsat 4、5、7、8 衛星影像，涵蓋了所有面積大于5 km2的陸地冰區，時間跨度為1985—2020年。該數據有120 m和240 m 兩種分辨率，本文使用由光學衛星圖像對生成的120 m 分辨率的速度數據（該數據獲取網站為https：//its-live.jpl.nasa.gov/）。

GoLIVE（Global Land Ice Velocity Extraction from Landsat 8）數據集［31］是將Landsat 8影像在空間域中進行歸一化互相關，進而提取的冰川運動速度，經過誤差調整后的數據精度可達到亞像素級別。該數據集覆蓋了82° N～82° S 之間的所有陸地冰川，空間分辨率為300 m，并提供多種時間跨度的數據（16 d、32 d、48 d、64 d、80 d和96 d）。本文利用多個時間跨度為32 d 的數據進行計算，獲取年尺度冰川運動速度（該數據集獲取網站為https：//nsidc.org/data/）。

Sentinel-1 冰川運動速度數據集［32］是通過強度偏移跟蹤技術，從Sentinel-1 SAR 影像中提取的冰川運動速度，該數據集涵蓋了極地冰蓋以外的12個主要冰川區域，并提供分辨率為200 m 的月度和年度速度數據。該數據集可從埃爾朗根-紐倫堡大學冰川入口網站獲?。╤ttp：//retreat.geographie.uni-erlangen.de.）。

實測數據是通過觀測冰面上布設的花桿所獲得，2015 年5 月至7 月在冰川表面海拔3 100～3 900 m區域布設了40根花桿（共布設40根花桿，其中10、12、39號桿丟失），2015年7月3日、9月28日、9月29日、9月30日分別對每根花桿采用2臺水平精度和垂直精度分別為5 mm+1 ppm 和10 mm+2 ppm的X20單頻靜態接收機進行GPS測量。實測中將一臺GPS 接收機作為基準站（位于科其喀爾巴西冰川大本營內），而另一臺GPS接收機作為流動站在冰川表面進行流動觀測。將觀測數據導入華測GPS數據后處理軟件中，進行求差解算坐標，從而獲取每根花桿的位置信息，即經緯度坐標，最后利用花桿坐標位移量以及測量時間間隔即可計算得到速度信息。

2.2 研究方法

提取冰川運動速度的流程主要包括兩部分，獲取冰川的位移量以及運動速度的計算。本研究使用COSI-Corr（Co-registration of Optically Sensed Images and Correlation）軟件包來獲取位移信息，該軟件包是加州理工學院基于IDL（Interactive Data Language）語言開發的，最初被應用于檢測同震形變［30］，后被廣泛用于處理Landsat、ASTER、SPOT 和Quickbird衛星影像。通過COSI-Corr軟件對影像進行配準和互相關，其結果精度可達到亞像元級別。在積雪區域，該方法仍能獲取精度較高的位移信息，已被證實是基于光學影像監測冰川運動速度的可靠方法［21，33-34］。

本研究使用的Landsat 8 L1T 產品是經過幾何精校正的數據，各景數據間的匹配精度較高，無需再次配準，可直接用于位移信息的提取。將參考影像和搜索影像進行傅里葉變換后，通過設置合適大小的窗口和步長，滑動計算兩景影像上的最大相關系數，進而獲得整個區域在影像成像期間的相對位移［10］。該方法相關參數設置如下，參考窗口設為128，搜索窗口設為32，步長設為8。此過程生成3個圖層：東西向位移（EW）、南北向位移（NS）、信噪比（SNR）。信噪比為從0到1的數值，值越大表示計算出的位移量越可信［2，35］。為保證結果的準確性，選擇SNR≥0.9的部分。受云和陰影的影響，部分值不能準確反映冰川表面的位移信息，對其進行剔除后，即可得到冰川表面位移量。

利用冰川表面位移量計算運動速度Vi：

式中：Di為每對影像之間測量的位移；Ti是兩景影像間的時間間隔。

2.3 精度評估

假設無冰區沒有運動，那么可將基于遙感圖像提取的無冰區運動速度看作誤差值［36］。誤差評估公式為：

式中：E為冰川運動速度誤差；MED為無冰區運動速度的平均值；SE為無冰區運動速度的標準誤差，其具體計算公式為：

式中：STDV為無冰區運動速度的標準差；Neff為去除空間自相關影響的有效像元個數，其計算公式為：

式中：Ntotal為無冰區像元的總數；R為像元分辨率；D為去除空間自相關距離，本研究使用Bolch 等［37］采用的保守估計，取空間去相關距離為像元分辨率的20倍。由此可得冰川運動速度誤差為0.007 m·d-1，遠小于冰川的實際運動速率。無冰區的誤差分布直方圖如圖2所示。

圖2 無冰區殘余運動頻率分布Fig.2 The frequency distribution of residual velocity in non-glacierized regions

3 結果分析

3.1 冰川運動速度空間分布特征

科其喀爾巴西冰川運動速度的空間分布如圖3所示，2014—2020 年科其喀爾巴西冰川運動速度的空間分布較為相似?？破淇柊臀鞅ㄟ\動速度沿中流線分布特征為中上游快、下游慢。冰川末端及下游（0～10 km）冰川運動較為緩慢（0.01 m·d-1）；冰川中上游（10～18 km）運動速度隨著海拔升高逐漸加快并達到最大值，最大速度為0.17～0.20 m·d-1；之后冰川運動速度隨著海拔升高逐漸減小至0.01 m·d-1。以上符合冰川運動速度的分布規律，即從積累區上游至物質平衡線附近冰川運動逐漸加快，物質平衡線附近冰川運動最快，從物質平衡線附近至末端逐漸減小。

圖3 2014—2020年科其喀爾巴西冰川運動速度空間分布Fig.3 Spatial distribution of flow velocity on the Koxkar Baxi Glacier from 2014 to 2020

科其喀爾巴西冰川運動速度沿橫向剖面線的分布情況如圖4（b）所示。冰川橫剖面X1 運動速度明顯呈現“中間快兩邊慢”的特點，中間最大運動速度為0.10～0.15 m·d-1。冰川橫剖面X2 和X3 特征雖不及X1 顯著，但也是中間快，中間最大運動速度分別為0.06～0.08 m·d-1和0.02～0.03 m·d-1，邊緣最小運動速度小于0.01 m·d-1。

圖4 冰川運動速度剖面圖Fig.4 Glacial velocity profile：the velocity along the centreline（a）；transects of surface velocity（b）

3.2 冰川運動速度季節性變化特征

為了探討冰川運動速度的季節性變化，本研究根據氣溫、降水的變化定義了暖季和冷季，將每年5—9 月作為暖季、9 月—次年5 月作為冷季。并提取該冰川沿中流線的運動速度對季節性變化特征進行分析，結果顯示（圖5），該冰川運動速度具有季節性波動特征，暖季運動速度均值為0.07 m·d-1左右，冷季運動速度均值為0.06 m·d-1左右，暖季比冷季快16.67%。

圖5 沿中流線季節運動速度箱線圖Fig.5 Box diagram of seasonal velocity along the center line

3.3 冰川運動速度年際變化特征

2014—2020 年科其喀爾巴西冰川年平均運動速度為0.04～0.05 m·d-1左右，運動速度的峰值均出現在距離冰川末端約18 km處（圖6），最大運動速度在0.17～0.20 m·d-1之間。對科其喀爾巴西冰川中流線上年均運動速度進行統計學分析，結果顯示（圖7），線性回歸模型R2為0.284，顯著性P值為0.218，2014—2020 年該冰川年際運動速度變化特征不顯著，但整體呈微弱減小趨勢，平均運動速度減小約0.01 m·d-1。其中2014—2016 年冰川運動速度呈微弱減慢趨勢（減小約21.3%），2016—2018 年冰川運動速度呈加快趨勢（加快約25.6%），2018—2020 年冰川運動速度呈減小趨勢（減小約27.6%）。

圖6 沿中流線冰川年運動速度Fig.6 Annual glacier velocity along the center line

圖7 2014—2020年均運動速度變化趨勢Fig.7 The variation trend of flow velocity from 2014 to 2020

4 討論

4.1 不確定性分析

由于實測數據有限，基于遙感方法的冰川運動速度誤差很難評估［10，38-39］。一般假設無冰區在研究時間內保持穩定且未發生運動，那么無冰川區域的測量值可以代表冰川運動速度的誤差。結果表明，冰川運動速度誤差值為0.007 m·d-1，誤差結果整體呈正態分布（圖2）。

為了進一步評估本研究的數據質量，以2015年冰川運動速度數據為例，將其與有限的花桿實測數據、ITS_LIVE、GoLIVE 以及Sentinel-1 冰川運動速度數據進行比較（數據獲取日期見表2）。結果表明（圖8），本研究利用COSI-corr 提取的冰川運動速度在總體空間分布上與其他數據產品具有較好的一致性。但由于數據源、數據時間間隔以及提取方法的不同，導致各結果之間存在差異。在冰川中上游（6～18 km）本研究結果與ITS_LIVE 數據最為接近，平均運動速度比ITS_LIVE 數據高1.07%，但在靠近冰川末端處（0～3 km）兩者結果相差較大，可能是因為ITS_LIVE 運動速度是一年內多期數據融合的結果；Sentinel-1 冰川運動速度數據雖然整體較低，但其不確定性為±0.08 m·d-1，在誤差允許范圍內，與本研究結果一致。GoLIVE 數據在冰川末段區域（0～8 km）與本研究的提取結果較為接近，但在冰川消融區中上游的差異較大，可能是由于選用間隔32天的GoLIVE 數據取平均，導致數據誤差疊加，使最終結果產生較大差異?；U實測數據在冰川末端區域（0～8 km）與本研究結果較為接近，但由于實測數據的獲取時間為夏季，較高的溫度和冰川融水導致冰的變形和滑動增大，從而使得末端的一些實測數據高于本研究結果。

圖8 不同方法的冰川運動速度比較Fig.8 Comparison of glacier flow velocity by different methods

表2 研究所使用的主要數據Table 2 List of the data used in this study

4.2 冰川運動速度季節變化特征影響因素

科其喀爾巴西冰川運動速度呈現顯著的季節變化，該變化趨勢與氣溫和降水變化趨勢一致（圖9）。冰川運動是冰的變形、底部滑動和底部巖屑層運動等共同作用的結果，與冰川物質平衡、冰溫、冰川的水力特征等密切相關［40-41］。冰川運動速度的季節性變化通常與表面融化引起的底部潤滑有關［42］，冰川在暖季雖會因消融變薄，但暖季（5—9 月）冰溫高，冰川表面消融強烈，大量冰川融水通過冰川裂隙流入冰川底部，減少冰川底部摩擦力，進而冰川運動速度增大［43-44］。降雨和融雪產生的水也可能與地表融水疊加，導致更大的融水輸入，促使冰川底部滑動，從而影響冰川速度［45-46］。此外，在底部滑動開始時，冰川與底部基巖之間的摩擦產生了大量的熱量，為冰川底部潤滑提供了額外的融水［47］，進一步削弱了冰川的抗滑力，導致冰川速度增加。

圖9 2014—2021年研究區月均氣溫及降水Fig.9 Monthly mean air temperature and precipitation in the study area from 2014 to 2021

4.3 冰川運動速度年際變化特征影響因素

科其喀爾巴西冰川年際運動速度呈現微弱下降趨勢（圖7），該現象可歸因于冰川厚度減薄。近幾十年來全球氣候變暖導致冰川呈現負物質平衡狀態，進而造成冰川運動呈現減速趨勢。聯合國政府間氣候變化專門委員會（IPCC）第六次評估報告顯示，近幾十年全球平均升溫速度為（0.2±0.1）℃·（10a）-1，山地冰川區域升溫速率達到（0.3±0.2）℃·（10a）-1［48］，而氣溫是控制冰川物質損失的主要因素。Brun 等［49］對亞洲高山區冰川物質平衡的研究結果表明，2000—2016 年間天山山區年平均物質平衡為（-0.28±0.20）m w.e·a-1，累計物質平衡始終為負，這與本研究中冰川年際運動特征呈減小趨勢的結論相一致。此外，液態降水也是加速冰川融化的重要驅動因素，當攜帶熱量的雨滴落在冰體時，雨滴可能會侵蝕冰體［42，50］，導致冰川進一步消融。氣溫和液態降水是導致冰川物質平衡改變的重要因素，當它們共同促進冰川物質虧損時，就會使得冰川運動速度減慢。

5 結論

本文利用2014—2021 年間多期Landsat 8 遙感影像，采用頻率域互相關方法提取了科其喀爾巴西冰川表面運動速度，并對冰川運動速度的空間分布、年際變化特征、季節性變化特征以及影響因素進行了分析，結果如下：

（1）科其喀爾巴西冰川運動速度由中流線向邊緣逐漸減小，在消融區隨著海拔的升高而加快，在積累區隨著海拔的升高而降低，物質平衡線附近運動最快（0.17～0.20 m·d-1）。

（2）科其喀爾巴西冰川具有季節性變化特征，較暖季節的冰川運動速度高于較冷季節，暖季運動速度比冷季快16.67%。

（3）2014—2020 年科其喀爾巴西冰川年平均運動速度顯示，該冰川年際變化特征并不顯著，但整體呈現微弱減小趨勢，平均運動速度減少約0.01 m·d-1。

（4）氣溫、降水對冰川運動速度的季節波動及年際變化具有一定的影響。