天地联合的区域山地管道地质灾害监测预警体系研究

邹永胜 李双琴 高建章 郑宇恒 李章青 沈忱 邹妍

[摘    要] 西南管道公司所辖山区管道占比达70%，由于穿越地形地貌复杂，管道地质灾害防控难度大，传统的监测手段很难实现管道全域覆盖，急需建立一套经济有效的区域管道地质灾害监测体系。文章提出了一种新的区域管道地质灾害监测预警体系，从GNSS和InSAR数据的融合、降雨-形变的归一化、天地联合的预警模型等方面详细阐述了该监测预警体系，并以某管道50km范围为试验区，采用Sentinel-1A、地表GNSS监测站和降雨量化监测站为数据源，对发现的82处潜在风险区进行了为期2年的跟踪监测预警，验证了该体系的有效性和可靠性。

[关键词] 天地联合;地表形变;InSAR;监测;管道;地质灾害;监测;预警

doi ： 10 . 3969 / j . issn . 1673 - 0194 . 2020. 15. 081

[中图分类号] TP311    [文献标识码]  A      [文章编号]  1673 - 0194（2020）15- 0192- 05

1      引    言

西南管道公司所轄管道作为保障国家四大油气战略通道之一，对改变高度依赖马六甲海峡的能源进口格局，实施油气进口多元化战略具有重大意义。西南管道公司所辖管道穿越七省一市，山区管道占比高达70%，沿途区域自然地理和地质环境复杂多样，管道受地质灾害威胁的安全风险高。特别是2008年汶川地震发生后，管道沿线发生了大量的滑坡、崩塌地质灾害，如K572滑坡、K0528+850滑坡、K0526+300山体变形。因此防治地质灾害是摆在管道安全管理面前的一项重要课题，而地质灾害监测预警是管道安全管理的一个重要环节。

目前的管道监测手段大多采用地面专业监测设备进行监测，通过持续采集地质灾害与管道实时的状态数据，从而判断地质灾害和管道的安全状态。这种传统的接触式监测以点为监测单元，虽然可以实现高频高精度的不间断监测，但无法应用于未发现的高位远程、高隐蔽性地灾隐患，很难实现对整个管廊区域开展大尺度全域监测。

目前，科学技术日新月异，随着遥感技术的不断进步，遥感技术在油气管道监测中也发挥着越来越重要的作用。近年来，遥感技术在油气管道选线、地质灾害监测等领域的应用已经得到高度重视。例如在国内已完成的6条管道选线工程中，陕甘宁气田至北京输气管道和内蒙二连浩特至山东日照输气管道均采用了遥感技术。王世洪  等基于多源遥感检测信息，以该段线路泥石流为例，采用模糊综合评判法分析了降雨作用下泥石流灾害发生概率，并对管道泥石流灾害危险性进行了评价，从而实现了管道泥石流灾害的快速预警。程涛利用布设人工角反射器（CR）PS-InSAR方法对西气东输管道工程靖边—临汾端子长县地区2004.01.10-2006.12.30的形变进行监测;徐小波也使用InSAR技术对子长县管道区域2009年5月至2010年9月的地表形变进行了监测。本文通过综合InSAR技术和地面专业监测技术，形成多源数据融合的天地联合的管道地质灾害监测预警体系，并以西南地区某管道晴隆段50km范围为试验区进行了验证，为InSAR技术在山区管道监测预警方面的应用提供参考。

2      天地联合的管道地质灾害监测预警体系

2.1   GPS和InSAR融合技术

目前针对地表形变的监测主要采用GPS测量和InSAR技术。这两种方法的精度、时间分辨率和空间分辨率各不相同，GPS连续观测测量技术时间分辨率高，InSAR技术空间分辨率较高。为了满足管道沿线地质形变监测的要求，需要对InSAR、GPS数据进行融合处理，充分利用各自数据和技术的优势，进行互补，获得高精度和高时空分辨率的数据成果，以便开展数据分析和综合预警。

将GPS接收机与角反射器安装在同一位置，确保地基监测的同时对SAR卫星进行修正。根据角反射器的高程值转换为相位值，并以其为约束，提高解缠精度，其公式为：

式中，Bh表示水平基线，By表示垂直基线，θ为视角，r1为主图像的斜距，λ为波长，h为GPS高程，Φ为GPS高程计算的相位值[1]。

对GPS形变结果进行空间域内插，GPS虽然可以提供高时间分辨率的地面形变资料，但受限于GPS监测点位的密度问题，其不可能反映出监测区域的全面沉降情况[2]。由于空间数据一般都具有相关性，即在一定的范围内，任何空间数据的变化一定会带动周边一起发生变化，其趋势和量级均具有一致性和渐变性的特征。

因此，结合InSAR的结果，通过对监测区域内的形变范围、趋势、量级等背景的GIS分析，将区域内离散GPS点的高时间分辨率特性按照InSARS数据的空间分辨率进行克里金插值，在此基础上，利用卡尔曼滤波对所有网格节点进行时间序列化处理，其方程式为：

L（tk）=H（tk）X（tk）+V（tk）（2）

式中，X（tk）为在tk时刻的形变量;H（tk）和F分别代表观测分配矩阵和状态转移矩阵;V（tk）表示观测噪声和;W（tk-1）表示状态噪声;Γ（tk-1）表示状态噪声分配矩阵。

2.2   潜在风险区识别

坡度、坡向以及山坡到管线的距离是产生山体滑坡从而危害石油管道的重要因素，故采用权重分析法综合评价三种因素对石油管道的影响程度，并后形成每个区域的关注等级，依据每个区域的权重因子划定形变量的风险阈值，根据InSAR处理结果的形变量，识别潜在风险区。同时，采用天地基数据迭代的方式，不断修正风险阈值。

（1）坡度计算重分类

通过DEM计算坡度，根据坡度等级进行重分类并赋予不同的权重。

坡度分类依据为：0°～5°为平坡，6°～15°为缓坡，16°～25°为斜坡，26°～35°为陡坡，36°～40°为急陡坡，41°～45°为急坡，46°以上为险坡。

（2）坡向计算

将坡向进行重分类并赋予不同的权重。由于石油管道两侧相同坡向发生滑坡时产生的危害程度不同，故将石油管道分为北侧和南侧，进行坡向重分类并赋予不同的权重。

（3）距离权重

以管道中线为中心进行欧氏距离计算，根据自然间断法将距离分级并赋予权重。

欧式距离是一个通常采用的距离定义，指在m维空间中两个点之间的真实距离，或者向量的自然长度（即该点到原点的距离）。在二维空间中的欧氏距离就是两点之间的实际距离，如图2所示。

二维空间的欧式距离计算公式如下：

式中，Z为（x1，y1），（x2，y2）两点之间的欧氏距离;|X|为点（x2，y2）到原点的欧氏距离。

自然间断法是基于数据中固有的自然分组，将对分类间隔加以识别，可对相似值进行最恰当地分组，并可使各个类之间的差异最大化。要素将被划分为多个类，对于这些类，会在数据值的差异相对较大的位置处设置其边界，不管采用的是InSAR分析的哪种方法，其结果数据都可以转换为GIS数据，因此在进行自然间断法分类时可以直接在ARCGIS软件中使用工具进行分组（如图3所示）。

（4）关注区等级划分

综合坡度、坡向以及距离权重结果进行综合权重计算，分为不同关注区等级，并根据形变速率等级划分为不同等级危险区域。

（5）潜在风险区识别

将标注出的沉降量或抬升量较大的区域叠加到三维地图中，结合实际坡度坡向以及到管线距离进行风险区域筛选。关注等级高的区域，其形变风险阈值相对较低（如图4所示）。

2.3   天地联合的区域地质灾害预警模型

对于区域地质灾害的预警分析，通过建立地面降雨自动化监测站获取实时的降雨信息，结合InSAR形变量和速率，通过改进切线角法建立降雨-形变归一化综合预警模型。

（1）地质灾害有效降雨量

受风险区内地形地貌、植被等因素影响，瞬时降雨量中只有部分对地质灾害发生起作用，因此在进行管道地质灾害预警时不能直接采用仪器测得的瞬时降雨量，而应该采用有效降雨量，其公式为：

Rc=R0+Ra=R0+kR1+k2R2+…+knRn（6）

式中，Rc是地质灾害发生的有效降雨量，R0是地质灾害发生当天的降雨量，Ra是地质灾害发生的前期有效降雨量，k为有效降雨系数，通常取0.9或0.8。Rn是地质灾害发生前第n天的降雨量，n通常取值为5天。

（2）降雨-形变归一化

InSAR数据和有效降雨量数据都是时序数据，因此采用时序归一化算法进行数据的融合，可以表示为：

式中，y（xm+n）表示第m+n时刻位移值，xm+n表示第m+n时刻前p个时刻的位移序列和有效降雨量，xi+p表示第i+p时刻前p个时刻的位移序列和有效降雨量，表示降雨影响天数，b为固定常数，a为有效降雨系数[3]。

（3）天地联合的区域地质灾害预警模型

滑坡变形-时间曲线的切线角在不同的阶段有不同的特点。初始变形阶段，切线角总体表现为由大减小，逐渐平稳;当进入等速变形阶段，切线角基本稳定在45°;在加速变形阶段，切线角增长速率急剧扩大，直至90°[4]，如图5所示。

根据图5的T-t曲线，可以得到改进的切线角αi的表达式：

式中，αi表示改进的切线角;ti表示某一监测时刻;Δt表示与计算S时对应的单位时间段（一般采用一个监测周期，如1天、1周等）;ΔT表示单位时间段内T（i）的变化量。

因此，本文在归一化值阈值预警基础上，综合考虑速率增量以及切线角指标，构建过程综合预警模型，如表1所示。

3      案例分析

3.1   研究区概况

研究区位于黔西南州普安县、晴隆县境内，属深切割岩溶侵蚀山区，区内大部分海拔在600～1 800 m之间，地势中部高而两端低。研究区属亚热带湿润季风气候区，具有高原性、季风性、湿润性的气候特点，年平均降水量1 438.9 mm，多集中在每年5～10月。区域内出露最老地层为石炭系（C），最新地层为第四系（Q）。研究区位于扬子准地台黔北台隆六盘水断陷普安旋钮构造变形区内。区域内褶皱、断层发育，岩体节理裂隙较为发育，受构造影响，岩体较破碎，岩层产状变化较大，倾角一般10°～40°之间，局部受断层的影响，岩层倾角较陡。

3.2   实验数据

本研究采用欧空局的Sentinel-1雷达卫星20 m分辨率数据，Sentinel-1数据监测时间为2016年8月6日～2018年8月8日，其基本信息如表2所示。

以查詢到的2016年8月～2018年8月8日数据为例，升轨及降轨遥感数据的覆盖范围见图6。由于所需数据量较大，监测持续时间长，选择下载Sentinel-1雷达卫星数据升轨影像数据44景和降轨影像数据42景。

研究区域内共布设了地面监测站9处，分别分布于沿线特点地质灾害点位上，与升降轨数据的位置叠加情况如图7所示，开展了为期半年到一年不等的专业监测。

通过地面监测站修正后的PS-InSAR数据共提取到361 028个PS点，形变速率的变化范围为[-24.59，24.40]，累计形变量的值域为[-72.34，71.51]。形变速率和累计形变量直方图大体呈现出左右对称分布的形态，绝大部分PS点的形变速率绝对值在4 mm/y以内（如图8～图10所示）。

通过可视区判定、坡度权重、坡向权重、欧式距离权重综合分析后，得到研究区內的潜在风险区，共82个，主要集中在管线首尾两端和中间部分区域（如图11所示）。

以2号风险区为具体分析对象，可以发现在2号风险区内有明显形变，采用本模型进行跟踪预警发现，该区域内的发育不稳定斜坡1处，从2016年8月6日开始监测以来，并无明显形变迹象，前15个月累计形变不足5 mm，从2017年12月开始发生形变，之后5个月内降雨-形变归一化平均值达到10.5 mm/d，引入改正的切线角模型，分析得出其切线角随时间的变化曲线，区间为2018年5月4日至6月9日，6月9日达到86.31°，发布了红色预警。根据后来现场调查，该斜坡于6月10发生局部滑塌（如图12～图13所示）。

4      结    语

本文提出的基于天地联合的区域山地管道地质灾害监测预警体系对于当前区域管道地质灾害监测预警是一种新的尝试，可以解决传统专业监测的空间局限性和InSAR监测时间与非地质分析局限性，充分发挥了二者的优势，形成了高时空精度的管道全域覆盖监测体系，并通过实验区的应用验证了该套体系的实际效果，对管道的日常管理和风险防控具有十分重要的科学意义。

主要参考文献

[1]刘冬，郑南山，王刚.基于GPS与D-INSAR 融合技术的矿区沉降监测[J].城市勘测，2011（4）：96-98，102.

[2]李更尔，周元华.InSAR、水准及GPS数据融合处理方法[J].测绘通报，2017（9）：78-82.

[3]俞俊平，陈志坚，武立军.基于多因素位移时序PSO-SVM的边坡变形预测[J].勘察科学技术，2015（1）：1-5.

[4]苏白燕.基于动态数据驱动技术的地质灾害监测预警研究[D].成都：成都理工大学，2018.

[5]亓星.突发型黄土滑坡监测预警研究——以甘肃黑方台黄土滑坡为例[D].成都：成都理工大学，2017.

[6]焦明连， 蒋廷臣. 基于GPS-InSAR集成技术地表形变的监测[J]. 测绘科学，2008（6）：58-60.

[7]焦明连. GPS与InSAR数据融合方法及其应用[J]. 全球定位系统，2010（3）：5-8.

[8]朱建军，李志伟，胡俊.InSAR变形监测方法与研究进展[J].测绘学报，2017，46（10）：1717-1733.

[9]张丽华.PS-InSAR技术在昆明长水机场地面形变的监测应用[J].昆明冶金高等专科学校学报，2018，34（1）：6-10.

[10]Hooper A， Segall P， Zebker H. Persistent Scatterer Interferometric Synthetic Aperture Radar for Crustal Deformation Analysis， with Application to Volcán Alcedo， Galápagos[J]. Journal of Geophysical Research， 2007，112.

[11]Ferretti A，Prati C，Rocca F.Permanent Scatterers in SAR Interferometry[J].IEEE Transactions on Geoscience and Remote Sensing，2001，39 （1）：8-20.

[12]Pritchard M E，Simons M. A Satellite Geodetic Survey of Large-scale Deformation of Volcanic Centres in the Central Andes[J].Nature，2002，418（ 6894）：167-171.

[13]陈强，刘国祥，李永树，等.干涉雷达永久散射体自动探测——算法与实验结果[J]. 测绘学报，2006，35（ 2） ： 112 - 117.

[14]李永生.高级时序InSAR地面形变监测及地震同震震后形变反演[D].哈尔滨：中国地震局工程力学研究所，2014.