[发明专利]一种基于多源数据集成的高铁沉降危害性评估方法

权利要求书

1.一种基于多源数据集成的高铁沉降危害性评估方法，其特征在于，包括以下步骤：

S1、利用InSAR技术采集地面形变数据，包括采集InSAR影像数据、光学遥感影像数据、水文地质资料数据、二级水准测量数据；

S2、对步骤S1中采集的数据进行预处理；

S3、采用PS-InSAR和SBAS-InSAR两种时序InSAR技术分析提取地面形变结果；

S4、针对高铁线路背景区域沉降分离；

S5、对影响地面形变的危害因子采用MIC算法赋权重；

S6、地面沉降危害性评估；

所述步骤S3中，针对有大量为稳定散射体的高铁沿线的处理方法为：进行PS-InSAR流程处理，在多时序SAR影像中选取时间空间基线小的影像为主影像，通过相干系数阈值和振幅离差阈值两种方法对影像进行选点处理，选取稳定的PS点；

针对有大量荒地或逐步建成的区域的处理方法为：进行SBAS-InSAR流程处理，除去必要的去噪之外，不必进行主影像选取，在影像基线优化处理时，引入最小二乘法优化干涉对组，在保障影像处理质量的同时，减少冗余干涉对，减轻计算负担；

进行PS-InSAR流程处理的具体方法如下：

将InSAR影像数据进行PS-InSAR处理，需要在N+1景影像中选取1景主影像m，其余N景InSAR影像将与主影像匹配，生成N幅干涉图，基于PS处理流程的特殊性，一般要求N≧20，且主影像的时间空间基线在N+1景影像中最小，在选取主影像时，要建立相关函数模型：

上述公式(1)、(2)中，是SAR影像N和M的空间基线，T^n,m为影像的时间基线，为多普勒基线，B_c，T_c，和f_c为临界基线，γ^mmax的SAR影像为主影像，则由式(1)可知，当空间基线和时间基线值最小时，γ^m值最大，因此主影像的时间空间基线应是多景SAR影像中最小的；

对于地形相位和平地相位，PS-InSAR流程中是采取加入DEM来减少影响，在干涉相位有5部分组成，在DEM校正中上述两种相位将被去除，其中高程上由于DEM选用的精度会存在误差，因此对于干涉图N的像元，其相位为：

式(3)中，为参考DEM校正中由于DEM精度引起的残余地形相位，为大气相位，为形变相位，为噪声，由式(3)可知，DEM精度越高则残余地形相位会偏小，则有：

将(4)，(5)中的带入(3)可知：

对于(6)中的大气相位和噪声非线性相位，在覆盖大量混凝土的城市、道路或高速铁路稳定地表的区域，二者的影像可视为线性形变；

经过上述步骤，会产生永久散射点即PS点，除了有形变信息的PS点外，还有许多对形变观测无用的点，所以需对PS点进行筛选；在本步骤中采取相干系数阈值法和振幅离差阈值法对PS点进行筛选，选取相干性高的PS点作为PS点，其中相干系数阈值法以一个局部矩阵进行像元的相干性计算，具体为：

式(7)中，γ是相干系数，M，S是两景SAR影像，m，n是局部矩阵的行列数，γ计算结果值为(0,1)，相干系数阈值是在进行运算之前给定阈值，大于阈值选为PS点，且相干系数值越接近于1表示其相干性越强；

相反的，振幅离差阈值是针对影像的单个像元，给定一个阈值范围，以像元振幅的稳定性来选取PS点，具体表现：

式(8)中，σ_φ和D_A表示相位标准差和振幅离差指数，μ_A为均值振幅，σ_A是振幅标准差，σ_nR为噪声实部标准差，σ_nl为噪声虚部标准差；

进行选点之后就可以进行解算，其中包括去除大气延迟、轨道误差，然后进行线性形变估算，最终得到PS点的形变速率，对其进行地理编码即可得到形变结果；

进行SBAS-InSAR流程处理的具体方法如下：

不同于PS-InSAR流程处理，SBAS-InSAR流程处理无需对N景SAR影像选取主影像，在用于生成干涉图的SAR影像时选取时间空间基线较小的SAR影像形成干涉组对，以减少失相干对结果的影响，N景SAR影像生成M幅干涉图，干涉图的相位表示：

式(9)中，x,r是干涉图坐标，t_A,t_B为生成干涉图的两景影像的获取时间，若t_At_B，式(9)中第一项即为形变相位，式中的d(t_B,x,r),d(t_A,x,r)为t_A,t_B时刻(x,r)相对于t0时刻视线方向的累积形变量；第二项即为残余高程相位；第三项为大气延迟相位；余式为噪声相位；若假设M景主影像和从影像的获取时间以序列IE＝[IE₁...IE_M]与IS＝[IS₁...IS_M],若对M幅干涉图先不考虑大气延迟，只求解高程误差和形变量，则：

且将形变相位转化为平均相位速度矢量[v]^T：

若速度矢量是线性变换的，即v＝Mp，则有矩阵方程：

DMp+CΔz＝Δφ (12)

式(12)中D是M×(N-1)的矩阵，p为速度模型的参数矢量，则参数个数i应小于等于(N-1)，M为(N-1)×i的矩阵；Δz为高程相位误差；对于式(12)采用最小二乘法求解，可得到形变相位和高程误差相位，并可以去除相位图中高程相位误差，则矩阵方程可修改为：

Dv+CΔz＝Δφ (13)

在SBAS处理中，由于考虑到空间失相干的问题，将影像分成若干基线对，导致式(13)中的D为秩亏矩阵，若基线对为L(L1)，则D的秩为N-L，使用SVD法对式(13)求解，对得到的速度矢量积分即可获得形变相位序列，需要注意的是此时结果中还有大气相位，采取高通滤波去除大气相位，得到研究区的形变序列,即研究区的地面形变信息；

所述步骤S4的具体方法如下：

S401、对步骤S3中采集的地面形变数据结果进行校准，针对高速铁路采用PS-InSAR流程处理，将结果与二等水准测量数据进行校正；

S402、对形变区域进行高频低频滤波的小波变换分析，用于反映出沉降区域的细节部分与趋势部分的差异；

其中，通过高频滤波组的用于表达高铁线路及周围沉降区域细节性的信息；通过低频滤波的用于反映整区域趋势性相关信息；

所述步骤S402中，小波变换分析的具体方法如下：

利用小波变换完成一次变换，将影像数据分成高频图组与低频图组，其中高频组会更好的反映锐化的部分；

对于一维的信号：

定义

其逆变换为：

上式中Ψ为傅里叶变换，C_Ψ取值为：

对于栅格影像的二维的小波变换，有一个前提条件是尺度可分离的，即信号函数Φ(x,y)＝Φ(x)Φ(y),其中Φ(x)Φ(y)分别对应的小波函数为Ψ(x)Ψ(y)其创立二维小波变换的三个基本小波可表达为：

ψ¹(x,y)＝φ(x)ψ(y),ψ²(x,y)＝φ(y)ψ(x),ψ³(x,y)＝ψ(x)ψ(y) (16)

式中构建的二维平方可积函数空间L²(R²)的正交归一，既有：

j≧0,l＝1,2,3,...且j,l,m,n均为整数；

对于一幅N×N的二维栅格影像f(x,y)开始，其中上标尺度N，N＝2n，对于j＝0，尺度2^j＝2⁰＝1，为源影像尺度；j每增大一次使尺度加倍，分辨率减半，影像分为三个高频组，一个低频组，低频组中包括大部分能量，还需进行分离，高频组则更好的反映锐化凸显线条的影像；对于第一个层次(j＝1)，则有：

通过小波变换分解后重构，将InSAR影响处理结果中有关高速铁路的部分提取出来，高频组中将高铁及沿线研究区的信息更加突出的表现；

所述步骤S5的具体方法如下：

S501、针对步骤S1中采集的数据，选择影响高铁行车地面形变的危害因子，包括地下水位、土壤压缩模量、交通流量、压缩层厚度；

S502、利用GIS软件工具，将数据进行统一的栅格化处理，对于数据点，首选克里金插值法进行插值，对于遥感地面沉降数据、光学遥感数据栅格数据，直接针对栅格灰度值进行栅格计算；

S503、将危害因子统一转为栅格数据后，根据MIC算法，对栅格灰度值之间进行计算，从而完成危害因子赋权重；

所述步骤S6的具体方法如下：

S601、基于步骤S5中得到的赋权重的危害因子与地面沉降图并结合随机森林算法建立决策树，并且将处理后的光学遥感影像数据进行分类，用以对高速铁路沿线不同地段的危害性进行评估；

S602、结合Boruta算法辅助修改各级决策树，Boruta算法通过对决策树建立副本并对比新建随机副本与原始变量的重要性对决策树进行修改，用于对高铁沿线行车安全地段进行提取标记，此方法可对赋值的栅格影像进行计算，根据权重提取特征地段，实现了对高铁沿线沉降危害性的分区与评估；

在执行步骤S601前，在建立决策树前，需对影像建立训练样本，通过构建训练样本以确定特征属性和目标变量的关系；采用CART算法构建决策树训练样本，CART算法通过基尼系数确定属性节点分割属性，其中选取属性基尼系数最小的进行二分化，最后构建一个二叉树的决策树：对于一个栅格数据集D，它的基尼系数为：

式中p_i是数据集D中任意数据被抽取成为训练样本的概率，对生成的二分数据集进行二分则有：

式中D₁，D₂为二分后的两个数据集，D₁为可进行二分的子数据集，并对D₁进行式(22)的计算再进行二分，并在节点处形成训练规则，直到所有二分子集都符合节点的规则时，结束算法，若有不符合的，可对其再次进行式(22)的计算，直到结束算法。

2.根据权利要求1所述的基于多源数据集成的高铁沉降危害性评估方法，其特征在于：所述步骤S2中，对采集的数据进行预处理的方法包括：

对光学遥感影像数据的辐射定标、大气校正、几何校正；

对水文地质资料数据的插值化、矢量化栅格化处理；

对InSAR影像数据的去噪、多视、选取主影像。