如果是学遥感的同学,应该都知道Landsat7的影像在2003年出现了问题,导致后期的数据都有条带。
最近我个人在使用Landsat7的数据,为了填补条带,就做了以下代码,大家如果有用得上的,就拿去用。选择好自己的研究区域,5秒钟内填补好landsat7的条带。该该方法使用的USGS提供的去条带原理,下面是5秒钟改善Landsat7条带步骤:
1.选择研究区;
2.选择需要填补的影像;
3.引入填补函数;
4.landsat填补,并下载 【GOOGLE|5秒钟用google earth engine(GEE)填补Landsat7条带问题】这个是代码,可以参考一下本人写的
https://code.earthengine.google.com/1da3330c6289b0dc75b0bb623c151da9
1.选择研究区: 研究区的选择,可以自己上传矢量文件,也可以像我这样直接在华北平原点上一个点。
文章图片
2.选择需要填补的影像: 我选择的是北京市2004年的landsat影像,作为需要填补的对象,只选择1-5波段并且去掉了landsat7的边缘缝隙,下面是具体代码。
//选择landsat7的影像
var image = landsat7
.select(['B3','B2','B1','B4','B5'], bands)
.filterBounds(roi)
.filterDate('2004-06-01', '2004-08-01')
.reduce(ee.Reducer.percentile([percentile]))
.rename(bands)
image = image.multiply(0.0001)
Map.addLayer(image, imageParams, "landsat7");
// 移除landsat7边缘的瑕疵
image = image.updateMask(image.select(0).mask().focal_min(90, 'circle', 'meters'))
有了需要填补的对象,我们接着将同时期的landsat5加载进来,用landsat5来填补landsat7的条带。
//使用landsat5填补条带
var fill = landsat5
.select(['B3','B2','B1','B4','B5'], bands)
.filterBounds(roi)
.filterDate('2004-06-01', '2010-08-01')
.reduce(ee.Reducer.percentile([percentile]))
.rename(bands)
fill = fill.multiply(0.0001)
Map.addLayer(fill, imageParams, "landsat5");
3.引入填补函数: 这个函数是usgs提供的,我们直接改写一下,放到GEE里面。
/* USGS提供的Landsat填补方法,可以自己下载下来看看
USGSLandsat填补方法链接:
https://landsat.usgs.gov/sites/default/files/documents/L7SLCGapFilledMethod.pdf */
//该函数进行landsat7条带的填补
var GapFill = function(src, fill, kernelSize) {//填补参数
var MIN_SCALE = 1;
var MAX_SCALE = 3;
var MIN_NEIGHBORS = 144;
var kernel = ee.Kernel.square(kernelSize * 30, "meters", false)//找到两张影像的相同的像素
var common = src.mask().and(fill.mask())
var fc = fill.updateMask(common)
var sc = src.updateMask(common)// 用回归法找出主要的比例因子。交错波段进行回归。这假设这些带具有相同的名称//
var regress = fc.addBands(sc)
regress = regress.select(regress.bandNames().sort())
var fit = regress.reduceNeighborhood(ee.Reducer.linearFit().forEach(src.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 src_stats = src.reduceNeighborhood(reducer, kernel, null, false)
var fill_stats = fill.reduceNeighborhood(reducer, kernel, null, false)
var scale2 = src_stats.select(".*stdDev").divide(fill_stats.select(".*stdDev"))
var offset2 = src_stats.select(".*mean").subtract(fill_stats.select(".*mean").multiply(scale2))
var invalid = scale.lt(MIN_SCALE).or(scale.gt(MAX_SCALE))
scale = scale.where(invalid, scale2)
offset = offset.where(invalid, offset2)// 当所有其他方法都失败时的处理方法
var invalid2 = scale.lt(MIN_SCALE).or(scale.gt(MAX_SCALE))
scale = scale.where(invalid2, 1)
offset = offset.where(invalid2, src_stats.select(".*mean").subtract(fill_stats.select(".*mean")))
var count = common.reduceNeighborhood(ee.Reducer.count(), kernel, null, true, "boxcar")
var scaled = fill.multiply(scale).add(offset)
.updateMask(count.gte(MIN_NEIGHBORS))return src.unmask(scaled, true)
}
4.landsat填补,并下载
填补就比较简单了,输入需要填补的图像,完成填补。并用export函数将填补好的的图像进行下载。
var result = GapFill(image, fill, 10);
//加载填补好的landsat7影像
Map.addLayer(result, imageParams, "landsat7 filled");
//导出已经处理好的landsat7影像
Export.image.toDrive({
image: result,
description: 'landsat7_result',
folder:'landsat7',
scale: 30,});
下面是填补之前与没有填补的对比:
文章图片
文章图片
可以看出,条带还是去掉了大部分的。有需要的同学可以直接运行我的代码试试效果
https://code.earthengine.google.com/1da3330c6289b0dc75b0bb623c151da9