[GIS] How to draw a point at an exact location with spplot

rshapefile

I imported a Shapefile holding a quite detailed world map from: http://www.diva-gis.org

library(sp)
library(rgdal)
world = readOGR('./countries_shp/', 'countries')

Using R with the sp-Packages, plotting is really easy:

spplot(world, 'UNREG2')
plot(world)

Now I want to

  • draw a customizable grid
  • add points derived from coordinates

For example importing coordinates of Vienna (48°12′32″N 16°22′21″E) can be done using char2dms of the sp-package:

vienna.N = char2dms("48d12'32\"N")
vienna.E = char2dms("16d22'21\"E")

But how can a draw a point at this exact location on my map?


Also adding a grid is straight forward:

plot(gridlines(world), add=TRUE)

assuming, that the last plot was created with "plot(world)".
How can this be achieved with the spplot output?

Best Answer

I figured out how to print a point for Vienna on the map.

Importing the coordinates from WikiPedia into R:

vienna.N = char2dms("48d12'32\"N")
vienna.E = char2dms("16d22'21\"E")

Building a sp-object

vienna.coord = data.frame(
        name='Vienna', 
        x=as.numeric(vienna.E), 
        y=as.numeric(vienna.N)
)
coordinates(vienna.coord) = c('x', 'y')

The projection is taken from http://toolserver.org/~geohack/geohack.php?pagename=Vienna&params=48_12_32_N_16_22_21_E_type:city_region:AT where the WikiPedia article links to:

proj4string(vienna.coord) = '+proj=latlong +datum=WGS84'

Printing this object gives the numeric representation which equals exactly the data from mentioned website:

> vienna.coord
         coordinates   name
1 (16.3725, 48.2089) Vienna

The resulting

> class(vienna.coord)[1]
[1] "SpatialPointsDataFrame"

can be used with base- and sp-graphics. First we transform the projection to the one used in the Diva worldmap:

vienna.coord = spTransform(vienna.coord, CRS(proj4string(world)))

Then a plot for Europe can be produced:

spplot(world, 'COUNTRY', colorkey=FALSE, 
    xlim=c(-24.5, 69.5), ylim=c(26.6, 71.2),
    sp.layout = c('sp.points', vienna.coord, col='red', pch=16),
    par.settings = list(panel.background=list(col='lightblue')),
    col.regions = rep('white', length(levels(world$COUNTRY)))
)

world map with red dot for Vienna

Related Question