0 引言

山地冰川是冰川的重要组成部分，是反映全球气候变化的敏感的记录器和指示器[1]。冰川运动是冰川的重要特征之一，能引起许多重大的自然灾害[2]。其中，冰川湖溃决及其引起的洪水灾害常给受影响区带来巨大的危害[3-4]。因此，监测山地冰川的运动状况对全球气候变化研究和冰川湖溃决危险性评估十分重要[5-6]。

合成孔径雷达(synthetic aperture Radar, SAR)可以全天时、全天候获取数据，能够长时间、大范围地监测山地冰川区地表环境[7]，相较受到时空限制的传统实地监测[8]和受到多云雨天气限制的光学遥感[9]更具优势，是冰川研究中的重要工具[10]。过去30多a以来，已有学者使用基于SAR影像的差分干涉测量技术(differential interferometry synthetic aperture Radar，DInSAR)、偏移量跟踪技术以及多孔径雷达干涉测量技术(multiple aperture InSAR，MAI)成功地监测了冰川表面运动情况[10-12]。但是，DInSAR技术和MAI技术主要利用SAR相位信息，对相干性要求较高，而冰川流速较快或SAR影像空间基线较大、时间基线较长时往往会使冰川区域失相干，限制了DInSAR和MAI技术的应用[9,13]。基于SAR影像强度信息的偏移量跟踪技术，在SAR影像之间不能保持良好相干性的情况下，仍可以直接获得冰川在方位向和距离向的位移[14]，在监测冰川运动方面比其他2种技术适用范围更大。

山地冰川大多分布于高海拔、地形陡峭的区域，传统的偏移量跟踪技术为了减小地形的影响，常采用空间基线较小的SAR影像数据来进行冰川运动速度的测量[15]，这大大限制了数据选取范围，不利于满足冰川研究需要。Yan等[16]针对地形崎岖地区、大空间基线SAR影像条件下地形引起的偏移量误差，提出一种基于外部数字高程模型(digital elevation model，DEM)的地形影响补偿算法，通过计算出一定基线条件下地形起伏引起的方位向和距离向偏移量，并从总的偏移量中予以去除，从而获得更为精确的结果。基于此算法使用L波段PALSAR数据在西昆仑山东部冰川区[16]、慕士塔格峰冰川区[17]以及塔吉克斯坦南部冰川区[18]获取了精确的冰川表面运动信息，证明了算法的可靠性。本文DEM辅助的偏移量跟踪技术利用DEM和SAR影像轨道信息建立起主从影像的配准查找表，将研究区地形引起的偏移量计算到查找表内，在空间基线较大时能更为精确地配准2幅影像，不仅在监测山地冰川表面运动时有效打破了基线长度的限制，也能应用于同一传感器不同波束模式的SAR影像，进一步扩大了监测冰川运动所使用的SAR数据范围。

本文使用2对不同空间基线的TerraSAR-X数据和一对不同入射角的COSMO-SkyMed数据，以卓莫拉日山系东段冰川分布区为研究区，对TerraSAR-X应用传统的以及DEM辅助的偏移量跟踪技术获取了研究区方位向和距离向的偏移量分布情况，通过对比2种结果在非冰川区的偏移量分布情况以及均方根误差，显示了本文方法的优势； 对COSMO-SkyMed进行数据处理获取了冰川表面运动信息，表明了本文方法可以进一步扩大可用的SAR影像来源。

1 技术方法

1.1 地形引起偏移量分析

偏移量跟踪技术利用归一化互相关算法来获取SAR影像方位向和距离向偏移量[8]，其原理是在主影像中选取一块特定的窗口，并将它和从影像中的搜索窗口进行匹配，并计算每一个位置的互相关系数，当互相关系数值最大时，匹配完成，从而获得匹配位置的总偏移量。

在地形起伏不大、基线较小，且不考虑大气扰动影响的情况下，总偏移量包括冰川位移引起的偏移量和轨道偏移量[19]。传统的偏移量跟踪技术在图像配准阶段基于图像亚像元互相关算法计算2幅影像的同名点，寻找影像的几何变换关系，一般使用最小二乘法计算主从影像的多项式变换模型，实现2幅影像的配准[20]。但在地形变化剧烈以及大空间基线条件下，地形引起的SAR影像之间的局部配准偏移量不可忽视[21]。传统的偏移量跟踪技术忽略了地形引起的偏移量，很难对受此影响的SAR影像进行精确配准，获得准确的运动监测结果。

地形在方位向和距离向引起的偏移量为[22]

(1)

(2)

式中： Atopo和Rtopo分别为地形引起的方位向偏移量和距离向偏移量； α为轨道交角； θ为入射角； Δh为高程差； a和r分别为方位向和距离向的分辨率； B⊥为垂直基线长度； R为斜距长。

根据以上公式可以看出，方位向偏移量与地形起伏、轨道交角以及入射角等有关； 在距离向上地形起伏导致的偏移量大小与垂直基线成正比。为了更直观地了解2分量与地形起伏之间的关系，参照实际情况设定TerraSAR-X数据θ为26°，r为0.9 m，a为2.0 m，R为560 km。则由式(1)可知，若主从影像对应的轨道交角α为0.025°，与平地条件相比，280 m的地形起伏将为同名地物在主从影像方位向上带来1/8个像元的额外偏移； 当然，若主从影像对应的轨道平行，地形起伏将不会引入相应偏移。由式(2)可知，100 m长度垂直基线情况下，当研究区地形起伏达到300 m时，将为同名地物在主从影像距离向上带来1/8个像元的偏移； 1 000 m长度垂直基线情况下，研究区地形起伏达到30 m时，将为同名地物在主从影像距离向上带来1/8个像元的偏移。因此，大空间基线条件下，地形起伏引起的距离向偏移量已经不可忽略，传统偏移量跟踪技术忽略了地形因素，大空间基线时无法对SAR影像精确配准，运动监测结果的精度会受到较大影响。

1.2 DEM辅助的偏移量跟踪技术

DEM辅助的偏移量跟踪技术，核心在于利用研究区的DEM和卫星轨道信息建立2幅SAR影像之间的映射关系，即配准查找表。因为该技术将地形引起的偏移量也计算在查找表内，所以可以获得更为精确的配准[23]。在地形变化剧烈及大空间基线条件下DEM辅助的偏移量跟踪处理流程如图1所示。主要分为3个步骤：

1)初始配准查找表的生成。地理编码获得基于主影像成像几何的SAR坐标系DEM及两者之间的映射关系； 然后建立起SAR坐标系DEM与从影像的映射关系； 结合2种映射关系，生成了记录主从影像偏移量的初始配准查找表，该查找表已将系统偏移与地形引起的偏移量同时计算在内。

2)初始配准查找表的精化。利用重采样获得具有相同成像几何的主、从影像，进一步使用图像互相关方法获取两者的配准多项式，据此精化配准查找表，降低可能存在的由轨道定位和DEM数据不精确带来的误差。

3)位移场的计算。依据精化配准查找表重采样从影像，按照传统偏移量跟踪技术流程计算出冰川方位向和距离向位移场。

此方法能够获得更精确的偏移量结果，对同一传感器不同波束模式的SAR影像也适用。

2 研究区概况与数据源

2.1 研究区概况

研究区位于西藏自治区康马县和浪卡子县之间喜马拉雅山脉东部、卓莫拉日山系东段，N28°07′～28°19′，E89°53′～90°22′之间，地处中国和不丹交界处，地形起伏较大，高程差超过2 000 m，平均海拔在5 000 m以上，包括大小近10条冰川。同时，研究区内有桑旺错和黄湖2个较大的冰川湖，其中桑旺错于1954年发生溃决，受灾人数达2万余人，死亡人数约400人，造成巨大损失[24]，现今仍有重大的成灾风险[25]。因此，监测研究区的冰川运动情况对冰川湖溃决危险评估和当地生产生活安全有重要意义。

2.2 数据来源与处理

TerraSAR-X 雷达卫星具有多极化、多入射角和精确的姿态和轨道控制能力,采用波长约3.2 cm的X波段SAR,可以全天时、全天候地对地观测[26]。COSMO-SkyMed对地观测系统是由4颗搭载了多个成像模式的X波段SAR的太阳同步轨道卫星组成，在轨道上相邻的不同卫星对同一地点的重访时间分别为1 d,3 d,4 d和8 d[27]。

选取2对不同长度空间基线、时间间隔为11 d的TerraSAR-X数据，主要成像参数为升轨右视、条带扫描模式，空间分辨率为3 m，其数据参数如表1所示。选取一对不同波束模式的COSMO-SkyMed影像，为降轨右视获取，空间分辨率为3 m，波束模式分别为HI01和HI02，2景影像在入射角和覆盖范围上均有差异，数据参数如表2所示。DEM为30 m空间分辨率的SRTM，绝对高程精度约为15 m[28]。研究区范围及数据覆盖范围如图2所示。

为了提高影像处理效率，将TerraSAR-X和COSMO-SkyMed影像范围适当裁剪。经过多次试验发现，偏移量跟踪搜索窗口过大时，不利于保存细节信息，而搜索窗口过小则会使冰川的位移情况过于细碎，不利于分析，故本文将偏移量跟踪的搜索窗口设为100像元×100像元，既可以较好地保留细节，又使得位移信息不太零散，同时也能防止遗漏偏移量较大的点。

3 结果分析

3.1 偏移量跟踪技术结果

对TerraSAR-X应用传统及DEM辅助的偏移量跟踪技术获取的非冰川区监测结果如图3和图4所示。图例中不同颜色表示偏移量值的变化，0值代表无偏移量，每一个色周代表3 m的偏移量，当研究区某处偏移量绝对值大于3 m时，偏移量结果图上将显示出色棒中颜色重复出现构成的颜色条纹。黑色区域为相关性阈值低于设定阈值的区域及SAR成像时的叠掩阴影区，即无结果区域。

(a) 传统偏移量跟踪技术方位向偏移量结果(b) DEM辅助偏移量跟踪技术方位向偏移量结果

(c) 传统偏移量跟踪技术距离向偏移量结果(d) DEM辅助偏移量跟踪技术距离向偏移量结果

图3 20160625—20160706像对提取结果

Fig.3 Offset tracking results from image pair 20160625-20160706

(a) 传统偏移量跟踪技术方位向偏移量结果 (b) DEM辅助偏移量跟踪技术方位向偏移量结果

(c) 传统偏移量跟踪技术距离向偏移量结果 (d) DEM辅助偏移量跟踪技术距离向偏移量结果

图4 20160808—20160819像对提取结果

Fig.4 Offset tracking results from image pair 20160808-20160819

从图上可以直观地看到DEM辅助的偏移量跟踪结果在非冰川区接近0值的颜色更加干净统一，表明这种方法所获得的冰川偏移量结果更加可靠。从传统方法在非冰川区方位向的偏移量监测结果来看，随着垂直基线的增加，非冰川区小部分区域误差虽略有增加，但整体仍趋于0值，说明方位向上地形起伏造成的偏移量误差受空间基线变化影响较小； 与之对比，DEM辅助的方法对监测结果有所改善。从传统方法在非冰川区距离向的偏移量监测结果来看，随着垂直基线的增大，非冰川区距离向的偏移量误差更加明显，甚至会超过一个像元大小，说明距离向上地形起伏造成的偏移量误差随着空间基线的增大而变大。与之对比可见，在空间基线较大时，DEM辅助的方法对距离向上非冰川区监测结果的改善效果非常明显。

对不同入射角、覆盖范围的COSMO-SkyMed应用DEM辅助的偏移量跟踪技术的监测结果如图5所示。由图5可知，监测结果在非冰川区接近0值，表明此方法同样适用于不同入射角、覆盖范围的SAR数据。

(a) DEM辅助偏移量跟踪技术方位向偏移量结果 (b) DEM辅助偏移量跟踪技术距离向偏移量结果

图5 20160630-20160731像对提取结果

Fig.5 Offset tracking results from image pair 20160630-20160731

3.2 误差分析

为了评价山地冰川流速估算的精度，对TerraSAR-X和COSMO-SkyMed结果非冰川区分别选取一块区域，计算区域内位移残差的均方根误差(root mean square error,RMSE)来衡量结果的精度。对TerraSAR-X数据处理结果选取的区域范围如图3(a)红色矩形所示，RMSE如表3所示； 对COSMO-SkyMed数据处理结果选取的区域范围如图5(a)红色矩形所示，RMSE如表4所示。

表3 TerraSAR-X数据传统和DEM辅助的

偏移量跟踪技术非冰川区RMSE

Tab.3 RMSE in ice-free region of Traditional and external DEM-assisted offset tracking techniques using TerraSAR-X data

(m)

表4 不同入射角COSMO-SkyMed数据

偏移量跟踪技术非冰川区均RMSE

Tab.4 RMSE in ice-free region of external DEM-
assisted offset tracking techniques using different incidence COSMO-SkyMed data

(m)

通过对比TerraSAR-X 2种方法的非冰川区结果RMSE可知，当SAR影像对垂直基线较小时，2种方法的RMSE都比较小； 并且，基于外部DEM辅助的方法RMSE更低，说明在小空间基线时，DEM辅助的方法也有助于提高结果的精度。当空间基线较大时，DEM辅助的方法在距离向上改善效果明显，方位向上也有提升。由公式(1)和(2)可知，主要是因为距离向上地形引起的偏移量误差受到空间基线的影响较大，方位向上受到空间基线的影响较小。依据表4可知，对不同入射角的COSMO-SkyMed数据应用本文方法时，监测结果在非冰川区的RMSE也很小，能获得精度较高的偏移量结果。以上分析表明，大空间基线下在地形起伏较大区域使用DEM辅助的方法能够有效降低地形引起的偏移量误差，减小空间基线长度的限制，提高偏移量跟踪结果的精度，扩大SAR数据的使用范围。

4 结论

1)偏移量跟踪技术不受冰川区域影像复相干性降低的影响，是一种监测山地冰川表面运动很有效的工具。

2)传统的偏移量跟踪技术忽略了地形起伏对结果的影响，在地形变化剧烈以及影像空间基线较大的情况下距离向结果的精度明显降低。

3)DEM辅助的偏移量跟踪技术在地形陡峭以及SAR数据对空间基线较大时能够获得可靠的结果，减小了地形以及空间基线长度带来的限制，有利于更好开展山地冰川表面运动监测。

4)不同波束模式的SAR影像应用DEM辅助的偏移量跟踪技术时仍能获得较好的冰川表面运动监测结果，进一步提高了SAR影像的利用率。

志谢： 衷心感谢GAMMA遥感公司Urs Wegmuller博士耐心的答疑和帮助。

参考文献(References)：

[1] Zhou J M,Li Z,Guo W Q.Estimation and analysis of the surface velocity field of mountain glaciers in Muztag Ata using satellite SAR data[J].Environmental Earth Sciences,2014,71(8):3581-3592.

[2] 谢自楚,刘潮海.冰川学导论[M].上海:上海科学普及出版社,2010:17-18.

Xie Z C,Liu C H.Introduction to Glacier[M].Shanghai:Shanghai Populao Science Press,2010:17-18.

[3] Bolch T,Buchroithner M F,Peters J,et al.Identification of glacier motion and potentially dangerous glacial lakes in the Mt.Everest region/Nepal using spaceborne imagery[J].Natural Hazards and Earth System Sciences,2008,8(6):1329-1340.

[4] 郭兆成,童立强,周成灿,等.基于遥感图像分析对金错冰川湖溃决泥石流事件的验证[J].国土资源遥感,2016,28(1):152-158.doi:10.6046/gtzyyg.2016.03.15.

Guo Z C,Tong L Q,Zhou C C,et al.Verification of glacial lake outburst debris flow events in Jincuo Lake of Tibet based on remote sensing image analysis[J].Remote Sensing for Land and Resources,2016,28(1):152-158.doi:10.6046/gtzyyg.2016.03.15.

[5] Kääb A.Combination of SRTM3 and repeat ASTER data for deriving alpine glacier flow velocities in the Bhutan Himalaya[J].Remote Sensing of Environment,2005,94(4):463-474.

[6] 邢 宇.青藏高原32年湿地对气候变化的空间响应[J].国土资源遥感,2015,27(3):99-107.doi:10.6046/gtzyyg.2015.03.17.

Xing Y.Spatial responses of wetland change to climate in 32 years in Qinghai-Tibet Plateau[J].Remote Sensing for Land and Resources,2015,27(3):99-107.doi:10.6046/gtzyyg.2015.03.17.

[7] 廖明生,王 腾.时间序列InSAR技术与应用[M].北京:科学出版社,2014:1-2.

Liao M S,Wang T.Time Series InSAR Technology and Application[M].Beijing:Science Press,2014:1-2.

[8] Huang L,Li Z.Comparison of SAR and optical data in deriving glacier velocity with feature tracking[J].International Journal of Remote Sensing,2011,32(10):2681-2698.

[9] Luckman A,Quincey D,Bevan S.The potential of satellite radar interferometry and feature tracking for monitoring flow rates of Himalayan glaciers[J].Remote Sensing of Environment,2007,111(2/3):172-181.

[10] Strozzi T,Luckman A,Murray T,et al.Glacier motion estimation using SAR offset-tracking procedures[J].IEEE Transactions on Geoscience and Remote Sensing,2002,40(11):2384-2391.

[11] Goldstein R M,Engelhardt H,Kamb B,et al.Satellite radar interferometry for monitoring ice sheet motion:Application to an antarctic ice stream[J].Science,1993,262(5139):1525-1530.

[12] Bechor N B D,Zebker H A.Measuring two-dimensional movements using a single InSAR pair[J].Geophysical Research Letters,2006,33(16):L16311.

[13] 王思胜,江利明,孙永玲,等.基于ALOS PALSAR数据的山地冰川流速估算方法比较——以喀喇昆仑地区斯克洋坎力冰川为例[J].国土资源遥感,2016,28(2):54-61.doi:10.6046/gtzyyg.2016.02.09.

Wang S S,Jiang L M,Sun Y L,et al.Evaluation of methods for deriving mountain glacier velocities with ALOS PALSAR images:A case study of Skyang glacier in central Karakoram[J].Remote Sensing for Land and Resources,2016,28(2):54-61.doi:10.6046/gtzyyg.2016.02.09.

[14] Riveros N,Euillades L,Euillades P,et al.Offset tracking procedure applied to high resolution SAR data on Viedma glacier,Patagonian Andes,Argentina[J].Advances in Geosciences,2013,35:7-13.

[15] Strozzi T,Kouraev A,Wiesmann A,et al.Estimation of Arctic glacier motion with satellite L-band SAR data[J].Remote Sensing of Environment,2008,112(3):636-645.

[16] Yan S Y,Liu G,Wang Y J,et al.Accurate determination of glacier surface velocity fields with a DEM-assisted pixel-tracking technique from SAR imagery[J].Remote Sensing,2015,7(8):10898-10916.

[17] Yan S Y,Guo H D,Liu G,et al.Mountain glacier displacement estimation using a DEM-assisted offset tracking method with ALOS/PALSAR data[J].Remote Sensing Letters,2013,4(5):494-503.

[18] Yan S Y,Liu G,Wang Y J,et al.Glacier surface motion pattern in the Eastern part of West Kunlun Shan estimation using pixel-tracking with PALSAR imagery[J].Environmental Earth Sciences,2015,74(3):1871-1881.

[19] 邓方慧,周春霞,王泽民,等.利用偏移量跟踪测定Amery冰架冰流汇合区的冰流速[J].武汉大学学报(信息科学版),2015,40(7):901-906.

Deng F H,Zhou C X,Wang Z M,et al.Ice-flow velocity derivation of the confluence zone of the Amery ice shelf using offset-tracking method[J].Geomatics and Information Science of Wuhan University,2015,40(7):901-906.

[20] 郑 茜,孙建宝,张 永.基于Landsat-8时间序列影像分析西昆仑山地区冰川滑移特征[J].大地测量与地球动力学,2016,36(7):604-608.

Zheng Q,Sun J B,Zhang Y.Fast and uniformly slipping Western-Kunlun glaciers from time-series deformation analysis using periodically captured Landsat-8 imagery[J].Journal of Geodesy and Geodynamics,2016,36(7):604-608.

[21] 闫世勇.山地冰川表面运动雷达遥感监测方法研究[D].北京:中国科学院大学,2013:53-55.

Yan S Y.Research on Extraction of Alpine Glacier Surface Movement by SAR Remote Sensing[D].Beijing:University of Chinese Academy of Sciences,2013:53-55.

[22] Sansosti E,Berardino P,Manunta M,et al.Geometrical SAR image registration[J].IEEE Transactions on Geoscience and Remote Sensing,2006,44(10):2861-2870.

[23] GAMMA Remote Sensing.Differential Interferometry and Geocoding Software-DIFF&GEO(Geocoding and Image Registration)[R].Switzerland:GAMMA Remote Sensing AG,2008.

[24] 戎丹雅,吴秀山,冉启华,等.冰湖溃决洪水演进模拟[J].水土保持,2014,2(3):29-36.

Rong D Y,Wu X S,Ran Q H,et al.Evolution simulation of glacial lake outburst flood[J].Open Journal of Soil and Water Conservation,2014,2(3):29-36.

[25] 刘春玲,童立强,祁生文,等.喜马拉雅山地区冰川湖溃决灾害隐患遥感调查及影响因素分析[J].国土资源遥感,2016,28(3):110-115.doi:10.6046/gtzyyg.2016.03.18.

Liu C L,Tong L Q,Qi S W,et al.Remote sensing investigation and influence factor analysis of glacier lake outburst potential in the Himalayas[J].Remote Sensing for Land and Resources,2016,28(3):110-115.doi:10.6046/gtzyyg.2016.03.18.

[26] 倪维平,边 辉,严卫东,等.TerraSAR-X雷达卫星的系统特性与应用分析[J].雷达科学与技术,2009,7(1):29-34,58.

Ni W P,Bian H,Yan W D,et al.System characteristics and application analysis of TerraSAR-X Radar satellite[J].Radar Science and Technology,2009,7(1):29-34,58.

[27] Italian Space Agency.COSMO-SkyMed Mission and Products Description[EB/OL].[2016-05-31].http://www.e-geos.it/cosmo-skymed.html.

[28] USGS Earth Resource Observation and Science(EROS) Center.Global Data Explorer(GDEx) Data Access User Guide[EB/OL].[2016-04].https://gdex.cr.usgs.gov/gdex/.