geom_raster_dap.R

#install.packages("rerddap")
#devtools::install_github("ropensci/plotdap")
library(rerddap)
library(plotdap)
library(ggplot2)
sstInfo <- info('jplMURSST41')

target_field <- 'analysed_sst'
murSST <- griddap(
   'jplMURSST41', latitude = c(22, 51), longitude = c(-140, -105),
   time = c('last', 'last'), fields = target_field
)
#p <- plotdap(crs = "+proj=robin")
#add_griddap(p, murSST, ~analysed_sst)
## the original data, in raster form
dres <- c(mean(diff(sort(unique(murSST$data$lon)))), mean(diff(sort(unique(murSST$data$lat)))))
prj <- "+init=epsg:4326"  ## manual as it doens't seem to be in the info or griddap object
prj_target <- "+proj=robin"
library(raster)
r_original <- 
   raster(extent(range(murSST$data$lon) + c(-1, 1) * dres[1]/2, 
                            range(murSST$data$lat) + c(-1, 1) * dres[2]/2), crs = prj, 
                     res = dres)
## project the abstraction of the raster
r_target <- projectExtent(r_original, prj_target)

r_original[] <- NA
r_original[cellFromXY(r_original, cbind(murSST$data$lon, murSST$data$lat))] <- murSST$data[[target_field]]

maxpixels <- 50000
## resize the abstraction of the target raster
dim(r_target) <- dim(r_target)[1:2] %/% sqrt(ceiling(ncell(r_target) / maxpixels))

## apply the reprojection to the raster data
r_target <- projectRaster(r_original, r_target)
names(r_target) <- target_field

ggplot() + 
  ## there's a way to use aes_ to have the variable name 
  geom_raster(data = as.data.frame(r_target, xy = TRUE), aes(x, y, fill = analysed_sst))  + 
  coord_sf(crs = prj_target)