地质灾害是自然或人为因素作用下形成的对人类生命财产、环境造成破坏或损坏的地质作用^[1]，主要包括斜坡变形与地形变形引起的地质灾害，如崩塌、滑坡、泥石流属于典型的斜坡变形灾害。通过多轮地质灾害风险普查调查与灾害隐患排查，以及建立全域覆盖的群测群防体系，发现了大量隐患点风险区。但近年来发生的多起重大地质灾害事件中有近70%并不在已发现的隐患库内^[2]。广东省2020年在册地质灾害隐患点4744处，其中，崩塌2646处、滑坡1654处、泥石流77处，崩滑流合计4377，占比92.3%，威胁24.95万人^[3]，全面识别排查潜在隐患仍是广东省当前和今后一段时期防灾减灾重要任务之一。

地质灾害隐患早期识别从以往光学遥感为主的“图谱测量”走向多种遥感手段的“图谱与几何测量”综合遥感应用^[4]。如吴森等^[5]提出的无人机倾斜摄影城市化建设地质灾害早期识别应用，贾虎军等^[6]在此基础上综合了无人机载LiDAR扫描技术，均表现出精度高投入少优势，但受限无人机载平台航时、空域资源等条件，大范围推广应用存在效率及成本等诸多难点。基于遥感影像纹理特征来定性或定量判释灾害位置范围是光学遥感^[7]主要的优势，但很难发现地质环境脆弱或稳定性较差区域，形变特征的缺失限制了光学遥感在地质灾害识别领域应用的深度。在地质灾害动态发展规律和特征基础上，许强等^[8]提出构建天－空－地一体化的“三查”体系进行重大地质灾害隐患的早期识别，并成功发现白格滑坡上下游多处地质灾害^[9]。特别是合成孔径雷达干涉测量技术（Interferometric synthetic aperture radar， InSAR）的深入应用^[10-11]，为时序准确发现地表形变特征提供支撑。图谱与几何测量正是结合了光学遥感（含激光雷达）与InSAR技术，综合获取了地质灾害潜在区域的位置以及形态、形变及形势的“三形”解译标识^[12-13]，成为贵州、甘肃、宁夏等地质灾害隐患早期识别试点主要方法^[14-16]。广东地处大陆东南沿海，地貌类型复杂、丘陵山地多、地势起伏大、地形切割强烈、地质构造较为复杂、地质环境脆弱、受人类工程活动影响大，地质灾害易发多发、点多面广、活动频繁、危害严重，是我国地质灾害多发省份之一。图谱与几何测量在广东沿海特殊地形地貌环境是否依然适用，特别是规模小、高位隐蔽的风险点能否发现是广东防灾减灾重要研究课题。

本文选取广州黄埔区为试验区域，采用2019年1月至2022年6月Sentinel-1A长时序雷达数据，综合利用高分辨率航空影像、基础地质调查普查等资料围绕崩滑流斜坡变形灾害开展综合遥感试验，论证图谱与几何测量在中小型、高位隐蔽潜在灾害隐患的有效性，为大范围推广提供可行技术参考。

黄埔区位于广州市东部，与白云、天河、海珠、增城和从化交界，总面积约484平方千米。区内古生界变质岩台地低丘地貌，主要由石英岩、片麻岩、斜长片麻岩等组成，分布在长洲岛。第三系中新统红色砂岩主要分布于南侧平原中，主要由凝灰质砾岩、砂岩、页岩组成。第四系冲积层分布于广深公路南侧及沙洲上，主要砂砾、砾石、砂质壤土等组成。地貌类型多样，叠加汛期强降雨和人类工程活动等影响引发崩塌、滑坡、泥石流斜坡变形为主。

1.1 雷达数据预处理

星载雷达数据是大范围获取地表形变主要数据源，其成像模式、分辨率、时间基线等参数对监测结果具有较大影响。L波数据穿透植被能力较强，适合用于多植被覆盖区；C波数据具备一定植被穿透能力，可以作为L波数据的补充；X波数据对地表微形变具有高灵敏性，可实现地表微小形变的监测，但植被区失相干较为明显^[17]，如表1所示。

表1 星载合成孔径干涉雷达参数及特性Table 1 Parameters and characteristics of spaceborne InSAR

本文主要采用欧空局免费发布的C波Sentinel-1A VH极化升轨数据，IW模式采用TOPS（Terrain Observation with Progressive Scans in azimuth）技术，一景SAR影像可以成像三个子条带及多个Burst，如图1所示，并根据研究区范围进行裁剪拼接，形成整幅单视复数影像。精密轨道数据来源于Sentinel-1A官方AUX_POEORB精确的轨道星历参数。Sentinel-1A使用精密轨道控制技术，以保证干涉对的空间基线在100 m范围内，采用成像21天之后发布的POD精轨数据，定位精度5 cm以内。

图1 IW成像模式形成的子条带及Burst Fig.1 Sub-bands and Burst formed by IW imaging mode

图2 雷达影像干涉处理流程Fig.2 Interference processing flow of radar images

数字高程模型DEM数据采用分辨率相对较高免费发布的AW3D30成果。对其中范围内部分无效数值区先进行探测，根据Kriging插值方法进行修正补充，对预处理后的DEM数据进行拼接及地理编码，统一到SAR影像统一坐标系。采用粗配准和精配准联合的方式处理，以消除多普勒质心偏移，以及由轨道不平行引起的系统误差及TOPS模式方位向特殊性的问题，干涉处理方法如图2所示。使用增强谱分集（ESD）方法进行高精度配准时，需要进行相位解缠，得到的解缠差分相位为[-π, π ]，相当于Sentinel-1A0.05个像元，因此粗配准精度应达到0.05个像元。

1.2 SBAS-InSAR小基线集时序干涉测量

小基线集方法（SBAS-InSAR）通过选取多主影像干涉，更好削弱由于时间基线过程而导致的时间失相干^[18]。假设共有N景影像，组成的小基线集含有M个干涉对，ϕ^T为各时间点上SAR影像中高相干点对应的相位（相对于参考点）组成的向量，即待求参数可表示为式（1）：

ϕ^T=[ϕ( t ₁ ),…,ϕ( t_N ) ]（1）

δϕ^T表示多时相干涉图解缠后相位组成的向量，δϕ_i为观测量，即：

δϕ^T =[ δϕ₁,…,δϕ_M]（2）

主影像（IE）和从影像（IS）对应的时间序列分别为：

IE=[IE₁，…,IE_M ] IS=[IS₁,…,IS_M ]（3）

假设主从影像是按照时间顺序排列的，即IE_j＞IS_j，其中j=1，…，M ，则第j幅干涉图对应的相位可表示为：

δϕ_j=ϕ( t_IEj ) -ϕ( t_ISj ) , j=1 ,…, M （4）

对于所有干涉图，将式（4）的线性模型表示为矩阵形式：

δϕ=Aϕ （5）

其中，系数矩阵A[M×N]每一行对应于一个干涉图，每一列对应于一景SAR图像。对于式（5），有M个方程，N个未知量。如果M≥N，且A的秩是N，则利用最小二乘法可得：

当矩阵A的秩小于N时，相应的法方程系数阵A^TA秩亏，因而根据最小二乘法得到的解不唯一。为解决系数阵关联和不同基线集之间的连接引起的法方程秩亏，采用奇异值分解（SVD）方法求解。

式（7）即为第j幅干涉图的相位值等于各时段速度在主从影像时间间隔上的积分，将其写成矩阵形式，即可得到一个新的矩阵方程：

δϕ=Bv（8）

和式（5）的矩阵A一样，B也是一个M×N矩阵。对第j行，位于主辅图像获取时间之间的列B(j,k）=t_k+1-t_k ，否则B（j, k）=0。将SVD分解应用于矩阵B，就可以得到速度矢量v的最小范数解。根据各个时间区间的沉降速度，对各时段速度在时间域上进行积分即可得到各个时间段的形变量。通过引入不同的线性相位贡献参数（DEM引起的误差，轨道误差等），可精确地估计各参数，使形变估计更准确。在线性模型的基础上，继续通过对残余相位进行时空滤波分离出大气相位和非线性形变相位，利用式（9）形变相位转换将地表形变相位场转换为形变距离场。

试验采用2019年1月至2022年4月合计94景Sentinel-1A长时序雷达数据通过SBAS-InSAR时序形变信息结合地形级实景三维模型解译风险斜坡，识别出研究区崩滑流斜坡变形潜在隐患风险84处，抽取15%比例开展外业调查，5处具有显著形变的崩滑流风险。地形级实景三维模型主要基于高分辨率航空正射影像叠加2×2m格网数字高程模型DEM，可以揭示斜坡体等致灾体和交通道路、房屋等承灾体形状、位置、地形地貌特征，利用人工判读InSAR异常形变区进一步揭示形变体时序形变趋势，如图3、图4、图5所示。

图3 疑似隐患图谱特征Fig.3 Map features of suspected hidden dangers

图4 疑似隐患形变几何特征Fig.4 Deformation geometric features of suspected hidden dangers

图5 疑似隐患现场调查图Fig.5 Site investigation diagram of suspected hidden danger

从2019年1月至2022年4月多年时序形变曲线分别反映了疑似隐患点长期稳定性趋势特征，如图3位于斜坡体上，坡度较陡，年平均速率-12.83 mm，累积形变超过40 mm，经现场调查验证，斜坡致灾体主要由切坡建房所致，以粘性土壤及以低矮植被为主。图5居民房屋致灾体位于多条乡村路之间，植被茂盛以灌木为主，坡度较陡，年平均速率-10.34 mm，累积形变超过35 mm，部分房屋出现明细裂缝。

顾及InSAR形变特征的斜坡变形灾害隐患早期识别遥感方法结合了光学遥感图谱特征，以及合成孔径雷达干涉测量形变几何特征，较全面的获取疑似隐患点形状位置和形变信息，对崩滑流等斜坡体形变灾害识别是可行的。在雷达升降轨综合应用、优化算法消除广东地区大气延迟相位等方面提出改进建议：

（1）优化雷达卫星观测几何，适当选择升降轨数据可发现不同视向斜坡变形，从灾害体的正面和背面两个方向获取雷达数据，选择可一定程度穿透植被L波雷达数据发现更多隐患区域，以便尽可能降低潜在隐患漏检概率。

（2）广东地处华南地区，大气延迟误差影响较大。利用大气延迟相位模型削除干涉对的大气延迟误差。利用PS-InSAR等方法可以通过获取永久散射体PS上长期形变来有效减弱空间时间基线等失相干问题。

（3）光学遥感提供的图谱信息和InSAR提供的几何形变信息有利于识别地质灾害隐患目标，结合密集激光点云数据可揭示地表精细化结构，为隐患稳定性提供更多解译标志。

（4）基于轻量级机器学习技术的疑似隐患图斑自动化提取有利于大范围推广应用。人工解译工作需要投入大量人力资源，需要作业员丰富经验。通过构建全灾种或主要灾种的样本库，建立轻量级机器学习模式。

致谢：感谢广东省自然资源厅广东省典型地质灾害隐患早期识别试点项目支持，感谢中南大学地球科学与信息物理学院在数据处理方面提供的支持。