融合北斗和InSAR数据的地表沉降监测方法及装置

摘要

本发明提供了一种融合北斗和InSAR数据的地表沉降监测方法及装置，涉及地球观测与导航的技术领域，包括：先获取监测点的北斗观测数据和PS点的SAR影像数据；然后分别对北斗观测数据和SAR影像数据进行沉降时序解算，得到北斗沉降时序和SAR沉降时序；将公共点的北斗沉降时序和SAR沉降时序进行融合，得到公共点的融合沉降时序；再构建监测区域的沉降场曲面方程；最后基于公共点的融合沉降时序和沉降场曲面方程，对其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以实现对监测区域的监测。本发明可以实现高精度，高时间分辨率，高空间分辨率，低成本，便于实现的有益效果。

权利要求

1.一种融合北斗和InSAR数据的地表沉降监测方法，其特征在于，包括：获取监测区域内监测点的北斗观测数据和监测区域内永久散射体PS点的合成孔径雷达SAR影像数据；

分别对所述北斗观测数据和所述SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序；

将所述PS点中与所述监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序；

基于所述公共点的融合沉降时序，构建监测区域的沉降场曲面方程；

基于所述公共点的融合沉降时序和所述沉降场曲面方程，对所述PS点中未与所述监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以基于所述公共点的融合沉降时序和所述其他PS点的融合沉降时序，实现对监测区域的监测。

2.根据权利要求1所述的方法，其特征在于，分别对所述北斗观测数据和所述SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序，包括：基于所述北斗观测数据，利用目标方法求取监测点的原始沉降时序，其中，所述目标方法包括：实时载波相位差分方法和精密单点定位方法；

对所述监测点的原始沉降时序进行第一预处理操作，得到所述监测点的北斗沉降时序；其中，所述第一预处理操作包括：粗差探测处理操作和剔除处理操作。

3.根据权利要求1所述的方法，其特征在于，分别对所述北斗观测数据和所述SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序，包括：在对所述SAR影像数据进行配准的过程中，按照预设稳定性条件选取PS点；

基于所述PS点对所述SAR影像数据进行第二预处理操作，得到PS点的形变相位信息；所述第二预处理操作包括：滤波操作，干涉操作和去地形处理操作；

基于所述PS点的形变相位信息，确定所述PS点的SAR沉降时序。

4.根据权利要求1所述的方法，其特征在于，将所述PS点中与所述监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序，包括：确定所述PS点中与所述监测点重合的公共点；

利用卡尔曼滤波算法对公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到所述公共点的融合沉降时序。

5.根据权利要求1所述的方法，其特征在于，基于所述公共点的融合沉降时序，构建监测区域的沉降场曲面方程，包括：基于所述公共点的融合沉降时序，提取所述公共点的属性值；其中，所述属性值包括平面坐标和沉降值；

基于所述公共点的平面坐标和沉降值，进行最小二乘法曲面拟合，得到监测区域的沉降场曲面方程。

6.根据权利要求1所述的方法，其特征在于，基于所述公共点的融合沉降时序和所述沉降场曲面方程，对所述PS点中未与所述监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，包括：从所述PS点中筛除与所述监测点重合的公共点，得到其他PS点；

从所述PS点的SAR沉降时序中获取所述其他PS点的SAR沉降时序；

基于所述其他PS点的SAR沉降时序和所述公共点的融合沉降时序，预测其他PS点在InSAR监测时刻的融合沉降时序；

基于所述沉降场曲面方程，预测所述其他PS点在非InSAR监测时刻的融合沉降时序；

将所述其他PS点在InSAR监测时刻的融合沉降时序和/或所述其他PS点在非InSAR监测时刻的融合沉降时序确定为其他PS点的融合沉降时序。

7.一种融合北斗和InSAR数据的地表沉降监测装置，其特征在于，包括：获取单元，用于获取监测区域内监测点的北斗观测数据和监测区域内永久散射体PS点的合成孔径雷达SAR影像数据；

解算单元，用于分别对所述北斗观测数据和所述SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序；

融合单元，用于将所述PS点中与所述监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序；

构建单元，用于基于所述公共点的融合沉降时序，构建监测区域的沉降场曲面方程；

计算单元，用于基于所述公共点的融合沉降时序和所述沉降场曲面方程，对所述PS点中未与所述监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以基于所述公共点的融合沉降时序和所述其他PS点的融合沉降时序，实现对监测区域的监测。

8.根据权利要求7所述的装置，其特征在于，计算单元，包括：筛除模块，用于从所述PS点中筛除与所述监测点重合的公共点，得到其他PS点；

获取模块，用于从所述PS点的SAR沉降时序中获取所述其他PS点的SAR沉降时序；

第一预测模块，用于基于所述其他PS点的SAR沉降时序和所述公共点的融合沉降时序，预测其他PS点在InSAR监测时刻的融合沉降时序；

第二预测模块，用于基于所述沉降场曲面方程，预测所述其他PS点在非InSAR监测时刻的融合沉降时序；

第一确定模块，用于将所述其他PS点在InSAR监测时刻的融合沉降时序和/或所述其他PS点在非InSAR监测时刻的融合沉降时序确定为其他PS点的融合沉降时序。

9.一种电子设备，包括存储器、处理器，所述存储器中存储有可在所述处理器上运行的计算机程序，其特征在于，处理器执行计算机程序时实现如权利要求1至6任一项所述的方法。

10.一种具有处理器可执行的非易失的程序代码的计算机可读介质，其特征在于，所述程序代码使所述处理器执行如权利要求1至6任一项所述的方法。

说明书

融合北斗和InSAR数据的地表沉降监测方法及装置

技术领域

[0001] 本发明涉及地球观测与导航的技术领域，尤其是涉及一种融合北斗和InSAR数据的地表沉降监测方法及装置。

背景技术

[0002] 我国地质灾害种类多样、分布广泛且活动频繁，是世界上地质灾害最为严重的国家之一。而位移形变是地质灾害发生的最直观的前兆信息，是地质灾害预警预报的主要依据，因此地表沉降监测是地质灾害监测最常用的技术。目前全球卫星导航系统（Global Navigation Satellite System ，GNSS）技术是连续实时获取地表三维矢量变形的主要技术手段。但是GNSS技术在大范围地表沉降监测中存在硬件成本高的缺陷，在复杂环境下存在监测精度受限的技术问题。

[0003] 针对地表沉降监测方法中的单一监测手段，其精度低，成本高，且时空分辨率低；另外，传统GNSS和InSAR结合的监测方法，因设计复杂，难以实际应用。

发明内容

[0004] 本发明的目的在于提供一种融合北斗和InSAR数据的地表沉降监测方法及装置，以缓解现有技术中存在的监测精度低，成本高，且时空分辨率低，设计复杂难以实际应用的技术问题。

[0005] 第一方面，本发明实施例提供了一种融合北斗和InSAR数据的地表沉降监测方法，其中，包括：获取监测区域内监测点的北斗观测数据和监测区域内永久散射体PS点的合成孔径雷达SAR影像数据；分别对所述北斗观测数据和所述SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序；将所述PS点中与所述监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序；基于所述公共点的融合沉降时序，构建监测区域的沉降场曲面方程；基于所述公共点的融合沉降时序和所述沉降场曲面方程，对所述PS点中未与所述监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以基于所述公共点的融合沉降时序和所述其他PS点的融合沉降时序，实现对监测区域的监测。

[0006] 进一步的，分别对所述北斗观测数据和所述SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序，包括：基于所述北斗观测数据，利用目标方法求取监测点的原始沉降时序，其中，所述目标方法包括：实时载波相位差分方法和精密单点定位方法；对所述监测点的原始沉降时序进行第一预处理操作，得到所述监测点的北斗沉降时序；其中，所述第一预处理操作包括：粗差探测处理操作和剔除处理操作。

[0007] 进一步的，分别对所述北斗观测数据和所述SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序，包括：在对所述SAR影像数据进行配准的过程中，按照预设稳定性条件选取PS点；基于所述PS点对所述SAR影像数据进行第二预处理操作，得到PS点的形变相位信息；所述第二预处理操作包括：滤波操作，干涉操作和去地形处理操作；基于所述PS点的形变相位信息，确定所述PS点的SAR沉降时序。

[0008] 进一步的，将所述PS点中与所述监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序，包括：确定所述PS点中与所述监测点重合的公共点；利用卡尔曼滤波算法对公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序。

[0009] 进一步的，基于所述公共点的融合沉降时序，构建监测区域的沉降场曲面方程，包括：基于所述公共点的融合沉降时序，提取所述公共点的属性值；其中，所述属性值包括平面坐标和沉降值；基于所述公共点的平面坐标和沉降值，进行最小二乘法曲面拟合，得到监测区域的沉降场曲面方程。

[0010] 进一步的，基于所述公共点的融合沉降时序和所述沉降场曲面方程，对所述PS点中未与所述监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，包括：从所述PS点中筛除与所述监测点重合的公共点，得到其他PS点；从所述PS点的SAR沉降时序中获取所述其他PS点的SAR沉降时序；基于所述其他PS点的SAR沉降时序和所述公共点的融合沉降时序，预测其他PS点在InSAR监测时刻的融合沉降时序；基于所述沉降场曲面方程，预测所述其他PS点在非InSAR监测时刻的融合沉降时序；将所述其他PS点在InSAR监测时刻的融合沉降时序和/或所述其他PS点在非InSAR监测时刻的融合沉降时序确定为其他PS点的融合沉降时序。

[0011] 第二方面，本发明实施例提供了一种融合北斗和InSAR数据的地表沉降监测装置，其中，包括：获取单元，用于获取监测区域内监测点的北斗观测数据和监测区域内永久散射体PS点的合成孔径雷达SAR影像数据；解算单元，用于分别对所述北斗观测数据和所述SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序；融合单元，用于将所述PS点中与所述监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序；构建单元，用于基于所述公共点的融合沉降时序，构建监测区域的沉降场曲面方程；计算单元，用于基于所述公共点的融合沉降时序和所述沉降场曲面方程，对所述PS点中未与所述监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以基于所述公共点的融合沉降时序和所述其他PS点的融合沉降时序，实现对监测区域的监测。

[0012] 进一步的，计算单元，包括：筛除模块，用于从所述PS点中筛除与所述监测点重合的公共点，得到其他PS点；获取模块，用于从所述PS点的SAR沉降时序中获取所述其他PS点的SAR沉降时序；第一预测模块，用于基于所述其他PS点的SAR沉降时序和所述公共点的融合沉降时序，预测其他PS点在InSAR监测时刻的融合沉降时序；第二预测模块，用于基于所述沉降场曲面方程，预测所述其他PS点在非InSAR监测时刻的融合沉降时序；第一确定模块，用于将所述其他PS点在InSAR监测时刻的融合沉降时序和/或所述其他PS点在非InSAR监测时刻的融合沉降时序确定为其他PS点的融合沉降时序。

[0013] 第三方面，本发明实施例提供了一种电子设备，包括存储器、处理器，所述存储器中存储有可在所述处理器上运行的计算机程序，其中，所述处理器执行所述计算机程序时实现如上述第一方面任一项所述的方法。

[0014] 第四方面，本发明实施例提供了一种具有处理器可执行的非易失的程序代码的计算机可读介质，其中，所述程序代码使所述处理器执行如上述第一方面任一项所述的方法。

[0015] 本发明提供的一种融合北斗和InSAR数据的地表沉降监测方法及装置，包括：先获取监测区域内监测点的北斗观测数据和监测区域内永久散射体PS点的合成孔径雷达SAR影像数据；然后分别对北斗观测数据和SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序；接着将PS点中与监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序；再基于公共点的融合沉降时序，构建监测区域的沉降场曲面方程；最后基于公共点的融合沉降时序和沉降场曲面方程，对PS点中未与监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以基于公共点的融合沉降时序和其他PS点的融合沉降时序，实现对监测区域的监测。本发明结合北斗导航卫星系统和InSAR技术，充分发挥北斗导航卫星系统的高精度、高时间分辨率的优点和InSAR技术高精度、高空间分辨率的优点，可以实现低成本、便于实现的有益效果。

[0016] 本发明的其他特征和优点将在随后的说明书中阐述，并且，部分地从说明书中变得显而易见，或者通过实施本发明而了解。本发明的目的和其他优点在说明书、权利要求书以及附图中所特别指出的结构来实现和获得。

[0017] 为使本发明的上述目的、特征和优点能更明显易懂，下文特举较佳实施例，并配合所附附图，作详细说明如下。

附图说明

[0018] 为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

[0019] 图1为本发明实施例提供的一种融合北斗和InSAR数据的地表沉降监测方法的流程图；

[0020] 图2为图1中步骤S102的一种流程图；

[0021] 图3为图1中步骤S102的另一种流程图；

[0022] 图4为图1中步骤S103的一种流程图；

[0023] 图5为图1中步骤S104的一种流程图；

[0024] 图6为图1中步骤S105的一种流程图；

[0025] 图7为本发明实施例提供的另一种融合北斗和InSAR数据的地表沉降监测方法的流程图；

[0026] 图8为本发明实施例提供的一种融合北斗和InSAR数据的地表沉降监测装置的结构示意图。

[0027] 图标：

[0028] 11-获取单元；12-解算单元；13-融合单元；14-构建单元；15-计算单元。

具体实施方式

[0029] 下面将结合实施例对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

[0030] 随着经济社会的快速发展，工程施工、基坑挖掘、矿业开发等人类活动不断加剧，极易引发重大地质灾害，对地质灾害进行监测预警显得迫在眉睫。而目前用于地表沉降监测的多种技术手段都在一定程度上存在缺陷与不足。全球卫星导航系统技术是连续实时获取地表三维矢量变形的主要技术手段。相比于其他卫星导航系统，北斗卫星导航系统（BeiDou Navigation Satellite System，BDS）在我国区域内具有卫星数目多、可视卫星高度角大等的特点，且其定位精度高，具备实时沉降监测的能力。但是在大范围地表沉降监测中，GNSS技术存在硬件成本高、复杂环境下监测精度受限等问题。多时相InSAR（Multi-Temporal InSAR，MT-InSAR）技术利用同一监测区域的多景SAR影像对时序上相位稳定的点（即PS点）进行精确分析，可以极大地降低由大气延迟等带来的测量误差，使得沉降监测精度达到厘米级甚至毫米级，可以对城市地表等进行大范围长时间高精度的连续监测。但是由于卫星的重访周期往往比较长，造成该技术时间分辨率低，不能够进行连续实时的监测预警，对沉降监测的实效性不强。目前针对地表沉降监测的单一的监测手段，但无法做到高精度、低成本、高时空分辨率的监测；而传统的GNSS和InSAR结合技术，也设计复杂，难以实际应用。

[0031] 基于上述现有技术的不足，本发明实施例提供了一种融合北斗和InSAR数据的地表沉降监测方法及装置，可以充分发挥北斗导航卫星系统的高精度、高时间分辨率的优点和InSAR技术高精度、高空间分辨率的优点，提高地表沉降监测的准确性。

[0032] 为便于对本实施例进行理解，首先对本发明实施例所公开的一种融合北斗和InSAR数据的地表沉降监测方法进行详细描述。

[0033] 实施例1：

[0034] 根据本发明实施例，提供了一种融合北斗和InSAR数据的地表沉降监测方法的实施例，需要说明的是，本发明实施例在附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。

[0035] 图1为本发明实施例提供的一种融合北斗和InSAR数据的地表沉降监测方法的流程图。如图1所示，该方法包括如下步骤：

[0036] 步骤S101，获取监测区域内监测点的北斗观测数据和监测区域内永久散射体PS点的合成孔径雷达SAR影像数据。

[0037] 在本发明实施例中，在监测区域内安置北斗监测接收机，该北斗监测接收机所在位置为监测点，因此该监测点可以理解为北斗监测点。北斗监测接收机用于获取所在位置的北斗观测数据，且获取方式可以为连续获取。获取到的SAR影像数据为监测区域内多景高分辨率的SAR影像。

[0038] 步骤S102，分别对北斗观测数据和SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序。

[0039] 在本发明实施例中，由于北斗观测数据和SAR影像数据具有不同的特点，因此对两者进行沉降时序解算的方式各不相同。

[0040] 步骤S103，将PS点中与监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序。

[0041] 在本发明实施例中，PS点中与监测点重合的点为公共点，也可以描述为：在监测点中与PS点重合的点为公共点，公共点既可以获取其北斗沉降时序，又可以获取其SAR沉降时序，进而可以得到公共点的融合沉降时序。

[0042] 步骤S104，基于公共点的融合沉降时序，构建监测区域的沉降场曲面方程。

[0043] 步骤S105，基于公共点的融合沉降时序和沉降场曲面方程，对PS点中未与监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序。

[0044] 在本发明实施例中，一方面在监测区域内进行空间插值，得到其他PS点在InSAR监测时刻的融合沉降时序。另一方面对其他PS点进行时间差值，得到其他PS点在非InSAR监测时刻的融合沉降时序，得到的上述两种时序都可以称为其他PS点的融合沉降时序。得到其他PS点的沉降时序的目的是，可以基于公共点的融合沉降时序和其他PS点的融合沉降时序，实现对监测区域的监测。

[0045] 为了能够对地表进行大范围、高密度、高精度且近实时的监测，弥补单一手段进行地表变形监测时的不足，本发明实施例将具有高时间分辨率特性的北斗定位技术与具有高空间分辨率特性的InSAR技术进行结合，利用时空序列数据处理方法和卡尔曼滤波技术，根据地表形变规律建立监测区域的曲面模型（沉降场曲面方程）等可以获取高时空地表形变，满足地表沉降监测的需求，为地质灾害预警预报提供依据。

[0046] 如图1和图7所示，本发明实施例提供的一种融合北斗和InSAR数据的地表沉降监测方法，包括：先获取监测区域内监测点的北斗观测数据和监测区域内永久散射体PS点的合成孔径雷达SAR影像数据；然后分别对北斗观测数据和SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序；接着将PS点中与监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序；再基于公共点的融合沉降时序，构建监测区域的沉降场曲面方程；最后基于公共点的融合沉降时序和沉降场曲面方程，对PS点中未与监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以基于公共点的融合沉降时序和其他PS点的融合沉降时序，实现对监测区域的监测。本发明结合北斗导航卫星系统和InSAR技术，充分发挥北斗导航卫星系统的高精度、高时间分辨率的优点和InSAR技术高精度、高空间分辨率的优点，可以实现低成本、便于实现的有益效果。因此本发明实施例提供的一种融合北斗和InSAR数据的地表沉降监测方法可以成为地表沉降监测的一种科学合理且高效的技术手段。

[0047] 在一个可选的实施例中，如图2所示，步骤S102，分别对北斗观测数据和SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序，包括以下步骤：

[0048] 步骤S201，基于北斗观测数据，利用目标方法求取监测点的原始沉降时序。

[0049] 在本发明实施例中，目标方法包括：实时载波相位差分方法（Real  Time Kinematic，RTK）和精密单点定位方法（Precise Point Position，PPP）。无论采用上述两种方法中的哪一种方法，均可以获取到监测点的原始沉降时序。

[0050] 步骤S202，对监测点的原始沉降时序进行第一预处理操作，得到监测点的北斗沉降时序。

[0051] 在本发明实施例中，第一预处理操作包括：粗差探测处理操作和剔除处理操作。北斗沉降时序也可以称之为北斗沉降监测时序。对原始沉降时序进行粗差探测与剔除，以获得高可靠性的北斗沉降时序。

[0052] 在一个可选的实施例中，如图3所示，步骤S102，分别对北斗观测数据和SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序，还包括以下步骤：

[0053] 步骤S301，在对SAR影像数据进行配准的过程中，按照预设稳定性条件选取PS点；

[0054] 步骤S302，基于PS点对SAR影像数据进行第二预处理操作，得到PS点的形变相位信息；第二预处理操作包括：滤波操作，干涉操作和去地形处理操作；

[0055] 在本发明实施例中，第二预处理操作包括：滤波操作，干涉操作和去地形处理操作。

[0056] 步骤S303，基于PS点的形变相位信息，确定PS点的SAR沉降时序。

[0057] 在本发明实施例中，先通过对SAR影像数据的配准、滤波、干涉等操作，然后根据差分干涉思想去除地形相位，得到高精度的PS点的形变相位信息，最后通过多时相InSAR技术对同一监测区域内的多景SAR影像进行时序稳定点（PS点）处理，进而获得该监测区域内PS点高精度的地表沉降时序（即PS点的SAR沉降时序）。

[0058] 在一个可选的实施例中，如图4所示，步骤S103，将PS点中与监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序，包括以下步骤：

[0059] 步骤S401，确定PS点中与监测点重合的公共点；

[0060] 步骤S402，利用卡尔曼滤波算法对公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序。

[0061] 在本发明实施例中，步骤S402，可以理解为是对公共点的北斗沉降时序和公共点的SAR沉降时序进行的沉降值的融合估计，也可以理解为对公共点的北斗沉降时序和公共点的SAR沉降时序进行的多传感器数据融合处理。将M个监测点在同一时刻的沉降量（沉降值）所组成的向量看作是状态向量，记为 ，步骤S402的具体步骤如下：

[0062] 预估计：

[0063] （1）

[0064] 其中， 为状态向量在当前历元的先验估计值， 为k-1时刻到k时刻的状态转移矩阵，一般为单位阵， 为在前一历元卡尔曼滤波器对系统状态 的最优
估计值。

[0065] 计算预估计协方差矩阵：

[0066] （2）

[0067] 其中， 为当前历元的先验协方差阵， 为前一历元的协方差阵， 为系统噪声的协方差阵。

[0068] 计算卡尔曼滤波增益矩阵：

[0069] （3）

[0070] 其中， 为卡尔曼滤波滤波增益， 为设计矩阵， 为量测噪声。

[0071] 更新估计：

[0072] （4）

[0073] 其中， 为k时刻公共点的SAR沉降时序中的沉降值。

[0074] 计算更新后的估计协方差矩阵：

[0075] （5）

[0076] 在以上公式（1）（5）中， 为待估参数，即公共点的融合沉降时序。在给定初值~和 的条件下，卡尔曼滤波器根据k时刻的观测值 进行递推计算，得到k时刻的状
态估计 。在估计过程中，利用公共点的SAR沉降时序及时对状态
参数进行修正，重复上述步骤进行递归计算，即可得到公共点的融合沉降时序。在本发明实施例中，公共点的融合沉降时序可以为其他PS点的融合沉降时序的计算作铺垫。

[0077] 在一个可选的实施例中，如图5所示，步骤S104，基于公共点的融合沉降时序，构建监测区域的沉降场曲面方程，包括以下步骤：

[0078] 步骤S501，基于公共点的融合沉降时序，提取公共点的属性值；其中，属性值包括平面坐标和沉降值；

[0079] 在本发明实施例中，从公共点的融合沉降时序中提取所有公共点的平面坐标和沉降值 。其中，公共点的融合沉降时序可以指经过地理编码后的干涉图产品。

[0080] 步骤S502，基于公共点的平面坐标和沉降值，进行最小二乘法曲面拟合，得到监测区域的沉降场曲面方程。

[0081] 在本发明实施例中，先获取公共点的融合沉降时序，然后执行上述步骤S501 步骤~S502，上述步骤S501 步骤S502即为三维曲面拟合沉降场的具体步骤，其中步骤S502，对公~
共点的属性值 进行最小二乘法拟合的方式具体如下：通常假设地表沉降场
是一个二次曲面，最小二乘法曲面拟合所用的曲面拟合方程为：

[0082] （6）

[0083] 其中， 为6个待求参数，且在不同InSAR监测时刻具有不同的参数值。

[0084] 为了求取待求参数 ，先将地表沉降场的第i个公共点的属性值 代入到公式（6）中，得到该公共点的拟合沉降值 ，然后求取拟合沉
降值 与真实沉降值 之间的误差 ，其中，

[0085] （7）

[0086] 由公式（7）可以得到曲面拟合方程的间接平差函数模型：

[0087] （8）

[0088] 其中，A为由 组成的矩阵， 为第i个公共点的真实沉降值， 为由 构成的矩阵。在沉降监测中，由于每个公共点的沉降值对二次曲面的影响是相同的，故可以取等权值。根据最小二乘原理进行求解，可得到曲面拟合方程的解：

[0089] （9）

[0090] 由上述步骤，可以获得曲面拟合方程的系数（ ）。在后续的计算中，只需要将检测区域内任一PS点的平面坐标（ ）代入公式（6），即可得到该PS点的融合沉降时序，在空间维度上实现了高分辨率。

[0091] 在一个可选的实施例中，如图6所示，步骤S105，基于公共点的融合沉降时序和沉降场曲面方程，对PS点中未与监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，包括以下步骤：

[0092] 步骤S601，从PS点中筛除与监测点重合的公共点，得到其他PS点；

[0093] 步骤S602，从PS点的SAR沉降时序中获取其他PS点的SAR沉降时序；

[0094] 步骤S603，基于其他PS点的SAR沉降时序和公共点的融合沉降时序，预测其他PS点在InSAR监测时刻的融合沉降时序；

[0095] 步骤S604，基于沉降场曲面方程，预测其他PS点在非InSAR监测时刻的融合沉降时序；

[0096] 步骤S605，将其他PS点在InSAR监测时刻的融合沉降时序和其他PS点在非InSAR监测时刻的融合沉降时序确定为其他PS点的融合沉降时序。

[0097] 在本发明实施例中，根据公共点的融合沉降时序和沉降场曲面方程即可实现对监测区域内其他PS点的时空插值。时空插值一般从两个维度上进行分析，一个是在时间维度上进行的插值，一个是在空间维度上进行的插值。其中时间维度上进行的插值，目的是得到其他PS点在InSAR监测时刻的沉降序列，而空间维度上进行的差值，目的是得到其他PS点在非InSAR监测时刻的沉降值。

[0098] 对时间维度上进行的插值进行如下分析：其他PS点的SAR沉降时序在时间上并不连续，因此时间分辨率较低。而本申请为了提高时间分辨率，先确定其他PS点的SAR沉降时序，然后基于其他PS点的SAR沉降时序和公共点的融合沉降时序，可以得到其他PS点在各个InSAR监测时刻的沉降时序。其中，公共点的SAR沉降时序为 ，其他PS点的SAR沉降时序为 ，公共点的融合沉降时序为 ，其他PS点在InSAR监测时刻的融合沉降时序为： 。对空间维度上进行的
插值进行如下分析：在非InSAR监测时刻，先使用公共点的融合沉降时序进行曲面拟合，确定曲面拟合方程的拟合系数（即待求参数），然后将其余PS点的坐标代入方程，即可得到其他PS点在非InSAR监测时刻的沉降值。

[0099] 通过步骤S601 步骤S604，可以获得地表监测区域内高时空分辨率的沉降时序，该~沉降时序融合了北斗技术的高时间分辨率和InSAR技术的高空间分辨率，包含了丰富的时空信息。对监测区域范围内的任意PS点进行提取，便可以得到该点连续的沉降值，有利于揭示其沉降演化规律。

[0100] 本实施例提供的融合北斗和InSAR数据的地表沉降监测方法融合了北斗定位技术与InSAR技术技术，可以充分发挥二者在点面监测、时空分辨率上的互补性，能够提高大范围沉降监测的效率，降低沉降监测成本，更加具有可行性和合理性。本发明实施例利用kalman滤波融合算法计算公共点在连续时间上的沉降值；然后以最小二乘曲面拟合方法确定InSAR监测时刻的沉降场曲面方程；再通过多景SAR影像确定的多个时刻的高相干点进行时空插值运算，得到高精度、高时空分辨率的地表形变数据序列，进而可以建立高精度的地质灾害形变场。本发明实施例结合了北斗定位技术与InSAR技术的优势，建立了两者之间的耦合关系，成本低、便于实现，能够更加精确细致地监测大范围地表沉降，可为后续的地质灾害预警预报提供科学依据。本实施例还可以用于GNSS其他系统（GPS、GLONASS、GALILEO）与InSAR融合数据处理。

[0101] 实施例2：

[0102] 本发明实施例还提供了一种融合北斗和InSAR数据的地表沉降监测装置，该融合北斗和InSAR数据的地表沉降监测装置主要用于执行本发明实施例上述内容所提供的融合北斗和InSAR数据的地表沉降监测方法，以下对本发明实施例提供的融合北斗和InSAR数据的地表沉降监测装置做具体介绍。

[0103] 图8为本发明实施例提供的一种融合北斗和InSAR数据的地表沉降监测装置的结构示意图。如图8所示，该融合北斗和InSAR数据的地表沉降监测装置主要包括获取单元11，解算单元12，融合单元13，构建单元14和计算单元15，其中：

[0104] 获取单元11，用于获取监测区域内监测点的北斗观测数据和监测区域内永久散射体PS点的合成孔径雷达SAR影像数据；

[0105] 解算单元12，用于分别对北斗观测数据和SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序；

[0106] 融合单元13，用于将PS点中与监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序；

[0107] 构建单元14，用于基于公共点的融合沉降时序，构建监测区域的沉降场曲面方程；

[0108] 计算单元15，用于基于公共点的融合沉降时序和沉降场曲面方程，对PS点中未与监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以基于公共点的融合沉降时序和其他PS点的融合沉降时序，实现对监测区域的监测。

[0109] 本发明实施例提供的一种融合北斗和InSAR数据的地表沉降监测装置，先利用获取单元11获取监测区域内监测点的北斗观测数据和监测区域内PS点的合成孔径雷达SAR影像数据；然后利用解算单元12分别对北斗观测数据和SAR影像数据进行沉降时序解算，得到监测点的北斗沉降时序和PS点的SAR沉降时序；接着利用融合单元13将PS点中与监测点重合的公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序；再利用构建单元14基于公共点的融合沉降时序，构建监测区域的沉降场曲面方程；最后利用计算单元15基于公共点的融合沉降时序和沉降场曲面方程，对PS点中未与监测点重合的其他PS点进行时空插值计算，得到其他PS点的融合沉降时序，以基于公共点的融合沉降时序和其他PS点的融合沉降时序，实现对监测区域的监测。本发明实施例结合北斗导航卫星系统和InSAR技术，充分发挥北斗导航卫星系统的高精度、高时间分辨率的优点和InSAR技术高精度、高空间分辨率的优点，可以实现低成本、便于实现的有益效果。

[0110] 可选地，计算单元，包括筛除模块，第一预测模块，获取模块，第二预测模块和第一确定模块，其中：

[0111] 筛除模块，用于从PS点中筛除与监测点重合的公共点，得到其他PS点；

[0112] 获取模块，用于从PS点的SAR沉降时序中获取其他PS点的SAR沉降时序；

[0113] 第一预测模块，用于基于其他PS点的SAR沉降时序和公共点的融合沉降时序，预测其他PS点在InSAR监测时刻的融合沉降时序；

[0114] 第二预测模块，基于沉降场曲面方程，预测其他PS点在非InSAR监测时刻的融合沉降时序；

[0115] 第一确定模块，用于将其他PS点在InSAR监测时刻的融合沉降时序和其他PS点在非InSAR监测时刻的融合沉降时序确定为其他PS点的融合沉降时序。

[0116] 可选地，解算单元，包括求取模块和第一预处理模块，其中：

[0117] 求取模块，用于基于北斗观测数据，利用目标方法求取监测点的原始沉降时序，其中，目标方法包括：实时载波相位差分方法和精密单点定位方法；

[0118] 第一预处理模块，用于对监测点的原始沉降时序进行第一预处理操作，得到监测点的北斗沉降时序；其中，第一预处理操作包括：粗差探测处理操作和剔除处理操作。

[0119] 可选地，解算单元，还包括选取模块，第二预处理模块和第二确定模块，其中：

[0120] 选取模块，用于在对SAR影像数据进行配准的过程中，按照预设稳定性条件选取PS点；

[0121] 第二预处理模块，用于基于PS点对SAR影像数据进行第二预处理操作，得到PS点的形变相位信息；第二预处理操作包括：滤波操作，干涉操作和去地形处理操作；

[0122] 第二确定模块，用于基于PS点的形变相位信息，确定PS点的SAR沉降时序。

[0123] 可选地，融合单元，可以包括：第三确定模块和融合模块，其中：

[0124] 第三确定模块，用于确定PS点中与监测点重合的公共点；

[0125] 融合模块，用于利用卡尔曼滤波算法对公共点的北斗沉降时序和公共点的SAR沉降时序进行融合，得到公共点的融合沉降时序。

[0126] 可选地，构建单元，包括提取模块和拟合模块，其中：

[0127] 提取模块，用于基于公共点的融合沉降时序，提取公共点的属性值；其中，属性值包括平面坐标和沉降值；

[0128] 拟合模块，用于基于公共点的平面坐标和沉降值，进行最小二乘法曲面拟合，得到监测区域的沉降场曲面方程。

[0129] 本发明实施例所提供的装置，其实现原理及产生的技术效果和前述方法实施例相同，为简要描述，装置实施例部分未提及之处，可参考前述方法实施例中相应内容。

[0130] 进一步地，本实施例还提供了一种电子设备，包括存储器、处理器，该存储器中存储有可在上述处理器上运行的计算机程序，其中，上述处理器执行计算机程序时执行前述方法实施例所提供的方法的步骤。

[0131] 进一步地，本实施例还提供了一种具有处理器可执行的非易失的程序代码的计算机可读介质，该程序代码使上述处理器执行前述方法实施例所提供的方法的步骤。

[0132] 本发明实施例所提供的融合北斗和InSAR数据的地表沉降监测方法及装置的计算机程序产品，包括存储了程序代码的计算机可读存储介质，所述程序代码包括的指令可用于执行前面方法实施例中所述的方法，具体实现可参见方法实施例，在此不再赘述。

[0133] 另外，在本发明实施例的描述中，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以具体情况理解上述术语在本发明中的具体含义。

[0134] 在本实施例的描述中，需要说明的是，术语“中心”、“上”、“下”、“左”、“右”、“竖直”、“水平”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本实施例的限制。此外，术语“第一”、“第二”、“第三”仅用于描述目的，而不能理解为指示或暗示相对重要性。

[0135] 所属领域的技术人员可以清楚地了解到，为描述的方便和简洁，上述描述的装置和单元的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。

[0136] 在本实施例所提供的几个实施例中，应该理解到，所揭露的装置和方法，可以通过其它的方式实现。以上所描述的装置实施例仅仅是示意性的，例如，所述单元的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，又例如，多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另一点，所显示或讨论的相互之间的耦合或直接耦合或通信连接可以是通过一些通信接口，装置或单元的间接耦合或通信连接，可以是电性，机械或其它的形式。

[0137] 所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本实施例方案的目的。

[0138] 另外，在本发明各个实施例中的各功能单元可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中。

[0139] 所述功能如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个处理器可执行的非易失的计算机可读取存储介质中。基于这样的理解，本实施例的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备（可以是个人计算机，服务器，或者网络设备等）执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：U盘、移动硬盘、只读存储器（ROM，Read-Only Memory）、随机存取存储器（RAM，Random Access Memory）、磁碟或者光盘等各种可以存储程序代码的介质。

[0140] 最后应说明的是：以上所述实施例，仅为本发明的具体实施方式，用以说明本发明的技术方案，而非对其限制，本发明的保护范围并不局限于此，尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，其依然可以对前述实施例所记载的技术方案进行修改或可轻易想到变化，或者对其中部分技术特征进行等同替换；而这些修改、变化或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围，都应涵盖在本发明的保护范围之内。