利用Cryosat-2数据确定格陵兰冰盖高程和体积变化

陈国栋，梁圣豪，孟子淇，朱家亨

（苏州科技大学 地理科学与测绘工程学院，江苏 苏州 215009）

随着人类社会的发展，人类活动对地球气候产生了巨大影响，其中，海平面上升是最受关注的环境问题之一，严重影响沿海地区的生产生活安全。格陵兰冰盖是世界第二大陆地冰盖，储存了全球十分之一的冰量，完全融化可以使全球海平面上升约6.5 m。进入21世纪以来，格陵兰冰盖消融的速度明显增加，对海平面上升的贡献已经远远超过南极冰盖[1]，因此，对格陵兰冰盖冰雪消融状况进行大范围、长时间的监测是十分必要的。

由于冰盖地理位置和自然环境特殊，一般的地面沉降手段[2]在冰雪覆盖的极地条件下难以实施，因此，极地冰盖高程变化的监测主要通过卫星测高技术实现，Envisat、ICESat等多颗测高卫星的数据都在格陵兰冰盖的监测中贡献了许多重要成果[3]。Cryosat-2是由欧空局（European Space Agency，ESA）研制的雷达测高卫星，于2010年4月8日发射，相比于过去的雷达测高卫星，Cryosat-2采用了新一代的雷达高度计SIRAL，不仅大大提高了沿轨分辨率，还具有三种不同的测量模式，以应对不同特点的地表类型[4]。因此，Cryosat-2卫星测高数据在两极海冰厚度监测、冰盖及冰川质量变化、重力场研究等方面都取得了令人瞩目的成果[5-7]。近年来，国内也有不少学者利用这颗卫星的观测数据开展了一系列的研究[8-11]。

笔者利用2011年1月至2017年12月间Cryosat-2卫星测高数据对格陵兰冰盖高程及体积变化趋势进行了确定，详细分析了冰盖各流域的变化情况，并将结果与2003—2009年间ICESat测高数据获得的冰盖变化情况进行了对比，结果表明格陵兰冰盖高程下降趋势正在加剧。

1 数据介绍及预处理

Cryosat-2卫星是ESA地球探测计划的组成任务之一，主要科学目标是为了监测地球海洋冰盖和大陆冰盖厚度的变化[3]。为了更好地应对不同地表形态，Cryosat-2搭载的新一代雷达高度计SIRAL采用了三种测量模式。其中在格陵兰冰盖表面采用了两种观测模式：在地形平坦的格陵兰岛内陆地区采用低分辨率模式（Low Resolution Mode，LRM）；在地形复杂的沿海区域则采用合成孔径雷达干涉测量模式（Interferometry Synthetic Aperture Radar，InSAR）。获得的数据经过不同类型的处理后，分为Level 1和Level 2两种级别，各级数据按照模式和处理方式及包含信息的不同又分为23类产品。笔者采用的数据为融合处理后的SIR_GDR测高数据，使用的数据时间跨度为2011年1月至2017年12月。该数据在处理过程中已经过了粗差剔除、波形重跟踪处理等各项误差改正，并将三种测量模式的结果以统一的形式进行存储，数据的坐标系统为WGS-84系统，时间系统为国际原子时。

为避免海洋、无冰覆盖的陆地等地表类型上的数据对格陵兰冰盖高程变化的计算产生影响，文中利用1 km分辨率格陵兰岛地表类型格网[12]对数据进行了筛选，只保留冰盖和冰川上的测高数据。此外，SIR_GDR产品中仍然存在一些粗差数据，为保证参与计算的观测值的质量，根据产品中给出的Height_error、SIN_x-track_angle_error和Calibration_Warning三项质量标记将相应的粗差数据进行剔除[13]。完成数据编辑后的格陵兰冰盖表面Cryosat-2测高数据分布如图1所示。

图1 LRM模式（左）和SARin模式（右）数据分布

2 确定冰盖表面高程和体积变化的方法

在利用卫星测高数据对冰盖变化的研究中，常用的方法包括交叉点法和重复轨道法，其中交叉点计算的高程变化精度较高，但交叉点分布稀疏，数据利用率低[14]，而重复轨道法虽然数据利用率高，但在计算时需要考虑地面坡度的影响，结果精度不高[15]。因此，为了能够更好地确定格陵兰冰盖表面高程的变化，文中采用了一种以地面格网为单位，在格网内平差确定地面高程变化的方法，将Cryosat-2数据中的大地坐标以切圆锥投影方式投影到平面，切纬度为71°N，中央经度为39°W，并将格陵兰区域分割成5 km×5 km分辨率的格网。图2为一个格网的示意图，其中的黑点表示Cryosat-2测高数据的平面位置。设格网中心点在起始时刻t0时的高程为h0，而每个Cryosat-2数据点高程hi与格网中心点起始高程h0之差包括两部分：一部分是由地面坡度及平面位置差异引起的空间上的高程差异；另一部分则是由于冰盖消融等原因造成的时间上的高程变化。因此，每个Cryosat-2测高值可表示为以下函数

图2 格网内测高点分布示意图

式中，Δx、Δy分别表示数据点与格网中心点在南北和东西方向上的坐标差，dh/dx和dh/dy分别表示由坡度引起的高程变化沿x轴和y轴方向上的分量，t表示数据采集的时间，dh/dt是高程随时间的线性变化率。

若每个格网内有n个测高点，则每个测高点都可以列出一个形如式（1）的方程，此时将方程改为误差方程，将每个误差方程联立成方程组

由图1可知，每个5 km×5 km格网内的测高数据个数约为数百个甚至上千个，显然式（2）可以通过最小二乘原理获得最或然值，其中dh/dt的估值即为该格网内冰盖表面高程的变化速率。

为了进一步避免观测值粗差对平差结果造成影响，文中对平差后获得的观测值改正数按照三倍中误差的限差进行检验，在剔除对应的观测值后重新进行平差计算，如此迭代，直到所有的改正值都在限差范围内，即可输出最终结果。获得每个格网的高程变化后，将其乘以格网投影面积，即可求得格网内的体积变化量，在一定范围内求和就可以获得相应区域内的冰盖体积变化。

3 结果及分析

文中利用Cryosat-2测高数据计算了2011年1月至2017年12月间格陵兰冰盖表面高程变化，并按照5 km×5 km分辨率进行格网化，结果如图3所示。参考国外相关研究[16]，采用自助抽样法对精度进行估计后，得到整个格陵兰地区冰盖、冰川的平均高程变化为（-12.3±0.38）cm·a－1，体积变化速率为（-220.38±6.0）km3·a－1。

格陵兰东南沿海地区冰川密集，是冰盖消融现象严重的区域[17－18]，但由于该区域地形起伏较大[19]，Cryosat-2数据质量不高，导致这一区域数据稀疏（如图1所示），因此，图3结果在该区域出现较多空白，无法体现当地冰川流失的情况，这一问题目前还没有较好的解决方案[19]，也是今后研究的重点之一。除东南沿海地区以外，文中在格陵兰冰盖的绝大多数区域都获得了有效的冰盖高程变化结果。由图3可知，格陵兰冰盖表面高程呈现上升或基本维持不变的区域主要集中在海拔较高的冰盖内陆地区；而体积流失现象则主要发生在冰盖边缘，特别是西北沿海、东南沿海以及雅各布港等冰川密集区域。

图3 格陵兰冰盖高程变化示意图

根据冰川流动情况，格陵兰冰盖可以划分为8个冰川流域系统，每个流域又分割为若干个子流域[20]（如图3中灰色实线所示），为更好地对格陵兰冰盖进行研究和分析，文中结合格陵兰冰盖冰川流域系统划分数据对各流域的变化情况进行了分析，其结果见表1。由表1可知，除去东南沿海数据缺失较多的DS3和DS4流域以外，在格陵兰岛的西部DS6、DS7、DS8流域是冰盖体积流失较多的地区，也是平均高程变化最大的地区，下降速度超过20 cm·a－1。而在北部地区，高程变化相较DS6、DS7、DS8这三个流域来说，冰盖消融较少，平均高程变化也比较小。

表1 2011—2017年格陵兰岛各流域冰盖体积变化表

文献[17]曾利用ICESat激光测高卫星对2003—2009年间格陵兰冰盖的体积变化进行过详细的分析，将其结果与文中获得的2011—2017年间的体积变化进行对比，结果如图4所示，由于Cryosat-2在东南沿海区域有较多数据缺失，故DS3、DS4、DS5三个流域未进行对比。数据表明，相比于2003—2009年的变化情况，2010年后的DS1流域中冰盖消融速度有所增长，但增长速度并不明显；而在DS6、DS7流域中，冰盖消融速度迅速增加，冰川流失情况更加严重。在DS2流域中，2003—2009年冰盖体积变化为体积增加，而2010年后，冰盖体积则在减少，即表明体积损失的区域进一步扩大，并延伸至气候寒冷的格陵兰北部。

图4 Cryosat-2与ICESat获得的格陵兰冰盖体积变化对比（单位：km3·a－1）

位于DS7.1流域的雅各布港冰川（图3中方框所示区域）是格陵兰体积损失最严重的冰川，因此，对该冰川所在区域进行进一步分析。图5所示为ICESat获得的2003—2009年以及Cryosat-2获得的2011—2017年雅各布港冰川地区高程变化速度，两者之差（Cryosat-2结果减去ICESat结果）如图6所示。由图5可知，雅各布港冰川地区呈明显的高程下降趋势，且下降速度较快，冰川中心区域超过了-5 m·a－1，对比两个时间段的结果可知，高程下降的区域和速度都有所增加。根据图6结果，其中变化趋势为正（即2011—2017年的消融速度相较于2003—2009年的速度有所减缓）的区域仅为3 700 km2，而变化趋势为负（即2011—2017年间消融速度相较于2003—2009年的速度有进一步增加）的区域面积则高达16 300 km2，由此可以明显看出该冰川及周边区域的体积损失在近10年中仍在进一步加剧，其范围也明显扩大。

图5 2003—2009年（左）与2011—2017年（右）雅各布港冰川高程变化（单位：m·a－1）

图6 2003—2009年与2011—2017年雅各布港冰川高程变化差异（单位：m·a－1）

4 结语

基于Cryosat-2卫星测高数据，以5 km×5 km分辨率确定了格陵兰冰盖在2011年1月至2017年12月间的高程和体积变化，对格陵兰冰盖流域进行了详细分析，并与ICESat激光测高卫星获得的2003—2009年间冰盖高程和体积变化进行了对比。结果表明，2011—2017年间，格陵兰冰盖仍处在明显的体积损失状态，其平均高程变化约为（-12.3±0.38）cm·a-1，对应的体积变化则为（-220.38±0.60）km3·a-1。进一步对各流域在不同时期的高程和体积变化速度的对比分析则表明，在格陵兰内部地区，冰川体积的积累速度已经有明显减缓，而在沿海地区，冰川流失现象则进一步加剧。

致谢：感谢欧空局（ESA）提供Cryosat-2数据产品、布里斯托大学Bamber教授提供格陵兰地区1 km分辨率地表类型格网数据及美国国家航空航天局Zwally教授提供格陵兰冰盖冰川流域系统划分。