Burn scar algorithm for Google Earth Engine, derived from Elvidge and Baugh (2014).

instructions.md

Set up for update

1. Go to fires download page - https://firms.modaps.eosdis.nasa.gov/download/
2. Zoom to Sumatra. It doesn't have to be perfect. We do some screening in EE to make sure there are no fires in Malaysia used.
3. Submit the download request (2013-03-30 to present) - csv is easiest - then wait for it to complete (usually < 30 minutes)
4. Upload the CSV file to Fusion Tables. Go here, then just click "new table" under the File menu and follow the instructions.
5. Get the docid from the url. For example, for docid=1SzJl73nj5IPVEOGqhGc8uv5Vkwb504uqK_YTnVGh just grab 1SzJl73nj5IPVEOGqhGc8uv5Vkwb504uqK_YTnVGh
6. Update dates in this script when you run it on EE - the second date in the POST variable may need to be extended.
7. Update the FIRES variable to equal to "ft:". You'll see how it is in the code.

Run + export

0. Paste code into EE playground.
1. Click the "run" button in EE. Then switch to the "Tasks" tab on the right.
2. Start the export process: Export settings: task name: nbr_indonesia_Elvidge mosaic name: NBR Indonesia Elvidge

Postprocessing

1. Check on NBR Indonesia Elvidge in data sources - make sure export completed successfully
2. add to NBR Indonesia Elvidge layer
3. remove existing update from layer
4. rename data source in layer as appropriate, eg. "burn scars Elvidge update 2014-09-22"
5. reprocess layer
6. get asset id, send to Blue Raster so they can update the site as necessary

nbr_elvidge.js

// updated February 2015

// Normalized Burn Ratio, Elvidge, production
// Hammer, Kraft, and Steele (Data Lab at WRI)
// GFW-Fires, prototype

// Reference
// Elvidge, Christopher D. and Kimberly Baugh. 2014. Burn scar mapping from Landsat 8.
// Presentation at APAN meeting in Bandung, Indonesia. January 20.
// URL: http://www.apan.net/meetings/Bandung2014/Sessions/EM/Elvidge_L8_burnscar_20140120.pdf

// Assets required:
//   LC8_L1T_TOA: Landsat 8 Top of Atmosphere, L1T product
//   Recent fires: ft:1UkooS2ir_Rx3Kbg7K2RWV9dTOCCtrjDWizZE_hL0

////////////////////////////////////////////////
// UPDATE THIS WITH NEW FIRES FUSION TABLE ID //
////////////////////////////////////////////////

var FIRES = 'ft:1m0tOYq4dTRhlOiFJxxnj66lS5gm10Z2AMM2q4x0n';

////////////////////////////////////////////////
// UPDATE END DATE
var INIT = ['2014-03-30', '2014-09-30'];
var POST = ['2014-10-01', '2015-01-31'];
//////////////////////////////////////////////// 
// BEGIN FUNCTIONS

var cloudDensity = function(img) {
  // Calculate the cloud score for the supplied L8 image

  var rescale = function(img, exp, thresholds) {
    // A helper to apply an expression and linearly rescale the output.
    return img.expression(exp, {img: img})
        .subtract(thresholds[0]).divide(thresholds[1] - thresholds[0]);
  };

  var score = ee.Image(1.0);
  // Clouds are reasonably bright in the blue band.
  score = score.min(rescale(img, 'img.B2', [0.1, 0.3]));

  // Clouds are reasonably bright in all visible bands.
  score = score.min(rescale(img, 'img.B4 + img.B3 + img.B2', [0.2, 0.8]));

  // Clouds are reasonably bright in all infrared bands.
  score = score.min(rescale(img, 'img.B5 + img.B6 + img.B7', [0.3, 0.8]));

  // Clouds are reasonably cool in temperature.
  score = score.min(rescale(img, 'img.B10', [300, 290]));

  // However, clouds are not snow.
  var ndsi = img.normalizedDifference(['B3', 'B6']);

  return(score.min(rescale(ndsi, 'img', [0.8, 0.6])));
};


function cloudMask(img, cloud_scores, thresh) {
  // Accepts an raw L8 image, a derived image of cloud scores corresponding 
  // to each pixel and a threshold to screen out excessively cloudy pixels.
  // Returns the original image, with clouds masked.
  var binary = cloud_scores.lte(thresh);
  return(img.mask(binary));
}


function makeCloudFree(img) {
  // Convenient wrapper `cloudMask` function that will only accept an L8
  // image and returns the cloud-masked image with the threshold set at 0.5
  var clouds = cloudDensity(img);
  return(cloudMask(img, clouds, 0.5));
}


function cloudScore(img) {
  // Calculates the cloud score for the supplied image, specifically the 
  // proportion of pixels in the image that are excessively cloudy.  If 
  // there are no clouds then the return value is null.
  var pix = cloudDensity(img);
  var binary = pix.gt(0.5);
  var num = binary.reduceRegion(ee.Reducer.mean());
  return(num.getInfo().constant);
}


function hsvpan(rgb, gray) {
  // Accepts the RGB and gray banded image for the same location.  Returns 
  // a pan-sharpened image.
    var huesat = rgb.rgbtohsv().select(['hue', 'saturation']);
    var upres = ee.Image.cat(huesat, gray).hsvtorgb();
    return(upres);
}


function createComposite(start, end) {
  // Accepts a start and end date, along with a bounding GEE polygon.  
  // Returns the L8 cloud composite. 
  var collection = ee.ImageCollection('LC8_L1T_TOA')
    .filterDate(start, end);
  var coll = collection.map(makeCloudFree);
  return(coll.min());
}


function generateBurnRatio(init_season, post_season) {
  // Accepts two tuples that delineate the beginning and end of the 
  // fire season for the previous and current seasons. Returns a binary 
  // image of the burn scars at the set parameter values.

  // create cloud free composites and mask out any remaining clouds 
  // from subsequent analysis
  var a = init_season[0]; var b = init_season[1];
  var c = post_season[0]; var d = post_season[1];
  var init = makeCloudFree(createComposite(a, b));
  var post = makeCloudFree(createComposite(c, d));

  // Add image layer as reference, visibility default set to false
  var rgb = post.select("B6","B5","B4");
  addToMap(rgb, {min:0.01, max:1.5, gamma:1.5}, 'post composite', false);

  // Generate the difference between the normalized burn ratios for 
  // the pre and post periods
  var nbr_init = init.normalizedDifference(['B7', 'B5']);
  var nbr_post = post.normalizedDifference(['B7', 'B5']);
  var nbr = nbr_post.subtract(nbr_init).gt(0.44);

  return(nbr.mask(nbr));
}


function loadFireBuffer(fire_location, meters) {
  // Accepts a fusion table location string and the distance of a buffer in 
  // meters and returns a feature of the unioned buffers.
  var fires = ee.FeatureCollection(fire_location);
  addToMap(fires, {color:'FF9900'}, 'fire points', false)
  var buffered = fires.map(function(f) { return f.buffer(meters); });
  return(buffered.union());
}


function burnRatio() {
  // Accepts the burn ratio image and the fire buffer features and returns
  // the burn scars clipped to the fire buffer extent.  
  var nbr = generateBurnRatio(INIT, POST);
  var fire_buffer = loadFireBuffer(FIRES, 5000);
  var img = nbr.clip(fire_buffer);

  return(img);
}


// Display the normalized burn ratio derived from L8 within 5km of fires.
centerMap(100.8, 1.7, 10);
addToMap(burnRatio(), {palette:'FF0000'}, 'NBR');



var burn = burnRatio();
//var vecs = burn.reduceToVectors({'bestEffort':true});
//addToMap(vecs);

var sumatra = [
          [
            [
              94.54833984375,
              6.206090498573885
            ],
            [
              95.91064453125,
              0.37353251022880474
            ],
            [
              97.93212890625,
              -1.4939713066293112
            ],
            [
              100.30517578125,
              -4.193029605360735
            ],
            [
              101.93115234375,
              -5.943899579425587
            ],
            [
              103.623046875,
              -5.900188795584159
            ],
            [
              104.3701171875,
              -6.380812331938274
            ],
            [
              105.3369140625,
              -6.1187077471908315
            ],
            [
              106.19384765625,
              -5.703447982149503
            ],
            [
              106.25976562499999,
              -3.3160183381615123
            ],
            [
              108.52294921875,
              -3.359889094873377
            ],
            [
              108.47900390625,
              -2.3943223575350774
            ],
            [
              106.72119140625,
              -2.174770572211874
            ],
            [
              106.14990234375,
              -1.1205340322500545
            ],
            [
              105.029296875,
              -1.3841426927920029
            ],
            [
              104.6337890625,
              -0.9008417889908996
            ],
            [
              105.1171875,
              0.13183582116662096
            ],
            [
              104.08447265624999,
              0.5053645409603004
            ],
            [
              103.07373046875,
              1.4500404973607948
            ],
            [
              100.26123046875,
              3.008869978848142
            ],
            [
              97.53662109375,
              5.572249801113911
            ],
            [
              95.82275390625,
              6.075011000682009
            ],
            [
              94.54833984375,
              6.206090498573885
            ]
          ]];

var burn_colored = burn.visualize({'palette':['FF0000'], 'forceRgbOutput':true});
exportImage(burn_colored, 'nbr_indonesia_Elvidge', {'maxPixels':2421986336, 'region':sumatra, 'scale':30});
// Display the normalized burn ratio derived from L8 within 5km of fires.
centerMap(100.8, 1.7, 10);