GEE处理哨兵数据并导出

博客提供了一个Google Earth Engine的代码链接https://code.earthengine.google.com/dc35e6894b8d75711ca0ec449edc0f73 ,可能与数据库、前端及JavaScript相关。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

https://code.earthengine.google.com/dc35e6894b8d75711ca0ec449edc0f73


// 定义研究区域和样本数据
var roi = ee.FeatureCollection('users/messi513112/xiamen2');
var samples = ee.FeatureCollection('users/messi513112/Expor_xiament_point'); // 替换为您的shapefile文件的实际路径

//导入sentinel-1和sentinel-2的图像集合
var s1 = ee.ImageCollection('COPERNICUS/S1_GRD');
var s2 = ee.ImageCollection('COPERNICUS/S2_SR');
// 获取NASA数据
var nasaDEM = ee.Image('NASA/NASADEM_HGT/001').clip(roi);


//重采样
var scale = 10;
// 重采样Sentinel-1影像集
var resampledS1 = s1.map(function(image) {
  return image.resample('bilinear').reproject({
    crs: image.select(0).projection().crs(),
    scale: scale
  });
});

// 重采样Sentinel-2影像集
var resampledS2 = s2.map(function(image) {
  return image.resample('bilinear').reproject({
    crs: image.select(0).projection().crs(),
    scale: scale
  });
});

// 重采样NASA DEM数据
var resampledDEM = nasaDEM.resample('bilinear').reproject({
  crs: nasaDEM.projection().crs(),
  scale: scale
});





// 定义时间范围
var start = ee.Date('2018-01-01');
var end = ee.Date('2018-12-31');

// 计算属性 'YSSZ_S' 中每个值的数量。
var valueCounts = samples.aggregate_histogram('YSSZ_S');
// 打印结果。
print(valueCounts);

// 使用简单的云和云影掩膜
function maskS2clouds(image) {
  var qa = image.select('QA60');

  // Bits 10和11是云和阴影的标识符
  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);
}


//纹理特征
var vh = resampledS1
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VV'))
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VH'))
  .filter(ee.Filter.eq('instrumentMode', 'IW'))
 .map(function(image) {
          var edge = image.lt(-30.0);
          var maskedImage = image.mask().and(edge.not());
          return image.updateMask(maskedImage).clip(roi);
        });
// 筛选不同观测角度的图像
var vhAscending = vh.filter(ee.Filter.eq('orbitProperties_pass', 'ASCENDING'));
var vhDescending = vh.filter(ee.Filter.eq('orbitProperties_pass', 'DESCENDING'));
var spring = ee.Filter.date('2018-01-01', '2018-12-31');
// 从不同极化和观测角度的平均值创建合成图像
var composite = ee.Image.cat([
  ee.ImageCollection(vhAscending.select('VH').merge(vhDescending.select('VH'))).filter(spring).mean(),
  ee.ImageCollection(vhAscending.select('VV').merge(vhDescending.select('VV'))).filter(spring).mean(),
]).focal_median();
print('composite',composite);
// 作为极化和后向散射特性的合成显示
Map.addLayer(composite, {min: -20, max: 0}, 'composite');
// 计算纹理特征
var glcm = composite.multiply(100).toInt().glcmTexture({
    'size': 1,
    'kernel': null
});

// 打印纹理特征
print('GLCM:', glcm);

glcm =glcm.select(['VH_asm',
                   'VH_contrast',
                   'VH_corr',
                   'VH_var',
                   'VH_idm',
                   'VH_savg',
                   'VH_sent',
                   'VH_ent',
                   'VH_dvar',
                   'VH_dent',
                   'VH_imcorr1',
                   'VH_imcorr2',
                   'VH_diss',
                   'VH_inertia',
                   'VH_shade',
                   'VH_prom',
                   //'VV_asm',
                   //'VV_contrast',
                  // 'VV_corr',
                   //'VV_var',
                   //'VV_idm',
                   //'VV_savg',
                   //'VV_sent',
                   //'VV_ent',
                   //'VV_dvar',
                   //'VV_dent',
                   //'VV_imcorr1',
                   //'VV_imcorr2',
                   //'VV_diss',
                   //'VV_inertia',
                   //'VV_shade',
                   //'VV_prom',
                    ]);



// 筛选Sentinel-2数据
var s2filtered = resampledS2
  .filterBounds(roi)
  .filterDate(start, end)
  .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 40));
// 应用去云函数到筛选后的影像集
var s2CloudRemoved = s2filtered.map(maskS2clouds);
// 打印筛选后的Sentinel-2数据集的元数据,以检查数据集
print('Sentinel-2 filtered metadata:', s2CloudRemoved);
//多时相合成
var s2composite = s2CloudRemoved
  .median() // 多时相合成,这里使用中位数
  .clip(roi); // 裁剪到研究区域
// 显示Sentinel-2多时相合成的结果
Map.addLayer(s2composite, {min: 0, max: 3000, bands: ['B4', 'B3', 'B2']}, 'Sentinel-2 Composite', false);
//提取Sentinel-2的光谱特征
var spectral = s2composite.select(['B4', 'B5', 'B6', 'B7', 'B8']);



// 定义一个函数来计算 NDVI
function addNDVI(image) {
  var ndvi = image.normalizedDifference(['B8', 'B4']).rename('NDVI');
  return ndvi;
}

// 定义一个函数来计算 EVI
function addEVI(image) {
  var evi = image.expression(
    '2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))', {
      'NIR': image.select('B8'),
      'RED': image.select('B4'),
      'BLUE': image.select('B2')
    }).rename('EVI');
  return evi;
}

// 定义一个函数来计算 NDWI
function addNDWI(image) {
  var ndwi = image.normalizedDifference(['B3', 'B8']).rename('NDWI');
  return ndwi;
}

// 定义一个函数来计算 NDRE
function addNDRE(image) {
  var ndre = image.normalizedDifference(['B8', 'B5']).rename('NDRE');
  return ndre;
}

// 定义一个函数来计算 SAVI
function addSAVI(image) {
  var savi = image.expression(
    '((NIR - RED) / (NIR + RED + L)) * (1 + L)', {
      'NIR': image.select('B8'),
      'RED': image.select('B4'),
      'L': 0.5
    }).rename('SAVI');
  return savi;
}

// 定义一个函数来计算 MSAVI
function addMSAVI(image) {
  var msavi = image.expression(
    '(2 * NIR + 1 - sqrt((2 * NIR + 1)**2 - 8 * (NIR - RED))) / 2', {
      'NIR': image.select('B8'),
      'RED': image.select('B4')
    }).rename('MSAVI');
  return msavi;
}

// 定义一个函数来计算 ARVI
function addARVI(image) {
  var arvi = image.expression(
    '(NIR - (2 * RED - BLUE)) / (NIR + (2 * RED - BLUE))', {
      'NIR': image.select('B8'),
      'RED': image.select('B4'),
      'BLUE': image.select('B2')
    }).rename('ARVI');
  return arvi;
}




// 应用函数并显示结果
var NDVI = addNDVI(s2composite);
var EVI = addEVI(s2composite);
var NDWI = addNDWI(s2composite);
var NDRE = addNDRE(s2composite);
var SAVI = addSAVI(s2composite);
var MSAVI = addMSAVI(s2composite);
var ARVI = addARVI(s2composite);

// 生长季开始(SOS)
var sos = NDVI.reduce(ee.Reducer.firstNonNull());

// 生长季结束(EOS)
var eos = NDVI.reduce(ee.Reducer.lastNonNull());

// 初花期
var startOfFlowering = NDVI.reduce(ee.Reducer.first().setOutputs(['StartOfFlowering']));

// 叶变色期
var endOfSenescence = NDVI.reduce(ee.Reducer.last().setOutputs(['EndOfSenescence']));

// 准备Sentinel-2数据
var ss2 = ee.ImageCollection('COPERNICUS/S2')
          .filterDate('2018-01-01', '2018-12-31') // 根据设置的时间范围过滤影像
          .filterBounds(roi) // 根据定义的区域过滤影像
          .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 40)) // 过滤云量小于40%的影像
          .map(function(image) {
            // 计算NDVI并将其作为一个波段添加
            var ndvi = image.normalizedDifference(['B8', 'B4']).rename('NDVI');
            // 将system:time_start属性转换为图像并重命名为'time'
            var time = ee.Image.constant(image.get('system:time_start')).toFloat().rename('time');
            // 返回包含NDVI和时间波段的图像
            return ndvi.addBands(time);
          });

// 对ImageCollection应用linearFit reducer
var linearFit = ss2.select(['NDVI', 'time']).reduce(ee.Reducer.linearFit());



// 计算地形特征
var terrain = ee.Algorithms.Terrain(resampledDEM);

// 提取地形特征波段
var elevation = terrain.select('elevation'); // 海拔
var slope = terrain.select('slope'); // 坡度
var aspect = terrain.select('aspect'); // 坡向
var hillshade = terrain.select('hillshade'); // 山体阴影



var img1 = ee.ImageCollection('COPERNICUS/S2')
    .filterDate('2018-2-13', '2018-2-14')
    .filterBounds(roi)
    .median();
  
var img6 = ee.ImageCollection('COPERNICUS/S2')
    .filterDate('2018-9-26', '2018-9-27')
    .filterBounds(roi)
    .median();

// 计算SVI
var svi = img1.addBands(img6)
            .rename(['B4_1', 'B8_1', 'B4_6', 'B8_6'])
            .expression(
              '((1/6) * timeDiff * (LB8 - LB4) * (B4_1 + 2 * B8_1 + 2 * B4_6 + B8_6))/100000', {
                'timeDiff': 225,
                //'centralWavelengths': centralWavelengths,
                'B4_1': img1.select('B4'),
                'B8_1': img1.select('B8'),
                'B4_6': img6.select('B4'),
                'B8_6': img6.select('B8'),
                'LB8':835.1,
                'LB4':664.5
              }
            ).rename('SVI');
            
// 计算B3与B8之间的SVI
var sviB3B8 = img1.addBands(img6)
  .rename(['B3_1', 'B8_1', 'B3_6', 'B8_6'])
  .expression(
    '(((1/6) * timeDiff * (835.1 - 560.0) * (B3_1 + 2 * B8_1 + 2 * B3_6 + B8_6))/10000)-6800', {
      'timeDiff': 225,
      'B3_1': img1.select('B3'),
      'B8_1': img1.select('B8'),
      'B3_6': img6.select('B3'),
      'B8_6': img6.select('B8')
    }
  ).rename('SVI_B3_B8');

// 计算B5与B8之间的SVI
var sviB5B8 = img1.addBands(img6)
  .rename(['B5_1', 'B8_1', 'B5_6', 'B8_6'])
  .expression(
    '((1/6) * timeDiff * (835.1 - 703.9) * (B5_1 + 2 * B8_1 + 2 * B5_6 + B8_6))/10000-2700', {
      'timeDiff': 225,
      'B5_1': img1.select('B5'),
      'B8_1': img1.select('B8'),
      'B5_6': img6.select('B5'),
      'B8_6': img6.select('B8')
    }
  ).rename('SVI_B5_B8');

// 计算B6与B8之间的SVI
var sviB6B8 = img1.addBands(img6)
  .rename(['B6_1', 'B8_1', 'B6_6', 'B8_6'])
  .expression(
    '((1/6) * timeDiff * (835.1 - 740.2) * (B6_1 + 2 * B8_1 + 2 * B6_6 + B8_6))/10000-2700', {
      'timeDiff': 225,
      'B6_1': img1.select('B6'),
      'B8_1': img1.select('B8'),
      'B6_6': img6.select('B6'),
      'B8_6': img6.select('B8')
    }
  ).rename('SVI_B6_B8');

// 计算B7与B8之间的SVI
var sviB7B8 = img1.addBands(img6)
  .rename(['B7_1', 'B8_1', 'B7_6', 'B8_6'])
  .expression(
    '((1/6) * timeDiff * (835.1 - 782.5) * (B7_1 + 2 * B8_1 + 2 * B7_6 + B8_6))/10000-1000', {
      'timeDiff': 225,
      'B7_1': img1.select('B7'),
      'B8_1': img1.select('B8'),
      'B7_6': img6.select('B7'),
      'B8_6': img6.select('B8')
    }
  ).rename('SVI_B7_B8');

// 计算B8A与B8之间的SVI
var sviB8AB8 = img1.addBands(img6)
  .rename(['B8A_1', 'B8_1', 'B8A_6', 'B8_6'])
  .expression(
    '((1/6) * timeDiff * (864.8 - 835.1) * (B8A_1 + 2 * B8_1 + 2 * B8A_6 + B8_6))/10000', {
      'timeDiff': 225,
      'B8A_1': img1.select('B8A'),
      'B8_1': img1.select('B8'),
      'B8A_6': img6.select('B8A'),
      'B8_6': img6.select('B8')
    }
  ).rename('SVI_B8A_B8');

// 特征融合 想导出的特征
var features = glcm//.addBands(spectral)
                   .addBands(NDVI)
                   .addBands(EVI)
                   //.addBands(NDWI)
                   //.addBands(NDRE)
                   //.addBands(SAVI)
                   //.addBands(MSAVI)
                   //.addBands(ARVI)
                   //.addBands(sos)
                   //.addBands(eos)
                   //.addBands(startOfFlowering)
                   //.addBands(endOfSenescence)
                   //.addBands(linearFit)
                   //.addBands(elevation)
                   //.addBands(slope)
                   //.addBands(aspect)
                   //.addBands(hillshade)
                   .addBands(svi)
                   //.addBands(sviB3B8)
                   //.addBands(sviB5B8)
                   //.addBands(sviB6B8)
                   //.addBands(sviB7B8)
                   //.addBands(sviB8AB8);
                   
                   
print(features);

// 假设'image'是您要导出的影像
var imageToExport = features.toFloat(); // 将所有波段转换为Float32类型

// 定义导出参数
var exportParams = {
  image: imageToExport,
  description: 'export_image',
  scale: 10,
  region: roi, // 确保已定义导出区域的几何形状
  maxPixels: 1e13
};

// 开始导出图像
Export.image.toDrive(exportParams);


// 定义您想要保留的样本数量。
var numSamples =200;
// 从特征集合中随机选择样本。
var randomSamples = samples.randomColumn('random').sort('random').limit(numSamples);
print(randomSamples);

// 训练数据准备
var training = features.sampleRegions({
  collection: randomSamples,
  properties: ['YSSZ_S'],
  scale: 10,
  tileScale:16
})


// 打印训练数据集,以检查样本和特征
print('Training data:', training);
// 随机化样本数据集并添加随机列
var samplesWithRandom = training.randomColumn('random');
// 设置训练和验证样本的比例
var split = 0.5; 
var trainingPartition = samplesWithRandom.filter(ee.Filter.lt('random', split));
var testingPartition = samplesWithRandom.filter(ee.Filter.gte('random', split));
// 训练随机森林分类器
var classifier = ee.Classifier.smileRandomForest(50).train({
  features: trainingPartition,
  classProperty: 'YSSZ_S',
  inputProperties: features.bandNames()
});
print(classifier.explain()); 
// 对整个区域进行分类
var classified = features.classify(classifier);

// 获取分类结果的特征集
var classifiedFeatures = ee.FeatureCollection(classified);

// 使用验证样本集评估模型性能
var test = testingPartition.classify(classifier);
// 计算混淆矩阵
var confusionMatrix = test.errorMatrix('YSSZ_S', 'classification');
// 打印混淆矩阵和精度评估指标
print('Confusion Matrix:', confusionMatrix);
print('Overall Accuracy:', confusionMatrix.accuracy());
print('Kappa Accuracy:', confusionMatrix.kappa());

// 定义一个颜色列表,每个类别一个颜色
var palette = [
  'ffd966', // 浅黄色,代表针叶林
  '8db360', // 橄榄绿,代表落叶阔叶林
  '006400', // 深绿色,代表常绿阔叶林
  'ffbb22', // 金色,代表混交林
  'f13c20', // 天蓝色,代表红树林
  '92cddc'  // 地红色,代表灌木丛
];
// 添加分类结果的图层到地图上
Map.addLayer(classified, {min: 1, max: 8, palette: palette}, 'Classification Result');

// 导出分类结果
Export.image.toDrive({
  image: classified,
  description: 'xiamen_classification',
  folder: 'xiamen_classification',
  scale: 10,
  region: roi,
  maxPixels: 1e13
});

### 使用 Google Earth Engine API 批量下载 Sentinel 卫星数据 #### 准备工作 在开始之前,确保已经安装配置好了 `earthengine-api` 和 `geemap` 库。这些库提供了访问 GEE 平台以及管理任务的功能。 对于 Python 用户来说,可以通过 pip 安装上述依赖项: ```bash pip install earthengine-api geemap ``` 初始化 Earth Engine: ```python import ee ee.Initialize() ``` #### 构建影像集合查询条件 定义感兴趣区域 (AOI),设置时间范围和其他筛选参数来构建特定的影像集合作为后续操作的基础[^1]。 ```python # 设置研究区 aoi = ee.Geometry.Polygon( [[[70.89, 26.54], [70.89, 27.30], [72.27, 27.30], [72.27, 26.54]]]) # 创建Sentinel-1 影像集合 sentinel1_collection = ( ee.ImageCollection('COPERNICUS/S1_GRD') .filterBounds(aoi) .filterDate('2023-01-01', '2023-12-31') # 时间区间可以根据需求调整 ) print(f'Found {sentinel1_collection.size().getInfo()} images.') ``` #### 下载单张图片作为测试 先尝试导出一张图片到 Google Drive 或本地文件系统验证流程是否正常工作[^2]。 ```python from geemap import ee_export_image image_to_download = sentinel1_collection.first() task_config = { 'region': aoi.getInfo()['coordinates'], 'scale': 10, 'folder': 'gee_downloads', } ee_export_image(image=image_to_download, filename='test_sentinel.tif', **task_config) ``` #### 实现批量下载功能 当确认单独图像可以成功下载之后,则可通过循环遍历整个影像集合完成批量化作业。这里展示了一个简单的例子,其中包含了基本错误处理逻辑以应对可能出现的问题。 ```python def batch_download(collection, folder_name): """Batch download all images from an image collection.""" try: size = collection.size().getInfo() if not size: print("No images found.") return for i in range(size): img = ee.Image(collection.toList(100).get(i)) task_params = dict( region=aoi.getInfo()["coordinates"], scale=10, maxPixels=1e13, fileFormat="GeoTIFF", description=f'sentinel_{i}', folder=folder_name ) export_task = ee.batch.Export.image.toDrive(**{ "image": img.select(['VV', 'VH']), **task_params }) export_task.start() print(f'Started exporting image {i}') except Exception as e: print(e) batch_download(sentinel1_collection, 'my_folder') ``` 此脚本会将符合条件的所有 Sentinel-1 GRD 类型的数据逐个上传至指定目录下保存成 GeoTiff 文件格式。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

燕南路GISer

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值