基于Landsat-8影像的2018—2021年音苏盖提冰川表面流速监测

王凯，孙永玲，孙世金，刘晓

（山东理工大学建筑工程学院，山东 淄博 255022）

0 引言

山地冰川是全球气候变化的敏感指示器［1-2］，20世纪90年代以来，全球气温升高导致冰川加速消融，冰川表面流速是研究冰川动力学的重要参数，因此监测冰川流速变化变得至关重要［3-4］。但山地冰川地处偏远、自然条件恶劣、人力测量成本昂贵等使实地监测变得困难，遥感技术的发展为冰川监测提供了低成本、大范围、高效率等有利监测条件［5］，特别是光学遥感监测具有视野广阔、获取信息量多、适应性强、可用于动态监测等突出优点，已广泛应用于全球范围内的冰川流速提取［6］。

冰川跃动是一种特殊的冰川运动现象，是冰川不稳定性的一种表现形式［7］。冰川在跃动期一般表现出冰川表面流速比平静期流速增大数十倍甚至数百倍［8］、冰川厚度在积累区减薄，消融区增厚和冰川末端前进等特征［9］。冰川跃动可能阻塞河流而引起水灾或冰川滑入阻塞湖引发冰湖溃决洪水，严重威胁下游区域人们的生命财产安全［10-11］。因此，研究跃动冰川运动特征对预警冰川跃动导致的冰川灾害具有重要的意义［12-13］。Copland等［14］确定音苏盖提冰川为跃动型冰川，并发现该冰川含有大量折叠状的冰碛物、冰川末端前进5 km，该冰川第一次跃 动发 生 在1973年 或 者更 早［15］。Sun等［16］根据Landsat ETM+数据获取音苏盖提冰川1996—2006年的年平均冰川流速，发现音苏盖提冰川在1996—2003年流速保持在20～80 m·a-1，2004—2006年冰川流速整体变大，流速最大约为100 m·a-1。蒋宗立等［17］利用ALOS PALSAR数据估计2007—2009年音苏盖提冰川流速，发现音苏盖提冰川主冰川年平均流速约为100 m·a-1，而南支流比较活跃，冰川流速较大（流速最大约为350 m·a-1）。Rankl等［18］利用ALOS PALSAR和TerraSAR-X数据提取喀喇昆仑区域流速，发现音苏盖提冰川在2003—2009年年平均流速相对较低（流速最大约为109 m·a-1），而2009—2011年年平均流速开始变大，该冰川在2011年发生跃动，最大冰川流速约为547 m·a-1。

本文以中国境内面积最大的音苏盖提冰川为研究区，选取了2018年1月1日—2021年11月16日20对Landsat-8影像，利用光学影像特征追踪方法提取音苏盖提冰川表面流速，评估流速不确定性，并分析音苏盖提冰川流速的时空变化特征，此外，本文还利用ITS_LIVE数据分析了音苏盖提冰川北支流（西）、南支流、北支流（东）的流速时间变化特征及最大流速与其高程分布。

1 研究区域和数据

1.1 研究区域

音苏盖提冰川（35°55′～36°15′N，75°55′～76°21′E），位于中国新疆喀喇昆仑山脉乔戈里峰的北坡（图1），总长约42 km，冰舌长约4 200 m，面积达390 km²，最高海拔约7 050 m，冰川末端海拔约4 000 m，是中国面积、储冰量最大的山谷冰川［19-20］。地形上呈NWW-SEE走向，该冰川是由4条主要支流和10余条小支流汇合而成的树枝状山谷冰川，该区域属于典型大陆型气候，气温较低，降水稀少［21］。

图1 音苏盖提冰川位置（黑色实线为音苏盖提冰川边界，白色实线为冰川中心线，白色虚线是冰川横剖面线，底图是2018年4月7日获取的Landsat-8影像）Fig.1 Location of Skamri Glacier（The black solid line is the border of Skamri Glacier.The white solid line is the central line of Skamri Glacier.The white dotted line is the glacier cross section.The background is Landsat-8 image from April 7，2018）

1.2 数据

本文选取L1T（Level 1T）完全覆盖音苏盖提冰川的Landsat-8 OLI全色波段影像数据，数据全部来源于美国地质调查局（USGS）（https：//earthexplorer.usgs.gov/）。OLI传感器的全色波段数据空间分辨率为15 m，卫星每16天可以实现一次全球覆盖［22］。本文选取云覆盖率小10%的影像，以消除云覆盖对冰川流速提取的影响，共计有20对Landsat-8影像用于提取2018年1月1日到2021年11月16日音苏盖提冰川表面流速（表1）。

表1 覆盖音苏盖提冰川的Landsat-8数据Table 1 Landsat-8 images covering Skamri Glacier

SRTM（Shuttle Radar Topography Mission）是美国奋进号航天飞机利用成像雷达于2000年2月11日至22日成功获取60°N至56°S间80%的陆地面积，对原数据处理后，获取平面精度±20 m，高程精度±16 m的全球数字高程模型（DEM）［23］。公开发布SRTM数字高程产品由美国国家航空航天局（NASA）（https：//earthexplorer.usgs.gov/）提供。由SRTM获取的DEM数据有30 m和90 m两种分辨率的产品［24］。本文选取的是SRTM 30 m DEM数据，用以获取音苏盖提冰川的表面高程。

2 方法

2.1 冰川边界提取

冰川边界是描述冰川位置和规模的重要特征，也是监测冰川变化信息的重要参数。目前冰川边界提取有目视解译和计算机解译两种方法。本文结合了全球冰雪中心（http：//www.glims.org/）发布的RGI 6.0（Randolph Glacier Inventory）［25］冰川边界数据为基础，结合音苏盖提冰川区域2018年1月1日—2021年11月16日的Landsat-8影像数据，采用目视解译方法对冰川边界进行调整，确定音苏盖提冰川边界。

2.2 冰川流速提取

本文利用光学影像特征追踪方法提取音苏盖提冰川表面流速，其方法的基本原理是对于不同时期冰川表面的遥感影像进行相关性分析。由于冰川表面特征点是随冰川一起运动的，计算两景影像特征点移动的位移，即可获取冰川表面流速［26-27］。图2为提取音苏盖提冰川年平均流速的流程图。由于本文选取的L1T数据已经进行过辐射校正和地面控制点几何校正，所以本文中的Landsat-8影像数据均不需要做正射校正，只需做影像配准和互相关计算。本文利用ENVI软件对获取时间较早的Landsat-8数据为参考影像，并以获取时间较晚的影像为匹配影像进行配准，并利用COSI-Corr（Co-registration of Optically Sensed Images and Correlation）软件采用频率域互相关算法对配准后的影像对进行互相关计算，获取音苏盖提冰川东西向和南北向偏移量和信噪比［28-29］。在互相关计算中，首先以初始窗口进行互相关计算，得到像素级的冰川偏移量，再以最终窗口进行互相关计算得到冰川亚像素级偏移量。本文初始搜索窗口设为128×128（像元），最终搜索窗口为64×64（像元），步长设为2（像元），最终获得30 m分辨率的东西向和南北向偏移量。根据光学影像特征追踪方法计算出偏移量的结果存在粗大误差，首先将结果中大于3倍标准差的异常值剔除，此外，利用COSI-Corr软件中的非局部均值滤波对冰川流速结果进行滤波处理，非局部均值滤波对图像中具有相同性质的区域进行分类并加权平均达到降噪效果，本文将搜索区域的窗口大小设为21×21像元，将邻域块大小设为7×7像元，噪声参数H设为1.6倍的标准差。最后根据滤波后的东西向和南北向偏移量以及影像对间隔时间获取音苏盖提冰川年平均流速。

图2 光学影像特征方法提取冰川表面流速流程图Fig.2 Flow chart showing how to derive the glacier surface flow velocity by using optical image feature tracking

由于音苏盖提冰川位于地形复杂的高海拔偏远山区，受到恶劣的自然条件限制，很难获得实测数据来对遥感影像所获取的表面流速进行验证。理论上，进行互相关计算时，冰川稳定区域（冰川内和周边基岩露头）不会产生偏移量，因此，可利用冰川稳定区域的偏移量来评估冰流结果的不确定性。其公式如下：

式中：eoff为冰川稳定区偏移量的误差；MED为冰川稳定区偏移量的平均值；SE为冰川稳定区偏移量的标准误差，其具体计算

式中：STDV为冰川稳定区偏移量平均值的标准差；Neff为空间非相关像元个数，其计算公式如下：

式中：Ntotal为冰川稳定区像元的总数；PS为像元分辨率；D为空间相关的最大距离，空间相关的最大距离一般是像元分辨率的20倍［30］。通过上述方法，本文提取的音苏盖提冰川表面流速不确定性约为12 m·a-1。

3 结果分析

本文利用获取的2018—2021年20对完全覆盖音苏盖提冰川遥感影像数据，采用频率域互相关算法，得到了年平均流速图（图3～5）。从图中看出，冰川流速从冰川上游到冰川下游逐渐减小，冰川中心流速比冰川边缘流速大，符合冰川运动规律。冰川末端由于气温高、海拔低和坡度小的原因导致末端冰消融强烈，几乎监测不到末端流速，所以冰川末端流速接近零［21］。如图所示，2018年1月到2019年6月，音苏盖提冰川南支流（图1的B-B'）年平均流速（最大达到441 m·a-1）大于北支流（西）（图1的AA'）年平均流速（最大达到138 m·a-1），而在2019年6月到2021年11月，北支流（西）年平均流速（最大达到1 008 m·a-1）大于南支流年平均流速（最大达到413 m·a-1）。音苏盖提冰川北支流（东）在2021年7月流速明显增大（流速最大达到280 m·a-1），可能发生跃动。此外，通过流速图发现，主冰川流速没有明显变化。

图3 2018年1月—2019年6月音苏盖提冰川年平均流速（m·a−1），底图为SRTM DEM山影效果图Fig.3 Annual mean velocity（m·a−1）of Skamri Glacier from January 2018 to June 2019 overlaid on the hillshade image of SRTM DEM

为了进一步分析冰川流速的时间变化特征，本文提取了音苏盖提冰川北支流（西）中心线流速随时间的变化（图6）。音苏盖提冰川北支流（西）流速在2018年1月—2019年6月相对比较稳定，流速最大约为138 m·a-1，而从2019年6月，北支流（西）冰川流速突然增大，流速最大达到333 m·a-1，流速持续增大到2020年7月，流速最大达到1 008 m·a-1［距北支流（西）冰川末端（图1的A'）4 500 m处］，是2018年1月份流速的7倍。随后流速开始减小，到2021年11月流速减小到300 m·a-1，但此时流速仍是2018年1月份流速的2倍。根据北支流（西）中心线流速变化，可推断音苏盖提冰川北支流（西）在2019年6月末开始发生跃动，至2021年11月仍处于跃动期。

图6 2018—2021年音苏盖提冰川北支流（西）中心线流速图（中心线位置见图1的A-A'剖面线）Fig.6 Annual velocity profiles along the centerline from terminus to up-glacier of the north（west）Skamri Glacier from 2018 to 2021（The location of the central line is marked in Fig.1 A-A'）

为了分析研究时间段内冰川末端变化，本文利用Landsat-8影像展示了音苏盖提冰川北支流（西）冰川末端变化（图7）。在2018年8月—2020年8月期间，北支流（西）的冰川末端没有明显变化，在此之后，该支冰川末端缓慢前进，并出现了支流末端变宽的现象，至2021年9月，北支流（西）冰川末端较2018年8月前进了约320 m。

图7 2018—2021年音苏盖提冰川北支流（西）冰川末端变化Fig.7 The changes of terminus in the north（west）of Skamri Glacier from 2018 to 2021

图8为音苏盖提冰川南支流的中心线流速随时间变化图。如图所示，2018年1—5月，南支流冰川流速明显增大，最大流速从342 m·a-1增大到441 m·a-1，随后流速开始减小，至2020年1月，流速减小到340 m·a-1，而从2020年1月—2021年2月，冰川流速持续增大，流速最大约为413 m·a-1，之后流速又开始减小，到2021年11月，流速约为299 m·a-1，音苏盖提冰川南支流在研究时段内流速一直比较活跃。

图8 2018—2021年音苏盖提冰川南支流中心线流速图（中心线位置见图1的B-B'剖面线）Fig.8 Annual velocity profiles along the centerline from terminus to up-glacier of the south Skamri Glacier from 2018 to 2021（The location of the central line is marked in Fig.1 B-B'）

为了分析各支流对主冰川流速的影响，本文提取了音苏盖提冰川主冰川横剖面线流速（图9）。通过T1、T2横剖面线流速变化发现，北支流（西）发生跃动之前，主冰川在北支流（西）汇入之前流速没有增大反而减小，最大流速从45 m·a-1减小到14 m·a-1［图9（a）］，而当北支流（西）发生跃动后，主冰川流速整体变大，并且横剖面T2流速（最大流速为100 m·a-1）大于T1的流速（最大流速为76 m·a-1）［图9（c）］，显然，北支流（西）的跃动影响了主冰川的流速；从图9（a）和9（b）中的T3、T4横剖面线流速变化发现，北支流（东）汇入前后对主冰川流速没有明显变化，但2021年10月，北支流（东）汇入主冰川之后，主冰川最大流速从53 m·a-1增大到73 m·a-1，可能跟北支流（东）发生跃动有关；通过T5、T6横剖面线流速变化发现，南支流汇入后主冰川流速增大，由此可推断，音苏盖提冰川主冰川流速受南支流影响，此外，北支流跃动会导致主冰川流速变大。

图9 音苏盖提冰川主冰川横剖面线图，流速为2018年11月—2019年3月（a）、2020年11月—2021年11月（b）、2021年10月—2021年11月（c），T1、T2、T3、T4、T5和T6横剖面线位置见图1Fig.9 The cross-sectional diagram of Skamri main glacier shows the annual average velocity from November 2018 to March 2019（a），from November 2020 to November 2021（b），and from October 2021 to November 2021（c）.The positions of T1，T2，T3，T4，T5 and T6 cross-sectional lines are shown in Fig.1

表2为音苏盖提冰川主冰川和各支冰川在研究时间段内的流速最大值及其高程分布。主冰川在2018年1—3月最大流速高程分布在4 980 m处，2018年3月—2020年1月最大流速高程上升到5 410 m处，随后最大流速高程又降低到4 980 m处，在2021年11月，最大流速高程分布约在5 410 m处；南支流在2018年1月—2021年11月最大流速的高程分布相对比较稳定，高程分布在4 860～4 960 m之间；北支流（西）在2018年1月—2020年7月最大流速的高程从5 058 m逐渐降低到5 001 m，随后又开始上升，在2021年11月，最大流速的高程约为5 106 m；北支流（东）在2018年1月—2021年11月，最大流速的高程分布在5 150～5 300 m之间。

表2 音苏盖提冰川主冰川和各支冰川流速最大值及其高程分布Table 2 Maximum velocity and elevation distribution of Skamri main glacier and its branches

图4 2019年6月—2021年2月音苏盖提冰川年平均流速（m·a−1），底图为SRTM DEM山影效果图Fig.4 Annual mean velocity（m·a−1）of Skamri Glacier from June 2019 to February 2021 overlaid on the hillshade image of SRTM DEM

图5 2021年2月—2021年11月年音苏盖提冰川年平均流速（m·a−1），底图为SRTM DEM山影效果图Fig.5 Annual mean velocity（m·a−1）of Skamri Glacier from February 2021 to November 2021 overlaid on the hillshade image of SRTM DEM

在空间上，主冰川和各支流冰川最大流速高程分布也存在差异。在2018年1—3月和2019年11月—2021年10月期间，北支流（东）最大流速的高程高于主冰川，南支流以及北支流（东）；而在2018年3月—2019年11月和2021年10—11月期间，主冰川最大流速的高程高于其他三条支流冰川。对于该区域流速最大值高程分布的时空差异的原因需要将来更多的研究去探讨。

4 讨论

ITS_LIVE数 据［31］（The inter-mission time series of land ice velocity and elevation）来源于NASA“研究环境制作地球系统数据记录2017”项目（https：//its-live.jpl.nasa.gov/），该数据包含高亚洲地区1985—2018年逐年冰川流速数据，本文选取2013—2018年流速数据来研究音苏盖提冰川表面运动。图10（a）～（f）是基于ITS_LIVE数据获取的2013—2018年的年均流速图，图10（g）～（i）是基于Landsat-8数据提取的2019—2021年的年平均流速图。根据ITS_LIVE流速结果和年平均流速结果发现2013—2016年北支流（西）流速相对稳定，最大流速介于55～74 m·a-1之间，从2017年流速缓慢增大，到2020年流速达到最大（流速约为811 m·a-1）。通过年平均流速变化可发现，北支流（西）在2019年3月—2021年11月期间发生跃动。刘星月等［32］利用2015—2016年Sentinel-1数据提取音苏盖提冰川表面流速，发现北支流（西）上部积累区在2016年5月至6月期间速度突增，最高流速达534 m·a-1，而本文发现音苏盖提冰川北支流（西）在2016年整体年平均流速都较小［图10（d）］，因此推断北支流（西）在2016年5月至6月可能是由于雪崩导致的流速增大。

图10 音苏盖提冰川2013—2018年ITS_LIVE流速（a）～（f），2019年3月—2021年11月年平均流速（g）～（i）Fig.10 ITS_LIVE Velocity（a）～（f）from 2013 to 2018，and mean annual velocity（g）～（i）from March 2019 to November 2021 of Skamri Glacier

此外，音苏盖提冰川主冰川流速在2013年较大，最大流速达到172 m·a-1，随后流速开始减小，一直减小到2015年（最大流速约为133 m·a-1），2015年之后流速没有明显变化，流速介于108～133 m·a-1之间。音苏盖提冰川北支流（东）在2013—2019年3月流速较稳定，最大流速介于33～63 m·a-1之间，从2019年3月流速开始增大，2021年11月流速达到之前的两倍（最大流速为121 m·a-1），因此可判断北支流（东）发生跃动。音苏盖提冰川南支流在2013—2021年期间流速比较活跃，但没有明显变化，最大流速介于349～358 m·a-1之间。

续表2

5 结论

本文基于频率域互相关算法提取了音苏盖提冰川2018年1月1日—2021年11月16日20期冰川表面流速，分析了音苏盖提冰川流速的时空变化特征，根据研究结果得出以下结论：

（1）音苏盖提冰川流速存在空间差异，在2018年1月—2019年6月期间，音苏盖提冰川南支流流速远大于其北支流（西）流速，而2019年6月北支流（西）发生跃动后，改变了空间流速分布，呈现北支流（西）流速大于南支流流速。

（2）根据冰川流速随时间的变化，本文首次通过冰川流速发现北支流（西）和北支流（东）发生跃动，北支流（西）冰川在2019年6月—2021年11月处于跃动期，且该支流冰川末端在2020年8月开始缓慢向主冰川推进，至2021年9月，约向主冰川推进320 m；而北支流（东）在2021年7—11月可能发生跃动。

（3）通过研究时间段内冰川横剖面线流速变化发现，音苏盖提冰川主冰川流速在南支流汇入后明显增大，此外，两条北支流的跃动也会导致主冰川流速增大。

（4）音苏盖提冰川主冰川和各支流冰川最大流速的高程分布在时间上和空间上均存在差异，将来需要更多的研究去探讨时空差异的原因。