Google Earth Engine(GEE)农作物莳植布局提取

[复制链接]
查看650 | 回复0 | 2023-8-23 11:44:20 | 显示全部楼层 |阅读模式
写在前面

前段时间由于考研的缘故因由不停没能更新,已经完成了农作物莳植布局的提取,现在给大家分享一下。
紧张也是联合前面写过的Google Earth Engine(GEE)利用土地利用数据(modis)上采样Landsat数据提取农田范围,将以上结果作为研究的基础,联合物候特性,光谱特性,地形特性,选择随机丛林算法进行农作物的提取。
不想看背面的同砚可以直接看代码:https://code.earthengine.google.com/3ed8a5303c610e063ae3a4517ca4e146
1.构建物候特性

在这里紧张通过不同农作物NDVI值的差异来选择符合的生长时间进行影像的选择。
通过构建典范地物的NDVI时序特性曲线,不同农作物随着时间的变革,NDVI值也在跟着变革,这种变革也反映了不同农作物在不同时期的光谱特性差异。因此,选择NDVI差异较大的越冬期和成熟期作为最佳辨认期进行农作物分类,可以有用进行农作物种类的区分。
  1. var winter = ee.ImageCollection('COPERNICUS/S2_SR')
  2.                   .filterDate('2021-11-30', '2022-01-01')
  3.                   .filterBounds(roi)
  4.                   .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 10))
  5.                   .map(maskS2clouds)
  6.                   .map(mndwi)
  7.                   .map(ndbi)
  8.                   .map(ndvi)
  9.                   .select(['B1','B2','B3','B4','B8','B12','QA60','MNDWI','NDBI','NDVI']);
  10. var winter1=winter.median()
  11.                  .addBands(elevation.rename("ELEVATION"))
  12.                  .addBands(slope.rename("SLOPE"))
  13.                  .clip(roi);
  14. Map.addLayer(winter1, {bands: ['B4', 'B3', 'B2'],min:0, max: 0.3}, 'winter1',false);
  15. var summer = ee.ImageCollection('COPERNICUS/S2_SR')
  16.                   .filterDate('2022-03-01', '2022-04-01')
  17.                   .filterBounds(roi)
  18.                   .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 10))
  19.                   .map(maskS2clouds)
  20.                   .map(mndwi)
  21.                   .map(ndbi)
  22.                   .map(ndvi)
  23.                   .select(['B1','B2','B3','B4','B8','B12','QA60','MNDWI','NDBI','NDVI']);
  24. var summer1=summer.median()
  25.                  .addBands(elevation.rename("sum_ELEVATION"))
  26.                  .addBands(slope.rename("sum_SLOPE"))
  27.                  .clip(roi);
复制代码
这里我把基础影像根据农田范围进行了裁剪,只留下农地步区,能有用的清除一些非农地步物范例。

2.构建光谱特性

也就是构建NDVI、NDWI、NDBI等光谱指数,来区分植被、修建,水体等。
  1. //构建光谱特征
  2. function mndwi(image){
  3.   return image.addBands(image.normalizedDifference(['B3', 'B8']).rename('MNDWI'))
  4. }
  5. function ndbi(image){
  6.   return image.addBands(image.normalizedDifference(['B12', 'B8']).rename('NDBI'))
  7. }
  8. function ndvi(image){
  9.   return image.addBands(image.normalizedDifference(['B8', 'B4']).rename('NDVI'))
  10. }
复制代码
3.将全部影像归并为一幅影像

将上面选择的越冬期和成熟期的多幅影像归并为一幅,在这里需要更改波段名,否则会出现波段名重复的情况。
  1. //更改波段名称,合成一幅影像
  2. var useimage = winter1.addBands([
  3.   summer1.select('B2').rename('sum_B2'),
  4.   summer1.select('B3').rename('sum_B3'),
  5.   summer1.select('B4').rename('sum_B4'),
  6.   summer1.select('B8').rename('sum_B8'),
  7.   summer1.select('B12').rename('sum_B12'),
  8.   summer1.select('MNDWI').rename('sum_MNDWI'),
  9.   summer1.select('NDBI').rename('sum_NDBI'),
  10.   summer1.select('NDVI').rename('sum_NDVI'),
  11.   summer1.select('sum_ELEVATION'),
  12.   summer1.select('sum_SLOPE')
  13.   ])
  14.   
  15. print(useimage,'useimage')
  16. //求交集,对影像进行掩膜
  17. var useimage1 = useimage.updateMask(esamask)
  18. // var useimage1 = useimage
  19. Map.addLayer(useimage1, {bands: ['sum_B4', 'sum_B3', 'sum_B2'],min:0, max: 0.3}, 'useimage1');
  20. var usebands = ['B2','B3','B4','B8','B12','MNDWI','NDBI','NDVI','ELEVATION','SLOPE',
  21.                 'sum_B2','sum_B3','sum_B4','sum_B8','sum_B12','sum_MNDWI','sum_NDBI','sum_NDVI','sum_ELEVATION','sum_SLOPE']
复制代码
4.构建随机丛林算法进行分类

构建算法也就是一些老套路了,无非就是创建样本,将样天职成训练样本和验证样本。
在这里想说一点,随机丛林的决定树棵树选择不能随便选择,需要进行不同棵数的实验,可以利用以下代码进行分析,选择精确率最高的颗数进行算法构建。
  1. //选取森林棵树
  2. var numTrees = ee.List.sequence(5, 50, 5);
  3. var accuracies = numTrees.map(function(t)
  4. {
  5.    var classifier = ee.Classifier.smileRandomForest(t)
  6.                      .train({
  7.                      features: trainingPartition,
  8.                      classProperty: 'landcover',
  9.                      inputProperties: usebands
  10.                      });
  11.    return testingPartition
  12.        .classify(classifier)
  13.        .errorMatrix('landcover', 'classification')
  14.        .accuracy();
  15. });
  16. print(ui.Chart.array.values({
  17.    array: ee.Array(accuracies),
  18.    axis: 0,
  19.    xLabels: numTrees
  20. }));
复制代码
得到的结果大概是如许,选择精确率最高的就可以了。

5.算法的存储

由于我们上面已经训练好了算法,以是这个算法是可以运用在其他年份的,我们将它生存在个人assets里,当缺少其他年份的样本时,可以利用该算法进行分类,也就是算法复用。
  1. var trees = ee.List(ee.Dictionary(testclassifier.explain()).get('trees'))
  2. var dummy = ee.Feature(roi.geometry())
  3. var col = ee.FeatureCollection(trees.map(function(x){return dummy.set('tree',x)}))
  4. print(col)
  5. Export.table.toAsset(col,'save_classifier','projects/dyb/assets/2022chuzhou')
复制代码
6.面积统计

上面已经做好了农作物的分类,接下来可以对农作物的面积做一个统计。
  1. var areaImage = ee.Image.pixelArea().addBands(classified)
  2. var clsdArea = areaImage.reduceRegion({
  3.     'reducer': ee.Reducer.sum().group({
  4.         'groupField': 1,
  5.         'groupName': 'class'
  6.     }),
  7.     'geometry': roi.geometry(),
  8.     'scale': 10,
  9.     'maxPixels': 1e13
  10. })
  11. print('Kerala landcover clsArea:', clsdArea.getInfo())
复制代码
最后把全部代码放在这里
https://code.earthengine.google.com/3ed8a5303c610e063ae3a4517ca4e146


来源:https://blog.csdn.net/qq_51118386/article/details/128572478
免责声明:如果侵犯了您的权益,请联系站长,我们会及时删除侵权内容,谢谢合作!

本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有账号?立即注册

x
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则