/*
This script uses Sentinel-2 imagery to monitor surface water in a specified 
area of interest over time. It calculates the NDWI for each image, applies 
a threshold to identify water pixels, and quantifies the total water-covered 
area for each date. The results are visualized as a time series graph and can be exported, 
along with NDWI and water mask images, to Google Drive.

author: marte.siebinga@acaciawater.com
date: Sept. 2025
*/

// Define your area of interest
// var aoi = ee.FeatureCollection('projects/ee-martesiebinga/assets/rubkona_county');
var aoi = geometry;
var aoi_name = 'test_AOI';

// Zoom to area of interest
Map.centerObject(aoi);
// Map.addLayer(aoi);

// Set dates and load Sentinel-2 data
var startDate = '2019-01-01';
var endDate = '2020-12-31';

// Set the date of the shown and exported layers
var exportDate = '2020-05-20';

// Returns all the images between the start date and the end date
// covering the area of interest and with less than 2% cloud cover.
var sentinelImageCollection = ee.ImageCollection("COPERNICUS/S2_SR_HARMONIZED")
                  .filterDate(startDate, endDate)
                  .filter(ee.Filter.lt("CLOUDY_PIXEL_PERCENTAGE", 2))
                  .filterBounds(aoi);
                  
print("Number of tiles = ", sentinelImageCollection.size());

// function to mosaic tiles by date if AOI is covering multiple tiles
function mosaicByDate(imcol){
  var imlist = imcol.toList(imcol.size())

  var unique_dates = imlist.map(function(im){
    return ee.Image(im).date().format("YYYY-MM-dd")
  }).distinct()

  var mosaic_imlist = unique_dates.map(function(d){
    d = ee.Date(d)

    var im = imcol
      .filterDate(d, d.advance(1, "day"))
      .mosaic()

    return im.set(
        "system:time_start", d.millis(), 
        "system:id", d.format("YYYY-MM-dd"))
  })

  return ee.ImageCollection(mosaic_imlist)
}

var sentinelMosaic = mosaicByDate(sentinelImageCollection); 

print("Number of mosaicked images = ", sentinelMosaic.size())

// Add True Color layer to map
Map.addLayer(
    sentinelMosaic.filterDate(exportDate).first(),
    {min: 0.0, max: 4000, bands: ['B4', 'B3', 'B2']},
    'RGB');


// Calculate the ndwi for each image in the image collection    
var ndwi = sentinelMosaic.map(function(image) {
  return image
    .normalizedDifference(['B3', 'B8'])  // (Green - NIR) / (Green + NIR)
    .rename('NDWI')
    .copyProperties(image, image.propertyNames()); 
});

// Define threshold and create NDWI mask 
var ndwiMasked = ndwi.map(function(image) {
  var waterMask = image
    .gte(0)   // Get pixels greater than or equal: this is the NDWI threshold 
                // which can be adjusted to fit local circumstances 
    .rename('mask');  
                            
  return image
    .addBands(waterMask)
    .updateMask(waterMask)        
    .copyProperties(image, image.propertyNames()); 
});

// Calculate the pixel area covered by water for each image
var waterAreas = ndwiMasked.map(function(image) {
  
  // Get the timestamp of the image and convert the date format
  var timestamp = image.get('system:time_start');
  var date = ee.Date(timestamp).format('YYYY-MM-dd');

  var areaDict = image.select('mask')
                      .multiply(ee.Image.pixelArea())
                      .reduceRegion({
                        reducer: ee.Reducer.sum(),
                        geometry: aoi,
                        scale: 10,
                        maxPixels: 1e10
                      });
  var areaKm2 = ee.Number(areaDict.get('mask')).divide(1e6); // convert m2 to km2
  
  return ee.Feature(null, {
    'date': date,
    'water_area_km2': areaKm2
  });
});

var waterAreaTable = ee.FeatureCollection(waterAreas);

// Plot results in graph
var chart = ui.Chart.feature.byFeature({
  features: waterAreaTable,
  xProperty: 'date',
  yProperties: ['water_area_km2']
  })
.setChartType('LineChart')           
.setOptions({
  title: 'Surface Water Area for ' + aoi_name,
  hAxis: {title: 'Date'},
  vAxis: {title: 'Water Area (km²)'},
  lineWidth: 2,
  pointSize: 4
});

print(chart);

// Add results to map
Map.addLayer(
    ndwi.filterDate(exportDate).first(),
    {palette: ['red', 'yellow', 'green', 'cyan', 'blue']},
    'NDWI');
    
Map.addLayer(ndwiMasked.select('mask').filterDate(exportDate).first(), {palette:['blue']}, 'NDWI mask first');


// Export the first image, specifying the CRS and region.
Export.image.toDrive({
  image: ndwiMasked.select('mask').filterDate(exportDate).first(),
  description: 'NDWI_mask',
  crs: 'EPSG:32636',
  scale: 10,
  region: aoi
});

// Export the first image, specifying the CRS and region.
Export.image.toDrive({
  image: ndwi.select('NDWI').filterDate(exportDate).first(),
  description: 'NDWI',
  crs: 'EPSG:32636',
  scale: 10,
  region: aoi
});


// Export table to CSV
Export.table.toDrive({
    collection: waterAreaTable,
    description: aoi_name,
    fileNamePrefix: aoi_name,
    fileFormat: 'CSV',
    selectors:['date', 'water_area_km2']
});