1.一种基于Sentinel‑1卫星SAR图像的大面积地表三维形变计算方法，其特征在于，该方法具体包括如下步骤：

步骤1：选定监测地点

从公开免费的对地遥感观测Sentinel卫星数据网站选定待监测地表三维形变的监测地区，配置数据参数，准备下载该待监测地区的SAR图像文件；

步骤2：定时下载SAR图像

SAR图像包括SLC模式和GRD模式的图像文件，每天定时从所述对地遥感观测Sentinel卫星数据网站获取该待监测地区的单视复数图像文件SLC和地距多视图像文件GRD，在后期数据处理中，使用差分干涉测量小基线集时序分析SBAS‑InSAR技术对SLC图像文件处理，获取视线向形变；使用偏移量追踪OT技术对GRD图像文件处理，获取方位向形变或距离向形变；

步骤3：判断轨道类型

在下载了待监测地区的SAR图像后，判断轨道类型，针对不同地区的轨道覆盖情况，分为四种情况：一个轨道覆盖类型的地区、两个交叉轨道覆盖的地区、两个平行轨道覆盖的地区和三个或更多轨道覆盖地区；

步骤4：对应解算方案

针对不同轨道类型的地区给出了不同的地表三维形变解算方案：a.针对一个轨道覆盖类型的地区，其地表三维形变解算方案为：S表示卫星，dLOS，ddistance，dazimuth分别代表视线向、距离向以及方位向上的形变，使用SBAS‑InSAR技术计算卫星S视线向的形变量dLOS，再使用OT技术计算出方位向形变dazimuth和距离向形变ddistance，然后再将视线向形变dLOS、方位向形变dazimuth和距离向形变ddistance分别投影到包含东西方向、南北方向和垂直方向的地表三维坐标系，便能够解算出地表的南北方向、东西方向以及垂直地面方向上的形变；

b.针对两个交叉轨道覆盖的地区，其地表三维形变解算方案为：S1，S2表示两颗轨道交叉的卫星，通过SBAS‑InSAR技术分别计算出卫星S1和卫星S2的视线向形变dLOS1和dLOS2，再使用OT技术计算出距离向形变ddistance，将dLOS1、dLOS2和ddistance分别投影到包含东西方向、南北方向和垂直方向的地表三维坐标系，便能够解算出地表的南北方向、东西方向以及垂直地面方向上的形变；

c.针对两个平行轨道覆盖的地区，其地表三维形变解算方案为：S1，S2表示两颗轨道平行的卫星，通过SBAS‑InSAR技术分别计算出卫星S1和卫星S2的视线向形变dLOS1和dLOS2，再使用OT技术计算出方位向形变dazimuth，将dLOS1、dLOS2和dazimuth分别投影到包含东西方向、南北方向和垂直方向的地表三维坐标系，便能够解算出地表的南北方向、东西方向以及垂直地面方向上的形变；

d.对于三个或更多轨道覆盖地区，其地表三维形变解算方案为：S1，S2表示两颗轨道平行的卫星，S3表示与S1，S2轨道交叉的卫星，通过SBAS‑InSAR技术计算出卫星S1、卫星S2和卫星S3各自的视线向形变dLOS1、dLOS2和dLOS3，将dLOS1、dLOS2和dLOS3分别投影到东西方向、南北方向和垂直方向便能够解算出地表的南北方向、东西方向以及垂直地面方向上的形变；

步骤5：解算形变序列

在不同轨道类型的地区的地表三维形变解算方案基础上，提出一种基于滑动窗口机制的处理方法，具体包括：

对于一个待监测地表三维形变的监测地区，设其包含n个成像区域，ar为成像区域的编号，且ar＝1，2，...，n， 表示成像区域ar的第j幅SAR图像， 表示成像区域ar的第j幅SAR图像的获取时间，采集的该待监测地表三维形变的监测地区的SAR图像序列P如下：为了计算整个图像序列P上的形变，采用滑动窗口机制，设置窗口大小k为1到n之间的整数，整体思想是：令k取遍1到n的所有值，当k＝n0时，将窗口从图像序列P起始点 开始向后移动，每次移动一幅SAR图像，直到 被纳入进窗口后结束，共移动n‑n0次，在移动过程中，每次移动均根据当前窗口内部包含的成像区域的轨道类型选用对应的地表三维形变解算方案，将计算出的形变按照时间顺序组织起来，覆盖上一次遍历的结果，逐步得到最终的地表三维形变序列；

先设窗口k＝1，代表每次只在该待监测地区的一个成像区域上计算，从 开始，认为该待监测地区为一个轨道覆盖类型的地区，对成像区域1的所有图像 中的任意两幅SAR图像使用一个轨道覆盖类型的地区的地表三维形变解算方案，得到成像区域1的所有图像的获取时间 上的形变序列L1，向后移动窗口，对成像区域2的任意两幅SAR图像使用一个轨道覆盖类型的地区的地表三维形变解算方案，得到上的形变序列L2，以此类推，直到该待监测地区的成像区域n计算完成；

再设窗口k＝2，代表每次在该待监测地区的两个成像区域上计算地表三维形变，从开始，窗口内部有 和 两幅SAR图像，如果窗口内部两幅SAR图像所属的两个成像区域属于同一轨道，则不进行操作；如果两个成像区域所属轨道相互平行，则使用两个平行轨道覆盖的地区的地表三维形变解算方案计算成像区域1和2上的地表三维形变，得到即公共时间段上的形变序列L1‑2；如果两个成像区域所属轨道交叉，则使用两个交叉轨道覆盖的地区的地表三维形变解算方案计算成像区域1和2上的地表三维形变，得到 即公共时间段上的形变序列L1‑2；再移动窗口，以此类推，直到该待监测地区的成像区域n计算完成，采用本次得到的公共时间段上的形变序列L1‑2覆盖k＝1时得到的上的形变序列；

再令k＝3，代表每次在该待监测地区的三个成像区域上计算地表三维形变，从 开始，窗口内部有 和 三幅SAR图像，如果窗口内部三幅SAR图像所属的三个成像区域对应三个轨道两两平行，则不进行操作；如果三个成像区域所属轨道存在交叉，则使用三个或更多轨道覆盖地区的地表三维形变解算方案计算成像区域1、2、3上的地表三维形变，得到即公共时间段上的形变序列L1‑3；再移动窗口，以此类推，直到该待监测地区的成像区域n计算完成，进一步采用本次得到的公共时间段上的形变更新k＝2时得到的 上的形

变序列；

再依次令k＝4，5......n，按照上述思路逐步计算，直到得到最终的地表三维形变序列L。

2.根据权利要求1所述的基于Sentinel‑1卫星SAR图像的大面积地表三维形变计算方法，其特征在于，所述步骤2中定时下载SAR图像具体包括：在所述对地遥感观测Sentinel卫星数据网站中，设置搜索待监测地区SAR图像的位置信息，点击“Filters”按钮，选择需要下载SAR图像的起止日期，“File Type”选择SLC和GRD类型分别用于SBAS‑InSAR和OT技术，“Beam Mode”选择IW模式，该模式主要用于陆地观测，“Direction”升降轨全选，“Subtype”表示需要获取的卫星图像来源，当前代表同时获取Sentinel‑1 A/B卫星图像，在点击搜索按钮后，检索出符合条件的SAR图像，点击“Queue”将SAR图像加入下载队列，在“Download”界面中获取到下载图像的python脚本，执行该脚本获取相关SAR图像。

3.根据权利要求2所述的基于Sentinel‑1卫星SAR图像的大面积地表三维形变计算方法，其特征在于，所述步骤3中判断轨道类型具体包括：在下载SAR图像前，先查看待监测地表三维形变的监测地区所有的SAR图像文件，在显示结果的地图中，点击成像区域后会展示出图像的Frame编号和Path编号，当两个成像区域Path编号相同而Frame编号不同则代表这两幅图像属于同一轨道的不同成像区域；当两个成像区域Path编号不同且两个成像区域的矩形框平行时，则代表这两个轨道平行；当Path编号不同且成像区域矩形框不平行时，则代表这两个轨道交叉；如果有三个成像区域Path编号不同，且对应的矩形框不两两平行，则说明这三个轨道交叉。

4.根据权利要求3所述的基于Sentinel‑1卫星SAR图像的大面积地表三维形变计算方法，其特征在于，所述步骤4中不同的地表三维形变解算方案具体为：采用SBAS‑InSAR技术降低相位噪声和误差，假定在时间t1至tM内获取同一地区的M幅SAR图像，然后根据干涉组合条件，在短基线距的条件下形成N幅干涉条纹图，且满足对tA和tB时刻两幅SAR图像生成的第nin幅干涉条纹图，在去除平地及地形相位影响后，第nin幅干涉条纹图中的第x个像素的干涉相位表示为：式中， 为tA～tB对应的视线向形变， 为tA～tB对应的地形相位误差， 为tA～tB对应的大气相位误差， 为tA～tB对应的基线轨道引起的相位误差； 为tA～tB对应的噪声误差，且t1＜tA＜tB＜tM，nin＝1，2，…，N；

对N辐干涉条纹图进行相位解缠、去除平地效应、去除大气相位、轨道矫正后便能够获取t1至tM的视线向形变dLOS；

OT技术则利用强度相关对两幅SAR图像逐像素配准，并从配准偏移量中解算地表位移，对于Sentinel‑1卫星而言，在进行轨道矫正之后认为配准偏移量即为地表偏移量，再将其分解到方位向和距离向，就得到了方位向形变dazimuth和距离向形变ddistance；

i i

令卫星的入射角为θ，其轨道倾角为α；θ、α分别第i个轨道对应卫星的入射角和轨道倾角，则有：

dU为垂直方向上的地表形变，dE为东西方向上的地表形变，dN为南北方向上的地表形变， 为使用SBAS‑InSAR技术计算出的第i个轨道的视线向形变，对于OT技术计算出的方位向形变dazimu和距离向形变ddistance有：i i

dazimuth＝dEcos(α‑3π/2)+dNsin(α‑3π/2)i i

ddistance＝dEsin(α‑3π/2)+dNcos(α‑3π/2)对于一个轨道覆盖类型的地区，解下列方程组便能够得到dU、dE、dN：dLOS＝dUcosθ‑dEsinθsin(α‑3π/2)‑dNsinθcos(α‑3π/2)dazimuth＝dEcos(α‑3π/2)+dNsin(α‑3π/2)ddistance＝dEsin(α‑3π/2)+dNcos(α‑3π/2)对于两个交叉轨道覆盖的地区，解下列方程组便能够得到dU、dE、dN：ddistance＝dEsin(α‑3π/2)+dNcos(α‑3π/2)对于两个平行轨道覆盖的地区，解下列方程组便能够得到dU、dE、dN：dazimuth＝dEcos(α‑3π/2)+dNsin(α‑3π/2)对于三个或更多轨道覆盖地区，以三轨道为例，解下列方程组便能够得到dU、dE、dN：

5.根据权利要求4所述的基于Sentinel‑1卫星SAR图像的大面积地表三维形变计算方法，其特征在于，所述步骤5中选择待监测地区中某四个成像区域A、B、C、D的情况，解算形变序列的具体内容包括：所述四个成像区域A、B、C、D中A、B成像区域轨道共线，A、B成像区域轨道与C和D成像区域所处轨道交叉，C、D成像区域轨道平行，并假设A、B、C、D四个成像区域的成像时间均匀分布，此时观测的图像序列如下：A1，B1，C1，D1，A2，B2，C2，D2其中，A1表示A的第1幅图像， 表示其获取时间；A2表示A的第2幅图像， 表示其获取时间；B1表示B的第1幅图像， 表示其获取时间；B2表示B的第2幅图像， 表示其获取时间；

C1表示C的第1幅图像， 表示其获取时间；C2表示C的第2幅图像， 表示其获取时间；D1表示D的第1幅图像， 表示其获取时间；D2表示D的第2幅图像， 表示其获取时间；最新图像为D2；设k代表滑动窗口大小，且1≤k≤4，k为整数；

令k＝1，从A1开始，对成像区域A使用一个轨道覆盖类型的地区的地表三维形变解算方案，即使用A1和A2获取 到 的地表三维形变d1，所述d1包含了东西方向、南北方向和垂直方向上的地表三维形变；再假定每次获取图像的时间间隔相同，得到 到 之间的形变为d1/4， 到 之间的形变为d1/2， 到 之间的形变为3d1/4；移动窗口至B1，获取 到 的地表三维形变d2， 到 之间的形变为d1/4+d2；移动窗口至C1，获取 到 的地表三维形变d3， 到 之间的形变为d1/2+d3；移动窗口至D1，获取 到 的地表三维形变d4， 到 之间的形变为3d1/4+d4，得到此时的地表三维形变序列L1；

然后，令k＝2，采用k＝2得到的公共时间段上的形变更新k＝1时得到的形变序列在该公共时间段上的形变；从A1，B1开始，由于A、B共线，仍相当于一个轨道覆盖类型的地区，因此不进行操作；移动窗口至B1，C1，B、C所属轨道交叉，相当于两个交叉轨道覆盖类型的地区，对B1，C1和B2，C2获取 到 公共时间段上的地表三维形变d5，则 到 之间的形变为d1/2+d5；

由于C、D轨道平行，相当于两个平行轨道覆盖类型的地区，对C1，D1和C2，D2获取 到 公共时间段上的地表三维形变d6，则 到 之间的形变为3d1/4+d6，得到此时的地表三维形变序列L2；

再令k＝3，采用k＝3得到的公共时间段上的形变进一步更新k＝2时得到的形变序列在该公共时间段上的形变；从A1，B1，C1开始，由于A、B共线，相当于两个交叉轨道覆盖类型的地区，不进行操作；移动窗口到B1，C1，D1，相当于多个轨道覆盖类型的地区，对B1，C1，D1和B2，C2，D2获取 到 的公共时间段上的地表三维形变d7，则 到 之间的形变为3d1/4+d7，得到此时的地表三维形变序列L3；

最后令k＝4，形变序列的更新规则和之前一样，由于此时仍然只包含三个不同方向的轨道，因此和k＝3时等同，得到最终的地表三维形变序列L。