## we can get that WKT by using Oblique Mercator params
prj <- "+proj=omerc +lonc=-122.6 +lat_0=40 +gamma=21 +x_0=1000000"
shp <- "/vsicurl/https://github.com/noaa-nwfsc/VMS-pipeline/raw/refs/heads/main/spatial_data/map_rotation/fivekm_grid_extent_rect_21.shp"
library(terra)
#> terra 1.8.52
v <- vect(shp)
wkt0 <- crs(prj)
## round trip with out proposed PROJ string form of the wkt
plot(project(project(v, "EPSG:4326"), wkt0))
plot(v, add = T, lwd = 2) ## it's the same
## get some raster data
dsn <- "/vsicurl/https://data.source.coop/alexgleith/gebco-2024/GEBCO_2024.tif"
target <- rast(ext(v), res = 2500, crs = prj) ## tweak resolution to be sensible (not too many pixels)
target
#> class : SpatRaster
#> size : 990, 826, 1 (nrow, ncol, nlyr)
#> resolution : 2499.146, 2500.01 (x, y)
#> extent : -231257, 1833038, -1252242, 1222769 (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=omerc +lat_0=40 +lonc=-122.6 +alpha=0 +gamma=21 +k=1 +x_0=1000000 +y_0=0 +datum=WGS84 +units=m +no_defs
topo <- project(rast(dsn), target, by_util = TRUE)
plot(topo); plot(v, add = TRUE)
Created on 2025-06-15 with reprex v2.0.2