GEE：对Sentinel-2遥感影像进行处理，水体提取与可视化

news2024/11/30 4:46:00

作者：CSDN @ _养乐多_

本文介绍了通过Google Earth Engine平台，并使用哨兵数据提取水体掩膜的方法和代码。通过裁剪和去除云等处理步骤，最终得到具有水体掩膜的影像，并进行可视化和导出。这种方法基于归一化水体指数（NDWI）和OTSU阈值计算技术，无需复杂的图像处理算法，适用于快速获取水体信息的需求。

图1 彩色合成

图2 水体掩膜

文章目录

- - 一、代码详解
  - 二、代码链接

一、代码详解

var roi = table          
Map.centerObject(roi,10)

// 可视化参数
var palette = ['FFFFFF', 'CE7E45', 'DF923D', 'F1B555', 'FCD163', '99B718',
               '74A901', '66A000', '529400', '3E8601', '207401', '056201',
               '004C00', '023B01', '012E01', '011D01', '011301'];
               
//可视化参数，按843波段合成                
var rgbVis = {
  min: 0.0,
  max: 0.35,
  bands: ['B8', 'B4', 'B3'],
};
               
//按矢量边界裁剪
function roiClip(image){
  return image.clip(roi)
}

//S2_SR去云
function remove_cloud(image) {
  var qa = image.select('QA60')
  var cloudBitMask = 1 << 10;
  var cirrusBitMask = 1 << 11;
  var mask = qa.bitwiseAnd(cloudBitMask).eq(0).and(qa.bitwiseAnd(cirrusBitMask).eq(0))
  return image.updateMask(mask).divide(10000).select("B.*").copyProperties(image, ["system:time_start"])
}

//时间范围
var startDate = ee.Date('2020-01-01');
var endDate = ee.Date('2020-02-01');

//数据集选择
var S2 = ee.ImageCollection("COPERNICUS/S2_SR")
           .filterDate(startDate, endDate)
           .select('B3', 'B4', 'B8', 'B11', 'QA60')
           .filterBounds(roi)
           .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 5))
           .map(remove_cloud)
           .map(roiClip)


//运行同一天影像
var newImgS2 = S2.max();
print(newImgS2)
Map.addLayer(newImgS2, rgbVis, 'newImg_RGB'); 

//计算NDWI
function calNDWI(image)
{
  var NDWI = image.normalizedDifference(['B3', 'B11']).rename('NDWI');
  NDWI = NDWI.updateMask(NDWI.gt(-1).and(NDWI.lt(1)));
  return image.addBands(NDWI);
}

//OTSU计算阈值
function otsu(histogram) {
  var counts = ee.Array(ee.Dictionary(histogram).get('histogram'));
  var means = ee.Array(ee.Dictionary(histogram).get('bucketMeans'));
  var size = means.length().get([0]);
  var total = counts.reduce(ee.Reducer.sum(), [0]).get([0]);
  var sum = means.multiply(counts).reduce(ee.Reducer.sum(), [0]).get([0]);
  var mean = sum.divide(total);
  var indices = ee.List.sequence(1, size);
  var bss = indices.map(function(i) {
    var aCounts = counts.slice(0, 0, i);
    var aCount = aCounts.reduce(ee.Reducer.sum(), [0]).get([0]);
    var aMeans = means.slice(0, 0, i);
    var aMean = aMeans.multiply(aCounts)
        .reduce(ee.Reducer.sum(), [0]).get([0])
        .divide(aCount);
    var bCount = total.subtract(aCount);
    var bMean = sum.subtract(aCount.multiply(aMean)).divide(bCount);
    return aCount.multiply(aMean.subtract(mean).pow(2)).add(
           bCount.multiply(bMean.subtract(mean).pow(2)));
  });
  return means.sort(bss).get([-1]);
}

//计算水体掩膜
function calWaterMask(image){
  var MNDWI_S2 = image.select("NDWI");
  var histogram = MNDWI_S2.reduceRegion({
    reducer: ee.Reducer.histogram(), 
    geometry: roi, 
    scale: 30,
    maxPixels: 1e13,
    tileScale: 8
  });
  var threshold_S2 = otsu(histogram.get("NDWI"));
  var mask = MNDWI_S2.gte(threshold_S2);
  return image.addBands(mask.rename('mask'));
}

//运行水体掩模函数
newImgS2 = calNDWI(newImgS2)
newImgS2 = calWaterMask(newImgS2).select('mask')
print(newImgS2)
Map.addLayer(newImgS2, {min:0, max:1, palette:['#DDDDDD', '#0099FF']}, 'S2_Water_body_mask');

//下载非水体概率影像
Export.image.toDrive({
  image:newImgS2,
  description: 'S2_water_mask',
  scale:30,
  region:roi,
  fileFormat: 'GeoTIFF',
  maxPixels:1e13,
});

这段代码的目的是对Sentinel-2遥感影像进行处理，提取水体掩膜，并将结果可视化和导出。

首先，定义了一个变量roi，表示感兴趣区域（Region of Interest），这是一个表格（table）对象。然后使用Map.centerObject(roi,10)将该区域设置为地图的中心，并指定缩放级别为10。
接下来，定义了一个调色板（palette）数组，用于可视化图像时的颜色映射。
然后，定义了一个可视化参数rgbVis，设置了最小值（min）和最大值（max），以及要显示的波段（bands）。
之后，定义了一个函数roiClip，用于将图像按照矢量边界进行裁剪。
接着是一个函数remove_cloud，用于去除图像中的云和雾。函数内部使用了位掩码（bitmask）的方式来选择需要保留的像素，然后将图像的像素值除以10000并选择特定波段，最后将处理后的图像的时间属性复制到结果中。
然后定义了起始日期和结束日期，用于选择要处理的影像时间范围。
接下来，使用ee.ImageCollection函数选择了名为"COPERNICUS/S2_SR"的Sentinel-2影像集合，并通过一系列操作对影像进行筛选、裁剪和去除云，最后将结果映射到ROI区域，并存储在变量S2中。
之后，使用S2.max()函数选取ROI区域内同一天的影像，并将结果存储在变量newImgS2中。
接着，定义了一个名为calNDWI的函数，用于计算归一化水体指数（NDWI）。
然后，定义了一个名为otsu的函数，用于计算阈值，该阈值用于将图像的像素分为水体和非水体。
接下来，定义了一个名为calWaterMask的函数，该函数通过计算NDWI的直方图，并使用OTSU方法计算阈值来生成水体掩膜。
然后，分别调用calNDWI和calWaterMask函数对newImgS2进行处理，并将结果可视化。
最后，使用Export.image.toDrive函数将非水体概率影像导出到Google Drive中。