ww3.md

fun stuff, compare two polygons from one dataset by plotting them together (local equal area projection ensures scale is right)

library(terra)
x <- terra::vect("/vsizip//vsicurl/https://github.com/wmgeolab/geoBoundaries/raw/main/releaseData/CGAZ/geoBoundariesCGAZ_ADM0.zip")

crs_local <- function(x) {
    cc <- terra::geom(terra::centroids(x))[, c("x", "y"), drop = FALSE]
    sprintf("+proj=laea +lon_0=%i +lat_0=%i", round(cc[,1]), round(cc[,2]))
    
}
crs_local(x)



## take two polygons
#i1 <- which(x$shapeName == "Brazil")
#i2 <- which(x$shapeName == "United States")


## repeat remaining code, or input your own two polygons
i1 <- sample(nrow(x), 1)
i2 <- sample(setdiff(1:nrow(x), i1), 1)
local1 <- project(x[i1, ], crs_local(x[i1, ]))
local2 <- project(x[i2, ], crs_local(x[i2, ]))
set.crs(local1, NA) ## make it official
set.crs(local2, NA)
plot(rbind(local1, local2), border = c("firebrick", "dodgerblue"))
expanse(x[c(i1, i2), ])/1e6