基于Sentinel-2(哨兵2)的杞县大蒜提取试验
一、前言
杞县大蒜已有2000多年的栽培历史是全国农产品地理标志之一,2009年9月14日,中华人民共和国农业部正式批准对“杞县大蒜”实施农产品地理标志登记保护。杞县已建成大蒜冷藏保鲜库600多座,年贮藏能力达到25万吨。建成蒜片、蒜米、蒜粒、蒜油等加工企业80多家,年加工能力5万多吨。已在全国20多个大中城市设立了大蒜销售网点,且远销到美国、欧洲、日本、韩国、东南亚、印尼等几十个国家和地区,每年直接间接出口大蒜制品25多万吨。
Sentinel-2(哨兵2)超大的幅宽,双星组合,较短的重访周期以及丰富的波谱信息都为农作物提取提供了宝贵的数据源。今天我们就来通过Sentinel-2来进行杞县大蒜提取试验。
二、数据选择
本次实验选取了三期的Sentinel-2 L1C数据,通过预处理(Sen2Cor大气校正)、SNAP采样输出到ENVI中进行 波段组合、裁剪等工作。
下面是本次实验的数据,分别为2019年3月14日、2018年4月8日、2019年4月28日,三期数据,下图是分别进行 B4B3B2 真彩色组合:
经过反复测试 B11B8B3 的假彩色组合显示最能突出地物间的差别,层次分明:
三、外业采样
本次选用手机APP:LocaSpace Viewer(LSV)软件来辅助外业采集样本,手机端打开GPS+WLAN+移动网络的定位模式,然后采集田块的角点作为坐标点。下面我们举例来说明:例如W1 W2 W3 W4为小麦样本的坐标点,G1 G2 G3 G4为大蒜样本的坐标点:
然后将样本点导入ENVI当中,蓝色点为大蒜地块的样本角点,红色点为小麦样本地块的样本角点:
采集的大蒜样本集:
采集的小麦样本集:
由上面的样本集图可见,大蒜的颜色在B11B8B3组合下呈现蓝紫色,小麦的颜色在B11B8B3组合下呈现亮绿色,下面我们分别来看两者的光谱曲线:
大蒜的光谱曲线:
小麦的光谱曲线:
对比两者的光谱曲线可以看到,在近红外波段小麦的反射率比大蒜要高,在绿波段和红波段大蒜的反射率比小麦要高,在短波红外波段,大蒜的反射率比小麦高。因此在B11 B8 B3的组合下,小麦要比大蒜更绿,大蒜要比小麦更蓝和更红。
四、2019年数据分类
我们先对2019年的数据进行分类,本次分类采用监督分类的方法,采用的是支持向量机的模型进行分类。因为数据时间是2019年3月14日,树木尚未返青,因此我们选取了大蒜、小麦、水域、居民地、空闲地五大类,并对样本进行了可分离度评价,结果如下:
由图可以看到各个样本之间的可分离度良好,基本都在1.9以上:
执行监督分类,下图依次是监督分类的结果和其中的大蒜结果:
五、2018年数据分类
5.1 监督分类
我们对2018年的数据进行监督分类,仍然采用支持向量机的模型进行分类。因为数据时间是2018年4月8日,树木已经返青,因此我们选取了大蒜、小麦、林地、水域、居民地、空闲地五大类,首先进行大蒜样本的选择:
大蒜在2018年4月8日的影像上依据B11 B8 B3组合显示下呈现蓝绿色:
接下来进行小麦样本的采集:
小麦在2018年4月8日的影像上依据B11 B8 B3组合显示下呈现鲜艳的亮绿色:
接着采集林地样本,林地多位于村庄的周边,呈现不规则形状:
后面我们依次采集水域、居民地、空闲地,并对样本进行了可分离度评价,结果如下:
各个样本之间可分离度良好,均在1.98以上:
执行SVM监督分类后的结果,其中蓝色为大蒜,黄色为小麦,红色为居民地,绿色是林地:
但是比对影像进行结果检查我们发现有一些田间的道路被错分到了大蒜里面:
5.2 基于规则的面向对象提取
下面我们换个方法进行分类,之前2019年3月14日的影像,因为时间还未处在植物的生长季节,因此其NDVI指数不高,但是本次2018年4月8日的影像,植物已经陆续拔尖返青,其NDVI指数处在较高的水平,因此我们可以考虑结合NDVI进行面向对象提取。
首先计算NDVI,采样B4,B8作为红和近红外进行计算,计算结果如下:
然后我们将之前监督分类选的样本加载进来进行统计,因为考虑到部分彩钢房的NDVI也比较高,基本跟植被的指数值相近,因此我们本次在上次的样本上增加了居民地2的样本:
统计后的样本NDVI如下:
由上图可见,大蒜的NDVI区间为0.5996-0.8434,小麦的NDVI区间为0.7938-0.9397,两者的NDVI区间存在一定的重叠,因此采用此NDVI作为指标进行提取可能会存在一定的误分,错分。
下面我们调整NDVI的计算参数,将B5(红边的第一波段)作为红,近红外依然为B8 进行NDVI计算:
跟之前的NDVI进行对比,两者进行对比发现将B5、B8作为参数进行计算的NDVI,图像地物的对比层次更加明显:
然后再次进行样本统计,统计后的样本NDVI如下:
由上图可见,大蒜的区间值为0.3888-0.6141,小麦的区间值为0.6655-0.7977,两者区间不重叠,可以作为分离指标,但是林地和居民地2的区间与大存在一定的重叠,因此接下来要找分离大蒜和林地、居民地2的区间。
下面我们统计样本在影像中的各个波段:
从图中可以看出大蒜和林地的分离可以采用B6( band5,因为本次实验剔除了60米波段,因此band5即为哨兵的B6波段),大蒜和居民地2的分离可以采用B12( band10,因为本次实验剔除了60米波段,因此band10即为哨兵的B12波段)。
下面我们进行基于规则的面向对象提取,计算NDVI我们选择band4和band7(因为本次实验剔除了60米波段,所以band4和band7即为B5、B8):
设置分割参数:
设置分离区间:
这个区间主要分离大蒜、小麦和居民地、水域、空闲地等,权重设为0.6。
这个区间主要分离大蒜和林地,权重设为0.2。
这个区间主要分离大蒜和居民地2,权重设为0.2。
输出结果:
但是我们放大后检查结果仍然会发现,田间小路的区分不是很好,仍会被错分到大蒜里:
这是因为一些田间道路的NDVI值正好处在大蒜的区间之内:
5.3 基于样本的面向对象提取
上面的监督分类和基于规则的面向对象提取均会存在将田间小路错分到大蒜的结果当中,为了改变这个结果,我们接下来尝试基于样本的面向对象提取。
为了充分的利用现有的数据,提高分类精度,我们先将B4、B8计算出来的NDVI与原影像进行波段组合,为影像增加一个可参考的图层:
然后我们打开基于样本的面向对象提取的时候,NDVI计算选择band4和band7(因为本次实验剔除了60米波段,所以band4和band7即为B5、B8),这样我们的数据里面就有了两个NDVI图层(B4、B8NDVI 和B5、B8 NDVI)。
接下来进行分割,为了分割出足够细小的地物,比如田间小路,我们将分割参数适当的设置小一些:
然后在分割出来的图像上选择样本,首先选择大蒜样本:
后面依次选择了小麦、林地、水域、居民地、空闲地样本,另外为了分离出田间道路,又增加了道路样本。然后在可选波段里面进行全选:
最后输出结果:
采用这种方法,道路基本被分离出来:
我们找一个田间道路比较多地方综合对比一下,由下图可见基于样本的面向对象提取结果的道路分割的更细致:
六、2017年数据分类
最后我们对2017年数据进行分类,为了提高运行速度,我们修改一下ENVI的配置参数,Cache Size缓存大小可以设置为物理内存的50%,Image Tile Size,假如内存为8G 64位系统,可以设置为50MB:
这次我们仍然采用基于样本的面向对象提取,但是这次为了提高效率,我们先选ROI采集样本,然后将样本导出为SHP格式,在ArcGIS里利用面转点的功能将ROI面转为点SHP:
接下来进行基于样本的面向对象提取操作,过程与上面2018年数据的操作一样,最后选择样本的时候,我们不在分割影像上选了,直接导入刚才的点SHP即可:
然后输出结果:
七、三期结果对比分析
对分类结果进行计算得出:
结合网络新闻,杞县大蒜常年稳定种植面积60万亩左右,本次三年的统计结果的平均值为63.3万亩,与之基本符合:
大家可能还记得2017年的“”蒜你狠“”的相关报道,2017年大蒜价格居高不下,这是跟种植面积少,产量低有关,2017年的高价大蒜直接导致了2018年的农民争相跟风种植:
这是一条2019年的报道:
报道当中的去年即为2018年,报道中为73.2万亩,基本与本次提取结果76万亩接近,报道中今年2019年62.17万亩,本次提取结果65万亩,也基本接近。
2018年农民跟风种植大蒜直接导致了蒜薹滞销,价格低迷,蒜农宁肯把蒜薹倒沟里:
好了,本次实验介绍就到这里。
- END -