基于PS‑InSAR技术的地铁沿线周边环境历史沉降风险评估方法与流程

本发明涉及一种基于ps-insar技术的地铁沿线周边环境历史沉降风险评估方法。

背景技术：

城市轨道交通工程诱发房屋开裂、管线断裂、道路坍塌等周边环境事故，除了与轨道交通建设引起的周边环境变形和位移有关外，还与周边环境自身前期历史变形和位移有较大关联。轨道交通工程建设周边环境事故多发的重要原因之一是对沿线道路、建(构)筑物前期历史沉降变形情况不明，并对其风险认识不足：1)对道路、房屋建筑、桥梁等既有变位及后续允许变位富余度缺乏认识；2)对导致地表沉降变形过大背后的不良地质条件调查不足。特别是，轨道交通建设前周边建(构)筑物、道路、管线已发生较大累计变形与位移时(导致变位富余度较小)，必将给轨道交通建设周边环境保护带来极大困难和风险。

常规变形监测技术包括采用经纬仪、水准仪、测距仪、全站仪等常规测量仪器测定点的变形值，其优点是：(1)能够提供变形体整体的变形状态；(2)适用于不同的监测精度要求、不同形式的变形体和不同的监测环境；(3)可以提供绝对变形信息。但外业工作量大，布点受地形条件影响，不易实现自动化监测。特殊测量手段包括应变测量、准直测量和倾斜测量，它具有测量过程简单、可监测变形体内部的变形、容易实现自动化监测等优点，但通常只能提供局部和相对的变形信息。在轨道交通发展建设节奏愈来愈快的今天，各地政府和轨道交通的建设单位越来越需要准确、及时的地面沉降基础资料进行轨道交通工程规划，同时，也越来越需要现代高新技术提供更全面、及时、科学的监测手段。

重复轨道雷达干涉测量(insar)是一门测量地表形变(塌陷、滑坡、地震、火山运动等引起的)的有效的技术。这种ps-insar的技术，通过sar影像，可以探测到地表形变的精度达到毫米级。该技术的基本原理是：不同时期、不同角度获得的影像的相位差与地形、获取期间的地表形变和大气变化有一定的关联。通过一些点来进行测量，这些点是散射特性比较稳定的点，称为ps点，通常分布在植被稀少区，对应着建筑、公路、水利等设施。

由于轨道交通在城市环境中呈线性走向分布，在其沉降变形监测过程中，具有的特点是：1)监测距离长，一条线少则十几公里，多则几十公里；2)监测项目多，包括基坑主体结构、沿线地面、周围建(构)筑物、管线、桥梁等。因此，采用常规的精密水准进行变形监测需要耗费大量的时问、人力、物力和财力，并且很难精确确定沿线的变形影响范围。同时，由于传统的精密水准观测和gps观测存在测点分布稀疏、作业周期长、劳动强度大的问题，难以适时、客观反映日益扩大的区域性地面沉降变化趋势。

技术实现要素：

为解决上述技术问题，本发明的目的是为石油库找到一种提高了地铁沉降的精度和准确性的基于ps-insar技术的地铁沿线周边环境历史沉降风险评估方法。

本发明基于ps-insar技术的地铁沿线周边环境历史沉降风险评估方法，包括：

以n幅sar图像为输入，通过图像配准的方法将所有sar图像配准到相同的网格内；

将所有获得的sar数据组合成若干个集合，原则是：集合内的sar图像基线距小，集合间的基线距大，并通过二轨法去除由地形起伏引入的干涉相位；

在图像中选择候选ps点(psc)，并利用psc的信息补偿由大气变化引入的误差相位和由轨道数据不精确引入的误差相位；

利用补偿后的相位信息对图像中所有像素点进行逐点分析，重新识别ps点，并估计其形变信息和高程误差信息。

进一步地，sar图像配准具体包括：在配准的过程中，选择n幅sar图像中的一幅sar图像为基准图像，将其他n-1幅sar图像都配准到基准sar图像的网格内；基准图像的选择需要综合考虑空间基线和时间基线两个指标，最佳的基准图像是到其他sar图像空间基线和时间基线的加权平均值最小的那幅sar图像；在处理过程中采用三级配准的方法：(1)基于卫星轨道数据的配准；(2)基于像素级的配准；(3)基于亚像素级的配准。

进一步地，二轨法处理具体包括：先利用主辅sar图像的卫星轨道数据和外部dem信息，计算出dem每个像素点的干涉相位，并将每个像素点投影到sar图像的坐标系中；此时，dem的像素点非均匀地分布在sar图像的网格中；然后，利用delaunay三角插值的方法对sar图像的均匀网格进行重采样，获取由地形信息模拟的干涉相位图；最后，再在真实的干涉相位中减去由外部dem模拟的干涉相位。

进一步地，参考ps点(psc)选择和提取psc处的干涉相位具体包括：psc的选择采用基于幅度统计特性的选择方法，此方法利用目标点的幅度离差信息来选择ps点，幅度离差的计算公式如下：

式中，σa表示目标点在输入的n幅sar图像中幅值的标准差，ma表示目标点在输入的n幅sar图像中幅值的均值；

先设定幅度离差门限dthreshold，然后将那些满足条件da＜dthreshold的像素点选为ps点，幅度离差门限dthreshold设置为0.3；输入的sar图像尽量多于25～30幅；

选出psc后，提取psc处的干涉相位，此时，提取的干涉相位为经过二轨法处理后的相位。

进一步地，三维空时相位解缠具体包括：雷达获取的相位数据是缠绕在区间[-π,π)内的数据；因此，为了恢复目标点的真实相位，需要对相位数据进行解缠绕处理。目标点真实相位和缠绕相位关系的数学表达式如式所示：

n为整数

相位解缠处理也就是估计未知整数n的过程；

在psinsar的处理过程中，需要对空时三维进行相位解缠，在空间二维的图像域，先根据psc的位置建立delaunay三角网格，然后再利用mcf算法获取空间二维的解缠结果。

进一步地，估计并补偿大气和轨道误差相位具体包括：大气相位和轨道误差相位是随空间缓变的，建模为以空间二维坐标为自变量的一阶函数或者建立为二阶或高阶函数)

式中，待估参数有三个a，b和c，而ε和η分别表示psc在sar图像中所对应的距离和方位二维坐标。

进一步地，ps点重新识别及形变反演和高程误差估计具体包括：识别一个点是否是ps点是看这个点是否与已知的形变和高程误差模型相比匹配，判断模型匹配的方法可以依据这个点的时间相关系数，其计算公式如式所示：

式中，表示目标点p在第i幅干涉图像中的干涉相位即补偿大气和轨道误差相位后，mi(p)表示由目标点p的形变运动模型和高程误差模型估计出的相位。最后通过设置相关系数门限，将大于门限的像素点选为最终的ps点，并估计这个ps点的形变量和高程误差。

与现有技术相比，本发明基于ps-insar技术的地铁沿线周边环境历史沉降风险评估方法具有以下优点：

雷达卫星干涉测量(ps-insar)技术具有非接触测量、不需要布设监测控制网、毫米级精度、效率高、成本低、覆盖范围广、不受天气影响、空间分辨率高等特点，可克服上述常规精密水准观测和gps观测的缺陷。

本发明采用ps-insar技术获得地铁工程施工前周边环境的沉降变形历史存档数据，通过反演计算获得高精度沉降变形成果，据以对地铁沿线条带缓冲区内的地表、建筑、道桥等设施进行沉降变化趋势分析，对补充勘查、管线及房屋调查、设计、施工等方面需要重点考虑的区段予以提示，可以协助建设单位全面做好地铁风险防范工作。

本发明通过分布在沿线地表不同位置的雷达干涉测量监测点，解算分析出各监测点的三维坐标，并根据历史沉降序列情况，通过数据分析分析区域内目标设施各监测点的沉降变化量、变化趋势，并结合其它资料对沿线整体的稳定性及沉降情况进行分析，解决传统监测技术无法监测地铁周边环境历史沉降风险的问题，提高地铁沉降的精度和准确性。

上述说明仅是本发明技术方案的概述，为了能够更清楚了解本发明的技术手段，并可依照说明书的内容予以实施，以下以本发明的较佳实施例并配合附图详细说明如后。

附图说明

图1是ps-insar形变测量的处理算法基本流程图；

图2是干涉图生成示意图比较；(a)传统的psinsar干涉图生成示意图；(b)新的psinsar干涉图生成示意图；

图3是德国柏林部分区域ps点的空间二维delaunay三角网格。

具体实施方式

下面结合附图和实施例，对本发明的具体实施方式作进一步详细描述。以下实施例用于说明本发明，但不用来限制本发明的范围。

如图1所示，本发明一种基于ps-insar技术的地铁沿线周边环境历史沉降风险评估方法的最佳实施例，包括：

以n幅sar图像为输入，先通过图像配准的方法将所有sar图像配准到相同的网格内。然后，选择时间和空间基线较小的sar图像对生成干涉图，并通过二轨法去除由地形起伏引入的干涉相位。同时，在图像中选择候选ps点(psc)，并利用psc的信息补偿由大气变化引入的误差相位和由轨道数据不精确引入的误差相位。最后，利用补偿后的相位信息对图像中所有像素点进行逐点分析，重新识别ps点，并估计其形变信息和高程误差信息。

sar图像配准

由于不同sar图像之间存在空间基线，其所对应的卫星位置存在差异，导致同一目标点在不同sar图像内的位置不同，这种偏差可能达到数百个像素点的程度。而insar形变测量技术的输入是同一目标点在不同sar图像之间的相位数据。因此，在数据处理的最初阶段，需要对输入的n幅sar图像进行配准处理。

sar图像配准

由于不同sar图像之间存在空间基线，其所对应的卫星位置存在差异，导致同一目标点在不同sar图像内的位置不同，这种偏差可能达到数百个像素点的程度。而insar形变测量技术的输入是同一目标点在不同sar图像之间的相位数据。因此，在数据处理的最初阶段，需要对输入的n幅sar图像进行配准处理。

在图像中，纵轴表示空间基线，横轴表示时间基线，每个点表示一幅sar图像，每条线表示用两幅sar图像生成的干涉图。

在干涉图生成的过程中，如果两幅sar图像的空间基线或时间基线过大，生成的干涉相位有很大概率会有较大误差。一般情况下，当形变监测的时间范围较大时，如果仍然采用传统的ps-insar干涉图生成方法，不可避免地会出现两幅sar图像空间基线或时间基线过大的情况，进而势必会降低最终的形变反演精度。

为了解决这个问题，如(b)所示，在新的ps-insar干涉图生成过程中，不再仅仅采用某幅sar图像作为主图像，而是选择空间基线和时间基线较小的sar图像对生成干涉图，并且在(b)中的时间基线和空间基线平面构建网格。由于生成了更多的干涉图，干涉数据会出现冗余，多幅干涉图在时间基线和空间基线的平面构成了“闭环”，通过对“闭环”内的干涉数据进行处理，能校准相位误差，最终提高相位精度。

二轨法处理

在二轨法的实际处理过程中，先利用主辅sar图像的卫星轨道数据和外部dem信息，计算出dem每个像素点的干涉相位，并将每个像素点投影到sar图像的坐标系中。此时，dem的像素点非均匀地分布在sar图像的网格中。然后，利用delaunay三角插值的方法对sar图像的均匀网格进行重采样，获取由地形信息模拟的干涉相位图。最后，再在真实的干涉相位中减去由外部dem模拟的干涉相位。

参考ps点(psc)选择和提取psc处的干涉相位

psc的选择主要包含两类方法：(1)基于幅度统计特性的选择方法；(2)基于相关系数的选择方法。一般而言，基于幅度统计特性的选择方法应用最为广泛。此方法主要利用目标点的幅度离差信息来选择ps点。幅度离差的计算公式如下：

式中，σa表示目标点在输入的n幅sar图像中幅值的标准差，ma表示目标点在输入的n幅sar图像中幅值的均值。在实际处理过程中，先设定幅度离差门限dthreshold，然后将那些满足条件da＜dthreshold的像素点选为ps点，幅度离差门限dthreshold可设置为0.3。同时，为了保证psc选择的质量，输入的sar图像尽量多于25～30幅，否则，需要利用其它信息，采用更加复杂的方法来选择psc。

选出psc后，可以提取psc处的干涉相位，此时，提取的干涉相位为经过二轨法处理后的相位。

三维空时相位解缠

在insar测量的过程中，雷达获取的相位数据是缠绕在区间[-π,π)内的数据。因此，为了恢复目标点的真实相位，需要对相位数据进行解缠绕处理。目标点真实相位和缠绕相位关系的数学表达式如式所示：

n为整数

相位解缠处理也可以理解是估计未知整数n的过程。

相位解缠是insar数据处理中最关键的步骤，其误差会在空间和时间范围内扩散，最终会影响整个形变反演结果。相位解缠算法一般分为三类，l⁰范数法、l¹范数法和l²范数法。

l⁰范数法中最基本的是枝切法，主要思路是寻找最优的解缠路径。这类方法虽然解缠效率较高，但在某些区域可能出现解缠失败。

l¹范数法是网格流法，代表方法是统计最小费用流法。它的基本思路是寻找整幅图像的最优解，由于网格的建立较为复杂，计算效率不高，但是它解缠的精度很好。

l²范数法和l¹范数法也是寻找最优解，但它的计算效率高于网格流法，在很多时候也有很大的作用。

一般而言，基于l¹范数法的统计最小费用流算法(minimumcostflow,mcf)是性能比较优越的算法，它在解缠精度上有很大的保证，且没有解缠失败的区域。

在psinsar的处理过程中，需要对空时三维进行相位解缠，在空间二维的图像域，先根据psc的位置建立delaunay三角网格，然后再利用mcf算法获取空间二维的解缠结果。其中，下图显示了德国柏林部分区域ps点的空间二维delaunay三角网格。在时间一维域，如图2(b)所示，由于在干涉图生成的过程中在时间维建立了网格(“闭环”)，也可以利用mcf算法实现相位解缠。

估计并补偿大气和轨道误差相位

相位数据解缠后，就能对大气相位和轨道误差相位进行估计和补偿。一般而言，大气相位和轨道误差相位是随空间缓变的，因此，它们可以建模为以空间二维坐标为自变量的一阶函数(如果需要的话，也可以建立为二阶或高阶函数)：

式中，待估参数有三个a，b和c，而ε和η分别表示psc在sar图像中所对应的距离和方位二维坐标。因此，轨道误差相位和大气相位估计的问题可以转化为三个待估参数的估计问题。根据参数估计理论，可以使用最小二乘方法，估计出轨道误差相位和大气相位，并最终将其补偿。

ps点重新识别及形变反演和高程误差估计

大气和轨道误差相位补偿后，就能对sar图像的每个像素点进行逐点分析，确认其是否是ps点，并估计其形变量和高程误差。一般来说，识别一个点是否是ps点的关键是看这个点是否与已知的形变和高程误差模型相比匹配。判断模型匹配的方法可以依据这个点的时间相关系数，其计算公式如式所示：

式中，表示目标点p在第i幅干涉图像中的干涉相位(补偿大气和轨道误差相位后)，mi(p)表示由目标点p的形变运动模型和高程误差模型估计出的相位。最后通过设置相关系数门限，将大于门限的像素点选为最终的ps点，并估计这个ps点的形变量和高程误差。

本发明：无需地面测站：由于雷达卫星干涉测量监测无需地面测站，因而可使监测时空范围的设计更为自由、方便。同时，可以避免地面控制点的限制，尤其是许多中间过渡点(采用常规大地测量方法进行变形监测时，为传递坐标经常要设立许多中间过渡点)，且不必建标，从而可节省大量的人力物力，大大提高监测效率。

主动发射微波：雷达卫星干涉测量由地面控制站根据监测任务安排，制定卫星数据获取计划，卫星根据编程指令，绕行地球通过制定区域时，向地面主动发射微波并接收回波完成测量。

观测点密度高：常规监测条件下，1平方公里内的监测点数量一般为1-100个，离散孤立的监测点，仅能近似地反映区域形变的情况。雷达干涉测量监测点数平均密度可达20000个/平方公里，高密度分布的观测点，为观测区域内不同目标的形变分析提供客观数据支持，进而实现区域内连续形变特征分析。

全天候观测：雷达卫星干涉测量不受气候条件的限制，在夜晚或是风雪雨雾条件下仍能进行有效观测。这一点对于汛期的崩塌、滑坡、泥石流等地质灾害监测是非常有利的。

全自动化观测：由于雷达卫星干涉测量的数据采集工作是自动进行的，同时卫星与接收站、接收站与用户之间通过数据链路进行联系，故用户可以较为方便地把雷达卫星干涉测量监测系统建成全自动化的监测系统。这种系统不但可保证长期连续运行，而且可大幅度降低变形监测成本，提高监测资料的可靠性。

历史存档数据：全国大中城市自2011年至今5年的历史存档雷达数据，为多时序分析历史形变分析提供支撑。

支持gis可视化分析接口：雷达卫星干涉测量监测结果具备支持地理信息系统(gis)平台接口，能够帮助用户将电子表格和数据库中无法体现数据之间的模式和变化趋势以图形的方式清晰直观地表现出来，并进行空间可视化分析，实现数据可视化管理、地理信息分析以及与相关业务工作应用有机集成，从而满足用户多元化的要求。

可以获得mm级精度：mm级的精度已可满足一般崩滑体变形监测的精度要求。需要更高的监测精度时应增加观测时间和时段数正因为gps定位技术具有上述优点，因而在滑坡、崩塌、泥石流等地质灾害的监测中得到了广泛的应用，成为一种新的有效的监测手段。

以上所述仅是本发明的优选实施方式，并不用于限制本发明，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和变型，这些改进和变型也应视为本发明的保护范围。