Monday, 2 October 2017

javascript - Smoothing/interpolating across images in an ImageCollection to remove missing data


I have some code that is working properly to calculate the monthly NDVI averages for the previous two years from a date but I sometimes have missing data in some of the months and I would like to get an average based on the mean of the previous two months. I suppose I could also interpolate from the closest geographic data point but it seems like it would be better to interpolate across time than distance given that sometimes the clouds cover reasonably large areas.


But I'm certainly open to other opinions, I just want to fill these holes as accurately as possible. Is there an efficient way to do this. Here is a link to my working code and here is the code itself:


// pick a landsat tile footprint to use as my geometry
var wrs2_descending = ee.FeatureCollection('ft:1_RZgjlcqixp-L9hyS6NYGqLaKOlnhSC35AB5M5Ll');

// use a manually defined point to pick the WRS2 tile
var wrs2_filtered = wrs2_descending.filterBounds(roi);
var layer1 = ui.Map.Layer(wrs2_filtered, {}, 'WRS2 filtered');
Map.layers().set(1, layer1);

var imageCollection = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR').filterBounds(wrs2_filtered);
var monthCount = ee.List.sequence(0, 11);

// Function to cloud mask from the pixel_qa band of Landsat 8 SR data.
function maskL8sr(image) {

// Bits 3 and 5 are cloud shadow and cloud, respectively.
var cloudShadowBitMask = 1 << 3;
var cloudsBitMask = 1 << 5;

// Get the pixel QA band.
var qa = image.select('pixel_qa');

// Both flags should be set to zero, indicating clear conditions.
var mask = qa.bitwiseAnd(cloudShadowBitMask).eq(0)
.and(qa.bitwiseAnd(cloudsBitMask).eq(0));


// Return the masked image, scaled to reflectance, without the QA bands.
return image.updateMask(mask).divide(10000)
.select("B[0-9]*")
.copyProperties(image, ["system:time_start"]);
}

// run through the image collection and generate monthly NDVI median images
var composites = ee.ImageCollection.fromImages(monthCount.map(function(m) {
var startMonth = 1; // thinking that I should always start from Jan so the series are similar

var startYear = ee.Number(2017-1); // number is one year before the current one

var month = ee.Date.fromYMD(startYear, startMonth, 1).advance(m,'month').get('month');
var year = ee.Date.fromYMD(startYear, startMonth, 1).advance(m,'month').get('year')

// filter by year and then filter by month to get monthly mosaics
var filtered = imageCollection.filter(ee.Filter.calendarRange({
start: year.subtract(1), // we want an average of the last two years
end: year,
field: 'year'

})).filter(ee.Filter.calendarRange({
start: month,
field: 'month'
}));
// mask for clouds and then take the median
var composite = filtered.map(maskL8sr).median();
return composite.normalizedDifference(['B5', 'B4']).rename('NDVI')
.set('month', ee.Date.fromYMD(startYear, startMonth, 1).advance(m,'month'));
}));
print(composites);


Map.addLayer(composites, {min: 0, max: 1}, 'check');

// stack the ImageCollection into a multi-band raster for downloading
var stackCollection = function(collection) {
// Create an initial image.
var first = ee.Image(collection.first()).select([]);

// Write a function that appends a band to an image.
var appendBands = function(image, previous) {

return ee.Image(previous).addBands(image);
};
return ee.Image(collection.iterate(appendBands, first));
};
var stacked = stackCollection(composites);
print('stacked image', stacked);

// Display the first band of the stacked image.
Map.addLayer(stacked.select(0).clip(wrs2_filtered), {min:0, max:1}, 'stacked');

Answer




You can replace the values using where().


// Replace masked pixels by the mean of the previous and next months 
// (otherwise, how to deal with the first images??)
var replacedVals = composites.map(function(image){
var currentDate = ee.Date(image.get('system:time_start'));
var meanImage = composites.filterDate(
currentDate.advance(-2, 'month'), currentDate.advance(2, 'month')).mean();
// replace all masked values:
return meanImage.where(image, image);
});


As an addition to your code, I added the property 'system:time_start', so you can use filterDate(). Furthermore, I made a time window of the two previous and next months to construct a meanImage and replace the values where the image is masked. See the link.


Link code


No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...