Converting USNG or MGRS to Latitude/Longitude in R

120 Views Asked by At

Is there a way to interpret (using sp, sf, etc.) grid locations expressed using the U.S. National Grid in R?

This is an example of a vector I would want to transform:

locations <- c("17T PH 8482 7766", "18T XL 6257 1639", "18T XL 6244 1636", "18T XL 6267 1673", "18T XL 6272 1671", "18T XL 6272 1674")

Ideally, I would want a way to get the Latitude/Longitude of the centroid of the grids, but even just being able to read them as a sf object would be half the work. Thank you!

I tried using the package mgrs (https://github.com/hrbrmstr/mgrs), but was unable to install it.

1

There are 1 best solutions below

0
margusl On BEST ANSWER

I can confirm that mgrs builds and installs just fine with remotes::install_github("hrbrmstr/mgrs"). At least on a Windows 10 machine with R 4.3.2 & RTools43 installed. However, if you still need an alternative, you could use JavaScript libraries through V8, e.g. https://github.com/proj4js/mgrs/ . It's kind of an overkill but perhaps an option if you can't build packages from source for whatever reason.

library(V8)
#> Using V8 engine 11.8.172.13
library(sf)
#> Linking to GEOS 3.11.2, GDAL 3.6.2, PROJ 9.2.0; sf_use_s2() is TRUE
locations <- c("17T PH 8482 7766", "18T XL 6257 1639", "18T XL 6244 1636", "18T XL 6267 1673", "18T XL 6272 1671", "18T XL 6272 1674")

# create a V8 environment for executing Javascript, 
# load mgrs from jsdelivr, assign js variable and
# apply mgrs.toPoint() on js array items
ct <- v8()
ct$source("https://cdn.jsdelivr.net/npm/[email protected]/dist/mgrs.min.js")
#> [1] "true"
ct$assign("locations", locations)
(lon_lat <- ct$get("locations.map(mgrs.toPoint)"))
#>           [,1]     [,2]
#> [1,] -78.72772 43.12917
#> [2,] -73.07329 40.78247
#> [3,] -73.07484 40.78222
#> [4,] -73.07202 40.78551
#> [5,] -73.07143 40.78532
#> [6,] -73.07143 40.78559

(points_sf <- sfheaders::sf_point(lon_lat) |> st_set_crs("WGS84"))
#> Simple feature collection with 6 features and 0 fields
#> Geometry type: POINT
#> Dimension:     XY
#> Bounding box:  xmin: -78.72772 ymin: 40.78222 xmax: -73.07143 ymax: 43.12917
#> Geodetic CRS:  WGS 84
#>                     geometry
#> 1 POINT (-78.72772 43.12917)
#> 2 POINT (-73.07329 40.78247)
#> 3 POINT (-73.07484 40.78222)
#> 4 POINT (-73.07202 40.78551)
#> 5 POINT (-73.07143 40.78532)
#> 6 POINT (-73.07143 40.78559)
mapview::mapview(points_sf)

Created on 2023-11-22 with reprex v2.0.2