基于 GEE 利用 Sentinel-1 双极化数据计算 SDWI 指数实现逐月提取水域面积

目录

1 SDWI 指数

2 研究方法及数据处理

3 完整代码

4 运行结果

1 SDWI 指数

Sentinel-1双极化数据SDWI水体提取指数公式：SDWI = ln ⁡(10×VV×VH)

Sentinel-1 Dual-Polarized Water Index (SDWI)水体信息提取方法对 Sentinel-1 双极化数据(VV 和 VH)之间水体信息提取的关系进行深入研究，达到增强水体特征的效果，水体区分较明显，同时消除土壤和植被对水体提取中造成的干扰。公式将 VV 和 VH 极化影像相乘，并且乘以 10，以此扩大水体与其他地物之间的差异，再以自然对数作为函数式。

2 研究方法及数据处理

• 水体识别方法

传统光学遥感因云雾遮挡往往难以实现连续监测，而本研究利用Sentinel-1雷达数据提出了一种综合方法，通过SDWI指数（Sentinel-1 Derived Water Index）实现高精度水体识别。

• SDWI指数计算

SDWI指数定义为VV和VH两个极化通道之比的对数值，SDWI的物理意义是利用水面对雷达波的镜面反射特性与植被和土地表面散射特性的差异，精确地区分水体与非水体。

• 掩膜过滤

为进一步提高精度，本研究引入了多个掩膜叠加条件：

（1）DEM坡度掩膜：利用SRTM数据计算坡度，排除坡度大于5°的地区，减少山区阴影的干扰。

（2）JRC全球水体掩膜：仅保留长期存在的水体区域（发生概率超过10%），减少短暂积水或湿地误判。

（3）VH极化通道筛选：设定VH阈值（< -15 dB）以进一步提高水体分类的纯度。

• 逐月水体面积动态监测

（1）Sentinel-1影像每月中值合成，降低单景图像噪声影响。

（2）应用综合掩膜提取水体区域。

（3）基于提取的水体掩膜，利用GEE进行逐月面积统计，并导出成表格。

3 完整代码

// 引入感兴趣区（ROI）
var roi_x = table;
Map.addLayer(roi_x);
Map.centerObject(roi_x, 8);

// 设置时间范围
var startYear = 2023;
var endYear = 2024;
var months = ee.List.sequence(1, 12);

// Sentinel - 1 数据集
var s1Collection = ee.ImageCollection('COPERNICUS/S1_GRD')
  .filterBounds(roi_x)
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VV'))
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VH'))
  .filter(ee.Filter.eq('instrumentMode', 'IW'));

// 引入 SRTM DEM 数据并计算坡度
var dem = ee.Image("USGS/SRTMGL1_003");
var slope = ee.Terrain.slope(dem);
var slopeThreshold = 5;
var landMask = slope.lt(slopeThreshold);

// **引入 JRC 全球水体数据**
var jrcWater = ee.Image("JRC/GSW1_3/GlobalSurfaceWater")
  .select("occurrence")
  .clip(roi_x);

// **JRC 水体阈值**
var jrcWaterThreshold = 10; // 仅保留出现概率 >50% 的水体区域
var jrcMask = jrcWater.gt(jrcWaterThreshold);

// 计算 Sentinel - 1 SDWI
var calculateSDWI = function (image) {
  var VV = image.select('VV');
  var VH = image.select('VH');

  // 计算 SDWI（Sentinel - 1 水体指数）
  var sdwi = VV.divide(VH).log10().rename('SDWI');

  return image.addBands(sdwi).addBands(VH.rename('VH')); // 添加 VH 以便后续筛选
};

// 计算逐月水体掩膜面积
var calculateMonthlyWaterArea = function (year, month) {
  var startDate = ee.Date.fromYMD(year, month, 1);
  var endDate = startDate.advance(1, 'month');

  // 过滤 Sentinel - 1 图像集
  var monthlyCollection = s1Collection
    .filterDate(startDate, endDate)
    .map(calculateSDWI);

  // 计算 SDWI 中值合成
  var sdwiComposite = monthlyCollection.median().clip(roi_x);
  var VH = sdwiComposite.select('VH');

  // 计算水体掩膜（结合 DEM 坡度、VH 极化、JRC 水体数据）
  var waterMask = sdwiComposite.select('SDWI').gt(-0.5)
    .and(sdwiComposite.select('SDWI').lt(0.3))
    .and(VH.lt(-15)) // 低 VH 值通常表示水体
    .and(landMask) // 过滤山体阴影
    .and(jrcMask) // 结合 JRC 长期水体数据
    .rename('WaterMask');

  // 计算水体面积（平方千米）
  var areaImage = waterMask.multiply(ee.Image.pixelArea()).divide(1e6);
  var waterArea = areaImage.reduceRegion({
    reducer: ee.Reducer.sum(),
    geometry: roi_x,
    scale: 10,
    maxPixels: 1e9
  }).get('WaterMask');

  return ee.Feature(null, {
    'Date': startDate.format('yyyy/M/d'),
    'Year': year,
    'Month': month,
    'WaterArea_km2': waterArea
  });
};

// 对每个月进行水体面积计算并转换为 FeatureCollection
var calculateWaterAreaByMonth = function (year) {
  return ee.FeatureCollection(months.map(function (month) {
    return calculateMonthlyWaterArea(year, month);
  }));
};

// 计算所有年份的水体面积（单位为平方千米）
var allYearsWaterArea = ee.FeatureCollection(
  ee.List.sequence(startYear, endYear).map(calculateWaterAreaByMonth)).flatten();

// 导出逐月水体面积到 CSV
Export.table.toDrive({
  collection: allYearsWaterArea,
  description: 'Refined_Monthly_Water_Area_' + startYear + '_' + endYear,
  fileFormat: 'CSV'
});

// 可视化参数
var visParams = {
  min: 0,
  max: 1,
  palette: ['FFFFFF', '0000FF'] // 白色表示非水体，蓝色表示水体
};

// 逐月显示优化后的水体掩膜，并导出影像
for (var year = startYear; year <= endYear; year++) {
  for (var month = 1; month <= 12; month++) {
    var startDate = ee.Date.fromYMD(year, month, 1);
    var endDate = startDate.advance(1, 'month');

    var monthlyCollection = s1Collection
      .filterDate(startDate, endDate)
      .map(calculateSDWI);

    var sdwiComposite = monthlyCollection.median().clip(roi_x);
    var VH = sdwiComposite.select('VH');

    // 计算水体掩膜
    var waterMask = sdwiComposite.select('SDWI').gt(-0.2)
      .and(sdwiComposite.select('SDWI').lt(0.2))
      .and(VH.lt(-18)) // 低 VH 表示水体
      .and(landMask) // 过滤山体阴影
      .and(jrcMask) // 结合 JRC 长期水体数据
      .rename('WaterMask');

    // 在地图上显示水体掩膜
    var label = 'Refined_WaterMask_' + year + '_' + month;
    Map.addLayer(waterMask, visParams, label);

    // 导出优化后的水体掩膜影像
    Export.image.toDrive({
      image: waterMask,
      description: 'Refined_WaterMask_' + year + '_' + month,
      folder: 'GEE/DTH/images',
      fileNamePrefix: 'Refined_WaterMask_' + year + '_' + month,
      region: roi_x,
      scale: 10,
      maxPixels: 1e13
    });
  }
}

4 运行结果

本文实现了2023年至2024年共24个月份的水体面积逐月变化监测，结果以.csv文件形式导出，便于进一步的数据分析、制图与展示。同时，导出的影像也可用于GIS软件进一步分析与可视化。