多源遥感技术的滑坡地灾隐患识别，以常澧地区为例

我国地理地质条件复杂,构造活动频繁,地质灾害频发,其中滑坡是造成人员伤亡和经济损失最主要的地质灾害类型之一,且往往数量多、规模大,对居民安全、区域经济发展等都会造成严重影响。特别的,相较于北方地理地质背景,南方丘陵-山地区地理地质环境更加复杂,浅层小型滑坡发育,具有点多、面广、零散、规模小、隐蔽、频率高等典型特征,防治难度更大。虽然当前我国地质灾害详查、排查方面的工作布置已经十分全面,但地质灾害一般发生在气候及地理地质环境条件复杂地区,现场调查较为困难,调查人员往往难以抵达,难以满足时效性要求,成本高、效率低。鉴于此,如何对滑坡地灾隐患进行早期识别和实时监测,以提前预防和统筹部署防治工作,仍然是当前地质灾害防治研究领域的难题。

基于航空摄影,现代遥感技术在20世纪60年代后逐渐发展起来,为灾害识别领域提供了良好的技术支持,具有规模大、时效性强、速度快等特点,且很大程度上节省了人力和物力。目视解译是最早的灾害信息提取手段,在21世纪初,3S技术在地质灾害领域得到了广泛应用,使滑坡识别研究得到充分开发,逐渐趋于多元化、精准化。特别是王冶华等学者提出的“数字滑坡”概念,根本性改变了传统滑坡遥感调查方法从获取到显示滑坡信息的方式方法。

近年来,计算机技术和测绘科学技术发展迅速,光学遥感、合成孔径雷达干涉测量(InSAR)、激光雷达探测等技术为地质灾害的识别和监测提供了海量的测绘数据。随着遥感技术的持续优化,信息量指数式增长、分辨率更高的遥感图像数据的不断加入,各种融合技术、遥感分类技术开始广泛应用于地质灾害领域当中,逐渐成为较目视解译识别能力更强、识别速度更快、精度更高、智能化程度更高、对人力和物质条件依赖性更低的解译方法。本研究基于InSAR形变速率数据、多光谱遥感图像和DEM数据对湖南常澧地区的滑坡地灾隐患进行实验性综合识别和提取,以期为滑坡灾害防治提供服务。

研究区位于湖南省常德市境内,111°20'~111°50'E,29°00'~29°40'N之间。所处区域属于洞庭湖四水(湘江、沅水、澧水、资江)流域范围,主要位于沅水、澧水流经区域(图1(a))。

研究区大地构造位于扬子板块,其东北部是第四纪形成的澧县坳陷,西南部是雪峰—武陵隆起,以燕山和喜山运动影响最明显^[12](图1)。区内地层主要有中元古界冷家溪群、新元古界板溪群以及震旦系至第四系发育,其中白垩-第四系地层发育最多,地表以第四系地层分布最为广泛。地层岩性主要由沉积岩和变质岩组成,其中沉积岩主要包括碎屑岩、碳酸盐岩、页岩和砂岩等; 变质岩主要为前寒武系板岩,在沅江流域出露明显。根据构造区划和地貌形态,研究区可划分为2部分: 约四分之三区域部分属于东部洞庭湖凹陷盆地,广泛分布中新生代沉积地层; 约四分之一西部区域以燕山期褶皱和断层为主,褶皱轴向和断层走向主要为NE和NNE,控制着山体和水系的延伸方向。据以往基础地质调查数据和前人研究成果,研究区滑坡地质灾害主要发生在较坚硬薄层-层块状碳酸盐岩和软硬相间薄层-层块状砂砾岩与泥(页)岩互层的地质环境,多为顺层滑坡,分析其工程力学原因,可能是这些岩组层面间夹泥质较多,遇水极易形成软弱夹层,岩土体力学指标急剧下降,易造成滑坡失稳滑动。地质构造对研究区滑坡的形成发育主要有2方面的影响作用: 一是控制着山地-丘陵地形地貌的形成发育,在构造抬升区往往山势险峻,沟谷深切,临空面发育,易于发生滑坡、泥石流等灾害; 二是破坏了岩土体结构的整体稳定性,尤其在褶皱轴部、转折端,断裂带及其两侧,风化层厚,岩石破碎,裂隙发育,易发生滑坡。

研究区主要为常德市中部的红色丘陵地貌,受沅水、澧水等地上水以及地下水、雨水等的冲刷侵蚀明显; 丘陵起伏,相对高度一般不超过200 m,坡度大多小于25°。中生代时期形成一些山间盆地,底部地层多为白垩-第三系红色岩层; 之后在晚第三纪受喜马拉雅运动影响,地壳有所回升,褶皱减轻,单斜构造普遍,侵蚀基准面相对降低,下切侵蚀加强,山间盆地底部的红色岩层遭到冲刷和侵蚀,使得原始的平坦地形呈现丘陵起伏地貌。研究区内还分布有红色丘岗盆地,以宽河谷和残丘为主,海拔多在200~500 m,丘陵的相对高度一般小于100 m,流域两岸仍可见河流冲积平原,居民点和耕地分布密集。地形地貌特征是影响地灾发育的主控因素之一。不同的地貌类型控制着地质灾害的发生类型: 高大的自由临空面是滑坡产生的有利地形,残坡积层较发育,软硬岩层相间的山区和丘陵区是滑坡多发区; 地形起伏大、沟谷发育、岩土体破碎、松散的山区也是滑坡的多发地区。由统计结果可知,滑坡与坡度的关系呈现一定的规律性。研究区滑坡、不稳定斜坡产生的坡度一般是 15°~25°,出现这一现象的主要原因是滑坡体发生较多的坡度范围正是松散堆积体滑动的有利角度范围: 一般来说,坡越高、坡度越大就越易发生滑坡,且规模也越大。

研究区土壤种类以潮土和潮泥土为主,土壤质地多为砂质,成土母质包括近代河流冲积物和湖泊沉积物。境内共有 49.6×10⁴ km²丘陵岗地,土种类型为红壤,质地粘稠,土层较厚; 山地面积共计 56.93×10⁴ km²,土壤缺乏养分,肥力较低; 有水域面积 31.27×10⁴ km²,其中河湖外滩达 66.67×10⁴ km²,土种主要为荒州湖潮土、潮土和河潮土,呈沙壤质地,土壤较厚且为中性。土壤的物理条件、结构条件和环境条件在某种条件达到临界状态时均有可能造成滑坡。其触发因素可以归纳为内部因素和外部因素2类, 内部因素包括土壤性质、地形地貌、植被和垂直节理等; 外部因素是指突然触发滑坡运动的事件,包括降雨、冻融以及人类活动等。

研究区位于中亚热带北部常绿阔叶林地带,山地植被垂直分布明显,大致可分为2个植被带: 从几十米到 800 m,是典型的中亚热带常绿阔叶林; 1 000~2 000 m属于北亚热带山地植被类型,是常绿阔叶树、落叶阔叶树、喜暖针叶树的大混交林,其中也有少量的落叶阔叶树的纯林,混交林遭人为破坏后形成落叶阔叶树^[13]。植被对滑坡稳定性的影响大致可分为力学效应和水文效应,具有正-负2方面的影响。力学作用主要与根系加固土体有关; 水文效应主要与植被根冠截留、增加入渗、 减少蒸腾等有关。

研究区属中亚热带向北亚热带过渡的湿润季风气候区,呈明显的大陆性和季风性气候。四季分明,冬凉夏热,年温差在23.2~24.0 ℃之间。年均总降水量为1 386.9 mm,时空分布不均匀,主要集中在4月初—7月初。全市年平均降雨量1 382.4 mm,是全国平均值的2.3倍,其中沅水流域1 445.5 mm,澧水流域1 335 mm。常德市头顶长江三口,腰跨两水,分别是阮水和澧水,脚踏东南洞庭湖,水系极为发达,河流众多,流域面积在10 km²、河长5 km以上的大小河流有432条^[14]。其中,流域面积200~3 000 km²河流有42条,流域面积3 000 km²以上主洪道有8条,洪道全长635 km,分别为长江三口六支290 km,以及沅水165 km,澧水180 km。湘、资、沅、澧4水在汇入洞庭湖前,蜿蜒于山地丘陵之间,其中上游强烈地切割了山地,形成了大量的又长又陡的谷坡,甚至是高大的悬崖陡壁。这些地形地貌条件,加上河水冲刷等其他因素的影响,导致这些区域的滑坡大量发生。对于降水导致滑坡地质灾害发生的情况分析,一般存在2种情况: 一是突然性的短时高强度降水,即为暴雨和特大暴雨等情况引发的地质灾害; 二是长时间的过程性降水,即为中雨、大雨等持续性降水引发的地质灾害。降水引发滑坡主要有 4 种破坏形式: ①雨水渗入下渗,增大了坡体重量,产生动静水压力; ②降水侵蚀坡脚破坏坡体,破坏了边坡原来的稳定性结构; ③雨水渗入,弱化岩体,粘土矿物的水化作用降低了粘着力,降低摩擦系数,甚至造成了粘着力的消失,降水改变边坡力学性能; ④滑动体的渐进性破坏和渗透力的作用^[15]。

本文以湘北常澧山地-丘陵地区为研究对象,首先采用2种决策树分类方法对多光谱数据进行了土地利用分类; 然后运用DEM数据提取了高程、坡度、坡向、起伏度和曲率等5项地形地貌因子; 再基于SBAS-InSAR技术进行了地表时序微形变测量; 最后在GIS系统内综合危险性评价结果和形变速率对研究区滑坡隐患进行圈定和识别。主要研究的技术路线如图2所示。

本次研究的光学遥感数据源为Sentinel-2A,于2015年6月23成功发射,和Sentinel-2B卫星在同一轨道上运行,2颗卫星协同工作,主要用于植被、土壤和水覆盖的土地观测、土地利用变化检测、气候变化检测、救灾支持等方面。数据为2017—2021年间的7幅影像,7幅影像的时序分别是2017年10月30日、2017年12月24日、2018年4月18日、2018年10月5日、2019年1月23日、2020年3月18日以及2021年3月28日。

Sentinel-2A多光谱影像通过几何精校正得到正射影像,同时完成辐射定标和大气校正等预处理过程; 另外,多光谱影像的覆盖范围远大于此次研究区域,为减少计算量,在正式处理影像之前有必要先对各幅影像按照研究区范围进行裁剪。

分类回归决策树算法(classification and regression tree, CART)是一种基于规则的面向对象分类方法^[17]。采用CART 决策树分类方法得到的面向对象的分类规则集,充分调用了多源分类依据,可选性好,更有利于找到不同地物类型的最佳分类特征及相应阈值,且CART 算法结构清晰,易于理解,没有外部参数的干扰。此次决策树分类操作选用的是2018年3月28日的Sentinel-2A影像,即形变速率监测研究的开始时间的影像,具体步骤如下:

1) 创建分类样本。根据在真彩色影像下观察到的研究区土地利用情况,在图像中选择良好的分类样本。此次分类工作创建和提取了6类样本: 城区、水体、林地、耕地、裸地和待耕种用地。实现方法是在待分类的影像下创建ROI文件,一个ROI文件代表一类地物,基于数量够多、分布均匀、代表性好等优选样本的标准,得到分类样本文件; 紧接着计算样本类别的可分性。经计算各类别的可分性指数均大于1.9(>1.8即代表可分性良好),符合分类样本的要求。

2) 准备各分类特征。为提高决策树的分类精度,应在CART模型自动生成决策树前尽可能多地考虑和提供帮助区分各地物类别的分类特征。此次准备的特征指数有: 归一化植被指数(normalized difference vegetation index,NDVI)、归一化水体指数(normalized difference water index,NDWI)、归一化建筑物指数(normalized difference building index,NDBI); 此外还提供了DEM高程数据以及多光谱影像本身的12个波段信息。3个归一化指数的提取方法则是运用波段运算工具(Band Math)输入: (float(B1)-float(B2))/(float(B1)-float(B2)),并从影像各波段中选择不同归一化指数各自对应的波段1(B1)和波段2(B2),输出即可。NDVI,NDWI及NDBI的具体要求见表1。

表1   归一化指数公式一览表

本次研究用到的DEM数据源为GDEM V3。1999年12月,高级星载热辐射反射辐射计(Advanced Space borne Thermal Emission and Reflection Radiometer,ASTER)搭载于NASA的EOS-AM1 (Terra) 平台被成功发射。ASTER传感器由3个子系统组成,其中VNIR子系统在第三波段(波长范围0.78~0.86 μm) 有一个沿轨迹方向 (along track)的后视波段,与其向下观测的镜头形成良好的立体观测,从而接收到良好的立体像对,帮助提取DEM数据。一般在保证良好控制点的精度、分布和数量的情况下,利用遥感数据源即可提取得到高质量DEM数据,数据分辨率为30 m。

滑坡的形成和演化往往不是由单一、明确的要素引起,而是由多种复杂要素耦合,受到内外动力的共同作用,处于不同时空条件下的滑坡行为往往具有与其特征相关的触发条件和特征要素。引起滑坡发育的地理条件主要包括地质构造条件、岩土条件、地下水条件、地形地貌条件、坡体覆被条件等,其中地形地貌条件对滑坡在空间格局上的分别起着至关重要的作用,是识别滑坡以及判断滑坡是否可能发生最直接最简捷的方法。本文以地形地貌作为基本的滑坡形成条件,基于DEM高程数据提取研究区的高程、坡度、坡向、起伏度、粗糙度、剖面曲率等特征参数,对研究区可能发生滑坡的区域中的潜在滑坡进行识别和提取。

1)坡度。斜坡坡度是两点的高程差与其水平距离的反正切值,坡度值的大小对斜坡的稳定性影响很大-坡度能够影响坡体内部应力场的分布,在重力作用下,当应力超过坡体岩土的最大承受能力时,就易导致坡体发生滑动或崩塌。本研究利用DEM数据提取了研究区的地形坡度值,得到坡度范围在0°~56°,根据研究区的分布特点作了5类分级,得到的坡度分布图如图5(a)所示。

2)坡向。坡向是指坡面法线在水平面上的投影的方向,相当于是由高到低的方向。坡向对坡体的日照时间和太阳辐射强度有明显的影响,引起坡体表面的湿度、覆被情况等有所差异,从而对滑坡发育产生影响。本文利用DEM数据对研究区提取了研究区的坡向分布情况,根据研究区的分布特点作了9类分级,得到的坡向分布图如图5(b)所示。

3)坡高。坡高即斜坡的高程,对斜坡有明显的影响。斜坡内各处的应力值均随坡高的增高而线性增大,斜坡坡高越高,其稳定性会越低,滑动范围或规模将逐渐增大^[21]。研究区内的坡高范围在6~595 m,根据地貌形态可将本研究区分为3个地貌单元: 平原区(<200 m)、丘陵区(200~500 m)、低山区(>500 m)。本文利用DEM提取研究区的坡高分布情况,根据研究区的分布特点作了5类分级,结果如图6(a)所示。

4)起伏度。地形起伏度是指在一定的小范围内,最高点与最低点之间海拔高度的差值。当地形起伏度增大,坡体稳定性会随之减小,坡体发生滑坡时的滑坡位移速率会更大。本文利用DEM数据提取研究区的起伏度分布情况,得到起伏度范围在0~189 m,根据研究区的分布特点作了6类分级,得到的起伏程度分布图如图6(b)所示。

5)坡面曲率。坡面曲率包括在水平方向上的平面曲率和垂直方向上的剖面曲率: 平面曲率是指地表任意一点,在经过其等高线时,等高线在该点处的曲率值,指示地表曲面沿水平方向的弯曲程度^[24]; 剖面曲率是高程变化的二次导数,指斜坡上任意一点坡度的变化率。曲率值为正时指示斜坡为凸形坡,为负时指示斜坡为凹形坡。斜坡凹凸不平反映了斜坡坡面粗糙,坡面面积大,更易渗水,且影响了岩土体应力的分布,从而较平面坡更易引起滑坡。本研究利用DEM 数据分布提取了研究区的平面曲率和剖面曲率的分布图,普通平缓山区的期望值普遍介于-0.4~0.4之间,因此以此作为类别划定的阈值,得到的曲率分布图如图7所示。

本次研究InSAR数据源为Sentinel-1A,Sentinel-1A 卫星在 2014 年4月3日发射升空,是欧空局发射的第一颗环境监测卫星。该卫星具有可观的全球覆盖和持续的数据获取能力,能够为连续的地表变化研究提供数据支撑。由于该卫星采用TOPS成像模式,传感器推扫式获取的多个burst条带能够提供十分可观的空间分辨率(4 m×14 m)和覆盖幅宽(250 km)。此外,对于同一地区,该卫星的重访周期仅有12 d,满足长时序高频次监测需求。工作波段为C波段,具有一定的抗植被干扰能力。本文共获取了2018/01/06—2021/05/14期间共96景IW成像模式下的Sentinel-1A影像数据,其Path为11,Frame为91,数据覆盖范围如图8(a)所示,其中黄色加粗框为本项目所采用的burst覆盖范围。根据短基线集数据处理要求,按照150 m空间基线阈值和36 d时间基线阈值进行SAR影像组对。为保证构网的时空完整性,少数干涉对的部分时间基线为48 d或96 d,同时,在部分影像受到的大气误差干扰较严重,最终采用的SAR影像数量为92景。时空基线图见图8(b),干涉测量所使用数据的空间-时间基线整体均匀,其中空间基线约为200 m,最长未超过300 m; 时间基线约为15 d,最长未超过45 d,时间序列采样密度高,数据选择合理,同时测量时间周期为2018—2021年,可方便监测分析。

SBAS-InSAR^{[13-14,25⇓⇓⇓⇓⇓⇓-32]}获取长时序地表形变实际上基于多期连续监测的D-InSAR结果进行平差处理所获取的,因此整理流程大致可分为3部分,分别为: 影像对组合与预处理、D-InSAR差分干涉、误差处理与时序形变序列解算,具体流程见图9。

1)影像组合与预处理。针对研究区域,获取精确覆盖相同区域的SAR数据并进行聚焦处理生成SLC数据。选取合适的影像作为配准主影像,结合精密轨道数据将所有影像精确配准到与主影像一致的空间下。并利用外部DEM模拟的SAR影像,对其进行地理编码。随后,对于所有配准好的SLC影像,结合监测区域的实际情况,设置合适的时间和空间基线进行组对。

2)差分干涉处理。对每一个干涉对,将主影像和从影像共轭相乘后得到干涉图,干涉图包括强度信息和相位信息。通过计算每个窗口内像元的相干系数得到相干图来评价2个像元的相关性。将外部DEM模拟的地形相位与干涉图作差,去除地形相位,得到差分干涉相位。随后采用自适应滤波对干涉相位和相干性进行滤波处理,并选取稳定区域采用最小费用流方法进行相位解缠。采用二次多项式拟合以及高程拟合去除轨道残余相位以及高程相关大气效应。

3)时序形变解算与误差处理。解缠的差分干涉相位是2景影像间的相对各干涉对主影像的变化量,需要采用SVD分解的方法,将这些相对变化量转换到单一参考时间的相对变化量,从而得到原始形变序列。此时的形变序列还包含大气误差、DEM残差和轨道误差等。由于原始形变序列中还包括了湍流大气、局部地形高程残差以及非线性形变信号的影响,还需要根据大气延迟误差和轨道误差在时空域上的不同信号特征,通过时间高通滤波和空间低通滤波将其分离开来。并将获取的线性和滤波得到的非线性形变进行相加,结合时间基线参数,重新估算形变时间序列以及形变速率。

研究区地理地质环境复杂,InSAR形变监测精度主要受植被覆盖、大气及失相干噪声3个方面的影响。由于区内暂无满足长时序监测需求的相关实测数据(如GNSS或水准观测等),因此,监测结果的定量精度难以验证。为保证形变测量结果的可靠性,本文从以下2个方面开展了相关工作。

1)SBAS-InSAR处理流程优化。地形、植被及失相干在InSAR观测中会造成像元干涉相位质量低下,而这一问题直接体现到干涉相干性指标上。为此,本文采取的SBAS处理流程从以下3个方面对其进行了相关优化。首先,较短时间间隔的时序多主干涉技术能够在一定程度上抑制由于植被生长及大气变化造成的时间去相干影响; 其次,采用了包括高低通滤波、非局部均值滤波等多种手段对干涉相位进行了优化,保证了干涉成果质量; 最后,采用单景(0.6)和平均(0.5)2个相干性阈值对变化较大的像元进行剔除,保证了最终结果的有效性。

2)形变测量精度定性评价。由于缺乏GNSS实测数据,本文的时序形变测量精度从定性层面展开了验证评价。在形变结果的时间变化层面上,对显著变化的像元时间序列进行提取,并发现不同地形条件下像元的时间序列变化与其对应的典型形变序列模式基本保持一致。在形变结果的空间分布层面,不同形变量级像元的空间分布同样存在差异。例如,较大形变的像元通常出现在坡度较陡的丘陵、山地等区域(图10(a)),符合孕灾要素特征。最后,根据InSAR时序形变结果,对重点形变区域进行光学遥感解译及实地调查验证,验证结果表明本文获取的重点形变区域内均存在明显形变表征及孕灾条件,对比结果一致性较好。

经过前述SBAS-InSAR解算以及相关误差优化等流程,本文获取了研究区全区范围内2018年1月—2021年5月的形变速率与形变时间序列结果。如图10(b)所示,获取的形变监测结果有效覆盖了常德-临澧城区及周边山地等绝大部分区域(研究区平地区域为非孕灾环境区,形变主要为地面沉降,非本此研究重点,未显示)。形变速率描述了地面目标在监测时段内发生形变的平均快慢,其随时间的变化是形变-时间曲线的一阶导数,有助于形变阶段研判,是评判地质灾害风险的有效指标。

从形变速率分布的空间范围而言,常德-临澧地区形变速率大部分均在0值附近,常德、临澧、澧县、石门等城市区域出现了一定的聚集性地面沉降,分析原因,主要可能是城市建设和地下抽水引起。山地地区最为显著的的聚集性形变主要表现在研究区中部、西北部、中东部和南东部山脉,形变速率在-1.5 ~ -3.5 cm/a之间,是重要的中大型滑坡地质灾害隐患分布区域; 丘陵区域形变一般面积小,零散,主要分布于研究区中部,是中小型滑坡等地灾隐患的重要分布区。从数值上看,可以看出常德-临澧地区在项目监测时段内基本保持稳定,这也与研究区地表形变现状基本吻合。监测范围内并未出现大片聚集性或空间不连续的形变特征,这表明本文所获取的结果较好地剔除了InSAR大气延迟等影响。因此,可以看出,研究区的形变基本表现出整体稳定、零星变化的分布特征。零散的形变信号主要集中在山地、丘陵、湿地、滩涂等区域,在空间上表现出一定的与地形特征/地物类型之间的关联。理论上来说,这些区域也往往更容易受到局部地形坡度/土壤含水量等因素的影响,表现出地表高程的不稳定变化。

不同时空条件下引发滑坡灾害的特征要素并不相同,因此如何合理地选择特征要素以及各要素的阈值选择是此次研究的关键。一些学者采用先基于形变速率信息和高分辨率光学遥感图像的目视解译结果综合圈定滑坡点,再综合分析这些滑坡点特征要素特点的方法,来对滑坡点附近的其他区域的潜在滑坡点进行排查和圈定^[13,25],该方法可行性高,效果良好。本文获得的研究区年形变率介于-1.5~-3.5 cm/a,研究区内大部分地方基本处于稳定状态,但在局部区域反映有明显的形变。鉴于此,本文首先在不加入其他参考数据的情况下,仅以形变率<-2 mm/a作为阈值条件,利用ArcGIS软件的栅格计算器功能,提取了研究区内满足条件的点位,并圈定为研究区的“滑坡危险区”,结果如图11的黄点标注位置所示。由图可知,得到的点位分散,主要分布在林地和耕作用地中。然而在平原地区,特别是人工活动较密集的城区和耕地区域,往往并不具有孕灾条件,不会导致相应的地质灾害隐患,这也说明在对一个大范围地区进行形变监测时仅凭单一的数据源并不能获得准确的结果,需要借助多源数据综合分析。

基于以上情况,本文将得到的“滑坡危险区”图像与CART决策树分类结果和7幅时间序列多光谱影像相结合,采用目视解译方法对明显的滑坡点位进行了初步提取。目视解译的主要判别标准是:

1)结合CART分类结果,排除掉位于城区和耕作用地的区域,在靠近林地附近的区域进行滑坡和滑坡隐患解译。

2)选择2018年3月—2021年3月的5幅多光谱影像对“滑坡危险区”进行多时相观察,以符合形变速率测量数据是介于2018年2月—2021年5月时间段内的实际情况。

3)比较各个区域在不同时相的植被覆盖变化情况以及非植被区域的土地利用类型,综合圈定滑坡地质灾害隐患点。

经对比观察,发现研究区内多数“滑坡危险区”在时间序列上发生过明显的植被覆盖度先降低后增加的现象。这种现象的成因既可能是雨水冲刷引起,属于滑坡的孕灾要素; 也可能是该区域植被在生长过程的自然变化,但无植被覆盖且地表为土质的情况仍符合滑坡孕灾条件。综上所述,此次目视解译滑坡隐患点的识别提取标志确定如下: 区域位于林地,随时间序列的变化区域的植被明显减少,且曾明显大面积地露出红褐色土质层,如图11。最终共圈定了16处滑坡隐患点。

本文将初步圈定的16处自然滑坡点作为样本点,将提取好的5项地形地貌特征参数以属性赋给各样本点,以此统计分析各个点的地形地貌要素的特征,对各要素的各个区间进行合理的赋值,并利用加权分析方法对研究区进行潜在滑坡点评价,帮助进行滑坡隐患点综合识别。

基于地形地貌要素的滑坡识别工作选用了研究区的坡度、起伏度、坡面曲率(剖面曲率、平面曲率)、坡高和坡向共计6个滑坡要素作为滑坡识别的评价因子,依据样本点的统计分析结果对各评价因子进行重分类并对各区间进行赋值,得到了如图12的评价标准。利用ArcGIS软件Spatial Analyst Tools-叠加分析-加权分析方法,对研究区进行滑坡危险性评价。得到的研究区滑坡危险性评价结果,如图12所示。

图12中6个评价因子的权重均为1,允许的分数范围为0~60分。滑坡危险性评价结果显示,评价结果分值较高的区域主要位于林地及其附近,且在林地区域的分布较为均匀,并无明显规律。鉴于此,本文拟将评价分值大于50的区域作为滑坡危险性评价阈值。

将滑坡危险性评价分值>50和形变率≥-0.02 m/a作为阈值条件,利用ArcGIS栅格计算器功能,对研究区进行了滑坡隐患点位的提取。经统计,此次提取的滑坡隐患的面积共计0.126 km²。

绝大多数点位分布在林地区域,其空间分布及局部图见图13,由于研究区植被覆盖严重,存在滑坡隐患的部分区域往往难以通过目视解译识别和提取,结合InSAR时序形变速率信息和滑坡致灾因子的危险性评价信息进行综合分析和提取,则有效降低了植被覆盖的影响,且多源信息的综合利用方法能够有效提高滑坡隐患提取的精度。

根据前文所建立的多源遥感技术集成的地灾隐患早期识别技术体系,圈定了常德-临澧应用示范区地灾隐患,并筛选了部分进行野外验证,以验证所建立技术方法体系的有效性和可靠性,验证结果见表2。从表2中可以清晰的看出,所圈定的地灾隐患点基本都发生了明显的形变迹象,一般地表地物对应着突变现象,部分严重地区已形成了明显的滑坡地质灾害,滑体与滑动结构面明显。证明本文基于InSAR、光学等多源遥感的滑坡地灾隐患的早期识别和提取技术在南方高植被覆盖区具有较高的准确度和精度,可进行推广示范应用。

表2   常德-临澧研究区部分地灾隐患点野外验证结果表

CART决策树分类结果图和研究区时序形变速率对比图如图14所示,图14(a)为研究区形变速率,其红色标志指示形变率<-0.01 m/a的区域; 图14(b)为CART决策树分类结果图,黑色标志指示城区,红色标志指示耕地,黄色标志指示未耕作用地。对比发现,城镇区形变速率一般≥-0.01 m/a,城郊区,即耕作用地区域(包括耕地和未耕作用地),出现了大量且明显的形变速率<-0.01 m/a的聚集区。

在进一步分析了研究区的水系分布情况以后,初步推测,除去圈定的滑坡点,此次研究范围内其他形变速率较高(<-0.01 mm/a)的区域主要有2种地理单元: 第一种位于沅水和澧水的干流流域,很可能受流水和雨水的影响明显,导致这些地方被冲刷剥蚀,但与滑坡孕灾条件并无明显相关性; 第二种为围绕研究区的洞庭凹陷盆地区域,其地面沉降形变可能是人类活动和地下水补排给引起。以上情况说明使用单一的形变速率数据进行滑坡隐患提取适合大范围初步筛选,综合地形地貌信息、土地利用分类信息等多源信息有助于提高了隐患的识别精度。

本文利用光学遥感数据、SBAS-InSAR时序形变速数据和DEM数据3种数据综合对湖南常澧地区的滑坡地灾隐患进行了研究。首先,基于2种决策树分类方法,利用光学遥感图像对研究区进行了土地利用分类; 然后运用DEM数据提取了高程、坡度、坡向、起伏度和曲率(剖面曲率和平面曲率)等6项地形地貌要素对研究区进行了基于滑坡危险性评价的滑坡识别; 最后综合危险性评价结果和时序形变速率在时间序列上对研究区滑坡隐患进行了提取和圈定,并基于CART决策树分类结果和研究区水系分布情况,对研究区滑坡隐患以外大于-0.01 m/a的区域进行了形变驱动因子分析。主要的结论如下:

1)相较传统决策树分类方法,CART决策树可充分利用多源数据用作分类规则,如NDVI,NDWI,NDBI等指数,以及DEM数据、坡度等,其分类结果精度更高,可较为有效地排除与滑坡地灾隐患不相关的土地利用区域,并成功识别出了与滑坡地灾隐患密切相关的林地用地类型。

2)本文采用光学遥感、InSAR形变测量、土地利用分类多源遥感技术实现了常澧地区滑坡隐患的综合遥感技术识别,建立了适用于该区的滑坡隐患识别的综合遥感解译标志,圈定了16处滑坡隐患高风险区,识别出了数处植被覆盖区滑坡隐患和4处裸露区与人工切坡相关的高隐蔽性滑坡隐患。证明了综合遥感技术在该区进行滑坡地灾隐患识别较单一方法更加有效,识别精度更高,为同类应用提供参考依据。