一种生产建设用地的高频次全自动遥感监测方法

1.本发明涉及遥感监测的技术领域，更具体地，涉及一种生产建设用地的高频次全自动遥感监测方法。


背景技术：

2.违法的生产建设用地侵占了大量优质农田，威胁着耕地红线和国家粮食安全，同时导致水土流失，对生态环境造成了极大破坏。生产建设用地开发在起始阶段通常只涉及土地平整，对环境造成的破坏较小且容易恢复。因此，开展高频次(每星期或每月)的生产建设用地开发监测，及时发现刚刚开始的违法开发行为，对永久基本农田管控、水土保持监管有着至关重要的作用。光学遥感影像可以提供丰富的地物光谱信息，是目前监测土地利用变化的主要数据源。然而，光学遥感监测周期一般较长(一年或多年)，难以及时发现并预防违法生产建设用地的开发，导致被侵占的土地难以恢复，对环境造成的破坏难以逆转。而传统的人工目视解译方法耗时耗人、自动化水平较低，导致监测成本较高、监测效率较低，很难适应高频次的监管任务。近年来，国内外越来越多的遥感卫星陆续发射，如我国高分(gf)系列卫星、欧空局sentinel-2以及美国landsat-8等，使遥感影像资源愈发丰富。结合多源光学遥感影像，可以显著提高对地观测频次，使生产建设项目开发实时监管成为可能。通过研究生产建设用地开发初期阶段土地平整导致的光谱变化特征发现，与其它地块相比，生产建设用地通常具有较高的土地平整强度。但是，农作物、自然植被、山体及建筑阴影的季节性变化也可导致土地平整强度较高，与生产建设用地混淆，给生产建设用地的高频次监测带来较多的误报。
3.2014年1月29日公开的中国专利申请cn103546726a提供了一种违章用地的自动发现方法，包括以下步骤：将不同时相的视频数据进行稳像处理；对稳像后的视频数据进行背景建模，并滤除视频数据中的流动前景，获得每一个预置点在不同时间节点下的背景图片；对背景图片的滤波处理，滤除背景图片中的噪声；将相同预置点不同时相的背景图片与历史背景图片匹配，获得差值图片，提取差值图片的变化特征并判断违章用地是否存在。本发明的技术方案利用数字图像、视频处理技术对目标监控区域的视频数据进行处理，建立标准监控目标背景历史数据库，通过分析比对当前目标背景和历史目标背景，提取差值图片的变化特征，然而农作物、自然植被、山体及建筑阴影的季节性变化也可导致相同的变化特征，发现的违章用地会存在大量误报，监测监督低。


技术实现要素：

4.本发明为克服上述现有技术生产建设用地检测精度低的缺陷，提供一种生产建设用地的高频次全自动遥感监测方法，可以消除农作物、自然植被、山体及建筑阴影季节性变化导致的误差，实现生产建设用地的高频次、全自动的监测，并显著提高监测精度。
5.为解决上述技术问题，本发明的技术方案如下：
6.本发明提供了一种生产建设用地的高频次全自动遥感监测方法，包括以下步骤：
7.s1：获取待检测区域三个不同时间的光学遥感影像，按照三个不同获取时间的先后顺序依次设置为矫正影像、参考影像和检测影像；
8.s2：计算检测影像和参考影像之间的土地平整强度，获得第一土地平整强度图像；计算检测影像和矫正影像之间的土地平整强度，获得第二土地平整强度图像；
9.s3：确定土地平整强度最佳阈值；
10.s4：根据土地平整强度最佳阈值，利用阈值法分别提取出第一土地平整强度图像上的生产建设用地和第二土地平整强度图像上的生产建设用地，获得第一检测结果图像和第二检测结果图像；
11.s5：根据第一检测结果图像和第二检测结果图像，获得最终生产建设用地的开发情况。
12.检测影像为当前时间待检测区域的光学遥感影像，参考影像的获取时间早于检测影像的获取时间，矫正影像的获取时间早于参考影像的获取时间；第一土地平整强度图像反映出参考影像的获取时间和检测影像的获取时间这段时间内，待检测区域的光谱信息；第二土地平整强度图像反映出矫正影像的获取时间和检测影像的获取时间这段时间内，待检测区域的光谱信息；利用阈值法，分别提取出第一土地平整强度图像上的生产建设用地和第二土地平整强度图像上的生产建设用地；由于检测影像和参考影像获取时间的时间间隔较短，无法避免农作物、自然植被、山体及建筑阴影季节性变化的影响，第一土地平整强度图像上的生产建设用地的区域较多，存在较多误报；矫正影像和检测影像获取时间的时间间隔较长，农作物、自然植被、山体及建筑阴影季节性变化不会被检测到，第二土地平整强度图像上的生产建设用地的区域误报明显减少；最后根据第一检测结果图像和第二检测结果图像上的生产建设用地的区域，确定最终生产建设用地。
13.优选地，所述步骤s1中，利用相同或不同遥感卫星获取待检测区域三个不同时间的光学遥感影像。利用相同或不同遥感卫星获取光学遥感影像，可以大大增加监测频次，减少天气因素对监测过程的影响。
14.优选地，所述步骤s1中，参考影像和检测影像的获取时间的时间间隔小于一年，检测影像和矫正影像的获取时间的时间间隔等于一年。
15.参考影像和检测影像的获取时间的时间间隔小于一年，根据高频次实时监测的需要，间隔为5-30天最优；由于农作物、自然植被、山体及建筑阴影虽然具有明显季节性变化，但是在不同年份的相同时间通常呈现相同的状态，将检测影像和矫正影像的获取时间的时间间隔设置为一年，最不易被矫正影像和检测影像上产生误报，精度更高；实际中可能无法获取时间间隔一年的矫正影像，应确保矫正影像的获取时间尽可能接近间隔一年的日期。
16.优选地，所述步骤s2中，计算检测影像和参考影像之间的土地平整强度的具体方法为：
[0017][0018]
式中，i
t12
表示检测影像和参考影像之间的土地平整强度，band(r)
t2
表示2时间获取的检测影像的红光波段，band(b)
t2
表示2时间获取的检测影像的蓝光波段；band(r)
t1
表示1时间获取的参考影像的红光波段，band(b)
t1
表示1时间获取的参考影像的蓝光波段；i
t12
值越大，表示t1—t2时间，进行土地平整的可能性越大，即进行生产建设用地开发可能性越
大。
[0019]
优选地，所述步骤s2中，计算检测影像和矫正影像之间的土地平整强度的具体方法为：
[0020][0021]
式中，i
t02
表示检测影像和矫正影像之间的土地平整强度，band(r)
t2
表示2时间获取的检测影像的红光波段，bαnd(b)
t2
表示2时间获取的检测影像的蓝光波段；band(r)
t0
表示0时间获取的矫正影像的红光波段，band(b)
t0
表示0时间获取的矫正影像的蓝光波段。i
t02
值越大，表示t0—t2时间，进行土地平整的可能性越大，即进行生产建设用地开发可能性越大。
[0022]
优选地，所述步骤s3中，确定土地平整强度最佳阈值的具体方法为：
[0023]
s3.1：获取多种遥感卫星拍摄的生产建设用地的参考图像和检测图像、非生产建设用地的参考图像和检测图像；
[0024]
s3.2：对多种遥感卫星拍摄的参考图像和检测图像进行组合，计算不同遥感卫星组合的土地平整强度；
[0025]
s3.3：利用决策树算法，根据不同遥感卫星组合的土地平整强度，获得土地平整强度阈值范围；
[0026]
s3.4：根据实际获取待检测区域的遥感卫星种类，在土地平整强度阈值范围内确定土地平整强度最佳阈值。
[0027]
优选地，所述步骤s4中，利用阈值法提取第一土地平整强度图像上的生产建设用地的方法为：
[0028]
比较第一土地平整强度图像与地平整强度最佳阈值的大小，i
t12
大于土地平整强度最佳阈值的区域，为第一土地平整强度图像上的生产建设用地，否则，为非生产建设用地。
[0029]
优选地，所述步骤s4中，利用阈值法提取第二土地平整强度图像上的生产建设用地的方法为：
[0030]
比较第二土地平整强度图像与地平整强度最佳阈值的大小，i
t02
大于土地平整强度最佳阈值的区域，为第二土地平整强度图像上的生产建设用地，否则，为非生产建设用地。
[0031]
优选地，所述步骤s5中，根据第一检测结果图像和第二检测结果图像，获得最终生产建设用地的具体方法为：
[0032]
对第一检测结果图像和第二检测结果图像取交集，即：
[0033]ct12
＝(i
t12
》t)&(i
t02
》t)
[0034]
式中，c
t12
表示最终生产建设用地，t表示土地平整强度最佳阈值，&表示取交集操作。c
t12
的取值为0或1，当第一检测结果图像上生产建设用地的区域与第二检测结果图像上生产建设用地的区域重叠时，所述区域作为最终生产建设用地。
[0035]
优选地，所述步骤s3.3中，土地平整强度阈值范围为[1.31，1.42]。
[0036]
与现有技术相比，本发明技术方案的有益效果是：
[0037]
本发明获取待检测区域三个不同时间的光学遥感影像，按照获取时间从早到晚依
次设置为矫正影像、参考影像和检测影像；通过计算检测影像和参考影像、检测影像和矫正影像之间的土地平整强度，获得第一土地平整强度图像和第二土地平整强度图像；利用阈值法，分别提取出第一土地平整强度图像上的生产建设用地和第二土地平整强度图像上的生产建设用地；由于检测影像和参考影像获取时间的时间间隔较短，无法避免农作物、自然植被、山体及建筑阴影季节性变化的影响，第一土地平整强度图像上的生产建设用地的区域较多，存在较多误报；矫正影像和检测影像获取时间的时间间隔较长，农作物、自然植被、山体及建筑阴影季节性变化不会被检测到，第二土地平整强度图像上的生产建设用地的区域误报明显减少，最后根据第一检测结果图像和第二检测结果图像上的生产建设用地的区域，确定最终生产建设用地，可以有效剔除农作物、自然植被、山体及建筑阴影季节性变化导致的误差，实现生产建设用地的高频次、全自动的检测，并显著提高监测精度。
附图说明
[0038]
图1为实施例所述的一种生产建设用地的高频次全自动遥感监测方法的流程图；
[0039]
图2为实施例所述的不同遥感卫星组合的土地平整强度阈值的示意图。
[0040]
图3为实施例所述的现有方法的检测结果示意图。
[0041]
图4为实施例所述的一种生产建设用地的高频次全自动遥感监测方法的示意图。
[0042]
图5为实施例所述的一种生产建设用地的高频次全自动遥感监测方法与现有方法的检测结果的精度对比图。
具体实施方式
[0043]
附图仅用于示例性说明，不能理解为对本专利的限制；
[0044]
为了更好说明本实施例，附图某些部件会有省略、放大或缩小，并不代表实际产品的尺寸；
[0045]
对于本领域技术人员来说，附图中某些公知结构及其说明可能省略是可以理解的。
[0046]
下面结合附图和实施例对本发明的技术方案做进一步的说明。
[0047]
实施例
[0048]
本实施例提供了一种生产建设用地的高频次全自动遥感监测方法，如图1所示，包括以下步骤：
[0049]
s1：获取待检测区域三个不同时间的光学遥感影像，按照三个不同获取时间的先后顺序依次设置为矫正影像、参考影像和检测影像；
[0050]
需要说明的可以利用相同或不同遥感卫星获取待检测区域三个不同时间的光学遥感影像，利用不同遥感卫星获取光学遥感影像，可以大大增加监测频次，减少天气因素对监测过程的影响；
[0051]
参考影像和检测影像的获取时间的时间间隔小于一年，检测影像和矫正影像的获取时间的时间间隔等于一年；
[0052]
参考影像和检测影像的获取时间的时间间隔小于一年，根据高频次实时监测的需要，间隔为5-30天最优；由于农作物、自然植被、山体及建筑阴影虽然具有明显季节性变化，但是在不同年份的相同时间通常呈现相同的状态，将检测影像和矫正影像的获取时间的时
间间隔设置为一年，最不易被矫正影像和检测影像上产生误报，精度更高；实际中可能无法获取时间间隔一年的矫正影像，应确保矫正影像的获取时间尽可能接近间隔一年的日期。
[0053]
s2：计算检测影像和参考影像之间的土地平整强度，获得第一土地平整强度图像；计算检测影像和矫正影像之间的土地平整强度，获得第二土地平整强度图像；
[0054][0055][0056]
式中，i
t12
表示检测影像和参考影像之间的土地平整强度，band(r)
t2
表示t2时间获取的检测影像的红光波段，band(b)
t2
表示t2时间获取的检测影像的蓝光波段；band(r)
t1
表示t1时间获取的参考影像的红光波段，band(b)
t1
表示1时间获取的参考影像的蓝光波段；i
t02
表示检测影像和矫正影像之间的土地平整强度，band(r)
t0
表示0时间获取的矫正影像的红光波段，band(b)
t0
表示0时间获取的矫正影像的蓝光波段；i
t12
值越大，表示t1—t2时间，进行土地平整的可能性越大，即进行生产建设用地开发可能性越大；i
t02
值越大，表示t0—t2时间，进行土地平整的可能性越大，即进行生产建设用地开发可能性越大。
[0057]
s3：确定土地平整强度最佳阈值；
[0058]
s3.1：获取多种遥感卫星拍摄的生产建设用地的参考图像和检测图像、非生产建设用地的参考图像和检测图像；
[0059]
s3.2：对多种遥感卫星拍摄的参考图像和检测图像进行组合，计算不同遥感卫星组合的土地平整强度；
[0060]
s3.3：利用决策树算法，根据不同遥感卫星组合的土地平整强度，获得土地平整强度阈值范围；
[0061]
s3.4：根据实际获取待检测区域的遥感卫星种类，在土地平整强度阈值范围内确定土地平整强度最佳阈值。
[0062]
s4：根据土地平整强度最佳阈值，利用阈值法分别提取出第一土地平整强度图像上的生产建设用地和第二土地平整强度图像上的生产建设用地，获得第一检测结果图像和第二检测结果图像；
[0063]
比较第一土地平整强度图像与地平整强度最佳阈值的大小，i
t12
大于土地平整强度最佳阈值的区域，为第一土地平整强度图像上的生产建设用地，否则，为非生产建设用地；
[0064]
比较第二土地平整强度图像与地平整强度最佳阈值的大小，i
t02
大于土地平整强度最佳阈值的区域，为第二土地平整强度图像上的生产建设用地，否则，为非生产建设用地。
[0065]
s5：根据第一检测结果图像和第二检测结果图像，获得最终生产建设用地的开发情况。
[0066]
检测影像为当前时间待检测区域的光学遥感影像，参考影像的获取时间早于检测影像的获取时间，矫正影像的获取时间早于参考影像的获取时间；第一土地平整强度图像反映出参考影像的获取时间和检测影像的获取时间这段时间内，待检测区域的光谱信息；第二土地平整强度图像反映出矫正影像的获取时间和检测影像的获取时间这段时间内，待
检测区域的光谱信息；利用阈值法，分别提取出第一土地平整强度图像上的生产建设用地和第二土地平整强度图像上的生产建设用地，由于检测影像和参考影像获取时间的时间间隔较短，无法避免农作物、自然植被、山体及建筑阴影季节性变化的影响，第一土地平整强度图像上的生产建设用地的区域较多，存在较多误报；矫正影像和检测影像获取时间的时间间隔较长，农作物、自然植被、山体及建筑阴影季节性变化不会被检测到，第二土地平整强度图像上的生产建设用地的区域误报明显减少；最后根据第一检测结果图像和第二检测结果图像上的生产建设用地的区域，确定最终生产建设用地。
[0067]
优选地，所述步骤s5中，根据第一检测结果图像和第二检测结果图像，获得最终生产建设用地的具体方法为：
[0068]
对第一检测结果图像和第二检测结果图像取交集，即：
[0069]ct12
＝(i
t12
》t)&(i
t02
》t)
[0070]
式中，c
t12
表示最终生产建设用地，t表示土地平整强度最佳阈值，&表示取交集操作。c
t12
的取值为0或1，当第一检测结果图像上生产建设用地的区域与第二检测结果图像上生产建设用地的区域重叠时，所述区域作为最终生产建设用地。
[0071]
在具体实施过程中，鉴于不同气候带的植被物候和类型不同，所以从全球多个气候带各选取一个城市，从这些城市及其周边地区选取生产建设用地和非生产建设用地的光学遥感影像；下表列出了不同遥感影像组合选取的生产建设用地和非生产建设用地的个数；
[0072][0073]
遥感影像组合中在前卫星型号表示获取参考影像的遥感卫星，在后卫星型号表示获取检测影像的遥感卫星；以landsat-5
–
landsat-7为例，表示获取参考影像为landsat-5遥感卫星，获取检测影像的为landsat-7遥感卫星。以上表的数据为样本集，利用r语言提供的决策树算法，确定土地平整强度最佳阈值；如图2所示，为不同遥感卫星组合的土地平整强度阈值的示意图，展示的是不同遥感影像组合计算的土地平整强度上样本数值分布和用决策树算法确定的土地平整强度阈值，landsat-5
–
landsat-5的土地平整强度最佳阈值为1.42，landsat-5
–
landsat-7的土地平整强度最佳阈值为1.39，landsat-5
–
landsat-7的土地平整强度最佳阈值为1.37，landsat-8
–
landsat-8的土地平整强度最佳阈值为1.31，landsat-8
–
sentinel-2的土地平整强度最佳阈值为1.39，sentinel-2
–
sentinel-2的土地平整强度最佳阈值为1.41，土地平整强度阈值范围为1.31到1.42；不同遥感影像组合对应
的土地平整强度最佳阈值略有不同，是由于不同遥感卫星的传感器之间波长的微小差距造成，例如，landsat-5和-7的红光波段为0.63
–
0.69μm，而landsat-8的红光波段为0.64
–
0.67μm；最后根据实际获取待检测区域的遥感卫星种类，在土地平整强度阈值范围内确定土地平整强度最佳阈值。土地平整强度最佳阈值确定后，便可在无需任何训练样本的情况下，实现对生产建设用地的全自动检测。
[0074]
下面以8个地区作为研究区域对本实施例提出的方法和现有的方法进行对比分析。
[0075]
测试区域
[0076]
8个地区分别选自全球8个气候带，分别为北京、柏林、开罗、芝加哥、广州、吉隆坡、孟买和罗马。
[0077]
测试数据
[0078]
本实施例测试数据有sentinel-2遥感卫星收集，sentinel-2是欧洲航天局哥白尼计划的组成部分之一，包含两颗相同卫星sentinel-2a和sentinel-2b，其获取的影像分辨率为10米，重访周期为5天。在每个城市，收集不同季节、不同时间间隔的生产建设用地的光学遥感影像。不同季节代表参考影像获取的季节，共包括四个季节，分别为季节1(1月至3月)、季节2(4月至6月)、季节3(7月至9月)和季节4(10月至12月)；不同时间间隔为参考影像和检测影像的时间间隔，共涉及5个不同时间间隔，分别为5天、1个月、3个月、6个月和1年。在每个城市共开展了20次生产建设用地检测，8个城市共开展160次检测。如下表所示，为每次检测时参考影像和检测影像获取的日期，其中r表示参考影响获取日期，d表示检测影像获取日期。
[0079]
[0080]
[0081]
[0082][0083]
如下表所示，为每次检测时验证样本数量，即实际生产建设用地数量及像素数，非生产建设用地及像素数；其中，ld表示生产建设用地数量，括号内为生产建设用地包含的像素数；ua表示非生产建设用地数量，括号内为非生产建设用地包含的像素数；na表示没有验证样本。
[0084]
[0085]
[0086][0087]
研究方法及效果
[0088]
利用现有方法进行，具体为计算参考影像和检测影像间的土地开发强度，提取建设用地开发地块，该方法没有考虑农作物、自然植被、山体及建筑物阴影季节性变化的影响。检测结果的精度评估主要基于三个统计量，分别为检测精度(da)、误报率(fr)和整体精度(oa)。检测精度表示正确检测出的生产建设用地的百分比。误报率表示错误标记的非生产建设用地的百分比。整体精度是正确标记的验证样本(即生产建设用地和非生产建设用地之和)的百分比。如图3所示，显示了现有方法的检测结果。图中，s1、s2、s3和s4代表不同季节；5d、1m、3m、6m、1y分别表示5天、1个月、3个月、6个月和1年的检测时间间隔；其中所有检测结果的检测精度均高于90.00％，但误报率变化很大，在短时间间隔(例如5天和1个月)的检测中，误报率通常较低，而在较长时间间隔(例如3和6个月)的检测中，误报率普遍较高。当检测时间间隔增加到1年时，误报率再次变低。在相同地区、相同季节开展检测时，5天时间间隔的检测得到的误报率均低于0.10％，而随着检测时间间隔的增加，误报率也开始变高，最高可达71.60％(如在孟买第三季度开展的时间间隔为3个月的检测)。经过分析误报原因来自三个方面，分别为农业活动、植被变化和阴影变化。季节性的农业生产活动，如作物的种植和收获，也需要进行土地平整，导致从植被到裸土的变化。区分作物生长导致的土地平整和生产建设用地开发导致的土地平整非常困难，因此很多农业生产活动导致的土地平整被错误地检测为生产建设用地开发，尤其是在北京、开罗、柏林、罗马和孟买等具有
大面积耕地的城市。在北京、柏林和芝加哥等地区，森林通常在温暖潮湿的夏季枝叶繁茂，而在冬季落叶。森林植被的季节性变化也会导致计算出的土地开发强度偏高，给生产建设用地开发检测带来大量误报。随着季节的变化，太阳方位角和天顶角也会发生变化，导致相同的建筑和山体阴影在不同季节差异明显，给生产建设用地开发检测带来误报。当时间间隔较短时(例如5天或1个月)，农作物、自然植被、山体和建筑物阴影变化较小，产生的误报也相对较少。随着时间间隔的增加，农作物、自然植被和阴影变化显著，导致在长时间间隔(例如3和6个月)检测中的误报较多。鉴于农作物、自然植被和阴影在不同年份相同时间的状态相似，因此在1年时间间隔的检测中导致的误报再次变少。
[0089]
针对误报率高于10.00％的检测，采用本实施例提出的方法进行检测。如图4所示，通过获取早于检测图像1年的矫正图像，消除检测图像与参考图像中季节性变化导致的误报。如图5所示，为本实施例提出的方法的检测结果与现有方法检测结果的精度对比图。图中，bj、bl、ca、mu和ro分别代表北京、柏林、开罗、孟买和罗马；s1、s2、s3和s4代表不同季节；3m、6m分别表示3个月和6个月的检测时间间隔；从图5中可以看出，本方法的在保持很高检测精度的同时，大幅度减少了误报率，最高可将误报率从71.60％较少到0.20％，说明本实施例可以有效减少农作物、自然植被、山体和建筑季节性变化导致的误报，显著提高生产建设用地的监测精度。本实时提出的方法可用于国土、水利部门，提高对生产建设用地的监测频次及监测精度，及时发现并预防违法用地和水土流失，确保耕地红线不被突破，避免难以逆转的环境问题和巨额的经济损失。
[0090]
显然，本发明的上述实施例仅仅是为清楚地说明本发明所作的举例，而并非是对本发明的实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式的变化或变动。这里无需也无法对所有的实施方式予以穷举。凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明权利要求的保护范围之内。