CN115032637A - 一种井工开采全生命周期地表沉陷监测方法 - Google Patents

Abstract

本发明公开了一种井工开采全生命周期地表沉陷监测方法，属于开采沉陷监测领域。利用不同时段多源SAR数据和时序InSAR技术生成矿区地表雷达视线向形变，并通过地理编码转换至相同坐标系下，结合SAR影像成像几何关系将视线向形变转为垂直向形变；当不同时段SAR影像存在时间重叠时，利用最小二乘方法对多源SAR数据监测结果进行平差；反之利用概率积分法对多源SAR数据监测结果进行融合；同时对获取的开采沉陷全生命周期形变场进行精度评估，若不满足精度要求，调整优化融合参数，重新进行数据融合处理，反之输出矿区全生命周期地表沉陷。其步骤简单，融合不同SAR影像获取开采全周期地表沉陷信息，克服了单一SAR影像源仅能获取部分时段地表沉陷的缺陷。

Description

一种井工开采全生命周期地表沉陷监测方法

技术领域

本发明涉及一种地表沉陷监测方法，尤其适用于种井工开采全生命周期地表沉陷监测方法，属于开采沉陷监测领域。

技术背景

虽然在煤炭供给侧改革和“碳中和”背景下，我国淘汰关闭煤矿数量剧增，但煤炭仍是我国的基础能源和重要原料，在我国能源体系中，煤炭一直起到压舱石和兜底保障作用。然而，随着煤炭资源的开发利用，在给社会带来经济效益的同时，也产生了一系列的环境地质灾害，如地面沉陷、地表塌陷、滑坡、地裂缝等，其沉陷灾害具有长期性、突陷性、隐蔽性和复杂性。因此，监测矿区全生命周期地表沉陷对研究地表沉陷时空演化规律、沉陷预测模型及灾害防控，促进煤矿区全生命周期安全绿色开发具有重要的理论和实际意义。

矿区全生命周期地表沉陷可分为三个阶段，分别为开采阶段、废弃(老)采空区阶段和矿井关闭阶段，地表沉陷持续时间可达十几年甚至上百年。目前矿区地表沉陷监测方法主要有常规地面变形监测和航天航空测量方法。常规地面变形监测方法主要包括水准测量、GNSS等，其测量精度高，但监测成本高、劳动强度大，不便于进行大区域、长时序地表变形监测。航空航天测量方法包括航空摄影测量、合成孔径雷达干涉测量(InterferometricSynthetic Aperture Radar,InSAR)等。航空摄影测量目前多采用无人机航空摄影测量系统，其优势是监测范围大、通过图像可再现变形体变形前状态，缺点是监测精度低，仅到分米级，不能满足开采沉陷监测要求。而InSAR技术具备全天时、全天候、大范围、高精度监测地表沉陷的能力，它利用SAR复数数据的相位信息，以毫米级或亚厘米级精度提取地表沉陷信息，同时利用存档SAR数据可追溯历史变形信息，是监测矿区全生命周期地表沉陷的首选方法。然而，由于SAR卫星平台寿命有限，导致单一平台SAR影像难以覆盖矿区全周期地表沉陷，仅能获取地表沉陷部分片段，无法全面反映全周期地表沉陷的时空演化规律。因此，融合不同时段多源SAR影像监测结果获取矿区全生命周期地表沉陷，可为揭示矿区地表沉陷时空演化规律、构建地表沉陷预测模型提供技术支持。

发明内容

针对上述技术问题，提供一种井工矿全生命周期地表沉陷监测方法，其利用不同时段多源SAR数据和时序InSAR技术获取矿区地表沉陷；并利用最小二乘法或概率积分法对存在或不存在时间重叠的多源SAR数据监测结果进行融合；同时对反演的全生命周期地表沉陷进行精度评估，获取最终的矿区全生命周期地表沉陷，有效地克服了单一SAR影像仅能获取部分地表沉陷的缺陷。

为实现上述技术目的，本发明的一种井工开采全生命周期地表沉陷监测方法，包括如下步骤：

步骤一、首先利用不同时段获取的多源SAR数据和时序InSAR技术生成矿区地表雷达视线向形变，然后对矿区地表雷达视线向形变进行地理编码转换至相同坐标系下，结合SAR影像成像几何关系将矿区地表雷达视线向形变转为垂直向形变；

步骤二、当不同时段SAR影像存在时间重叠时，利用最小二乘方法对多源SAR数据监测结果进行平差；若不存在时间重叠时，则利用概率积分法对多源SAR数据监测结果进行融合；

步骤三、同时对获取的开采沉陷全生命周期形变场进行精度评估，若不满足精度要求，调整最小二乘法收敛阈值或概率积分法预测参数，重新利用最小二乘法或概率积分法对多源SAR数据监测结果进行融合，其中收敛阈值越小精度越高，概率积分法预测参数通过调整求参点数量进行优化；若满足精度要求，则输出矿区全生命周期地表沉陷，实现了融合不同时段多源SAR影像获取井工开采全生命周期地表沉陷。

进一步，步骤一具体为：利用不同平台卫星获取被测因井工开采导致地表沉陷的矿区不同时段SAR影像，对多源SAR影像分别利用时序InSAR技术获取被测矿区不同时段地表视线向形变，同时对多源SAR影像监测结果分别进行地理编码，将监测结果统一至相同的标准坐标系下，并结合各SAR影像成像几何关系，将被测矿区不同时段视线向形变转换至垂直向形变，在不同时段SAR影像监测结果中选取同名的相干点，即选取相同位置上的相干点进行处理，以获取时间上完整的地表形变。

进一步，将不同时段SAR影像监测结果获取的地表沉陷归化至全周期地表沉陷部分具体步骤为：

设共获取到K个不同时段SAR影像监测结果，因每一个SAR影像集均以第一景SAR影像为参考影像获取地表时序形变，即：

其中

表示第i个SAR影像集中第一景影像时间，

为第i个SAR影像集中起始时间的形变量；因此，在每组SAR影像监测结果间均存在偏移量：Λⁱ(i＝1,2,…,K)，因每个相干点的形变速率不同，故不同相干点的偏移量不同；

为了获取全周期地表沉陷，将每组SAR影像中的SAR数据按照获取时间进行排序，并将获取时间最早的SAR数据定义为第一组SAR影像，将第一组SAR影像作为全局参考影像，即Λ¹≡0；然后通过偏移量Λⁱ(i＝1,2,…,K)，利用下式将不同时段SAR影像获取的地表形变起始时间归化至第一组SAR影像：

Disp_full(tⁱ)＝Disp(tⁱ)+Λⁱ(i＝1,2,…,K)，

式中Disp_full(tⁱ)表示第i个SAR影像集获取的地表沉陷归化至全周期地表沉陷部分，Disp(tⁱ)表示第i个SAR影像集获取的地表沉陷。

进一步，当判断不同时段SAR影像存在两个相邻时间段重叠时，利用最小二乘法获取偏移量Λⁱ(i＝1,2,…,K)；

设两个相邻时间段的SAR影像集为第A组与第B组，第A组SAR影像集重叠部分为T_A＝(t_n-x，t_n-x+1，t_n-x+2，···，t_n)，式中的t是影像帧对应的时间，其中x为重叠部分SAR影像个数，n为第A组SAR影像集总数；第B组SAR影像集重叠部分为T_B＝(t₁，t₂，t₃，···，t_y)，其中y≤m为重叠部分SAR影像个数，m为第B组SAR影像集总数；

则对于第A组与第B组两组数据中同名，即代表相同位置的相干点P的偏移量

通过下式进行估计：

式中，

表示A组SAR数据中相应时间段(t_n-t_n-1)的形变速率，

表示B组SAR数据中相应时间段(t_y-t_y-1)的形变速率。

进一步，当判断不同时段SAR影像不存在两个相邻时间段重叠时，利用概率积分模型获取偏移量Λⁱ(i＝1,2,…,K)；

设两个任意相邻时间段的SAR影像集为第A组与第B组，A、B两组SAR数据不存在时间重叠，通过A组SAR数据获取的地表沉陷数据结合概率积分模型反演形变参数：

式中，W₀＝mqcosα为地表最大下沉值，r＝H/tanβ为主要影响半径；(x，y)为地表任一点坐标，l，L分别为工作面走向和倾向长度，m为采厚，H为采深，α为煤层倾角，为已知量可由地质采矿条件确定；q为动态下沉系数，tanβ为主要影响角正切，为待求参数，可根据地表沉陷进行反演；

获取动态下沉系数q和主要影响角正切tanβ后，结合《建筑物、水体、铁路及主要井巷煤柱留设与压煤开采规程》对动态下沉系数q和主要影响角正切tanβ进行修正，

利用修正后的下沉系数q和主要影响角正切tanβ结合概率积分法对地表任一点最大形变量进行预测，结合Knothe时间函数，则对于A、B两组数据中同名，即表示相同位置的相干点P的偏移量

利用下式估计：

其中W_P为同名相干点P的最大下沉值，通过概率积分法预计；

分别为A、B组数据最后一景和第一景SAR影像时间；c为模型参数，通过矿区地质采矿条件确定。

进一步，步骤三中的精度要求根据实际应用场景进行设定，

对于应用于开采沉陷影响范围圈定时，因开采沉陷以下沉-10mm为界，因此设定监测精度为10mm；

对融合监测结果获取的全周期地表形变进行精度评估，当监测误差小于10mm时，输出当前结果；

监测误差大于10mm时，则调整优化步骤二中最小二乘法或概率积分法的求解参数，迭代运行步骤二直至监测精度满足精度标准。

有益效果：本方法针对传统变形监测方法获取大范围、长时序地表沉陷的不足，单一SAR影像仅能获取部分地表沉陷的缺陷的问题，充分利用不同时段多源SAR影像的优势，面向矿区沉陷的三个阶段(开采阶段、废弃(老)采空区阶段和矿井关闭阶段，获取矿区全生命周期地表沉陷，扩展了多源SAR数据监测开采沉陷的应用范围，填补了矿区全生命周期地表沉陷监测方法的空白；技术流程结构清晰，具有监测精范围广、精度高、易实现、省时、省力、经济实惠的优势，对揭示矿区开采沉陷时空演化规律、构建沉陷预测模型、预防控制沉陷灾害及促进煤矿区全生命周期安全绿色开发等具有重要的理论和实际价值。

附图说明

图1为本发明井工开采全生命周期地表沉陷监测方法流程图；

图2(a)为本发明井工开采全生命周期地表沉陷监测方法的ENVISAT干涉对时空基线图；

图2(b)为本发明井工开采全生命周期地表沉陷监测方法的ALOS-1PALSAR干涉对时空基线图；

图2(c)为本发明井工开采全生命周期地表沉陷监测方法的Sentinel-1A干涉对时空基线图；

图3(a)为本发明井工开采全生命周期地表沉陷监测方法的ENVISAT数据沉陷监测结果；

图3(b)为本发明井工开采全生命周期地表沉陷监测方法的ALOS-1PALSAR数据沉陷监测结果；

图3(c)本发明井工开采全生命周期地表沉陷监测方法的Sentinel-1A数据沉陷监测结果；

图3(d)为本发明井工开采全生命周期地表沉陷监测方法的同名相干点；

图4为本发明井工开采全生命周期地表沉陷监测方法的同名相干点全生命周期地表时序累计形变；

图5为本发明井工开采全生命周期地表沉陷监测方法的地表沉陷精度评估。

具体实施方式

以下将结合图和具体实施过程对本发明做进一步详细说明：

如图1所示，本发明的井工开采全生命周期地表沉陷监测方法步骤如下：

步骤1：利用不同时段多源SAR影像和InSAR技术获取井工开采地表沉陷：

收集或购买被测矿区不同时段多源SAR影像，依据相干性最大原则，设置多源SAR影像时间、空间阈值，选取SAR影像干涉对，并结合外部DEM数据进行差分干涉；利用时序InSAR技术获取被测矿区不同时段地表视线向时序形变，并将监测结果地理编码至相同坐标系下；根据不同时段多源SAR影像所选高相干点的坐标信息，选取同名相干点；根据多源SAR影像的成像几何关系，将视线向形变转换至垂直向：

式中，D_sub为垂直向形变，单位米；D_los为雷达视线向形变，单位米；α为雷达入射角，单位度。

步骤2：利用最小二乘法或概率积分法融合不同时段多源SAR监测结果：

不同时段多源SAR监测结果融合可分为两种情况，多源SAR影像集间存在和不存在时间重叠。

假设共获取K组不同时段SAR影像，每组影像通过步骤1处理均可获取一段时间内地表沉陷信息，共可获取K个不同时段地表沉陷。因每组SAR影像通过InSAR技术处理时均以第一期影像为参考获取地表形变，即：

式中，

为第i组SAR影像第一景数据获取时间；

为第i组SAR影像第一景数据对应的形变量。

因此，在获取的K组地表沉陷中前一组的最后一期和后一组的第一期之间存在形变偏移Λⁱ(i＝1,2,3,…K)，无法正确地反映出全生命周期地表沉陷的时空演化规律。同时，因每个相干点的形变速率不同，故不同相干点的形变偏移量不同。为了获取矿区全周期地表沉陷，将第一组SAR影像作为全局参考影像，即以第一组SAR影像的第一期影像为时间基点，令：

Λ¹≡0 (3)

而后面K-1组不同时段SAR影像监测结果可通过估计形变偏移量Λⁱ(i＝2,3,…K)，将参考时间规划至第一组SAR影像监测结果，从而实现矿区全周期地表沉陷监测，即：

Disp_full(tⁱ)＝Disp(tⁱ)+Λⁱ(i＝2,3,…K) (4)

式中，tⁱ表示第i组SAR影像获取时间；Disp(tⁱ)表示第i组SAR影像获取的地表沉陷，单位米；Disp_full(tⁱ)表示第i组SAR影像获取的地表沉陷归化至全周期地表沉陷部分，单位米。

由于不同时段多源SAR影像获取时间可能不连续，存在获取时间重叠和不重叠两种情况。

对于两个相邻时间段的SAR影像获取时间存在重叠时，形变偏移量Λⁱ(i＝2,3,…K)可利用最小二乘法进行估计。

假设A、B两组SAR影像集间存在时间重叠，第A组SAR影像获取时间在第B组之前，其中第A组SAR影像重叠部分的影像获取时间T_A为：

T_A＝(t_n-x,t_n-x+1,t_n-x+2,…,t_n) (5)

式中，x为重叠部分SAR影像个数；n为第A组SAR影像总数；t_n-x为SAR影像获取时间。

第B组SAR影像重叠部分的影像获取时间为：

T_B＝(t₁,t₂,t₃,…,t_y) (6)

式中，y≤m为重叠部分SAR影像个数；m为第B组SAR影像总数；t_y为SAR影像获取时间。

对于A、B两组SAR影像监测结果中同名相干点P的形变偏移量

可通过最小二乘法进行估计：

式中，

表示A组SAR数据中相应时间段(t_n-t_n-1)的形变速率；

表示B组SAR数据中相应时间段(t_y-t_y-1)的形变速率，形变速率可通过步骤1获取。在获得每个同名点的形变偏移量后，利用公式(4)可将B组SAR影像监测结果的时间基准归化至A组SAR影像监测结果，同理，其它组存在时间重叠的SAR影像监测结果也可利用公式(7)计算同名点的形变偏移量，再利用公式(4)统一时间基准。

对于两个相邻时间段的SAR影像获取时间不存在重叠时，形变偏移量Λⁱ(i＝2,3,…K)可利用概率积分法进行估计。

假设A、B两组SAR影像不存在时间重叠，则形变偏移量Λ^A-B无法通过SAR影像实测地表沉陷数据利用最小二乘法进行估计，而此时可利用矿区沉陷预测模型概率积分法进行估计。首先利用A组SAR影像获取的地表沉陷数据结合概率积分法反演模型参数，即：

式中，W_SAR(x,y)为A组SAR影像获取的地面任一点(x,y)的地表沉陷，单位米；W₀＝mqcosα为地表最大下沉值，单位米；r＝H/tanβ为主要影响半径，单位米；l，L分别为工作面走向和倾向长度，单位米；m为采厚，单位米；H为采深，单位米；α为煤层倾角，单位度；q为下沉系数；tanβ为主要影响角正切。

其中，W_SAR(x,y)、l、L、m、H、α为已知量可由SAR监测结果和地质采矿条件确定，q、tanβ为待求参数，可根据SAR数据获取的地表沉陷利用式(8)进行反演。

由于此时利用A组SAR影像监测数据获得动态下沉系数q和主要影响角正切tanβ并不是充分采动后的下沉系数和主要影响角正切，因此，需根据《建筑物、水体、铁路及主要井巷煤柱留设与压煤开采规程》将动态下沉系数q和主要影响角正切tanβ修正为充分采动条件下的下沉系数和主要影响角正切，即：

式中，下沉系数q和主要影响角正切tanβ可利用SAR监测结果反演；L_goal为工作面走向长度或倾向宽度；H为采深；k₁、k₂为修正系数可根据《建筑物、水体、铁路及主要井巷煤柱留设与压煤开采规程》进行确定；q_充、tanβ_充为修正后充分采动条件下的下沉系数和主要影响角正切。

利用修正后的下沉系数q_充和主要影响角正切tanβ_充结合概率积分法(式(8))预测地表任一点充分采动时地表沉陷，同时结合Knothe时间函数可获取地表任一点g在任意时间点地表沉陷，即：

式中，

为地表点g的最大沉陷值，单位米；c为模型参数，可通过矿区地质采矿条件确定；W(t)为任意时间点地表沉陷，单位米。

因此，根据式(8)和式(10)，不存在时间重叠的A、B两组数据中同名相干点P的偏移量

可通过下式估计：

式中，

同名相干点P的最大下沉值，单位米，可通过概率积分法预计；

分别为A、B组SAR影像最后一景和第一景数据时间；c为模型参数，可通过矿区地质采矿条件确定。

在获得每个同名点的形变偏移量后，利用公式(4)可将B组SAR影像监测结果的时间基准归化至A组SAR影像监测结果，同理，其它组不存在时间重叠的SAR影像监测结果也可利用公式(11)计算同名点的形变偏移量，再利用公式(4)统一时间基准。

步骤3：评估获取的全周期地表沉陷精度：

利用水准数据或GNSS数据评估获取的全生命周期地表沉陷精度，利用标准偏差对精度验证，当精度满足要求时，则输出结果，当精度不满足要求时，则调整优化步骤2中最小二乘法或概率积分法参数，重复步骤2，直至输出结果满足精度。

实施例一：

步骤1：利用不同时段多源SAR影像和InSAR技术获取井工开采地表沉陷：

本实例选取被测井工开采矿区，采用不同卫星平台获取的被测矿区14景C波段3.92m×7.80m(方位向×距离向)分辨率的ENVISAT影像，影像获取时间为2006年6月-2008年2月；9景L波段3.17m×4.68m(方位向×距离向)分辨率的ALOS-1PALSAR影像，数据获取时间为2007年2月至2011年1月；54景C波段13.94m×2.33m(方位向×距离向)分辨率的Sentinel-1A影像，数据获取时间为2016年10月至2018年8月。分别选取ENVISAT、ALOS-1PALSAR和Sentinel-1A数据集中2007年2月16日、2009年1月10日和2017年12月22日的SAR影像为主影像，其余影像均与其进行配准并重采样。依据相干性最大原则，设置时空基线，三组数据分别选取了19、17和143个干涉对，本发明实施例形成的干涉对组成情况结果如图2(a)，图2(b)，图2(c)所示。利用DInSAR技术获取差分干涉图，设置相干性阈值为0.4选取高相干点，采用时序InSAR技术基于不同时段ENVISAT、ALOS-1PALSAR及Sentinel-1A影像分别获取矿区地表视线向形变速率，并将监测结果地理编码至相同坐标系下，本发明实施例ENVISAT、ALOS-1PALSAR及Sentinel-1A影像获取的不同时段地表沉陷分别如图3(a)、图3(b)和图3(c)所示。同时基于SAR影像的成像几何关系，将视线向形变转换至垂直向形变；根据不同时段多源SAR影像所选高相干点的坐标信息，选取同名相干点，本发明实施例同名相干点如图3(d)所示。

步骤2：基于最小二乘法或概率积分法融合监测结果

由步骤1中采用的多源SAR数据获取时间可知，ENVISAT数据和ALOS-1PALSAR数据之间存在时间重叠，而ALOS-1PALSAR数据和Sentinel-1A数据之间不存在时间重叠，因此，对于存在时间重叠的ENVISAT和ALOS-1PALSAR数据可利用公式(7)，采用最小二乘法估计同名相干点的形变偏移量；对于不存在时间重叠的ALOS-1PALSAR和Sentinel-1A数据可利用公式(8)、(10)和(11)，采用概率积分法预计同名相干点的形变偏移量。然后，结合公式(4)对三组不同时段SAR影像监测结果进行融合就可获取矿区全生命周期地表沉陷，图4上下两部分为本发明实施例获取的被测矿井地表A、B两个同名相干点的全生命周期地表沉陷，同名相干点A、B位置见图3(d)。

步骤3：评估获取的全周期地表沉陷精度：

根据步骤2中获取的同名相干点全生命周期地表沉陷结果，结合地表实测水准数据，评估监测结果的精度，二者的对比结果如图5所示，由图5可知，利用最小二乘法和概率积分法和不同时段多源SAR影像获取的矿区全生命周期地表沉陷与实测数据吻合较好，形变趋势一致，二者的标准偏差为8.5mm，满足开采沉陷监测精度要求。

从本发明实施例获取的全生命周期地表沉陷和实施流程可以看出，本发明的井工开采全生命周期地表沉陷监测方法可精准地获取矿区全生命周期地表沉陷，扩展了多源SAR技术在矿区开采沉陷监测中的应用范围，为揭示开采沉陷时空演变规律、构建沉陷预测模型提供了新的思路和技术。

Claims (6)

1.一种井工开采全生命周期地表沉陷监测方法，其特征在于，具体步骤如下：

步骤一、首先利用不同时段获取的多源SAR数据和时序InSAR技术生成矿区地表雷达视线向形变，然后对矿区地表雷达视线向形变进行地理编码转换至相同坐标系下，结合SAR影像成像几何关系将矿区地表雷达视线向形变转为垂直向形变；

步骤二、当不同时段SAR影像存在时间重叠时，利用最小二乘方法对多源SAR数据监测结果进行平差；若不存在时间重叠时，则利用概率积分法对多源SAR数据监测结果进行融合；

步骤三、同时对获取的开采沉陷全生命周期形变场进行精度评估，若不满足精度要求，调整最小二乘法收敛阈值或概率积分法预测参数，重新利用最小二乘法或概率积分法对多源SAR数据监测结果进行融合，其中收敛阈值越小精度越高，概率积分法预测参数通过调整求参点数量进行优化；若满足精度要求，则输出矿区全生命周期地表沉陷，实现了融合不同时段多源SAR影像获取井工开采全生命周期地表沉陷。

2.根据权利要求1所述井工开采全生命周期地表沉陷监测方法，其特征在于步骤一具体为：利用不同平台卫星获取被测因井工开采导致地表沉陷的矿区不同时段SAR影像，对多源SAR影像分别利用时序InSAR技术获取被测矿区不同时段地表视线向形变，同时对多源SAR影像监测结果分别进行地理编码，将监测结果统一至相同的标准坐标系下，并结合各SAR影像成像几何关系，将被测矿区不同时段视线向形变转换至垂直向形变，在不同时段SAR影像监测结果中选取同名的相干点，即选取相同位置上的相干点进行处理，以获取时间上完整的地表形变。

3.根据权利要求1所述的井工开采全生命周期地表沉陷监测方法，其特征在于将不同时段SAR影像监测结果获取的地表沉陷归化至全周期地表沉陷部分具体步骤为：

设共获取到K个不同时段SAR影像监测结果，因每一个SAR影像集均以第一景SAR影像为参考影像获取地表时序形变，即：

其中

表示第i个SAR影像集中第一景影像时间，

为第i个SAR影像集中起始时间的形变量；因此，在每组SAR影像监测结果间均存在偏移量：Λⁱ(i＝1,2,…,K)，因每个相干点的形变速率不同，故不同相干点的偏移量不同；

为了获取全周期地表沉陷，将每组SAR影像中的SAR数据按照获取时间进行排序，并将获取时间最早的SAR数据定义为第一组SAR影像，将第一组SAR影像作为全局参考影像，即Λ¹≡0；然后通过偏移量Λⁱ(i＝1,2,…,K)，利用下式将不同时段SAR影像获取的地表形变起始时间归化至第一组SAR影像：

Disp_full(tⁱ)＝Disp(tⁱ)+Λⁱ(i＝1,2,…,K)，

式中Disp_full(tⁱ)表示第i个SAR影像集获取的地表沉陷归化至全周期地表沉陷部分，Disp(tⁱ)表示第i个SAR影像集获取的地表沉陷。

4.根据权利要求3所述的井工开采全生命周期地表沉陷监测方法，其特征在于当判断不同时段SAR影像存在两个相邻时间段重叠时，利用最小二乘法获取偏移量Λⁱ(i＝1,2,…,K)；

设两个相邻时间段的SAR影像集为第A组与第B组，第A组SAR影像集重叠部分为T_A＝(t_n-x，t_n-x+1，t_n-x+2，···，t_n)，式中的t是影像帧对应的时间，其中x为重叠部分SAR影像个数，n为第A组SAR影像集总数；第B组SAR影像集重叠部分为T_B＝(t₁，t₂，t₃，···，t_y)，其中y≤m为重叠部分SAR影像个数，m为第B组SAR影像集总数；

则对于第A组与第B组两组数据中同名，即代表相同位置的相干点P的偏移量

通过下式进行估计：

式中，

表示A组SAR数据中相应时间段(t_n-t_n-1)的形变速率，

表示B组SAR数据中相应时间段(t_y-t_y-1)的形变速率。

5.根据权利要求3所述的井工开采全生命周期地表沉陷监测方法，其特征在于当判断不同时段SAR影像不存在两个相邻时间段重叠时，利用概率积分模型获取偏移量Λⁱ(i＝1,2,…,K)；

设两个任意相邻时间段的SAR影像集为第A组与第B组，A、B两组SAR数据不存在时间重叠，通过A组SAR数据获取的地表沉陷数据结合概率积分模型反演形变参数：

式中，W₀＝mqcosα为地表最大下沉值，r＝H/tanβ为主要影响半径；(x，y)为地表任一点坐标，l，L分别为工作面走向和倾向长度，m为采厚，H为采深，α为煤层倾角，为已知量可由地质采矿条件确定；q为动态下沉系数，tanβ为主要影响角正切，为待求参数，可根据地表沉陷进行反演；

获取动态下沉系数q和主要影响角正切tanβ后，结合《建筑物、水体、铁路及主要井巷煤柱留设与压煤开采规程》对动态下沉系数q和主要影响角正切tanβ进行修正，

利用修正后的下沉系数q和主要影响角正切tanβ结合概率积分法对地表任一点最大形变量进行预测，结合Knothe时间函数，则对于A、B两组数据中同名，即表示相同位置的相干点P的偏移量

利用下式估计：

其中W_P为同名相干点P的最大下沉值，通过概率积分法预计；

分别为A、B组数据最后一景和第一景SAR影像时间；c为模型参数，通过矿区地质采矿条件确定。

6.根据权利要求1所述的井工开采全生命周期地表沉陷监测方法，其特征在于：步骤三中的精度要求根据实际应用场景进行设定，

对于应用于开采沉陷影响范围圈定时，因开采沉陷以下沉-10mm为界，因此设定监测精度为10mm；

对融合监测结果获取的全周期地表形变进行精度评估，当监测误差小于10mm时，输出当前结果；

监测误差大于10mm时，则调整优化步骤二中最小二乘法或概率积分法的求解参数，迭代运行步骤二直至监测精度满足精度标准。

Images