本发明涉及采煤沉陷积水区水下地形测量技术领域，尤其涉及一种利用insar监测数据的采煤沉陷积水区水下地形反演方法。
背景技术：
两淮地区是我国重要的煤炭基地，该地区地下潜水位高、松散层厚、下沉系数大，煤炭资源的开采引起地表下沉形成盆地，加上大气降水与地下水补给等因素的共同影响，沉陷盆地内形成大面积积水，导致沉陷区水土流失严重，村庄、耕地被积水淹没，原有的陆地生态系统转化为水陆复合生态系统，对当地生态环境、社会经济造成了严重的影响。
当前针对高潜水位沉陷积水区的沉陷积水区的研究主要集中在水面的范围，对沉陷积水区水下地形的观测主要通过水深测量杆、水下地形测量船等传统方式进行。由于沉陷积水区范围大、水下地形坡度小，且积水区多处于农田或无人居住区域，交通不便，较大型的测量船只难以进入积水区进行测量作业，常规的观测手段成本高、耗时长、且有一定的危险性。除了实测手段获取水下地形外，也有通过数学方法进行水下地形插值，主要有：基于经验参数的概率积分法预计、线性插值等方法，这些方法由于缺少实测数据，所以预计结果误差较大。
针对沉陷区积水水下地形测量困难，数学方法预测精度不高等问题，故提出一种利用insar监测数据的采煤沉陷积水区水下地形反演方法。
技术实现要素：
本发明的目的在于克服传统技术中存在的上述问题，提供一种利用insar监测数据的采煤沉陷积水区水下地形反演方法。
为实现上述技术目的，达到上述技术效果，本发明是通过以下技术方案实现：
一种利用insar监测数据的采煤沉陷积水区水下地形反演方法，包括如下步骤：
s1.获取采煤工作面上方地表沉陷前和沉陷后的两副insar数据，将沉陷前insar数据作为主影像、沉陷后insar数据作为辅影像，利用insar的处理技术，经过对两副insar影像配准、干涉处理、去平地、解缠相位获得采煤沉陷区未积水区dem，用h1(x,y)未积水表示，其中(x,y)未积水为沉陷区未积水区点的坐标，h1为沉陷后点的高程；
s2.生成采煤沉陷区沉陷前dem；
s3.在地面建立具有wgs84和国家坐标系同名点，获取点位坐标并进行七参数解算得到坐标转换参数，运用坐标转换参数将雷达坐标系下的采煤沉陷区未积水区dem转换到国家坐标系下的dem；
s4.将采煤沉陷区沉陷前dem与获得的采煤沉陷区未积水区dem的h0(x,y)未积水叠加、差值计算，获得沉陷区边缘未积水区地表下沉值
w(x,y)未积水＝h0(x,y)-h1(x,y)未积水；
s5.利用沉陷区部分下沉监测数据即获得的沉陷区边缘未积水区数据，反求出概率积分法模型中参数，概率积分法模型为：
w0＝mqcosα，
r＝h/tanβ，
其中，w(x,y)为地表下沉值，w0为地表最大下沉值；m为煤层采厚，q为下沉系数；r为主要影响半径，h为采深，tanβ是主要影响角正切；d是地下开采范围；
s6.根据采煤工作面参数，利用概率积分法模型以及s3.中已根据部分实测沉陷数据求取的参数，对沉陷区积水区域进行插值计算，获得积水区域下沉值w(x,y)；
s7.与采煤沉陷区沉陷前dem叠加差值计算，获取沉陷区积水区水下地形
h1(x,y)积水＝h0(x,y)-w(x,y)积水。
进一步地，步骤s1.中，insar的处理技术包括对两副insar影像配准、干涉处理、去平地和解缠相位。
进一步地，步骤s5.中，概率积分法模型中参数包括下沉系数、主要主要影响角正切、开采影响传播角和拐点偏移距。
进一步地，步骤s6.中，采煤工作面参数包括工作面煤层倾角、工作面走向长、工作面倾向长和工作面采深。
本发明的一个优选方案为：步骤s2.中，采煤沉陷区沉陷前dem的生成方法为：收集大比例地形图，提取高程点和等高线信息，运用空间分析功能生成数字高程三角模型tin，并建立高精度的dem。
本发明的另一个优选方案为：步骤s2.中，采煤沉陷区沉陷前dem的生成方法为：下载全球数字高程数据srtm，对覆盖大的研究区需要镶嵌并根据区域范围进行矢量裁剪。
本发明的有益效果是：
1、基于d-insar获取采煤沉陷区未积水区地表沉降信息，无需工作人员现场测量，不受天气等因素影响，且sar数据获取便利、成本低、有历史数据，d-insar获取的沉降数据精度高，点位多。
2、基于d-insar半实测数据的概率积分法参数较概率积分法经验参数，更适合指定采煤沉陷区地表移动盆地的沉陷预计。基于概率积分法模型的水下地形是基于半实测数据，故预测精度较常用的概率积分法预计、常规插值法更高，大大减少人力实测成本。
当然，实施本发明的任一产品并不一定需要同时达到以上的所有优点。
附图说明
为了更清楚地说明本发明实施例的技术方案，下面将对实施例描述所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明的流程框图；
图2为本发明利用d-insar技术获取沉陷区未积水区地表下沉值的示意图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其它实施例，都属于本发明保护的范围。
实施例一
请参阅图1-2所示，本实施例为一种利用insar监测数据的采煤沉陷积水区水下地形反演方法，包括如下步骤：
s1.获取采煤工作面上方地表沉陷前和沉陷后的两副insar数据，将沉陷前insar数据作为主影像、沉陷后insar数据作为辅影像，利用insar的处理技术，经过对两副insar影像配准、干涉处理、去平地、解缠相位获得采煤沉陷区未积水区dem，用h1(x,y)未积水表示，其中(x,y)未积水为沉陷区未积水区点的坐标，h1为沉陷后点的高程；
s2.生成采煤沉陷区沉陷前dem；采煤沉陷区沉陷前dem的生成方法为：收集大比例地形图，提取高程点和等高线信息，运用arcgis空间分析功能生成数字高程三角模型tin，并建立高精度的dem；
s3.在地面建立具有wgs84和国家坐标系同名点，获取点位坐标并进行七参数解算得到坐标转换参数，运用坐标转换参数将雷达坐标系下的采煤沉陷区未积水区dem转换到国家坐标系下的dem；雷达坐标系是wgs84坐标系统，需要统一到国家坐标系下，国家坐标系为cgcs2000、西安80或者北京54；
s4.将采煤沉陷区沉陷前dem与获得的采煤沉陷区未积水区dem的h0(x,y)未积水叠加、差值计算，获得沉陷区边缘未积水区地表下沉值
w(x,y)未积水＝h0(x,y)-h1(x,y)未积水；
s5.利用沉陷区部分下沉监测数据即获得的沉陷区边缘未积水区数据，反求出概率积分法模型中参数，概率积分法模型为：
w0＝mqcosα，
r＝h/tanβ，
其中，w(x,y)为地表下沉值，w0为地表最大下沉值；m为煤层采厚，q为下沉系数；r为主要影响半径，h为采深，tanβ是主要影响角正切；d是地下开采范围；
s6.根据采煤工作面参数，利用概率积分法模型以及s3.中已根据部分实测沉陷数据求取的参数，对沉陷区积水区域进行插值计算，获得积水区域下沉值w(x,y)；
s7.与采煤沉陷区沉陷前dem叠加差值计算，获取沉陷区积水区水下地形
h1(x,y)积水＝h0(x,y)-w(x,y)积水。
本实施例步骤s5.中，概率积分法模型中参数包括下沉系数、主要主要影响角正切、开采影响传播角和拐点偏移距。
本实施例步骤s6.中，采煤工作面参数包括工作面煤层倾角、工作面走向长、工作面倾向长和工作面采深。
实施例二
请参阅图1-2所示，本实施例为一种利用insar监测数据的采煤沉陷积水区水下地形反演方法，包括如下步骤：
s1.获取采煤工作面上方地表沉陷前和沉陷后的两副insar数据，将沉陷前insar数据作为主影像、沉陷后insar数据作为辅影像，利用insar的处理技术，经过对两副insar影像配准、干涉处理、去平地、解缠相位获得采煤沉陷区未积水区dem，用h1(x,y)未积水表示，其中(x,y)未积水为沉陷区未积水区点的坐标，h1为沉陷后点的高程；
s2.生成采煤沉陷区沉陷前dem；采煤沉陷区沉陷前dem的生成方法为：在nasa官网下载全球数字高程数据srtm，对覆盖大的研究区需要镶嵌并根据区域范围进行矢量裁剪；
s3.在地面建立具有wgs84和国家坐标系同名点，获取点位坐标并进行七参数解算得到坐标转换参数，运用坐标转换参数将雷达坐标系下的采煤沉陷区未积水区dem转换到国家坐标系下的dem；雷达坐标系是wgs84坐标系统，需要统一到国家坐标系下，国家坐标系为cgcs2000、西安80或者北京54；
s4.将采煤沉陷区沉陷前dem与获得的采煤沉陷区未积水区dem的h0(x,y)未积水叠加、差值计算，获得沉陷区边缘未积水区地表下沉值
w(x,y)未积水＝h0(x,y)-h1(x,y)未积水；
s5.利用沉陷区部分下沉监测数据即获得的沉陷区边缘未积水区数据，反求出概率积分法模型中参数，概率积分法模型为：
w0＝mqcosα，
r＝h/tanβ，
其中，w(x,y)为地表下沉值，w0为地表最大下沉值；m为煤层采厚，q为下沉系数；r为主要影响半径，h为采深，tanβ是主要影响角正切；d是地下开采范围；
s6.根据采煤工作面参数，利用概率积分法模型以及s3.中已根据部分实测沉陷数据求取的参数，对沉陷区积水区域进行插值计算，获得积水区域下沉值w(x,y)；
s7.与采煤沉陷区沉陷前dem叠加差值计算，获取沉陷区积水区水下地形
h1(x,y)积水＝h0(x,y)-w(x,y)积水。
本实施例步骤s5.中，概率积分法模型中参数包括下沉系数、主要主要影响角正切、开采影响传播角和拐点偏移距。
本实施例步骤s6.中，采煤工作面参数包括工作面煤层倾角、工作面走向长、工作面倾向长和工作面采深。
该方法中insar数据获取成本低，d-insar技术获取的地表下沉数据精度高、数据量大，参数反演基于半实测数据，参数反演结果适用于研究区，该方法中基于概率积分法模型的水下地形插值，可以反演出工作面采动过程中积水区域水下地形，指导沉陷区水陆复合生态系统的建立，为沉陷区地表的土地利用规划、土地复垦、生态恢复等提供重要参考。
以上公开的本发明优选实施例只是用于帮助阐述本发明。优选实施例并没有详尽叙述所有的细节，也不限制该发明仅为具体实施方式。显然，根据本说明书的内容，可作很多的修改和变化。本说明书选取并具体描述这些实施例，是为了更好地解释本发明的原理和实际应用，从而使所属技术领域技术人员能很好地理解和利用本发明。本发明仅受权利要求书及其全部范围和等效物的限制。
1.一种利用insar监测数据的采煤沉陷积水区水下地形反演方法，其特征在于，包括如下步骤：
s1.获取采煤工作面上方地表沉陷前和沉陷后的两副insar数据，将沉陷前insar数据作为主影像、沉陷后insar数据作为辅影像，利用insar的处理技术获得采煤沉陷区未积水区dem，用h1(x,y)未积水表示，其中(x,y)未积水为沉陷区未积水区点的坐标，h1为沉陷后点的高程；
s2.生成采煤沉陷区沉陷前dem；
s3.在地面建立具有wgs84和国家坐标系同名点，获取点位坐标并进行七参数解算得到坐标转换参数，运用坐标转换参数将雷达坐标系下的采煤沉陷区未积水区dem转换到国家坐标系下的dem；
s4.将采煤沉陷区沉陷前dem与获得的采煤沉陷区未积水区dem的h0(x,y)未积水叠加、差值计算，获得沉陷区边缘未积水区地表下沉值
w(x,y)未积水＝h0(x,y)-h1(x,y)未积水；
s5.利用沉陷区部分下沉监测数据即获得的沉陷区边缘未积水区数据，反求出概率积分法模型中参数，概率积分法模型为：
w0＝mqcosα，
r＝h/tanβ，
其中，w(x,y)为地表下沉值，w0为地表最大下沉值；m为煤层采厚，q为下沉系数；r为主要影响半径，h为采深，tanβ是主要影响角正切；d是地下开采范围；
s6.根据采煤工作面参数，利用概率积分法模型以及s3.中已根据部分实测沉陷数据求取的参数，对沉陷区积水区域进行插值计算，获得积水区域下沉值w(x,y)；
s7.与采煤沉陷区沉陷前dem叠加差值计算，获取沉陷区积水区水下地形
h1(x,y)积水＝h0(x,y)-w(x,y)积水。
2.根据权利要求1所述的利用insar监测数据的采煤沉陷积水区水下地形反演方法，其特征在于：步骤s1.中，insar的处理技术包括对两副insar影像配准、干涉处理、去平地和解缠相位。
3.根据权利要求1所述的利用insar监测数据的采煤沉陷积水区水下地形反演方法，其特征在于：步骤s5.中，概率积分法模型中参数包括下沉系数、主要主要影响角正切、开采影响传播角和拐点偏移距。
4.根据权利要求1所述的利用insar监测数据的采煤沉陷积水区水下地形反演方法，其特征在于：步骤s6.中，采煤工作面参数包括工作面煤层倾角、工作面走向长、工作面倾向长和工作面采深。
5.根据权利要求1-4任一项所述的利用insar监测数据的采煤沉陷积水区水下地形反演方法，其特征在于：步骤s2.中，采煤沉陷区沉陷前dem的生成方法为：收集大比例地形图，提取高程点和等高线信息，运用空间分析功能生成数字高程三角模型tin，并建立高精度的dem。
6.根据权利要求1-4任一项所述的利用insar监测数据的采煤沉陷积水区水下地形反演方法，其特征在于：步骤s2.中，采煤沉陷区沉陷前dem的生成方法为：下载全球数字高程数据srtm，对覆盖大的研究区需要镶嵌并根据区域范围进行矢量裁剪。
技术总结