mikejohnson51 / climateR

An R 📦 for getting point and gridded climate data by AOI
https://mikejohnson51.github.io/climateR/
MIT License
168 stars 40 forks source link

Extract climate values to point locations #4

Closed meaghanregina closed 5 years ago

meaghanregina commented 5 years ago

Hello-- I am looking for a way to extract daily climate data (from the MACA climate dataset) from 5 GCMs to ten latitude/longitude point locations, instead of an AOI extent. Thanks!

mikejohnson51 commented 5 years ago

Hi Meaghan,

With the new update, point can be found using a lat/long pair.

A single point can be passed at a time. Certainly, a for loop could pass through many points but for more speed, I would suggest the 'foreach package. This example should retrieve 3 parmeters for 10 points in Colorado over 1 year across 5 random GCMs (54,750 values):

# Get Bounding Box of Colorado
AOI = getAOI(state = 'CO') %>% AOI::bbox_st()

# Generate Random Points
n = 10
pts = data.frame(lat = runif(n, AOI$ymin, AOI$ymax),
                 lon = runif(n, AOI$xmin, AOI$xmax))

system.time({
  no_cores  <- parallel::detectCores() - 1
  doParallel::registerDoParallel(no_cores)

  output = foreach::foreach(i = 1:NROW(pts), .packages = 'climateR') %dopar% {
    getMACA(
      AOI = pts[i,],
      param = c('tmax', 'prcp', 'shum'),
      model = 5,
      startDate = "2018-01-01",
      endDate   = "2018-12-31"
    )

  }
})

  user  system elapsed 
  2.511   2.519  14.482 

The output will look something like this, a list of data.frames with (# of models) x (# of parameters) COLUMNS and (# number of days) ROWS.

str(output)

List of 10
 $ :'data.frame':   365 obs. of  19 variables:
  ..$ source             : chr [1:365] "maca" "maca" "maca" "maca" ...
  ..$ lat                : num [1:365] 37.4 37.4 37.4 37.4 37.4 ...
  ..$ lon                : num [1:365] -105 -105 -105 -105 -105 ...
  ..$ date               : Date[1:365], format: "2018-01-01" "2018-01-02" ...
  ..$ bcc-csm1-1_tmax    : num [1:365] 274 274 274 274 277 ...
  ..$ bcc-csm1-1_prcp    : num [1:365] 0 0 0.915 0 0 ...
  ..$ bcc-csm1-1_shum    : num [1:365] 0.00126 0.00178 0.00179 0.00189 0.00198 ...
  ..$ CNRM-CM5_shum      : num [1:365] 0.00261 0.0022 0.0013 0.00263 0.00228 ...
  ..$ CNRM-CM5_prcp      : num [1:365] 3.668 0.352 1.67 23.092 0.376 ...
  ..$ CNRM-CM5_tmax      : num [1:365] 276 274 275 274 275 ...
  ..$ HadGEM2-ES365_shum : num [1:365] 0.00191 0.00239 0.00292 0.0033 0.00271 ...
  ..$ HadGEM2-ES365_prcp : num [1:365] 0 0 1.55 5.47 1.77 ...
  ..$ HadGEM2-ES365_tmax : num [1:365] 277 278 278 275 275 ...
  ..$ IPSL-CM5A-MR_tmax  : num [1:365] 273 270 275 275 275 ...
  ..$ IPSL-CM5A-MR_shum  : num [1:365] 0.00169 0.00168 0.00171 0.00189 0.0021 ...
  ..$ IPSL-CM5A-MR_prcp  : num [1:365] 4.3 4.49 1.36 0 0 ...
  ..$ MIROC-ESM-CHEM_prcp: num [1:365] 0 0 0 0 0.291 ...
  ..$ MIROC-ESM-CHEM_tmax: num [1:365] 275 277 281 278 276 ...
  ..$ MIROC-ESM-CHEM_shum: num [1:365] 0.00127 0.00236 0.00277 0.00181 0.00179 ...

...

plot(x = var[[1]]$date, 
       y = var[[1]]$`bcc-csm1-1_tmax`, 
       type = 'l', main = "TMAX Point 1", xlab = "Date", ylab = "TMAX (K)" )
image

Hope its usefull!

Mike

meaghanregina commented 5 years ago

Thank you!!