基于InSAR的气田潜在地质灾害隐患识别技术研究

张 苏，徐兴雨

(1.中国石油化工集团有限公司安全监管部，北京 100728 2.中国石化胜利油田技术检测中心，山东东营 257000)

0 前言

天然气资源是关系国民经济发展的重要战略资源，近年来，随着东部陆相盆地油气勘探开发难度的日益增大、国际油价的节节攀升以及我国油气对外依存度屡创新高，以西部海相天然气为代表的油气资源在我国一次能源中的地位正变得越来越重要。普光气田作为我国最大规模海相整装高含硫气田，在促进我国能源结构调整、缓解沿海发达地区用气短缺、减小大气污染等方面发挥了重要作用，是中国石化乃至我国天然气资源增储上产的“主阵地”、“压舱石”。普光气田开发区块主体位于我国川东北地区，近年来受区域地震活动增强、极端降水天气增多等因素影响，滑坡、崩塌等各类地质灾害发生频率、数量明显上升，因地质灾害造成停工停产事件时有发生，集输管道等设施一旦受损导致高含硫化氢原料气泄漏，极可能引发严重安全事故，造成重大的社会影响。

当前气田地质灾害隐患辨识主要依靠定期踏勘，周期较长且无法及时捕捉微小变形，依靠水准仪、GNSS等技术开展精确测绘，则存在工作量大、通视要求高等不足。随着卫星测量技术的发展，基于合成孔径雷达干涉技术(InSAR)在建构筑物变形、山体滑坡形变监测等领域的应用日渐广泛，其优点在于利用卫星的定期巡访，不需要人工投入和部署任何设施即可实现不受地形、天气影响的精确成像，获得稳定点毫米级的变形数据。近年来部分学者[1～3]利用InSAR技术对山区地质灾害前兆信息监测开展了相应研究，但对于油气设施特别是管道等线性工程的相关研究较少，特别是PS-InSAR、SBAS-InSAR等技术在气田地质灾害中的应用有待进一步研究。本文基于时序InSAR(SBAS-InSAR)技术并结合高分辨率遥感和地面精细核查对山区气田开展潜在地质灾害隐患辨识研究，为摸清气田范围内地质灾害隐患底数，指导地面重点调查与风险分级管控提供依据。

1 技术原理概述

合成孔径雷达(SAR)是一种全天候、全天时对地观测的主动式高分辨率雷达成像系统。传统InSAR技术通过差分干涉测量(D-InSAR)，对同一地区不同时间两幅雷达影像进行处理，进而获取监测点位变形信息，这种方法受卫星轨道、大气效应、失相干等影响，常存在较大误差，为了消除干涉过程中的失相干与大气延迟带来的测量误差，意大利学者费雷蒂首先提出基于永久散射体进行差分干涉(PS-InSAR)的技术思路。该项技术通过分析不受时空变化影响的点目标来提取形变信息，进而对多幅影像进行处理，削弱或部分消除上述误差源的影响，提高形变信息获取的精度。SBAS-InSAR技术则采用设置相应的时空基线阈值，使影像自由组合构成干涉对，可以进一步削弱空间失相关的影响，处理后数据可达到毫米级分辨率。

本次研究为了有效地降低时空去相干的影响，综合考虑研究区域地貌地形特征、气候条件以及SAR数据量等因素，以遥感地质学、工程地质学、水文地质学为指导，以空间分辨率优于15 m的分辨率遥感数据为数据源，以工作区地质环境条件及地灾分布情况为研究对象，采用SBAS技术对观测数据进行处理解译，获得毫米级变形数据，进而分析工作区内地表形变规律，解译变形区域，确定典型灾害范围，并结合遥感、地面调查等结果进行分析验证。项目总体实施的技术路线如图1所示。

图1 总体技术路线

2 数据选取与处理

2.1 数据来源

研究区位于四川省达州市宣汉县普光镇，属超深、高含硫、高压、复杂山地气田[4]，其整体区块位于四川盆地东北大巴山南麓，地处渝川陕鄂结合部。研究选用欧空局回访周期12天的sentinel-1A合成孔径雷达卫星数据作为基础数据源，轨道为升轨，入射角值约39.456°，采集2019—2021年91景升轨数据。以日本航天局发布12.5 m分辨率外部DEM作为地理参考，选取稳定基岩和人工基础设施等点作为GCP点，为轨道精炼提供依据。采用Goldstein滤波算法进行地表植被的滤波处理，用最小费用流算法进行奇异值分解，最终获得研究区的形变范围[5]。

2.2 生成连接图

通过对研究区升轨数据进行时序干涉处理，共生成455对相干像对，时间基线和空间基线连接图见图2。所有的像对用连接线表示，每个线段表示所对应的干涉像对，连线两端表示不同时间的影像数据，其中黄色点即为超级主影像，其余绿色点均为从影像。研究表明，较为理想的连接图中一幅影像应至少和5幅影像数据发生干涉。通过时空基线和后期质量较差的部分影像数据剔除，基于SAR图像基线，并结合短时空基线原则进行阈值确定，生成连接图。

图2 SBAS-InSAR时空基线

2.3 空间干涉

在短基线子集干涉测量中，相干性阈值是表征雷达影像数据相干性质量的重要参数，相干性大小决定了干涉图中相位信噪比情况。相干性的高低还受雷达成像时天气、地形及地物类别影响。采用二分法在[0，1.0]范围内取依次取值0.8，0.4，0.2，0.1等值作为相干性阈值，基于对比分析可知，研究区管线跨越范围内，相干性阈值为0.2时，干涉强度图亮度较高，相干性良好。本文在升轨数据中采用的阈值为0.2，多视视数Azimuth Looks设置为1，Range Looks设置为4，基于研究区丘陵地区相干性较低的特征，采用最小费用流解缠方法有效抑制噪声。如图3所示，为了削弱大气、植被等引起的噪声，选择Goldstein算法进行滤波处理。在剔除此类相干性不良、解缠效果差的影像像对后，保留优质干涉像对用以进一步估算地表变形速率[6]。

图3 影像相干性图与滤波后的差分干涉、解缠结果

2.4 轨道精炼与形变反演

受地形起伏和地球曲率影响，干涉像对解缠后还存在一定的残余相位和相位坡道，通过确定研究区内未变形的点并将其作为控制点，从而确定变形区域的形变情况。基于野外踏勘和干涉结果，选择稳定非变形、干涉图上选择没有形变条纹的部分道路、建筑、岩石部位作为GCP点用以轨道精炼。

两次反演在整个形变提取过程中占有极为重要的位置。第一次反演主要是估算地表的形变速率和残余地形，用来对干涉图进行去平地效应。考虑到研究区整体相干性较低的特点，首先采用奇异值分解法估算第一次形变速率和残余地形，同时进行二次解缠优化输入的干涉图，为第二次反演提供较优的干涉图，第二次反演中，通过滤波估算和去除大气相位，得到最优形变时间序列结果(图4)。

图4 地形相位与大气相位估计示意

3 研究区地质灾害识别结果分析

根据上述InSAR数据处理流程，本次解译形变点共58处，为了对本次解译成果准确性进行对比验证，选取解译形变量较大典型地区，结合高分辨率遥感和地面现场调查等方法，进行现场交叉核查验证。

3.1 毛坝1井-1#阀室管线西侧

如图5所示，本点通过变形解译分析并结合气田地形地貌数据，初步判断该处年平均变形速率约0.038 m/a，根据其变形情况，圈绘疑似滑坡点位置，最大变形累计量约0.106 m。野外利用现场无人机调查和实地勘察等方法，结合现场的土质及工程建设情况，发现该处道路有明显变形迹象，表现为道路出现裂缝，局部有掉土掉块的现象，其变形范围约0.01 km2，地势较陡，平均坡度约20°～30°，地表主要以第四系覆盖层为主，未发现基岩出露，判定为滑坡灾害。

图5 毛坝1井-1#阀室管线西侧地质灾害InSAR解译及现场示意

3.2 P104集气站-17#阀室段

该点位于P104集气站附近，基于InSAR干涉结果和地形地貌数据，判断该处明显变形，年平均变形速率约0.024 m/a，根据其变形情况，勾绘变形范围面积约0.012 km2。在滑坡区域选取形变速率较大剖面绘制历史形变累计结果，发现该点在3年期间形变量达0.071 m，且变形体中部变形大于变形体前缘和后缘(图6)。实地勘察表明地形坡度在15°～25°之间，主要为第四系覆盖层，结合现场土质及工程建设情况，其变形迹象主要为坡体局部开裂下沉，坡表管线桩倾斜、步梯开裂，初步判定该处为一不稳定斜坡。

图6 P104集气站-17#阀室段地质灾害InSAR解译及现场示意

3.3 总结分析

综上所述，本研究通过InSAR变形解译分析共计识别出疑似地质灾害隐患点58处，通过高分遥感对比、现场验证对典型变形点确认、核对，现场验证，剔除因耕种削土、采矿沉降等因素造成的变形点合计20处，确定为地质灾害隐患点共计38处，其中距各类油气较远不威胁设施安全地质灾害点7处，其余25处地表形变点纳入普光气田地质灾害隐患台账。

4 地质灾害分布特征

4.1 高程分析

山体高程决定了地质灾害体初始重力势能，进而影响灾害运动范围、破坏能量等关键因素，此外从空间上，人类活动、植被分布、降雨、冻融基本都发生在一定高程范围内，故高程变化对斜坡的稳定性影响不可忽视。

研究区以低山地形为主，东北端和西南侧分别为中山地形和丘陵地形。通过提取地表数字高程信息分析表明，研究区海拔范围在105～1 300 m之间，地质灾害主要发育高程在105～740 m范围内，高程在450 m以下、450～580 m、580～740 m、740～960 m、960 m以上范围内地质灾害占比分别为37.9%，39.6%，17.2%，5.1%，0。这说明滑坡在高程较低的区域较为集中，而随着高程增加，灾害分布数量呈减小趋势(图7)。

图7 地质灾害分布特征分析

4.2 坡度分析

边坡的坡度影响着岩土体应力分布且与灾害体有效临空面密切相关。在应力方面，随着坡度的增大，灾害体坡脚处应力集中越明显，失稳概率相应增大。故坡度对斜坡类灾害的发育起着至关重要作用。如图7所示，研究区整体坡度在0°～89°之间，存在多级平原以及夷平面，其中地质灾害多发育在17°～37°的坡度中，该范围内坡体在饱水条件下极易发生破坏，其中滑坡主要发育于26°～37°范围内，不稳定斜坡主要发育于17°～26°范围内。

4.3 坡向分析

坡向对坡体的湿度与温度产生影响，不同方向斜坡受太阳照射强度与时长的差异导致温度和湿度、植被不同，进而影响外界的物理、生物风化作用。通过统计各个坡向滑坡的分布情况，针对滑坡、崩塌、不稳定斜坡、水毁进行坡向分析，据图7(c)可知，滑坡、不稳定斜坡主要发育于N、NE与W方向，这是由于工作区的构造和工程岩组对其影响的结果。其中，N、NE方向斜坡以不稳定斜坡为主(N向5处，NE向4处)，而E/W方向斜坡以滑坡发育数量较多(E向3处，W向4处)，其余坡向灾害发育规律不明显。

5 结论

通过获取最新的Sentinel-1A高分辨率雷达数据，利用InSAR技术对普光气田2019—2021年的地表形变时序进行反演，通过对结果进行分析研究得到以下结论。

a) Sentinel-1A影像数据时间基线短、影像数量较多，为InSAR时序分析提供了有利条件，本研究采用Goldstein滤波算法进行地表植被的滤波处理，用最小费用流算法进行奇异值分解，获得了高质量的差分干涉图，相比传统地质灾害调查手段，InSAR技术具有重访周期短，探测面积广、识别精确、自动化水平高等优势。

b) 基于对2019—2021年监测时间段地表变形速率分析，在研究区内识别出58处地表形变速率异常区域，结合高分遥感与地面详细调查核实，剔除因耕种削土、植被生长、采矿沉降等因素造成的变形点合计20处，确定为地质灾害隐患点共计38处，准确率65.5%，有效验证了InSAR在大范围地质灾害隐患筛查中的作用，相关成果可为后期地质灾害精细化地面调查提供有效指导。

c) 通过提取研究区内地表数字高程信息分析表明，地质灾害多发育于105～740 m高程和17°～37°的坡度中，其中滑坡主要发育于26°～37°范围内，不稳定斜坡主要发育于17°～26°范围内，随高处和坡度增加，灾害数量呈下降态势。