In paleolimbot/crs2crs: Generic Coordinate System Transforms

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.path = "man/figures/README-",
  out.width = "100%"
)

crs2crs

The goal of crs2crs is to provide a flexible framework for reproducible coordinate transforms. With the last few releases of PROJ it can be difficult to predict which transform will be applied when converting between spatial reference systems. The crs2crs package provides tools to make these more explicit or to define explicit transforms between reference systems for future reproducibility. It is capable of using PROJ in several forms but doesn't require it.

Installation

You can install the development version from GitHub with:

# install.packages("devtools")
devtools::install_github("paleolimbot/crs2crs")

If you can load the package, you're good to go!

library(crs2crs)

Example

Use the PROJ command-line engine to transform some coordinates:

crs_set_engine(crs_engine_proj_cmd())
pos <- wk::xy(-64, 45, crs = "OGC:CRS84")
crs_transform(pos, "EPSG:3857")

You can also define an engine that uses R functions if you only need a few transforms:

engine <- crs_engine_fun()
engine <- crs_engine_fun_define(engine, "EPSG:3857", "OGC:CRS84", function(coords) {
  r <- 6378137
  coords$x <- coords$x * pi / 180 * r
  coords$y <- log(tan(pi / 4 + coords$y * pi / 180 / 2)) * r
  coords
})
crs_set_engine(engine)
crs_transform(pos, "EPSG:3857")

The crs_transform() function is built on wk::wk_transform(), so you can pass it any object that defines a wk::wk_handle() method (like sf objects!).

library(sf)
nc <- read_sf(system.file("shape/nc.shp", package = "sf"))
nc %>% 
  st_transform("OGC:CRS84") %>% 
  crs_transform("EPSG:3857")

You can also combine engines to, for example, force a custom pipeline to be used for a given conversion:

proj_engine <- crs_engine_proj_cmd()
engine <- crs_engine_fun_define(engine, "EPSG:3857", "OGC:CRS84", function(coords) {
  pipe <- paste(
    "+proj=pipeline",
    "+step +proj=unitconvert +xy_in=deg +xy_out=rad",
    "+step +proj=webmerc +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84"
  )
  xy_out <- crs_transform_pipeline(wk::as_xy(coords), pipe, engine = proj_engine)
  coords[c("x", "y")] <- as.data.frame(xy_out)[c("x", "y")]
  coords
})

crs_set_engine(engine)
nc %>%
  st_transform("OGC:CRS84") %>% 
  crs_transform("EPSG:3857")

If you're running PROJ >= 7.1 (sorry Ubuntu 20.04 without ubuntugis!) you can also use sf's interface to PROJ which is much faster (especially for sf objects):

# for non-ballpark datum transforms
sf_proj_network(TRUE, "https://cdn.proj.org/")

crs_set_engine(crs_engine_proj_sf())
nc %>%
  crs_transform_pipeline("+proj=axisswap +order=2,1") %>% 
  crs_set("NAD27") %>% 
  crs_transform("OGC:CRS84")

You can also use crs_engine_sf(), which just wraps sf::st_transform():

crs_set_engine(crs_engine_sf())
crs_transform(nc, "OGC:CRS84")

paleolimbot/crs2crs documentation built on Jan. 8, 2022, 6:25 a.m.