Plotting country borders with maptools - R

Question

Following the excelent post from Joe Wheatley (http://joewheatley.net/ncep-global-forecast-system/) I managed to produce a temperature global map. But, instead of only plotting coastline I've tried to use maptools package to plot country borders. The problem comes when only eastern hemisphere country borders are plotted. I should be missing something I can't figure out, still looking for on stackoverflow and google. Hope you can help.

Here is the code I'm using (most coming from Joe's post)

loc=file.path("ftp://ftp.ncep.noaa.gov/pub/data/nccf/com/gfs/prod/gfs.2013052100/gfs.t00z.sfluxgrbf03.grib2")
download.file(loc,"temp.grb",mode="wb")

system("wgrib2 -s temp.grb | grep :TMP: | wgrib2 -i temp.grb -netcdf TMP.nc",intern=T)
system("wgrib2 -s temp.grb | grep :LAND: | wgrib2 -i temp.grb -netcdf temp.nc",intern=T)

library(ncdf)
landFrac <-open.ncdf("LAND.nc")
lon <- get.var.ncdf(landFrac,"longitude")
lat <- get.var.ncdf(landFrac,"latitude")

temp=open.ncdf("TMP.nc")
t2m.mean <- get.var.ncdf(temp,"TMP_2maboveground")


library("fields")
library("sp", lib.loc="/usr/lib/R/site-library")
library("maptools", lib.loc="/usr/lib/R/site-library")

day="DIA"

png(filename="gfs.png",width=1215,height=607,bg="white")

rgb.palette <- colorRampPalette(c("snow1","snow2","snow3","seagreen","orange","firebrick"), space = "rgb")#colors
image.plot(lon,lat,t2m.mean,col=rgb.palette(200),main=as.expression(paste("GFS 24hr Average 2M Temperature",day,"00 UTC",sep="")),axes=T,legend.lab="o C")
data(wrld_simpl)
plot(wrld_simpl, add = TRUE)

dev.off()

and this is the image produced

This is a global map, should I use xlim and ylim in image.plot to extract a region (i.e. Europe)

EDIT: Added url for temp.nc file

http://ubuntuone.com/29DKAeRjUCiCzLblgfSLc9

Any help would be appreciated, thanks

Answer 1

The data set wrld_simpl is aligned on longitudes [-180, 180] but your grid data is clearly [0,360]. Since you have maptools loaded try this to add a modified copy to your plot:

plot(elide(wrld_simpl, shift = c(360, 0)), add = TRUE)

There are other tools to crop and move/combine and recentre data like this, including ?rotate in the package raster. It really depends on what you need.

Another graphics-alternative is to just use "world2" in the maps package

library(maps)
map("world2", add = TRUE)

Either of those will work to finish your plot started above.

Here's another discussion related to this: Fixing maps library data for Pacific centred (0°-360° longitude) display

Answer 2

You can do this:

library(raster)
r <- raster("temp.nc")
r <- rotate(r)
plot(r)