hypertidy/ceramic

Cannot use `sf` object as `loc` argument

vronizor opened this issue · 4 comments

Hi and thanks for the great package.

I used to use sf objects as loc arguments inside the cc_location method, but it does not seem to work this time around. The method still works when feeding a point as argument, though. Below is a reprex of the issue. Not sure what is going on with the proj4::ptransform error, any clues?

library(sf)
#> Linking to GEOS 3.8.1, GDAL 3.2.1, PROJ 7.2.1
library(ceramic)
library(tmap)

#### Load NC shapefile and NYC long/lat ####

nc = st_read(system.file("shape/nc.shp", package="sf"))
#> Reading layer `nc' from data source 
#>   `/Library/Frameworks/R.framework/Versions/4.1/Resources/library/sf/shape/nc.shp' 
#>   using driver `ESRI Shapefile'
#> Simple feature collection with 100 features and 14 fields
#> Geometry type: MULTIPOLYGON
#> Dimension:     XY
#> Bounding box:  xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
#> Geodetic CRS:  NAD27

nyc = cbind(-73.941130, 40.714414)

#### Download with Ceramic ####

map.nc = cc_location(loc = nc, buffer = 1e8,
                      base_url = "https://api.mapbox.com/styles/v1/mapbox/light-v10/tiles/{zoom}/{x}/{y}")
#> Warning in to_proj(source): not a proj-like string

#> Warning in to_proj(source): not a proj-like string
#> Error in proj4::ptransform(x, source, target, ...): no arguments in initialization list

map.nyc = cc_location(loc = nyc, buffer = 1e4,
                     base_url = "https://api.mapbox.com/styles/v1/mapbox/light-v10/tiles/{zoom}/{x}/{y}")
#> Preparing to download: 12 tiles at zoom = 12 from 
#> https://api.mapbox.com/styles/v1/mapbox/light-v10/tiles/
#> Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj =
#> prefer_proj): Discarded ellps unknown in Proj4 definition: +proj=merc +lon_0=0
#> +k=1 +x_0=0 +y_0=0 +R=6378137 +units=m +no_defs +type=crs
#> Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj =
#> prefer_proj): Discarded datum unknown in Proj4 definition

#### Map ####

tm_shape(map.nc, raster.downsample = F) +
  tm_rgb()
#> Error in as.list.environment(environment()): object 'map.nc' not found

tm_shape(map.nyc, raster.downsample = F) +
  tm_rgb()

Created on 2021-06-28 by the reprex package (v2.0.0)

yeh, I'm afraid it's very difficult to get a basic extent and crs for use for transforming with these big player packages ... the core problem is the raster::projection(spex::spex(nc)) now gives only the string "NAD27", which is some combination of sp and raster ... spex was my approach to a "extent with crs" back in the day, and I'll need to have a think about how that gets updated. Thanks for the note.

I'd recommend just using a raster extent in longlat if you have trouble with sf objects, something like

raster::extent(sf::st_transform(x, 4326))

that will be pretty reliable (we don't need sub-metre accuracy for a basic extent to get a rough area ...

Concretely:

map.nc = cc_location(loc = raster::extent(sf::st_transform(nc, 4326)), buffer = 1e8,
                     base_url = "https://api.mapbox.com/styles/v1/mapbox/light-v10/tiles/{zoom}/{x}/{y}")

that redundantly attempts transformation "longlat" (different datum but no one cares at this scale), but for general sf objects will also work because ceramic can guess well enough from longlat data. Unless you have tricky polar regions or very weird projections you should be fine in most places.

HTH, and sorry there's no obvious fix - I'll have a think about where this goes but there is no fundamental support for these basic shared things for dev sad to say.

Hey @mdsumner, thanks so much for taking the time to answer so thoroughly!

Thanks as well for the fix, looks way better than my own terrible hack (getting the center of the bbox and eyeballing the buffer...).

I guess this is related: I just suppressed warnings in the slopes package, which auto gets elevations thanks to this great package. In case it's of use/interesting, and in case I'm doing something wrong! ropensci/slopes#37

I think I've fixed this, raster::projection now returns those weird strings that only work in some situations - in this case it's "NAD27" so I catch it.