2025-0412学习记录14——土地覆盖分类代码解析

最近不是没学习，最近在忙着尽快把分类结果做出来，一天都当成两天过的在努力，毕竟都中期了，结果还没做出来，是有点着急了，论文也要尽快开始写了。

不过最近代码写的挺好，也可以稍微做个记录和展示了。

// 导入 Sentinel-2 影像和云概率数据集
var s2 = ee.ImageCollection('COPERNICUS/S2_SR_HARMONIZED');
var s2_cloud = ee.ImageCollection('COPERNICUS/S2_CLOUD_PROBABILITY');

// ========== 参数设置 ==========
var point = table; // 研究区域
Map.centerObject(point, 8);

var seasonalPeriods = [
  {name: 'winter', start: '2024-02-01', end: '2024-03-31'},
  {name: 'spring', start: '2024-04-01', end: '2024-05-30'},
  {name: 'summer', start: '2024-07-01', end: '2024-08-31'},
  {name: 'autumn', start: '2024-10-01', end: '2024-10-31'}
];

// ========== 云去除和归一化函数 ==========
function rmCloudByProbability(image, threshold) {
    return image.updateMask(image.select('probability').lte(threshold));
}
function scaleImage(image) {
    return image.divide(10000).set('system:time_start', image.get('system:time_start'));
}
function getMergeImages(primary, secondary) {
    return ee.ImageCollection(ee.Join.inner().apply(primary, secondary, 
        ee.Filter.equals({ leftField: 'system:index', rightField: 'system:index' })
    )).map(function(image) {
        return ee.Image(image.get('primary')).addBands(image.get('secondary'));
    });
}
function fillMissingValues(image, radius, iterations) {
    return image.unmask(image.focal_mean({
        kernel: ee.Kernel.square({ radius: radius, units: 'pixels' }),
        iterations: iterations
    }));
}

// ========== DEM 特征 ==========
var dem_10m = ee.ImageCollection("COPERNICUS/DEM/GLO30")
             .filterBounds(qhh_border)
             .limit(10)
             .map(function(image) {
               return image.select('DEM')
                           .clip(qhh_border)
                           .resample('bilinear')
                           .reproject({crs: 'EPSG:4326', scale: 10});
             })
             .mosaic();

Map.addLayer(dem_10m, {min: 0, max: 5000, palette: ['blue', 'green', 'red']}, 'DEM');

var demFeatures = ee.Image.constant(1).rename('constant')
    .addBands([
      dem_10m.lt(3700).rename('is_urban_elev'),
      dem_10m.gte(3800).rename('is_desert_elev'),
      dem_10m.gte(3200).and(dem_10m.lte(4100)).rename('is_beach_elev'),
      dem_10m.gte(3197).rename('is_river_elev')
  ]);

// ========== 构建季节性特征 ==========
function extractSeasonalFeatures(period) {
    // Sentinel-2
    var s2 = ee.ImageCollection('COPERNICUS/S2_SR_HARMONIZED')
        .filterDate(period.start, period.end)
        .filterBounds(point)
        .map(scaleImage);
    var s2_cloud = ee.ImageCollection('COPERNICUS/S2_CLOUD_PROBABILITY')
        .filterDate(period.start, period.end)
        .filterBounds(point);
    var merged = getMergeImages(s2, s2_cloud).map(function(image) {
        return rmCloudByProbability(image, 20);
    });
    var medianS2 = merged.median().clip(point).toFloat();
    var filled = fillMissingValues(medianS2, 2, 3);
    var clipped = filled.clipToCollection(qhh_border);

    // Sentinel-2波段重命名
    var bandNames = ['B2','B3','B4','B5','B6','B7','B8', 'B11'];
    var renamedBands = bandNames.map(function(band) {
        return band + '_' + period.name;
    });
    clipped = clipped.select(bandNames).rename(renamedBands);

    // 指数计算
    var indices = {
        NDVI: clipped.normalizedDifference(['B8_' + period.name, 'B4_' + period.name]).rename('NDVI_' + period.name),
        NDBI: clipped.normalizedDifference(['B11_' + period.name, 'B8_' + period.name]).rename('NDBI_' + period.name),  // 注意：这行用的是原始波段名
        MNDWI: clipped.normalizedDifference(['B3_' + period.name, 'B11_' + period.name]).rename('MNDWI_' + period.name),
        BSI: clipped.normalizedDifference(['B11_' + period.name, 'B4_' + period.name]).rename('BSI_' + period.name)
    };

    // Sentinel-1
    var s1 = ee.ImageCollection('COPERNICUS/S1_GRD')
        .filterDate(period.start, period.end)
        .filterBounds(point)
        .filter(ee.Filter.eq('instrumentMode', 'IW'))
        .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VV'))
        .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VH'))
        .filter(ee.Filter.eq('orbitProperties_pass', 'DESCENDING'))
        .select(['VV', 'VH']);

    var s1_median = s1.median().clip(qhh_border)
        .rename(['VV_' + period.name, 'VH_' + period.name]);

    // 拼接所有波段并返回
    var result = clipped
        .addBands(indices.NDVI)
        .addBands(indices.NDBI)
        .addBands(indices.MNDWI)
        .addBands(indices.BSI)
        .addBands(s1_median);

    return result;
}


// ========== 构建时间序列合成影像 ==========
var seasonalImages = seasonalPeriods.map(function(period) {
  return extractSeasonalFeatures(period);
});

var timeSeriesImage = ee.ImageCollection(seasonalImages).toBands().float();

// 添加 DEM 特征
var image = timeSeriesImage
    .addBands(dem_10m.rename('Elevation'))
    .addBands(demFeatures);

// ========== 样本集合 ==========
var classNames = Temperate_Steppe_11_0
    .merge(Alpine_Steppe_12_1)
    .merge(Alpine_Meadow_13_2)
    .merge(Alpine_Shrubland_20_3)
    .merge(Lake_31_4)
    .merge(River_32_5)
    .merge(Alpine_Desert_41_6)
    .merge(Sand_Dunes_42_7)
    .merge(River_Beach_43_8)
    .merge(Alpine_Wetland_50_9)
    .merge(Arable_Land_60_10)
    .merge(Urban_Land_70_11)
    .merge(Snow_and_Ice_80_12);

// 特征名自动收集
var bands = image.bandNames();

var training = image.select(bands).sampleRegions({
    collection: classNames,
    properties: ['landcover'],
    scale: 10
});

var split = 0.7;
var withRandom = training.randomColumn('random');
var trainingPartition = withRandom.filter(ee.Filter.lt('random', split));
var testingPartition = withRandom.filter(ee.Filter.gte('random', split));

// 分类器训练
var classifier = ee.Classifier.libsvm().train({
    features: trainingPartition,
    classProperty: 'landcover',
    inputProperties: bands
});

var classified = image.select(bands).classify(classifier);

// 精度评估
var test = testingPartition.classify(classifier);
var confusionMatrix = test.errorMatrix('landcover', 'classification');
print('Confusion Matrix:', confusionMatrix);
print('Overall Accuracy:', confusionMatrix.accuracy());
print('Kappa Accuracy:', confusionMatrix.kappa());

// 添加图层显示
Map.addLayer(classified, {
    min: 0, max: 12,
    palette: [
      'b3ca1f', '52f132', '015c14', '20a315', '193bd6',
      '16b7ff', 'ffdfca', 'fff518', '2bffe9', '821299',
      'ffc82d', 'ff0000', 'fcfff3']
}, 'SVM Classified (Time Series)');

print('Training sample count per class:', training.aggregate_histogram('landcover'));

// ========== UI 可视化 ==========

// 定义季节与对应影像字典
var seasonalDict = {
  'Winter': extractSeasonalFeatures(seasonalPeriods[0]),
  'Spring': extractSeasonalFeatures(seasonalPeriods[1]),
  'Summer': extractSeasonalFeatures(seasonalPeriods[2]),
  'Autumn': extractSeasonalFeatures(seasonalPeriods[3])
};

// 可视化参数
var visParams_S2