以下文章来源于GEE学习室 ,作者GEEStudyRoom
光学影像由于受到天气因素(云、雨和雾等)影响,导致单张影像数据存在大量坏死像元。此处,坏死像元指由于受到云遮挡等导致下垫面地物覆盖不能准确被卫星信息捕捉从而不能正常用于实际应用的像素(云识别等研究除外,因为这类研究就是需要有云像素)。坏死像元的存在造成实际应用中数据需求难以得到满足,因此有必要考虑时序影像合成等技术来补充/弥补影像。
时序影像合成,从粗到细,可以分为年合成、半年合成、季度合成、月度合成和半月合成等,在不同学科研究中都有广泛应用,例如半月合成在降雨侵蚀(B因子计算等)中被广泛采用。按照学科和研究需求不同,可以自主选择适合自己研究的合成方式开展研究。
下面就分别以北京主城区为例,使用Landsat-8影像数据集分别实现年合成、半年合成、季度合成、月度合成和半月合成等其中合成方式选择median(大家也可以选择mean、max、min等其他合成方式)。
具体代码:https://code.earthengine.google.com/7e7cb1e0f819ff8bcd5401c59a686a77
var Beijing = ee.Geometry.Polygon([[[116.29730339279362, 39.99599891062294],[116.29730339279362, 39.84222494267003],[116.48681755294987, 39.84222494267003],[116.48681755294987, 39.99599891062294]]], null, false);
var roi = Beijing;
Map.addLayer(roi, {'color':'grey'}, 'studyArea');
Map.centerObject(roi);var year = 2020;// reomove cloud for Landsat-8
function rmL8Cloud(image) { var cloudShadowBitMask = (1 << 3); var cloudsBitMask = (1 << 5); var qa = image.select('pixel_qa'); var mask = qa.bitwiseAnd(cloudShadowBitMask).eq(0) .and(qa.bitwiseAnd(cloudsBitMask).eq(0));return image.updateMask(mask).toFloat().divide(10000).copyProperties(image).copyProperties(image, ['system:time_start','system:time_end','system:footprint']);
}var l8 = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR').filterBounds(roi).filterDate(ee.Date.fromYMD(year,1,1),ee.Date.fromYMD(year+1,1,1)).map(rmL8Cloud).select(['B2','B3','B4','B5','B6']);// year composite
var l8_year = l8.median().clip(roi);
print("l8_year",l8_year);// half-year composite
var l8_halfyear = ee.List.sequence(0, 1*11,6).map(function(n) {var start = ee.Date.fromYMD(year,1,1).advance(n, 'month');var end = start.advance(6, 'month');var tmpMedian = l8.filterDate(start, end).median().clip(roi).set("system:time_start", start.millis());return tmpMedian;
}).flatten();
var l8_halfyear = ee.ImageCollection.fromImages(l8_halfyear);
print("l8_halfyear",l8_halfyear);// season composite
var l8_season = ee.List.sequence(0, 1*11,3).map(function(n) {var start = ee.Date.fromYMD(year,1,1).advance(n, 'month');var end = start.advance(3, 'month');var tmpMedian = l8.filterDate(start, end).median().clip(roi).set("system:time_start", start.millis());return tmpMedian;
}).flatten();
var l8_season = ee.ImageCollection.fromImages(l8_season);
print("l8_season",l8_season);// month composite
var l8_month = ee.List.sequence(0, 1*11,1).map(function(n) {var start = ee.Date.fromYMD(year,1,1).advance(n, 'month');var end = start.advance(1, 'month');var tmpMedian = l8.filterDate(start, end).median().clip(roi).set("system:time_start", start.millis());return tmpMedian;
}).flatten();
var l8_month = ee.ImageCollection.fromImages(l8_month);
print("l8_month",l8_month);// l8_halfmonth
var l8_halfmonth = ee.List.sequence(1,12).map(function(item){var start_date = ee.Date.fromYMD(year, item, 1);var mid_date = start_date.advance(15,'day');var end_date = start_date.advance(1,'month');var pre_month_img = l8.filterDate(start_date,mid_date).mean().set('system:time_start',start_date.millis()).set('time',start_date.format('YYYY-MM-dd'));var pro_month_img = l8.filterDate(mid_date,end_date).mean().set('system:time_start',start_date.millis()).set('time',mid_date.format('YYYY-MM-dd')); return ee.ImageCollection([pre_month_img,pro_month_img]);
});// define iterate function
function accumulate(imageCol, list) {return ee.ImageCollection(list).merge(imageCol);
}
// flatten each element in the l8_halfmonth collection
var l8_halfmonth = l8_halfmonth.iterate(accumulate,ee.ImageCollection([]));
var l8_halfmonth= ee.ImageCollection(l8_halfmonth);
print("l8_halfmonth",l8_halfmonth);// show the image
var visParam = {bands:['B6','B5','B4'],min:0,max:0.3};
Map.addLayer(l8_year, visParam, 'l8_year');
代码解析
年合成相对来说简单,只需要使用mean函数就可以。具体代码见31-32行。
2、半年合成比年合成相对复杂一些,但是也不复杂,主要是增加List和Date类别里面有关函数的使用,其中List涉及到sequence使用,主要用于生成列表序列用于表征2个月。具体代码见34-42行。
3、季度合成和半年合成类似,只需要修改月度List里面的生成间隔即可,即将ee.List.sequence(0, 111,6)变成ee.List.sequence(0, 111,3)(Note,map函数里面也需要对应修改)。当然,这里的按季合成只是简单地将1-3月当成了一个季度,更合理的方式可能是按照“春-夏-秋-冬”对应月份合成。具体代码见44-52行。
4、月度合成参考前面的半年合成和季度合成,只是再次修改ee.List.sequence参数就行。具体代码见54-62行。
5、半月合成和之前的合成方式有些许不同,在于一次循环需要得到2张影像。即一个月需要分成上旬和下旬来计算,这样返回值就需要返回两张影像(但是只能有一个返回值)。所以一开始得到的是12个对象,每个对象含有两张影像数据。为了顺利得到24张影像数据,这里还需要巧妙地借助iterate函数,将影像数据都抽出来组成一个24个影像,具体代码见64-85行。这块是重点内容,还请您多多留意~
具体结果图如下所示。