My goal is to create a world map (using ggplot2) with a Pacific perspective.
I have managed to re-center the map, by following this solution. Code below:
# world country layer sppolys <- rnaturalearth::ne_countries(scale = "medium", returnclass = "sp") wgs84 <- CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0") # ========================================================================= # Split sPolygonsDF object by "split line" # ========================================================================= # shift central/prime meridian towards west - positive values only shift <- 152 # create "split line" to split worldmap (split at Prime Meridian) split.line <- sp::SpatialLines(list(Lines(list(Line(cbind(180 - shift, c(-90, 90)))), ID="line")), proj4string = wgs84) sppolys <- gBuffer(sppolys, byid=TRUE, width=0) # intersecting line with world polys line.gInt <- gIntersection(split.line, sppolys) # create a very thin polygon (buffer) out of the intersecting "split line" bf <- gBuffer(line.gInt, byid=TRUE, width=0.000001) # split world polys using intersecting thin polygon (buffer) sppolys.split <- gDifference(sppolys, bf, byid=TRUE) #========================================================== ## plot #========================================================== countries <- data.table(map_data(as(sppolys.split, "SpatialPolygonsDataFrame"))) # Shift coordinates to fall correctly on shifted map countries$long.shift <- countries$long + shift countries$long.shift <- ifelse(countries$long.shift > 180, countries$long.shift - 360, countries$long.shift) # plot shifted map ggplot() + geom_polygon(data=countries, aes(x=long.shift, y=lat, group=group), colour="black", fill="gray80", size = 0.25) + coord_equal()
Which results in the map below, providing a Pacific perspective, as desired.
However, polygons along the anti-meridian (now the map center) are still split. I suppose these could be unioned together somehow, which I tried to do with the following code (basically the inverse of the above splitting).
# create SpatialLine where along anti-meridian old.split.line <- sp::SpatialLines(list(Lines(list(Line(cbind(180, c(-90, 90)))), ID="line")), proj4string = wgs84) # ID intersection btwn this line and polygons layer old.line.gInt <- gIntersection(old.split.line, sppolys.split) # buffer it to ensure intersection old.bf <- gBuffer(old.line.gInt, byid=TRUE, width=0.000001) # union polygons which intersect the line sppolys.union <- gUnion(sppolys.split, old.bf)
The problem with this solution (as far as I can tell), is that the re-centering is done after converting to a data.table object, and therefore there is no actual intersection between the polygon parts. The only solution I can think if is to reconstruct a polygons object (
sf) after the shifting is done, but this is a huge hassle just to perform the unioning.
So in a nutshell, I'm looking for a way to shift/re-center the polygons THEN combine polygons along the old anti-meridian.
This can either be before or after converting to a data.table, whatever is most efficient. (Simple Features solutions are also welcome)