2014-2018年中巴经济走廊地表变形数据集

白艳萍，康建芳，李萌，韩守富，赵宝强，艾鸣浩，马金辉＊

1. 兰州大学资源环境学院，兰州 730000

2. 中国科学院西北生态环境资源研究院，兰州 730000

数据集（库）基本信息简介表

引 言

中巴经济走廊区域地势北高南低，区域构造活跃、地势险峻、冰川活跃，地质、水文等工程地质条件极为恶劣，公路沿线雪崩、泥石流、坍方、滑坡、岩崩等地质灾害频发，严重影响公路的安全运营[1]。2005年克什米尔7.6级地震导致公路洪扎河（Huza River）以北段的岩体结构严重破坏，这种历时短、规模大的灾害已导致中巴公路的路面与路基严重损毁[2]。类似地质灾害问题在中巴经济走廊区域普遍存在，而且随着未来工程建设的进行和气候环境的变化，地质灾害有可能进一步频发。

地质灾害识别、活动性调查与分析以及危险性评价是走廊内公路、铁路、通信、油气管道等重大工程建设、安全运行与科学减灾的基础。基于此，学术界和工程界采用多种手段开展了该区域地质灾害研究。但由于该区域自然环境复杂、野外实际勘测难度很大，基于先进的卫星遥感技术和时序InSAR技术进行地表变形的反演和地质灾害的早期识别、监测预警，已成为值得重视的研究手段。星载合成孔径雷达干涉测量（Interferometric synthetic aperture radar，InSAR）因其高精度、高分辨率、全天候等优点已迅速成为常用的大地测量技术之一[3-4]，旨在通过对同一区域两次或多次过境的合成孔径雷达（Synthetic Aperture Radar，SAR）影像的复共轭相乘，来提取地物目标的地形或者形变信息[5]。如易靖松通过时序 InSAR技术进行地表形变监测，分析潜在滑坡的时空分布特征[6]。朱建军提出 SAR数据的全球无缝覆盖和高密度重返都会为潜在灾害识别和地质监测提供重要的数据源[5]。Hilleyge利用永久散射体星载合成孔径雷达干涉测量（Persistent scatterer InSAR，PSInSAR）技术对美国Berkley地区的滑坡展开监测，结果较好地反应了滑坡形变的历史演变进程[7]。均可以说明利用时序InSAR变形监测技术，识别并监测潜在地质灾害的可能性。潜在灾害的识别和监测是地质灾害成灾机制研究风险评估和工程防治的基础，也是科学预测自然灾害发展、演化趋势，合理评估潜在自然灾害影响的关键科学问题。

本数据集为中巴经济走廊全部区域时序地表形变数据，通过PSInSAR技术对2014-2018年间覆盖中巴经济走廊全部区域的Sentinel 1A数据进行地表形变监测而获取，为中巴经济走廊的地质灾害识别和经济建设提供参考数据。

1 数据采集和处理方法

1.1 原始数据说明

1.1.1 DEM数据

DEM 数据是进行地形因子计算的基础，本数据集选择 SRTM 90 m分辨率的 DEM 数据（https://dds.cr.usgs.gov/srtm/version2_1/SRTM3/Eurasia/），基础影像共19景，将其进行拼接、裁剪和空值填补后建立了覆盖走廊范围的完整 DEM，高程范围为 986～7770 m，空间参考椭球体为WGS_1984。

1.1.2 InSAR数据

本研究根据SAR数据的分辨率、波长、观测角度、性价比和可获得性等，针对研究区地形和变形特征选用欧空局合成孔径雷达（C波段）Sentinel 1A数据（https://scihub.copernicus.eu/dhus/#/home。）。所涉及轨道为相对轨道107（下降轨道）、27（上升轨道）。覆盖范围如图1所示。时间范围从2014年10月到2018年6月。其中107轨道共400景、覆盖74天，27轨道共340景、覆盖85天，数据基本特征如表1所示。空间覆盖范围如图1所示，时间覆盖范围如图2所示。107轨道和27轨道影像的采样时间分布如图2。虽然107轨道缺少20161014、20161026、20161107、20161119、20161201、20161213、20161225等日期的数据，但是采样时间分布相对均匀，可以满足图像高精度配准要求。

图1 数据集覆盖范围

图2 影像数据时间分布图

表1 Sentinel 1A 数据基本信息

1.2 数据集处理方法及过程

1.2.1 采用原理和技术路线

本文采用PSInSAR技术获得地表形变量，其基本原理是利用多景同一地区的SAR影像，通过统计分析时间序列上幅度和相位信息的稳定性，探测不受时间、空间基线去相关影响的稳定点目标[8]。PSInSAR技术通过对主辅影像进行干涉处理，将获取的强度与相位信息进行统计分析，选取在时间上散射特性相对稳定、回波信号较强的永久散射体（Persistent Scatterer，PS）目标点作为观测对象。这些目标点可能是人工建筑物、裸露的岩石、人工布设的角反射器等。由于它们在时间序列SAR影像中几乎不受斑点噪声的影响，经过长时间间隔仍然保持稳定的散射特性，所以被称作PS。构建PS点网络，网络中任意相邻PS点间经过差分处理的方式分离出形变信息，进而求出相邻点间的高程差值和形变速度差值[9]。以高程差和形变速度差作为原始数据，利用相位解缠算法求解出所有PS点的形变、高程参数，再通过最小二乘法等求得研究区内PS点的形变量和形变速度[10]。

本文采用SARProZ软件——SAR影像处理系统，主要面向雷达影像数据PSInSAR处理技术及相关应用，基于并行多线程处理技术，能够快速处理分析海量雷达影像数据，可识别区分建筑物和地面上的PS点，区分建筑物形变和地表形变；精确反演影像采集时间的大气相位图，剔除大气相位残值；对于城市区域或强反射硬目标形变速率监测精度可达3～5 mm。其数据处理流程如图3。

1.2.2 中巴经济走廊PSInSAR的数据处理过程

（1）中巴经济走廊计算范围设计

计算范围选取107轨道和27轨道覆盖的中国境内的经济走廊区域。具体计算时，采用分幅计算的方式，各分幅计算范围如图4所示。本文仅展示境内重点区域的计算结果。

境内重点区域的覆盖范围为图4中紫色区域所示，涵盖了慕士塔格山和塔什库尔干塔吉克自治县两个典型沉降区，其四角点经纬度分别为：（38°42′12″N，75°03′41″E）、（38°34′16″N，75°48′54″E）、（37°26′35″N，75°29′34″E）、（37°31′22″N，74°46′01″E）。面积约为 8000 km2。

图3 PSInSAR数据处理流程图

图4 样本数据覆盖范围

（2）InSAR主影像选择和影像连接拓扑图构建

针对研究区Sentinel 1A数据集，获取所有影像对的时间、空间基线，将时间、空间基线之和为最小的影像作为主影像，构建时间和空间影像连接拓扑图（图5），其中主影像日期为2016年12月31日。

图5 PSInSAR影像连接拓扑图

（3）影像配准

配准过程中先采用精确轨道参数进行粗配准，为精配准提供初始偏移量，以加快精配准速度。然后参考90 m分辨率的SRTM数据，选择整景影像上均匀分布的同名像素点建立多项式方程进行精配准，配准精度小于1/8个像元。

（4）干涉相位计算

采用Goodstein mode滤波器，滤波窗口为15×15。对已经过滤波的主辅影像像元对生成干涉相位值，逐像元计算干涉相位，且采用DEM去除地形相位，生成时间序列干涉图集。干涉图用于评价InSAR的计算结果，干涉环越明显，形变效果越好。根据图6所示干涉图（共85幅干涉图，本文选择部分干涉图进行展示），可以看出残留地形相位不明显，且存在明显干涉环，表明各时期影像与主影像相干性较好，干涉结果可以用于后继PSInSAR计算。

图6 干涉处理结果

（5）干涉相位解缠

由于干涉图整体相干性较高时，采用枝切法进行相位解缠，解缠效果良好。

（6）候选PS点选择、空间网络构建、大气效应剔除

在进行候选PS点计算时，采用Amplitude Stable. index进行PS点选择，阈值设置为0.65，线性趋势估值设置为-20和20，高程的估值参数设置为-50和50，多项式项数设置为1，计算过程中将相干性（Coherence）作为权重，以此估算大气效应，并反演得到地表形变速度和各时期形变量，最终PS点的相干性大于0.7，符合国际PS点选择的惯例。

2 数据样本描述

本数据集采用PSInSAR技术获取中巴经济走廊地表形变数据，本文选择两个典型的境内重点区域对数据集样本进行描述。第一个典型重点区域位于塔什库尔干塔吉克自治县，以东经 75°17′，北纬37°39′为区域中心，沿轨道方向长为98公里、宽42公里的长方形区域。图7所示为27轨道数据计算结果，图8所示为107轨道数据计算结果。对比两条轨道数据的计算结果可以看出，库孜滚村的东南方向2.2公里的位置均发生明显的下沉，下沉位置完全吻合。其余地方没有明显变形。

图7 PSInSAR计算结果样本数据（27轨道）

图8 PSInSAR计算结果样本数据（107轨道）

第二个典型重点区域位于克孜勒苏柯尔克孜自治州阿克陶县，有着世界冰山之父之名的慕士塔格山，以东经75°08′，北纬38°13′为区域中心。本文选用27轨道的计算结果，图9所示的慕士塔格山的冰川下缘，沿着冰川触水径流方向有明显的下沉现象，其在山脚下西南方向也有局部明显沉降。

3 数据质量控制

本数据集基础数据选用欧空局合成孔径雷达（C波段）Sentinel 1A数据，它具有良好的基线控制和短暂的重返周期，能够在短暂时间积累大量影像数据，有利于InSAR时序分析，并且提升雷达图像参数反演的精度。在PSInSAR数据处理过程中，通过影像配准、干涉成像、相位解缠、高程解算和地理编码等关键步骤，有效地去除地平、高程相位以及大气效应，获得了理论精度毫米级的形变信息。并且文中采用了107轨道和27轨道两个轨道，分别对研究区进行计算，不同轨道计算结果吻合较好。

4 数据使用方法和建议

使用数据时，可利用ArcGIS地理信息系统软件对shape文件进行查询、编辑等操作，其shape中的属性表可作为分析与统计地表形变的基本数据。用Google Earth可对KML文件进行查询、分析等操作。本数据集可靠性高，代表性较强，可以结合地质岩性结构数据，用于中巴经济走廊地质变化研究，满足用户对中巴经济走廊自然灾害科学数据的需求。