CODEWITHSUNDEEP
rgisrgdal
andybega
andybega

Reputation: 1437

Distortions when reprojecting cshapes world map

I am trying to plot a reprojected world map using data from the cshapes package and sp::spTransform, but the projection leads to distorted plots. How can I correctly reproject and plot a cshapes map?

Here is an example that shows that the map plots fine by itself (code adapted from this blog post):

library("cshapes")
library("ggplot2")
library("rgdal")

wmap <- cshp(date=as.Date("2012-06-30"))
wmap_df <- fortify(wmap)

ggplot(wmap_df, aes(long,lat, group=group)) + 
  geom_polygon() + 
  labs(title="World map (longlat)") + 
  coord_equal()

ggsave("~/Desktop/map1.png", height=4, width=7)

Map without distortion

And here is the distorted version when I reproject to Robinson:

wmap_robin <- spTransform(wmap, CRS("+proj=robin"))
wmap_df_robin <- fortify(wmap_robin)
ggplot(wmap_df_robin, aes(long,lat, group=group)) + 
  geom_polygon() + 
  labs(title="World map (robinson)") + 
  coord_equal()

ggsave("~/Desktop/map2.png", height=4, width=7)

Distorted map

Some additional info:

  1. I know there are other data sources for country borders, but I need maps that reflect changes in country borders, which cshapes does.
  2. My guess is the the problem is related to issues with the underlying map polygons, but I have no idea where to begin looking and it's probably better to ask what I want to get in the end, not how to fix a hunch.
  3. The problem is not with ggplot2, plotting the map with base graphics shows the same distortions (plot(wmap_robin)).

Upvotes: 3

Views: 747

Answers (2)

andybega
andybega

Reputation: 1437

Update for 2018, using a solution with sf:

library("cshapes")
library("sf")

cshp(as.Date("2015-01-01")) %>%
  st_as_sf() %>%
  st_crop(ymin = -90, ymax = 90, xmin=-180, xmax=180) %>%
  st_transform(crs = "+proj=robin") %>%
  `[`(1) %>% plot()

The output is:

enter image description here

Upvotes: 0

Robert Hijmans
Robert Hijmans

Reputation: 47146

You can use raster::crop to remove nodes that are just smaller than -180 or larger than 180 

library(cshapes)
library(raster)

wmap <- cshp(date=as.Date("2012-06-30"))
w <- crop(wmap, extent(-180, 180,-90,90))
w_robin <- spTransform(w, CRS("+proj=robin"))
plot(w_robin)

Upvotes: 1

Related Questions