Try a couple more ways to streamline map rotation in R

[brooke-l-hawkins-noaa]

I'm making a new issue, since I didn't get to this as part of the spatial_data reorganization ticket #45.

This ticket is about rotating maps, and getting a nice mapping workflow as part of the VMS pipeline. I worked on this some and got a workaround, but I'd like something more streamlined.

Here are my thoughts for outputs:

• Create a notebook (Rmd or Qmd) with a couple methods explored for rotating maps
  • Is there a way to use rot = function(a) matrix(c(cos(a), sin(a), -sin(a), cos(a)), 2, 2) from the sf documentation that doesn’t remove geometry from the sf object, so the map can still be made with ggplot? Or is there some other workaround?
  • Is there a way to programmatically adjust the CRS using WKT (example file), rather than needing to load an adjusted shapefile for every degree of rotation (example file)?
• Optionally follow up on sf issue that Brooke started and never replied to…
  • Having supporting code and a non-confidential example would help the discussion.
  • Edzer responded by suggesting I use st_transform() and the new projection, and that’s what I’m already doing, but I think it shouldn’t be necessary. What I actually want is to be able to rotate the output map without needing to upload a new CRS, just by adding a rotation argument in sf_transform or geom_sf. (Is that ignorant spatially of me? I understand that a rotation is technically a projection, but it also seems a bit silly to not be able to just rotate some pixels…)

This is low priority. If we use something interactive like leaflet in an R Shiny app for an interactive visualization tool, then the rotation is less important. This would be relevant for static figures in slides or a paper, though.

[mdsumner]

I thought it might be helpful to generate that rotated projection (Oblique Mercator) via PROJ string, here's an example that gets the shapefile and plots it with topography.

https://gist.github.com/mdsumner/e5cf8d10574bf61f99c6a84f5dc1378e

(your crs has a false easting x_0 1000000 where the defaults are 0, this is entirely arbitrary and is just a shift applied after reprojection to move where 0,0 is)

happy to follow up and help if I can, it can be painful to determine an extent with the right aspect ratio and so forth but you can use project on an extent fwiw if you provide the source projection as well:

prj <- "+proj=omerc +lonc=-122.6 +lat_0=40 +gamma=21 +x_0=1000000"
project(ext(c(xmin = -132, xmax = -116, ymin = 30, ymax = 49)), to = prj, from = "EPSG:4326")


prj <- "+proj=omerc +lonc=-122.6 +lat_0=40 +gamma=15 +x_0=1000000"
project(ext(c(xmin = -132, xmax = -116, ymin = 30, ymax = 49)), to = prj, from = "EPSG:4326")

https://proj.org/en/stable/operations/projections/omerc.html#central-point-and-azimuth-method

seeing your question in recent weeks I tried a few things but wasn't sure the context you were needing so I didn't pursue further, here that is in case it's of interest:

https://gist.github.com/mdsumner/274f191f8718ecb9a07bca0a41e572c8

definitely omerc is a good way to "rotate" an arbitrary location, and if you have a crs for that and an extent you can drive any other data into that coordinate system.

[befeist]

@mdsumner, thanks so much for looking into this, you're too kind! Rotating maps that use any CRS is a trivial exercise using conventional GIS software, but it seems that that simplicity is not readily available in R. So I had to find a CRS where the rotation could be achieved using the parameters, and the Oblique Mercator projection (aka, the rectified skew orthomorphic) was the only option I could find. I think we have settled on a workflow that reads the .wkt file and then the user sets the rotation to the desired angle. I would imagine that this workflow could also use the proj string instead of the .wkt, but I've heard buzz about proj strings potentially becoming deprecated, so just went with .wkt to play it safe.