项目简介
该项目使用Google Earth Engine (GEE)平台，对Landsat 5和Landsat 7卫星影像进行预处理与影像填补操作。主要功能包括影像的选取、波段处理、缺失影像的填补以及最终影像的导出。代码中的核心功能是通过空间回归方法对Landsat 7和Landsat 5影像进行时序配准与插值填补。
主要步骤
1. 定义感兴趣区域 (ROI)：
- 使用四个坐标点绘制矩形多边形，确定需要处理的地理区域。
2. 选择影像数据集：
- 使用Landsat 5和Landsat 7数据集 (`ee.ImageCollection`)，分别选择时间范围为2010年5月1日至2010年5月31日之间的影像。
3. 影像预处理：
- 定义影像预处理函数，将选择的波段进行裁剪、波段重命名，并转换像元值范围，以便与其他数据进行处理和分析。
4. 影像填补：
- 使用空间回归方法，计算两个影像集之间的关系，通过滑动窗口对像元进行回归分析并填补缺失区域。主要步骤包括：
- 选择匹配的填补影像；
- 计算回归模型，生成比例和偏移量；
- 处理填补失败的情况，使用默认值填补。
5. 显示与导出结果:
- 显示原始Landsat 7、Landsat 5影像和填补后的Landsat 5影像，并将结果导出为影像文件，便于后续分析。
核心代码说明
- 预处理函数 `preprocessImage`:
该函数处理Landsat影像集中的波段并应用裁剪区域。通过乘以系数与增加偏移量将影像值缩放至0-1范围。
- 影像填补函数 `GapFill`:
这是实现影像插值填补的主要函数，输入Landsat 7影像集和Landsat 5影像集，并对它们进行空间回归分析，最终生成无缺失的填补影像。
- 导出函数 `Export.image.toDrive`:
使用GEE的导出功能将填补后的影像导出为文件，以便本地存储和处理。
如何运行
1. 打开Google Earth Engine的代码编辑器，并将上述代码粘贴到编辑器中。
2. 修改`startDate`和`endDate`参数以更改需要处理的影像时间范围。
3. 运行代码，预处理并显示影像，同时将填补后的影像导出为本地文件。
var roi =
/* color: #b50000 */
/* shown: false */
/* displayProperties: [
{
"type": "rectangle"
}
] */
ee.Geometry.Polygon(
[[[117.29699355756334, 40.12436832412577],
[117.29699355756334, 39.890323465395774],
[117.65610915814928, 39.890323465395774],
[117.65610915814928, 40.12436832412577]]], null, false);
// 定义数据集
var landsat7 = ee.ImageCollection('LANDSAT/LE07/C02/T1_L2');
var landsat5 = ee.ImageCollection('LANDSAT/LT05/C02/T1_L2');
var landsat8 = ee.ImageCollection('LANDSAT/LC08/C02/T1_L2');
// 定义波段
var bands = ['red', 'green', 'blue', 'nir', 'swir'];
var percentile = 30;
var imageParams = {min: 0.0, max: 0.3, bands: ['red', 'green', 'blue']};
// 输入时间范围
var startDate = '2010-05-01';
var endDate = '2010-05-31';
// 影像预处理函数
var preprocessImage = function(image) {
var percentileImage = image
.select(['SR_B3', 'SR_B2', 'SR_B1', 'SR_B4', 'SR_B5'])
.multiply(2.75e-05).add(-0.2)
.clip(roi)
.rename(bands); // 重命名波段，使其保持与原始影像一致
return percentileImage
//.updateMask(image.select(0).mask().focal_min(90, 'circle', 'meters'))
.copyProperties(image, ['system:time_start']); // 复制时间等属性;
};
// 选择 Landsat 7 影像
var landsat7Collection = landsat7
//.select(['B3', 'B2', 'B1', 'B4', 'B5'], bands)
.filterDate(startDate, endDate)
.map(preprocessImage);
print('landsat7Collection:',landsat7Collection.limit(10))
// 选择 Landsat 5 影像
var landsat5Collection = landsat5
//.select(['SR_B3', 'SR_B2', 'SR_B1', 'SR_B4', 'SR_B5'], bands)
.filterDate(startDate, endDate)
//.reduce(ee.Reducer.percentile([percentile]))
.map(preprocessImage);
print('landsat5Collection:',landsat5Collection.limit(10))
/*// 选择 Landsat 8 影像
var landsat8Collection = landsat8
.select(['SR_B4', 'SR_B3', 'SR_B2', 'SR_B5', 'SR_B6'], bands)
.filterBounds(roi)
.filterDate(startDate, endDate)
.reduce(ee.Reducer.percentile([percentile]))
.rename(bands)
.map(preprocessImage);
*/
// 填补函数
var GapFill = function(srcCollection, fillCollection, kernelSize) {
var MIN_SCALE = 1/3;
var MAX_SCALE = 3;
var MIN_NEIGHBORS = 144;
var kernel = ee.Kernel.square(kernelSize * 30, "meters", false);
var fillFunction = function(image) {
var date = ee.Date(image.get('system:time_start'));
// 筛选匹配的填补影像
var fillImage = fillCollection
.filterBounds(roi)
.filterDate(date.advance(-10, 'day'), date.advance(10, 'day'))
.median()
.select(bands);
var srcImage = srcCollection
.filterBounds(roi)
.filterDate(date.advance(-10, 'day'), date.advance(10, 'day'))
.median()
.select(bands);
// 计算共同掩膜
var commonMask = srcImage.mask().and(fillImage.mask());
var fillMasked = fillImage.updateMask(commonMask);
var srcMasked = srcImage.updateMask(commonMask);
// 进行回归计算
var regression = fillMasked.addBands(srcMasked);
regression = regression.select(regression.bandNames().sort());
var fit = regression.reduceNeighborhood(ee.Reducer.linearFit().forEach(srcImage.bandNames()), kernel, null, false);
var offset = fit.select(".*_offset");
var scale = fit.select(".*_scale");
// 计算二级比例因子
var reducer = ee.Reducer.mean().combine(ee.Reducer.stdDev(), null, true);
var srcStats = srcImage.reduceNeighborhood(reducer, kernel, null, false);
var fillStats = fillImage.reduceNeighborhood(reducer, kernel, null, false);
var scaleFactor = srcStats.select(".*stdDev").divide(fillStats.select(".*stdDev"));
var offsetFactor = srcStats.select(".*mean").subtract(fillStats.select(".*mean").multiply(scaleFactor));
var invalidScale = scale.lt(MIN_SCALE).or(scale.gt(MAX_SCALE));
scale = scale.where(invalidScale, scaleFactor);
offset = offset.where(invalidScale, offsetFactor);
// 处理所有其他方法都失败的情况
var fallbackScale = scale.lt(MIN_SCALE).or(scale.gt(MAX_SCALE));
scale = scale.where(fallbackScale, 1);
offset = offset.where(fallbackScale, srcStats.select(".*mean").subtract(fillStats.select(".*mean")));
var neighborCount = commonMask.reduceNeighborhood(ee.Reducer.count(), kernel, null, true, "boxcar");
var scaledImage = fillImage.multiply(scale).add(offset).updateMask(neighborCount.gte(MIN_NEIGHBORS));
return image.unmask(scaledImage, true).copyProperties(image, ['system:time_start']);
};
return srcCollection.map(fillFunction);
};
// 应用填补函数
var filledLandsat5 = GapFill(landsat7Collection, landsat5Collection, 10);
print('filledLandsat5:',filledLandsat5.limit(10))
//var filledLandsat8 = GapFill(landsat7Collection, landsat8Collection, 10);
// 组合填补结果
//var resultCollection = ee.ImageCollection.fromImages([filledLandsat5, filledLandsat8]);
//var finalResult = resultCollection.mosaic();
// 显示结果
Map.centerObject(roi, 12); // 调整地图中心和缩放级别
Map.addLayer(landsat7Collection.mean(), imageParams, "Landsat 7");
Map.addLayer(landsat5Collection.mean(), imageParams, "Landsat 5");
Map.addLayer(filledLandsat5.mean(), imageParams, "Landsat 5 Filled");
//Map.addLayer(filledLandsat8, imageParams, "Landsat 8 Filled");
//Map.addLayer(finalResult, imageParams, "Combined Filled");
// 导出处理后的影像
Export.image.toDrive({
image: filledLandsat5.mean(), // 选择Landsat 5填补结果进行导出
description: 'landsat7_filled_landsat5',
folder: 'landsat7',
scale: 30,
});
