Sentinel-1卫星是由Sentinel-1A和Sentinel­1B两颗位于同一轨道平面卫星组成的双星座卫星,两颗卫星在轨道平面内的相位差180°^[7]。其中Sentinel-1A在地震形变监测和灾害应急响应中具有重要作用。

Sentinel-1A的所有数据由欧空局免费向所有用户开放,该卫星的Level-1级数据产品包括 SLC(Single Look Complex)和GRDH(Ground Range Detected)数据 ,因其具有地理参照系数､时间参考､辐射系数等而适用于绝大多数用户。 Sentinel-1A卫星的C波段成像系统具有四种成像模式:宽幅模式､波模式､条带模式和超宽幅模式,其中干涉宽幅模式的影像可按沿轨干涉获取方位向形变^[6]。

Sentinel-1A为近极地太阳同步轨道卫星,轨道倾角为98.18°,轨道平面近似南北向。在卫星运动过程中,卫星轨道平面与地球赤道面有两个交点,分别称为升交点和降交点(图1)。由于地球的自转运动和卫星绕地运动,使得在卫星成像范围内,地面上任意区域都会存在卫星两种相向的飞行轨迹,由此形成了卫星的两种成像数据:降轨数据和升轨数据。升轨数据即卫星在由南向北飞行时对地观测所形成的影像,反之成像结果为降轨数据。

图1 卫星运行轨道
Fig.1 Satellite orbit

本文选取了如表1所示的台湾地震前后的 6景SAR影像,包括 3景降轨数据和3景升轨数据,每景SAR影像数据均覆盖以震中(22.94°N,120.54°E)为圆心的 20 km半径内区域,1景震前数据和1景震后数据组成1个干涉对,共得到 4个干涉对 ,用于研究升降轨D-InSAR的差异及其对结果的影响,获取台湾地震地表形变分布。其中 ,T069为升轨数据,T105为降轨数据(如图2(a)所示)。

表1 实验采用的 Sentinel-1A数据信息表
Table 1 Sentinel -1A data information table used in the experiment

台湾省位于中国大陆架的东南延伸处,处在菲律宾海板块与欧亚板块之间的挤压地带 ,由长久的造山运动挤出海面所形成。在台湾南部,欧亚板块正向菲律宾海板块不断靠近 ,两个板块约以每年80 mm的速度相互聚合 ,板块与板块之间的挤压作用使得台湾地区地形南北狭长,并且在岛内的中央山脉､东部的海岸山脉､西部的麓山带以及西部海岸平原区域形成了大量的呈现南北方向的褶皱断层,这些断层大都具有很高的活动性 ,因此造成了台湾地区地震多发的历史事实^[5]。

根据台湾地区的特殊地质构造情况,科学界目前把台湾划分为三个不同的地震构造带 ,分别是西部､东部及东北部地震构造带。有专家指出,本次高雄地震发生在台西部地震构造带,震中位于左镇断层､旗山断层之间的地震密集区(图2(b)),地震密集带和人口密集带叠加,加剧了人员和经济损失^[8]。而旗山断层恰好位于台湾西部地震构造带,其构造类型属于活动断层,此次高雄地震由旗山断层的错动所引发^[9]。

图2 研究区域范围及地质概况
Fig.2 Study area and general geology

InSAR利用雷达回波信号的相位信息来获取地表的高程信息^[10]。InSAR的工作方式主要有三种:交轨干涉测量(Cross-Track Interferometry,XTI)､顺轨干涉测量(Along-Track Interferometry,ATI)､重复轨道干涉测量 (Repeat-Track Interferometry,RTI)^[11]。重复轨道干涉测量又称为双轨干涉测量,即同一颗卫星在不同时刻,两次经过几乎相同的轨道时对地表目标成像,获得两幅具有一定相干性的影像 ,进行差分处理,即影像对应像素复数值共轭相乘,得到干涉图。在实际处理中 ,通常结合已知的 DEM数据 ,消除干涉图中地形起伏因素对相位的影响,从而检测出地表的微小形变。利用双轨D-InSAR进行差分干涉测量处理获取地表形变的流程如下图3所示:

图3 双轨D-InSAR差分干涉处理流程图
Fig.3 Flow chart of dual-track D-InSAR differential interference processing

按照图3所示流程,对 4个干涉对分别进行升降轨 D-InSAR处理 ,主要包括差分干涉图生成､差分干涉图滤波､相位解缠､轨道精炼和重去平､相高转换等。

(1)SAR影像裁剪

实验所用的Sentinel -1A数据是在干涉宽幅成像模式下得到的,数据幅宽 250 km,一景数据覆盖范围很大,为避免失相干区域以及噪声较多的非形变区对处理带来较大误差,需要对影像数据进行裁剪 ,保留震中附近区域 ,便于处理 ,并得到较好的效果。

(2)干涉图生成

将各干涉对主､辅两景影像的同名像素做复共轭相乘 ,并利用外部 DEM数据去除地形相位的影响 ,得到如图4所示的差分干涉图。从初始干涉图中均可以看出地表有一定的形变发生 ,升轨数据中干涉纹图更清晰 ,降轨数据噪声较多且形变不明显。

图4 差分干涉图
Fig.4 Differential interferogram

(3)干涉图滤波噪声。

初始干涉图包含有较多的由大气延迟､轨对滤波后的干涉图进行地理编码 ,得到地道误差等产生的噪声 ,因此要对干涉图进行滤理坐标系下的滤波后差分干涉图5。从图中可以波处理,实验采用Goilden滤波方法进行2次滤看出,升降轨数据均能反应出地表形变 ,在降波处理,在最大限度保留有用信息情况下滤除轨数据中,雷达视线向有两个相邻的形变区。

图5 滤波后差分干涉图(经地理编码后)
Fig.5 Differential interferogram after filtering (after geocoding)

(4)相位解缠

相位解缠是解决相位模糊问题,即在差分干涉图相位差主值的基础上确定整周期的个数,以获得地表形变的真实相位。目前主要的解缠方法有:最小费用流法､区域增长法等。当存在大范围的低相干或是有其他限制增长的原因而致使解缠困难时 ,最小费用流法可以取得较好的结果,因此本实验采用自适应最小费用流法进行相位解缠。

(5)轨道精炼和重去平､相高转换

轨道精炼及重去平即轨道精炼和相位偏移的计算 ,轨道精炼方法采用健壮性好､在基线微小情况也可以使用的多项式优化的方法 ,以此解决多项式的卫星轨道相位偏移的纠正问题,之后将相位值转化为形变值,得到升轨数据和降轨数据在雷达视线向的形变分布图6。

图6 研究区域地表形变分布
Fig.6 Distribution of surface deformation in the study area

(1)从图6(a)和(b)可以看出 ,两个升轨干涉对的 D-InSAR在震中的左侧都监测出了地表抬升 ,其抬升中心一致 ,在 22°56′00″N,120°24′00″E附近,但由于两个干涉对的时间和空间基线的差异,二者监测到的最大抬升量和地表抬升的范围并不相同,第一个干涉对监测到的最大抬升量为16.5 cm,抬升范围为5~16.5 cm,第二个干涉对监测到的最大抬升量为14 cm,抬升范围为4~14 cm。

(2)从图6(c)和(d)可以看出 ,两个降轨干涉对的D-InSAR在震中的左､右两侧分别监测出了地表抬升和下沉 ,抬升中心在 22°57′00″N,120°25′00″E附近 ,下沉中心在 22°57′00″N,120°37′00″E附近,同样,由于两个干涉对的时间和空间基线不同,二者监测到的最大抬升､最大下沉量､抬升或下沉的程度和范围等也不太相同。

(3)在地震发生后 ,地表出现了不同程度的变形,升降轨数据处理得到的结果也存在一定的差异。升轨数据监测的地表形变主要表现为地表的抬升。降轨数据监测的地表形变表现为震中一侧抬升 ,一侧下沉,抬升区地理位置和抬升量与升轨数据较为一致,下沉量最大达到-11 cm。

(4)Sentinel-1A卫星是右侧视几何成像模式(图7),这种模式会导致升降轨数据获取的同一地表影像存在一定的差异。台湾中部为东北-西南走向的中央山脉,平均海拔3000 m,西部地区较为平坦,主要为平原､丘陵地区 ,最大的嘉南平原平均海拔低于100 m,地形起伏较小。卫星升轨时右侧视成像 ,距离向上依次为平原地区､高山地区 ,平原地区起伏小,卫星成像时能够获取视线向地表抬升形变,而沉陷区在抬升区东侧,抬升区对其有一定遮挡,卫星未能探测到 ,因此,在升轨数据获取的形变分布图上,只有明显的抬升形变。卫星降轨成像时,距离向上依次为高山地区､平原地区 ,海拔较高的中央山脉对侧视成像有一定遮挡 ,在山脉西侧存在雷达阴影 ,成像质量较差,因此在降轨数据中,地表形变不明显。沉陷区在抬升区东侧,卫星能够探测到 ,形成两相邻的形变区(图8)。

图7 右侧视几何成像模式
Fig.7 Geometric imaging mode of right view

图8 升降轨飞行示意图
Fig.8 Schematic diagram of flying on ascending and descending orbit