display tiles - flood

display.R

viz_patches <- function(pid) {
  metadata_list <- read.csv("https://storage.googleapis.com/flood_dataset/flood-training-metadata.csv")
  date <- metadata_list[metadata_list$chip_id %in% pid, 6][1]
  
  url <- sprintf("https://storage.googleapis.com/flood_dataset/%s.tif", pid)  
  local_tif <- tempfile(fileext = ".tif")
  download.file(url, local_tif)
  
  # BOX
  ee_extent <- local_tif %>% 
    raster() %>%
    projectRaster(crs = "EPSG:4326") %>% 
    extent() %>%
    as.matrix() %>%
    as.numeric() %>% 
    ee$Geometry$Rectangle()
  bbox_map <- Map$addLayer(ee_extent, name = "box", shown = FALSE)
  
  
  Map$centerObject(ee_extent)
  

  # # -----------------------------------------------------------------------
  # JRC
  # # -----------------------------------------------------------------------
  surface_water <- ee$Image("JRC/GSW1_3/GlobalSurfaceWater")
  
  ocurrence_vizparams <- list(bands = "occurrence", min = 0, max = 100,
                    palette = list('ffffff', 'ffbbbb', '0000ff'))
  ocurrence_map <- Map$addLayer(surface_water, ocurrence_vizparams, 'Occurrence', shown = FALSE)
  
  recurrence_vizparams <- list(bands = "recurrence", min = 0, max = 100,
                              palette = list('ffffff', 'ffbbbb', '0000ff'))
  recurrence_map <- Map$addLayer(surface_water, recurrence_vizparams, 'Recurrence', shown = FALSE)
  
  seasonality_vizparams <- list(bands = "seasonality", min = 0, max = 12,
                               palette = list('ffffff', 'ffbbbb', '0000ff'))
  seasonality_map <- Map$addLayer(surface_water, seasonality_vizparams, 'Seasonality', shown = FALSE)

  max_extent_vizparams <- list(bands = "max_extent", min = 0, max = 1,
                                palette = list('ffffff', 'ffbbbb', '0000ff'))
  max_extent_map <- Map$addLayer(surface_water, max_extent_vizparams, 'Max_extent', shown = FALSE)
  
  # # -----------------------------------------------------------------------
  # Sentinel-1
  # # -----------------------------------------------------------------------
  init_date <- as.Date(date) %>% as.character()
  last_date <- (as.Date(date) + 2) %>% as.character()

  rescale_values <- function(image) {
    
    minMax <- image$reduceRegion(
      reducer = ee$Reducer$minMax(),
      geometry = ee_extent,
      scale = 10
    )
    
    unit_scale <- ee_utils_pyfunc(function(name) {
      name <- ee$String(name)
      band <- image$select(name)
      band$unitScale(ee$Number(minMax$get(name$cat('_min'))), ee$Number(minMax$get(name$cat('_max'))))
    })
    
    unitScale <- ee$ImageCollection$fromImages(
      image$bandNames()$map(unit_scale))$toBands()$rename(image$bandNames())
    
    vv <- unitScale$select("VV")
    vh <- unitScale$select("VH")
    vv_div_vh <- ee$Image$divide(vv, vh)
    ee$Image(list(vv, vh, vv_div_vh))$rename(c("R","G","B"))
  }
  
    
  sar_img <- ee$ImageCollection('COPERNICUS/S1_GRD')$
    filter(ee$Filter$listContains('transmitterReceiverPolarisation', 'VV'))$
    filterDate(init_date, last_date)$
    filterBounds(ee_extent)$
    filter(ee$Filter$eq('instrumentMode', 'IW'))$
    select(c('VV', 'VH'))
  dates <- ee_get_date_ic(sar_img)$time_start %>% as.Date()
  sar_map <- Map$addLayer(sar_img$map(rescale_values)$first(), name = paste0("SAR-", dates), shown = FALSE)
  
  to_display <- recurrence_map + ocurrence_map + 
    seasonality_map + max_extent_map + 
    sar_map + bbox_map
  
  # # -----------------------------------------------------------------------
  # Sentinel-2
  # # -----------------------------------------------------------------------
  dataset <- ee$ImageCollection('COPERNICUS/S2_SR')$
    filterDate(init_date, last_date)$
    filterBounds(ee_extent)
  
  if (dataset$size()$getInfo() != 0) {
    dataset <- dataset$first()$divide(10000) 
    visualization <- list(min = 0, max = 0.3, bands =  c('B4', 'B3', 'B2'))
    rgb <- Map$addLayer(dataset, visualization, name = "sentinel2", shown = FALSE)
    to_display <- to_display + rgb
  }
  
  # # -----------------------------------------------------------------------
  # Landsat-8
  # # -----------------------------------------------------------------------
  dataset <- ee$ImageCollection$Dataset$LANDSAT_LC08_C01_T1_SR$
    filterDate(init_date, last_date)$
    filterBounds(ee_extent)
  
  if (dataset$size()$getInfo() != 0) {
    dataset <- dataset$first()$divide(10000)
    visualization <- list(min = 0, max = 0.3, bands =  c('B4', 'B3', 'B2'))
    rgb <- Map$addLayer(dataset, visualization, name = "landsat", shown = FALSE)
    to_display <- to_display + rgb
  }
  
  # # -----------------------------------------------------------------------
  # DEM
  # # -----------------------------------------------------------------------
  dem <- ee$Image("NASA/NASADEM_HGT/001")
  
  rscale2 <- function(image) {
    image <- image$select("elevation")
    minMax <- image$reduceRegion(
      reducer = ee$Reducer$minMax(),
      geometry = ee_extent,
      scale = 10
    )
    unit_scale <- ee_utils_pyfunc(function(name) {
      name <- ee$String(name)
      band <- image$select(name)
      band$unitScale(ee$Number(minMax$get(name$cat('_min'))), ee$Number(minMax$get(name$cat('_max'))))
    })
    
    ee$ImageCollection$fromImages(
      image$bandNames()$map(unit_scale))$toBands()$rename(image$bandNames())  
  }
  
  demr <- rscale2(image = dem)
  palette <- c("#04D7FF", "#F1BEBE", "#E19D9D", "#D07D7D", "#BC5D5D", "#A94141", "#962929", "#861919", "#750B0B", "#640000")
  map_dem <- Map$addLayer(demr, list(max = 0.5, palette=palette), name = "dem", shown = FALSE)
  to_display <- map_dem + to_display
  target <- raster(local_tif)
  target_map <- suppressWarnings(mapview(target, to_display, legend = FALSE))
  target_map
}