Google Earth Engine（GEE）使用土地利用数据（modis）上采样Landsat数据提取农田范围

前言

最近因为大创项目，一直在研究GEE平台的一些操作，主要在GEE中利用海量的遥感影像，来实现关于滁州市的农作物种植结构的长时间序列分析。本篇主要利用GEE中发布的土地利用数据（modis）500m分辨率进行上采样到Landsat的30m影像中，提取农田范围，为后续农作物分类提供帮助。

1.导入研究区域

首先在Assets中导入研究区域边界的shp文件，导入到代码中。

2.选择影像数据

选择MCD12Q1全球土地覆盖产品和Landsat7影像，并按时间进行过滤筛选。

var modisLandcover =  ee.ImageCollection('MODIS/006/MCD12Q1')
    .filterDate('2019-01-01', '2019-12-31')
    .first()
    .select('LC_Type1')
    .subtract(1)       //将波段排序从0开始排列
    .clip(geometry);   //切割出研究区范围

// 设置IGBP分类系统配色方案
var landcoverPalette = '05450a,086a10,54a708,78d203,009900,c6b044,dcd159,' +
 'dade48,fbff13,b6ff05,27ff87,c24f44,a5a5a5,ff6d4c,69fff8,f9ffa4,1c0dff';
//设置土地分类可视化参数
var landcoverVisualization = {palette: landcoverPalette, min: 0, max: 16, format: 'png'};

var l7 = ee.ImageCollection('LANDSAT/LE07/C01/T1')
    .filterBounds(geometry)
    .filterDate('2019-01-01', '2020-01-01');
// 最小云量影像合成
var landsatComposite = ee.Algorithms.Landsat.simpleComposite({
  collection: l7,
  asFloat: true
});

MCD12Q1土地利用数据中LC_Type1的第12个波段表示农田，所以我们后续重点提取bands12属性。

3.选取研究波段

根据modis数据属性可知，bands12表示农田数据，所以通过编写函数，进行指定波段提取。

//选取农田波段
var mask=function(image,my_type){
  var mask=image.eq(my_type);
  var masked_image=image.updateMask(mask);
  return masked_image;
}
var modismask = mask(modisLandcover,11)  //调用mask函数，筛选出bands12农田数据

4.通过modis上采样到Landsat影像中

首先要创建数据集，在MCD12Q1中随机采样1000个随机点作为数据集，并将数据集通过smileCart方法进行训练，将训练器运用到Landsat影像中。

// 创建训练数据集，采样1000个随机采样点
var training = modisLandcover.addBands(landsatComposite).sample({
  region: geometry,
  scale: 30,
  numPixels: 1000
});

// 利用训练数据集进行训练，选取IGBP分类系统
var classifier = ee.Classifier.smileCart().train({
  features: training,
  classProperty: 'LC_Type1',
});

// 将训练器运用到Landsat影像中
var upsampled = landsatComposite.classify(classifier);

5.输出结果分析

这是modis数据通过波段筛选，筛选出的农田范围。

这是通过modis数据进行采样，进行训练后的landsat影像中的农田范围。

结果可以看到通过modis数据500m分辨率进行上采样到landsat影像30m分辨率上，精度有着明显的提高，可以分辨出更细小的农田范围。

6.输出影像存入Assets中

Export.image.toAsset({
  image: result1,
  scale:10,
  maxPixels:5e8,
  region: geometry
});

总结

提取出了滁州市的农田范围数据，下一步将进行农作物分类，通过modis数据上采样到landsat数据中，提高了后续农作物分类的分辨率问题，只选择了表示农田的bands12，避免了其他数据对农作物分类的干扰。
完整代码：

//选择2001年数据的IGBP分类系统
var modisLandcover =  ee.ImageCollection('MODIS/006/MCD12Q1')
    .filterDate('2019-01-01', '2019-12-31')
    .first()
    .select('LC_Type1')
    .subtract(1)
    .clip(geometry);
//选取农作物波段
var mask=function(image,my_type){
  var mask=image.eq(my_type);
  var masked_image=image.updateMask(mask);
  return masked_image;
}
var modismask = mask(modisLandcover,11)
// 设置IGBP分类系统配色方案
var landcoverPalette = '05450a,086a10,54a708,78d203,009900,c6b044,dcd159,' +
 'dade48,fbff13,b6ff05,27ff87,c24f44,a5a5a5,ff6d4c,69fff8,f9ffa4,1c0dff';
//设置土地分类可视化参数
var landcoverVisualization = {palette: landcoverPalette, min: 0, max: 16, format: 'png'};
//显示影像
Map.centerObject(geometry, 9);
Map.addLayer(modismask, landcoverVisualization, 'MODIS landcover');

// 筛选Landsat影像
var l7 = ee.ImageCollection('LANDSAT/LE07/C01/T1')
    .filterBounds(geometry)
    .filterDate('2019-01-01', '2020-01-01');
// 最小云量影像合成
var landsatComposite = ee.Algorithms.Landsat.simpleComposite({
  collection: l7,
  asFloat: true
});
//影像显示
Map.addLayer(landsatComposite, {min: 0, max: 0.3, bands: ['B3','B2','B1']}, 'Landsat composite');

// 创建训练数据集，采样1000个随机采样点包含MCD12Q1与Landsat属性
var training = modisLandcover.addBands(landsatComposite).sample({
  region: geometry,
  scale: 30,
  numPixels: 1000
});

// 利用训练数据集进行训练，选取IGBP分类系统
var classifier = ee.Classifier.smileCart().train({
  features: training,
  classProperty: 'LC_Type1',
});

// 将训练器运用到Landsat影像中
var upsampled = landsatComposite.classify(classifier);
var result = mask(upsampled,11)
var result1 = result.clip(geometry);
// 可视化分类结果
Map.addLayer(result1, landcoverVisualization, 'Upsampled landcover');
// 显示分类区域
Map.addLayer(geometry, {}, 'Training region', false);

Export.image.toAsset({
  image: result1,
  scale:10,
  maxPixels:5e8,
  region: geometry
});