89_Simple
89_Simple

Reputation: 3805

extract values of multiple lat lon for image collection in earth engine using rgee

I have a dataframe of lat-lon for which I want to extract the daily rainfall data from the "NASA/NEX-GDDP" dataset which is daily temporal resolution and 0.25 degree spatial resolution. My thought process is to

  1. download daily raster on my local drive. So for a given year, there will 365 rasters
  2. extract the values for each daily raster for the my lats lons
  3. remove the daily raster to save space on disk

I need to do this for multiple years. The above steps are very time consuming esp the step 1. Is there any way I can extract the daily data directly instead of having to download rasters. Below is how I am doing currently

  library(rgee)
  ee_Initialize(email = '[email protected]')
  
  # define a bounding box within which my lat lon lies 
  polygon <- ee$Geometry$Polygon(
    list(
      c(-76.12, -36.69),
      c(-15.42,-36.69),
      c(-15.42, 6.71),
      c(-76.12, 6.71),
      c(-76.12, -36.69)
    ))
  
  latlon <- structure(list(lat = c(-11.125, -11.125, -10.875, -10.875, -10.875, -10.875, -10.875, -10.875), 
                           lon = c(-37.375, -37.125, -70.625, -70.375, -70.125, -69.875, -69.625, -69.375)), 
                      class = "data.frame", row.names = c(NA, -8L))
  
  model_ref <- 'ACCESS1-0' 
  rcp_ref <- 'historical'
  year_ref <- 2005
  
  Imagecollection_pr <-
    ee$ImageCollection("NASA/NEX-GDDP")$ # read the data
    select('pr')$ # select rainfall
    filter(ee$Filter$calendarRange(year_ref, year_ref, "year"))$ # filter the year
    filter(ee$Filter$eq("model", model_ref))$ # filter the model
    filter(ee$Filter$eq("scenario", rcp_ref)) # filter the scenario

  # Downloading as raster from Earth Engine to Local. This downloads 365 images for 2005
  ee_imagecollection_to_local(
    ic = Imagecollection_pr,
    region = polygon,
    dsn = paste0('pr'),
    via = 'getInfo')
  
  # extract values for the given latlon
  library(raster)
  library(rgdal)
  
  ll <- latlon
  
  coordinates(ll) <- ~ lon + lat
  proj4string(ll) <- "+proj=longlat +datum=WGS84 +no_defs"
  
  raster_ls <- list.files(file.path(getwd()), pattern = '.tif')
  
  temp_ls <- list()
  
  for(rr in seq_along(raster_ls)){

    raster_temp <- raster_ls[rr]
    raster_doy <- raster(raster_temp)
    
    crs(raster_doy) <- "+proj=longlat +datum=WGS84 +no_defs"
    
    ymd <- 
      as.numeric(gsub(model_ref,"",raster_temp) %>% 
                   gsub('pr',"",.) %>% 
                   gsub(rcp_ref,"",.) %>%
                   gsub(".tif","",.) %>%
                   gsub("___","",.)) 
    
    latlon$val <- extract(raster_doy, ll)
    temp_ls[[rr]] <- latlon %>% dplyr::mutate(ymd = ymd)
  }  
  
  file.remove(raster_ls)
  results <- do.call('rbind', temp_ls)
  head(results)      
  
  lat     lon          val      ymd
  1 -11.125 -37.375 1.677004e-05 20050101
  2 -11.125 -37.125 1.582177e-05 20050101
  3 -10.875 -70.625 4.219083e-05 20050101
  4 -10.875 -70.375 3.825104e-05 20050101
  5 -10.875 -70.125 3.438132e-05 20050101
  6 -10.875 -69.875 2.946513e-05 20050101

Upvotes: 5

Views: 455

Answers (0)

Related Questions