Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Plot Map of Pacific with filled countries

Tags:

r

map

I'm trying to plot a map of the Pacific using World2Hires in R's mapproj library but there's an odd glitch when I try to fill the countries. How do I fix this?

library(maps)
library(mapproj)
library(mapdata)
map("world2Hires", 
    xlim=c(120, 260), 
    ylim=c(-60, 40), 
    boundary=TRUE, 
    interior=TRUE,
    fill=TRUE,
    col="gray30",
)
map.axes()

Here's the output:

Broken Map Image

like image 294
Puzzled79 Avatar asked Sep 25 '13 00:09

Puzzled79


2 Answers

The issue seems to be with a small subset of areas which cause wrapping. From some trial and error saving the original map call like mapnames <- map(...) and then passing subsets of this list to the regions= argument in a new call, I could avoid the wrapping around of the fills. E.g.:

library(maps)
library(mapproj)
library(mapdata)
map("world2Hires", regions=mapnames$names[c(1:7,14:641)],
    xlim=c(120, 260), 
    ylim=c(-60, 40), 
    boundary=TRUE, 
    interior=TRUE,
    fill=TRUE
)
map.axes()

enter image description here

As to a more thorough or sensible solution to prevent this happening, I am stumped. Playing with the wrap= option does nothing helpful, and likewise for the other options. As a side-note, this issue does not appear using the "world" database, but only pops up for "world2" and "world2Hires".

like image 75
thelatemail Avatar answered Oct 24 '22 02:10

thelatemail


The answer by @thelatemail is the best and simplest solution I've seen to this problem. To make it more universal, though, it is better to remove the polygons by name. This is because depending on the limits you give your first call to map(), the indices of the polygon names can be different.

library(maps)
library(mapproj)
library(mapdata)

mapnames <- map("world2Hires", xlim=c(120, 260), ylim=c(-60, 40),
                fill=TRUE, plot=FALSE)

mapnames2 <- map("world2Hires", xlim=c(100, 200), ylim=c(-20, 60),
                 fill=TRUE, plot=FALSE)

mapnames$names[10]
[1] "Mali"
mapnames2$names[10]
[1] "Thailand"

There are 8 countries that are intersected by the Prime Meridian: United Kingdom, France, Spain, Algeria, Mali, Burkina Faso, Ghana, and Togo. By matching these country names with mapnames$names, you can remove the polygons regardless of your original extent:

remove <- c("UK:Great Britain", "France", "Spain", "Algeria", "Mali",
            "Burkina Faso", "Ghana", "Togo")

map("world2Hires", regions=mapnames$names[!(mapnames$names %in% remove)],
    xlim=c(120, 260), 
    ylim=c(-60, 40), 
    boundary=TRUE, 
    interior=TRUE,
    fill=TRUE
)
map.axes()

You could also use grepl() but because polygons are named heirarchically, you may remove some sub-polygons of the nations in question. For example, mapnames$names[grepl("UK", mapnames$names)] returns 34 matches.

I would have suggested this as an edit, but I don't have privileges yet.

like image 44
Cotton.Rockwood Avatar answered Oct 24 '22 00:10

Cotton.Rockwood