2

When I download a large map within R (using the OpenStreetMap package) and I reproject the map (ideally to +proj=longlat), the map will always turn out distorted, making country labels etc. near the poles illegible.

Q1: Is there a way to obtain such a map (ideally still in R) in such a way that you can still read the labels? Secondly, if I keep the original map (in mercator projection) but "reproject" the points, and I use the map in a (ggplot2) plot, the labels will still look strange - probably because the size of the map image changed between downloading and plotting.

Q2: Is there a way to either download the perfectly sized tiles or to use the original tile size in the plot, so that the labels don't look blurred?

My code (for Q1, but Q2 is based on this as well) looks similar to this:

library(OpenStreetMap)
map <- openmap(c(85,-179.99999),c(-60,179.99999), zoom=2, type = "osm")
# - type="nps" is about the only type that looks okay because it has no labels
map <- openproj(map)
plot(map)

...and the map looks like this (if you look closely, you will see that "Brazil" is clearly visible, but "Russian Federation" (in Russian) not so much): enter image description here

2 Answers 2

2

If you reproject a raster with labels, you will obviously get squeezed labels.

The only way to avoid this is to render the raster from vector data directly into the desired projection. You might want to look into mapnik, tilemill or maperitive to do this from Openstreetmap raw data (which is vector data).

The R openstreetmap package only offers raster data from tiles, which is not suitable for this task. The osmar package might do a better job, but you will not get lucky with it on a worldwide scale. It's just too much data to process.

2

The following code might seem a little long-winded as it represents a manual approach rather than relying on OpenStreetMap, but maybe it's of any help to you anyway. I took the country boundaries and the referring country labels from the wrld_simpl dataset (class 'SpatialPolygonsDataFrame', projected in EPSG:4326) that comes with maptools. The shapefile data has been cropped by the extent you mentioned above.

## sample data
library(maptools)
data("wrld_simpl")

## remove antarctica
library(raster)
ext <- extent(c(-179.99999, 179.99999, # xmin, xmax
                -60, 85))              # ymin, ymax
spy_world <- crop(wrld_simpl, ext)

With regard to label placement, only countries with an area larger than 250,000 (* 1,000 hectares) are labeled at the moment, but you might wish to adjust this to your personal preferences.

## countries larger than 250,000 (* 1,000 ha)
int_id_large <- spy_world@data$AREA > 250000
spy_world@data$LABEL <- " "
spy_world@data$LABEL[int_id_large] <- as.character(spy_world@data$NAME[int_id_large])

I wrote a small code chunk that identified the largest sub-polygon per country and calculated the center coorinates of it based on rgeos::gCentroid. Otherwise, some labels would have been poorly positioned (e.g. US with Alaska, Hawaii, etc. that quite distorted the center coordinates as compared to the mainland area).

## centroid coordinates of largest polygon per country
library(rgeos)

ls_spt_cntr <- lapply(spy_world@polygons, function(i) {

  # areas of sub-polygons per country
  num_area <- sapply(i@Polygons, function(j) {
    attr(j, "area")
  })

  # polygons to spatial polygons
  py_largest <- i@Polygons[which.max(num_area)]
  py_largest <- Polygons(py_largest, ID = i@ID)

  spy_largest <- SpatialPolygons(list(py_largest), 
                                 proj4string = CRS(proj4string(spy_world)))

  # center coordinates
  gCentroid(spy_largest)
})
spy_cntr <- do.call("rbind", ls_spt_cntr)

Finally, I visualized the polygons using spplot along with lattice-powered label placement via ltext.

## visualize
library(latticeExtra)

png("worldmap.png", width = 16, height = 12, units = "cm", res = 300)
p <- spplot(spy_world, "NAME", col.regions = "cornsilk", colorkey = FALSE,
            par.settings = list(panel.background = list(col = "cadetblue")), 
            col = "grey75", scales = list(draw = TRUE))

p + 
  layer(ltext(spy_cntr@coords, labels = spy_world@data$LABEL, 
                font = "bold", cex = .8))
dev.off()

worldmap

1
  • +1 for your effort, and gCentroid. I didn't know that. I don't want to use the maptools maps though.
    – maj
    Commented Aug 15, 2015 at 9:07

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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