2

I would like to limit a map to a particular range of values (i.e. northern hemisphere) under an azimuthal or orthographic projection in ggplot2. When I set ylim, I receive an error message that relates to geom_path() and certain portions of the map are missing (note the triangle cut out of the United States). Any ideas on how to plot a particular region using a non-cartesian projection without losing portions of the shapefile?

world<-map_data('world') 

p<-ggplot(data=world,aes(x=long,y=lat,group=group)) +
ylim(30,90) +
geom_polygon(fill='skyblue') +
geom_path(color='gray60',aes(x=long,y=lat),size=0.1) +
coord_map("orthographic")
p

Warning message:
Removed 276 rows containing missing values (geom_path). 

example

2

What I suggest below is to generate a polygon with the area of interest and clip your polygons.

library(ggplot2)
library(sp)
library(rgeos)
library(maptools)

data(wrld_simpl) # use this since it is already a SpatialPolygons object

# make a bounding box for the northern hemisphere
g <- gridlines(wrld_simpl, easts=0, norths=30, ndiscr = 200)
x <- c(coordinates(g)[[1]][[1]][,1], rev(coordinates(g)[[1]][[1]][,1]), coordinates(g)[[1]][[1]][1])
y <- c(rep(90, 200), rep(30, 200), 90)
bnds <- cbind(x=x, y=y)
SP <- SpatialPolygons(list(Polygons(list(Polygon(bnds)), "1")), proj4string=CRS(proj4string(wrld_simpl)))

# the code snippet below is courtesy of Roger Bivand https://stat.ethz.ch/pipermail/r-sig-geo/2012-June/015340.html

gI <- gIntersects(wrld_simpl, SP, byid=TRUE)
out <- vector(mode="list", length=length(which(gI)))
ii <- 1
for (i in seq(along=gI)) if (gI[i]) {
  out[[ii]] <- gIntersection(wrld_simpl[i,], SP)
  row.names(out[[ii]]) <- row.names(wrld_simpl)[i]; ii <- ii+1
}
out1 <- do.call("rbind", out)

# fortify SP object to plot wiht ggplot
mymap <- fortify(out1, region='id') 
mymap <- mymap[order(mymap$order), ]

ggplot(data=mymap,aes(x=long,y=lat,group=group)) +
  geom_polygon(fill='skyblue') +
  geom_path(color='gray60',aes(x=long,y=lat),size=0.1) +
  coord_map("orthographic")

enter image description here

2

I'm not sure if ggplot depends on proj.4, but I solved similar problems in QGIS by defining the projection on a sphere instead an ellipsoid:

Manipulating Azimuthal Equidistant Projections in QGIS

For orthographic projections, it might be better to crop the polygons to the visible hemisphere:

Where did the polygons go after projecting a map in QGIS?

See also: http://docs.ggplot2.org/current/coord_map.html (the last remark).

And here you may find some hints on how to use proj strings: use proj4 to specify Robinson projection with R ggmap and ggplot2 packages?

Here are some more tips:

https://stackoverflow.com/questions/11490493/preventing-partial-mapping-of-coastlines

https://stackoverflow.com/questions/20143851/ggmap-is-unable-to-plot-specific-region

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.