CODEWITHSUNDEEP
rggmapr-sf
geotheory
geotheory

Reputation: 23630

Cannot plot sf layer on top of ggmap

Subject of some previous questions - e.g. this. I cannot get the sf layer to render over the ggmap layer.

require(ggmap)
require(raster)
require(sf)

target_bbox <- c(left = -0.65, bottom = 51.2, right = 0.45, top = 51.8)
map = get_stamenmap(target_bbox, zoom = 9, maptype = 'toner-2010')

# can test this with: ggmap(map)

e = extent(list(x = c(-0.65, 0.45), y = c(51.2, 51.8)))
r = raster(volcano)
crs(r) = CRS("+proj=longlat +datum=WGS84")
p = raster::rasterToPolygons(r)
sf = st_as_sf(p) %>% st_transform(3857)

# test with: ggplot() + geom_sf(data = sf, aes(fill = layer), col = NA)

ggmap(map) + 
   geom_sf(data = sf, aes(fill = layer), col = NA, alpha = 0.4, inherit.aes = FALSE)
Coordinate system already present. Adding new coordinate system, which will replace the existing one.

enter image description here

Coordinate system clash but no errors, yet no layer is visible. Any suggestions much appreciated.

This github solution - extracting the bounding box from the ggmap object with custom function ggmap_bbox() - only succeeds in adding another warning to that above:

Warning message:
In st_is_longlat(x) :
  bounding box has potentially an invalid value range for longlat data

Upvotes: 2

Views: 1395

Answers (1)

mrhellmann
mrhellmann

Reputation: 5489

It looks like you're trying to plot the sf volcano object over London, but you've defined it to be at lon/lat 0,0. Other than that, your call to ggmap and then geom_sf works fine. The problem was just that the sf object wasn't in the bounding box of the ggmap.

I had to change it to a 'Spatial' object to move it using maptools::elide, and then back to an 'sf' object as I couldn't find an easy way to shift sf polygons. The shift is an eyeball approximation.

sf_moved <- sf %>% 
  st_transform(3857) %>% # <- may be unnecessary
  as('Spatial') %>% 
  maptools::elide(.,shift = c(-80000, 6650000)) %>% 
  as('sf') %>% 
  st_set_crs(3857) %>% 
  st_transform(4326) 

ggmap(map) + 
  geom_sf(data = sf_moved, inherit.aes = FALSE, 
          col = NA, aes(fill = layer), alpha = .4) + 
  scale_fill_viridis_c()

enter image description here

Upvotes: 1

Related Questions