app_func

// released under Open Source GPL License Copyright © 2023 Ka Hei Chow

// This script stores functionalities of all available app

// import module
var palettes = require('users/gena/packages:palettes');
var helper = require('users/pinkychow1010/WB_IntraUrban:helper');
var analysis = require('users/pinkychow1010/WB_IntraUrban:analysis_utils');
var basemap = require('users/pinkychow1010/WB_IntraUrban:basemap_resources');
var intraUrban = require('users/pinkychow1010/WB_IntraUrban:app_interface');

// app functions

exports.set_up = function set_up() {
  // config base map
  basemap.addCustomBasemap('MutedMonotone');
  Map.setCenter(70, 30, 4.5);
  
  // turn off drawing tool bar
  Map.setControlVisibility({drawingToolsControl:false});
};

exports.pop_count = function pop_app() {
  var geometry = ee.FeatureCollection("projects/ee-pinkychow1010/assets/WB_GEE/karachi_union_council_subset");
  analysis.zonal_pop_sum(geometry);
};

exports.fvc_lst_app = function fvc_lst_app() {
  intraUrban.app(fvc_lst_impact, "Land Surface Temperature based on Vegetation Change");
};

exports.landuse_ratio_app = function landuse_ratio_app() {
  intraUrban.app(calculate_landuse_ratio, "Land Use Cover Ratio (2020 ESA World Cover)");
};

//##

exports.map_lst_app = function landuse_lst_app() {
  intraUrban.app(map_lst, "Land surface Temperature (2022 Summer MODIS-based)");
};

function map_lst(aoi) {
  
  var bbox = aoi.geometry();
  var centroid = bbox.centroid({'maxError': 1});
  var listCoords = ee.Array.cat(centroid.coordinates(), 0);
  var yCoords = ee.List(listCoords).getInfo()[1];
  
  // south/ north
  var south = ee.Number(yCoords).lt(0);
  var shift = south.multiply(6).getInfo();
  
  // evaluation year
  var start_year = ee.Date("2022-01-01");
  var end_year = ee.Date("2022-01-01");
  
  // filter dataset to evaluation year
  var yearFilter = ee.Filter.date(
    start_year.advance(-1, "year"), 
    end_year.advance(1, "year")
    );
  
  // get summer months depends on location
  var summerFilter = ee.Filter.calendarRange(5+shift, 9-shift, "month");
  
  var dataset = ee.ImageCollection('MODIS/061/MOD11A1')
                  .filter(yearFilter)
                  .filter(summerFilter)
                  .mean().multiply(0.02).add(-273.15)
                  .clip(bbox);
  
  var dayTemp = dataset.select('LST_Day_1km');
  var nightTemp = dataset.select('LST_Night_1km');
  var tempDiff = dayTemp.subtract(nightTemp);
  
  var palette = ['blue','white','yellow','orange','red'];
  var vis = {min: 20, max:40, palette: palette};
  Map.addLayer(dayTemp, vis, "Day Time Temperature");
  Map.addLayer(nightTemp, {min: 10, max:30, palette: palette}, "Night Time Temperature");
  Map.addLayer(tempDiff, {min: 0, max:15, palette: palette}, "Daily Temperature Difference");
  
}


//##

exports.landuse_lst_app = function landuse_lst_app() {
  intraUrban.app(calculate_landuse_lst, "Temperature by Land Use Cover (2022 Summer)");
};

function calculate_landuse_lst(aoi) {
  
  var bbox = aoi.geometry();
  var centroid = bbox.centroid({'maxError': 1});
  var listCoords = ee.Array.cat(centroid.coordinates(), 0);
  var yCoords = ee.List(listCoords).getInfo()[1];
  
  // south/ north
  var south = ee.Number(yCoords).lt(0);
  var shift = south.multiply(6).getInfo();
  
  // evaluation year
  var start_year = ee.Date("2022-01-01");
  var end_year = ee.Date("2022-01-01");
  
  // filter dataset to evaluation year
  var yearFilter = ee.Filter.date(
    start_year.advance(-1, "year"), 
    end_year.advance(1, "year")
    );
  
  // get summer months depends on location
  var summerFilter = ee.Filter.calendarRange(5+shift, 9-shift, "month");
  
  var dataset = ee.ImageCollection('MODIS/061/MOD11A1')
                  .filter(yearFilter)
                  .filter(summerFilter)
                  .mean().multiply(0.02).add(-273.15)
                  .clip(bbox);
  
  var dayTemp = dataset.select('LST_Day_1km');
  
  var vis = {min: 20, max:40, palette: ['blue','white','yellow','orange','red']};
  Map.addLayer(dayTemp, vis, "MODIS-based Land Surface Temperature");
  
  var chart = reduce_by_lulc(dayTemp, bbox, "LST by LULC");
  
  var lulc = ee.ImageCollection('ESA/WorldCover/v100').first().clip(bbox);
  
  var visualization = {
    bands: ['Map'],
  };
  
  Map.addLayer(lulc, visualization, 'ESA Landcover (2020)');
  
  return chart;
}


   
    
function calc_ls_lst(landsat, bbox) {
    // #################################
    // # calculating LST from landsat ##
    // #################################
    // https://gis.stackexchange.com/questions/314374/calculating-lst-from-landsat-8-in-google-earth-engine

    
    // ndvi    
    var ndvi = landsat.normalizedDifference(['B5', 'B4']).rename('NDVI');

    //select thermal band 10(with brightness tempereature), no calculation 
    var thermal= landsat.select('B10').multiply(0.1);
    
    // find the min and max of NDVI
    var min_ndvi = ee.Number(ndvi.reduceRegion({
      reducer: ee.Reducer.min(),
      geometry: bbox,
      scale: 30,
      maxPixels: 1e9
      }
      ).values().get(0));
    
    var max_ndvi = ee.Number(ndvi.reduceRegion({
      reducer: ee.Reducer.max(),
      geometry: bbox,
      scale: 30,
      maxPixels: 1e9
      }
      ).values().get(0));
    
    // Fractional Vegetation
    var fv = (ndvi.subtract(min_ndvi).divide(max_ndvi.subtract(min_ndvi)))
      .pow(ee.Number(2))
      .rename('FV'); 
    
    // Emissivity
    var EM = fv.multiply(ee.Number(0.004)).add(ee.Number(0.986)).rename('EMM');
    
    // LST in Celsius Degree bring -273.15
    var LST = thermal.expression(
      '(Tb/(1 + (0.00115* (Tb / 1.438))*log(Ep)))-273.15', 
      {
        'Tb': thermal.select('B10'),
        'Ep': EM.select('EMM')
      }
      ).rename('LST');
    
    return LST;
}


function reduce_by_lulc(image, bbox, title) {
    // #################################
    // # plot band by land use class  ##
    // #################################
    
    var lulc = ee.ImageCollection("ESA/WorldCover/v200").first().clip(bbox);
    var combine_bands = image.addBands(lulc);
    
    // Grouped a mean reducer: LST by land cover category.
    var means = combine_bands.reduceRegion({
      reducer: ee.Reducer.mean().group({
        groupField: 1,
        groupName: 'landUseClass',
      }),
      geometry: bbox,
      scale: 1000,
      maxPixels: 1e8
    });
    
    // Helper function: Convert number to string
    function toString(number) {
      return ee.Number(number).format('%d')
    }
    
    // Create a dict for look up class names
    var lookup_names = ee.Dictionary.fromLists(
      ee.List(lulc.get('Map_class_values')).map(toString),
      lulc.get('Map_class_names')
    );
    
    // Create a dict for look up class colors
    var lookup_palette = ee.Dictionary.fromLists(
      ee.List(lulc.get('Map_class_values')).map(toString),
      lulc.get('Map_class_palette')
    );
    
    function createChartSliceDictionary(fc) {
      return ee.List(fc.aggregate_array("landcover_class_palette"))
        .map(function(p) {
          return {
            'color': p
          };
        }).getInfo();
    }
    
    // Create feature without geometry
    function createFeature(stats) {
      var stats_dict = ee.Dictionary(stats);
      var class_number = stats_dict.get('landUseClass');
      var result = {
        landcover_class_number: class_number,
        landcover_class_name: lookup_names.get(class_number),
        landcover_class_palette: lookup_palette.get(class_number),
        mean: stats_dict.get('mean')
      };
      return ee.Feature(null, result); // Creates a feature without a geometry.
    }
    
    // Object to list of group info
    var stats = ee.List(means.get('groups'));
    
    // Convert list to feature collection
    var by_lulc = ee.FeatureCollection(stats.map(createFeature));
    
    // Create bar chart
    var landcover_summary_chart = ui.Chart.feature.byFeature({
        features: by_lulc,
        xProperty: 'landcover_class_name',
        yProperties: ['mean']
      })
      .setChartType('BarChart')
      .setOptions({
        hAxis: {viewWindow: {min: 10, max: 50}},
        title: title,
        slices: createChartSliceDictionary(by_lulc),
        sliceVisibilityThreshold: 0 // Don't group small slices.
      });
    
    return landcover_summary_chart;
}


function maskL8sr(col) {
  // #################################
  // ###  Landsat 8 cloud mask   #####
  // #################################
  
  // 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 = col.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 col.updateMask(mask);
}


//##

function calculate_landuse_ratio(aoi) {
  
  var bbox = aoi.geometry();

  var lulc = ee.ImageCollection('ESA/WorldCover/v200')
  .first()
  .clip(bbox);
  
  var visualization = {
    bands: ['Map'],
  };
  
  Map.addLayer(lulc, visualization, "Landcover");
  
  // Generate an image in which the value of each pixel is the area of that pixel in square meters
  var area_image = ee.Image.pixelArea().addBands(lulc);
  
  // Group image of target area with class values
  var reduction_results = area_image.reduceRegion({
    reducer: ee.Reducer.sum().group({
      groupField: 1,
      groupName: 'Map_class_value',
    }),
    geometry: bbox,
    scale: 30,
    bestEffort: true,
  });
  
  // Get stats and create feature collection
  var roi_stats = ee.List(reduction_results.get('groups'));
  var landcover_fc = ee.FeatureCollection(roi_stats.map(createFeature));
  
  // Summarize land use with a pie chart
  var landcover_summary_chart = ui.Chart.feature.byFeature({
      features: landcover_fc,
      xProperty: 'landcover_class_name',
      yProperties: ['area_m2', 'landcover_class_number']
    })
    .setChartType('PieChart')
    .setOptions({
      title: 'Summary of landcover class areas',
      slices: createPieChartSliceDictionary(landcover_fc),
      sliceVisibilityThreshold: 0 // Don't group small slices.
    });
  
  return landcover_summary_chart;
}


// Create a dict for ESA land use cover look up class names
var worldcover_ds = ee.ImageCollection('ESA/WorldCover/v200').first();

var lookup_names = ee.Dictionary.fromLists(
  ee.List(worldcover_ds.get('Map_class_values')).map(toString),
  worldcover_ds.get('Map_class_names')
);

// Create a dict for look up class colors
var lookup_palette = ee.Dictionary.fromLists(
  ee.List(worldcover_ds.get('Map_class_values')).map(toString),
  worldcover_ds.get('Map_class_palette')
);

// Create a feature for a landcover class that includes the area covered.
function createFeature(landcover_class_stats) {
  landcover_class_stats = ee.Dictionary(landcover_class_stats);
  var class_number = landcover_class_stats.get('Map_class_value');
  var result = {
    landcover_class_number: class_number,
    landcover_class_name: lookup_names.get(class_number),
    landcover_class_palette: lookup_palette.get(class_number),
    area_m2: landcover_class_stats.get('sum')
  };
  
  return ee.Feature(null, result); // Creates a feature without a geometry.
}

// Create a JSON dictionary that defines piechart colors based on the
// landcover class palette.
// https://developers.google.com/chart/interactive/docs/gallery/piechart
function createPieChartSliceDictionary(fc) {
  return ee.List(fc.aggregate_array("landcover_class_palette"))
    .map(function(p) {
      return {
        'color': p
      };
    }).getInfo();
}

// Helper function: Convert number to string
function toString(number) {
  return ee.Number(number).format('%d')
}


function fvc_lst_impact(aoi){
  
  function rollingMean(collection) {
    function doIt(idx) {
      var d1 = date0.advance(ee.Number(idx).add(2), "year");  // last year will be excluded
      var d0 = d1.advance(-3, "year");
      return collection.filterDate(d0, d1)
      .reduce(ee.Reducer.median())
      .set("system:time_start", d0.advance(1, "year"));
    }
    return doIt;
  }
  
  // check location
  var bbox = aoi.geometry();
  var centroid = bbox.centroid({'maxError': 1});
  var listCoords = ee.Array.cat(centroid.coordinates(), 0);
  var yCoords = ee.List(listCoords).getInfo()[1];
  
  // south/ north
  var south = ee.Number(yCoords).lt(0);
  var shift = south.multiply(6).getInfo();
  
  // evaluation year
  var start_year = ee.Date("2000-01-01");
  var end_year = ee.Date("2021-01-01");
  
  // filter dataset to evaluation year
  var yearFilter = ee.Filter.date(
    start_year.advance(-1, "year"), 
    end_year.advance(1, "year")
    );
  
  // get summer months depends on location
  var summerFilter = ee.Filter.calendarRange(5+shift, 9-shift, "month");
  
  // modis land surface temperature (day and night)
  var modisLST = ee.ImageCollection("MODIS/061/MOD11A1")
    .filter(summerFilter).filter(yearFilter) // filter images
    .map(function(img) {
      return img.multiply(0.02).subtract(273.15) // apply scale factor and offset
      .copyProperties(img, ["system:time_start"]);
  });
  
  // extract Landsat time series and merge to single collection
  var ndviLandsat5 = ee.ImageCollection("LANDSAT/LT05/C01/T1_32DAY_NDVI").filter(yearFilter);
  var ndviLandsat7 = ee.ImageCollection("LANDSAT/LE07/C01/T1_32DAY_NDVI").filter(yearFilter);
  var ndviLandsat8 = ee.ImageCollection("LANDSAT/LC08/C01/T1_32DAY_NDVI").filter(yearFilter);
  var ndvi = ndviLandsat5.merge(ndviLandsat7.merge(ndviLandsat8)).filterBounds(aoi);
  
  var fvc = ndvi.map(function(img) {
    var f = ee.Image(
      (img.subtract(0.1))
      .divide(ee.Number(0.5).subtract(ee.Number(0.1)))
      )
      .multiply(100);
    return f.where(f.gt(100), 100).where(f.lt(0), 0).copyProperties(img, ["system:time_start"]);
  });
  
  // check out images with before and after scanerios
  var before = fvc.filter(ee.Filter.date('2001-01-01', '2001-12-01')).median().clip(aoi);
  var after = fvc.filter(ee.Filter.date('2020-01-01', '2020-12-01')).median().clip(aoi);
  

  // ################################################ //
  // ########## VEGETATION LOSS GROUP ############### //
  // ################################################ //
  
   // get urban areas
  var urban_mask = ee.ImageCollection('ESA/WorldCover/v200').first().clip(bbox).eq(50);
  
  // calculate urban density
  var reduceNeighborhood = urban_mask.reduceNeighborhood(
    ee.Reducer.mean(), 
    ee.Kernel.circle(5000, 'meters')
  ).multiply(100).toByte();

  // focus on dense areas
  var urban_threshold = reduceNeighborhood.gt(2);
  
  // Map.addLayer(before, {min: 0, max: 20, palette: ['white', 'green']}, "FVC 2001");
  // Map.addLayer(after, {min: 0, max: 20, palette: ['white', 'green']}, "FVC 2020");
  Map.addLayer(after.subtract(before).updateMask(urban_threshold), {min: -10, max: 10, palette: ['red', 'white', 'green']}, "FVC Difference");
  
  var threshold = (after.subtract(before)).lt(-20).and(urban_threshold);
  
  // FVC at different temporal aggregations (create sequence)
  var date0 = start_year.advance(1, "year");
  var n_years = end_year.difference(date0, "year").toInt().subtract(2);
  var idxs = ee.List.sequence(0, n_years);
  
  // Fractional Vegetation Cover, LST Day Time, LST Night Time
  // create smoothed temperol sequence
  var lstDay = ee.ImageCollection.fromImages(idxs.map(rollingMean(modisLST.select("LST_Day_1km"))));
  var lstNight = ee.ImageCollection.fromImages(idxs.map(rollingMean(modisLST.select("LST_Night_1km"))));
  var fvc_ds = ee.ImageCollection.fromImages(idxs.map(rollingMean(fvc)));
  
  // generate line chart for visualization
  
  // FRACTIONAL VEGETATION CHANGE 
  var fvcStyle = {
    title: 'Fractional Vegetation Cover Time Series (Vegetation Loss)',
            hAxis: {title: 'Year', titleTextStyle: {italic: false, bold: true}},
            vAxis: {
              title: 'FVC',
              titleTextStyle: {italic: false, bold: true}
            },
            lineWidth: 2,
            colors: ['green']
  };
  
  var fvcChart1 = ui.Chart.image.series({
      imageCollection: fvc_ds.map(function (image) {
        return image.updateMask(threshold)
      }),
      region: aoi,
      reducer: ee.Reducer.mean(),
      scale: 500
    }).setOptions(fvcStyle);
  
  // DAY TEMP
  var dayStyle = {
            title: 'MODIS Day LST Time Series (Vegetation Loss)',
            hAxis: {title: 'Year', titleTextStyle: {italic: false, bold: true}},
            vAxis: {
              title: 'Day Time Temp',
              titleTextStyle: {italic: false, bold: true}
            },
            lineWidth: 2,
            colors: ['red']
          };
  
  var dayChart1 = ui.Chart.image.series({
      imageCollection: lstDay.map(function (image) {
        return image.updateMask(threshold)
      }),
      region: aoi,
      reducer: ee.Reducer.median(),
      scale: 500
    }).setOptions(dayStyle);
  
  // NIGHT TEMP
  var nightStyle = {
    title: 'MODIS Night LST Time Series (Vegetation Loss)',
            hAxis: {title: 'Year', titleTextStyle: {italic: false, bold: true}},
            vAxis: {
              title: 'Night Time Temp',
              titleTextStyle: {italic: false, bold: true}
            },
            lineWidth: 2,
            colors: ['blue']
  };
  
  var nightChart1 = ui.Chart.image.series({
      imageCollection: lstNight.map(function (image) {
        return image.updateMask(threshold)
      }),
      region: aoi,
      reducer: ee.Reducer.median(),
      scale: 500
    }).setOptions(nightStyle);
  
  // ################################################ //
  // ########## VEGETATION GAIN GROUP ############### //
  // ################################################ //
  
  var threshold = (after.subtract(before)).gt(20).and(urban_threshold);
  
  // VEGETATION FRACTIONAL COVER
  var fvcStyle = {
    title: 'Fractional Vegetation Cover Time Series (Vegetation Gain)',
            hAxis: {title: 'Year', titleTextStyle: {italic: false, bold: true}},
            vAxis: {
              title: 'FVC',
              titleTextStyle: {italic: false, bold: true}
            },
            lineWidth: 2,
            colors: ['green']
  };
  
  var fvcChart2 = ui.Chart.image.series({
      imageCollection: fvc_ds.map(function (image) {
        return image.updateMask(threshold)
      }),
      region: aoi,
      reducer: ee.Reducer.mean(),
      scale: 500
    }).setOptions(fvcStyle);
  
  // DAY TEMP
  var dayStyle = {
            title: 'MODIS Day LST Time Series (Vegetation Gain)',
            hAxis: {title: 'Year', titleTextStyle: {italic: false, bold: true}},
            vAxis: {
              title: 'Day Time Temp',
              titleTextStyle: {italic: false, bold: true}
            },
            lineWidth: 2,
            colors: ['red']
          };
  
  var dayChart2 = ui.Chart.image.series({
      imageCollection: lstDay.map(function (image) {
        return image.updateMask(threshold)
      }),
      region: aoi,
      reducer: ee.Reducer.median(),
      scale: 500
    }).setOptions(dayStyle);
  
  // NIGHT TEMP
  var nightStyle = {
    title: 'MODIS Night LST Time Series (Vegetation Gain)',
            hAxis: {title: 'Year', titleTextStyle: {italic: false, bold: true}},
            vAxis: {
              title: 'Night Time Temp',
              titleTextStyle: {italic: false, bold: true}
            },
            lineWidth: 2,
            colors: ['blue']
  };
  
  var nightChart2 = ui.Chart.image.series({
      imageCollection: lstNight.map(function (image) {
        return image.updateMask(threshold)
      }),
      region: aoi,
      reducer: ee.Reducer.median(),
      scale: 500
    }).setOptions(nightStyle);
    
  return [
    fvcChart1, dayChart1, nightChart1,
    fvcChart2, dayChart2, nightChart2,
    ];
};