First, load the stationaRy and oce packages:

library(stationaRy)
library(oce)
## Loading required package: gsw

We can get the information on all NCDC stations (including the lon/lat) by doing:

info <- get_ncdc_station_info()

Let’s make some maps. First, load the coastlineWorld dataset from oce, and the use mapPlot() to make a map using the default mollweide projection, and then add the station locations, as semi-transparent points:

data(coastlineWorld)
mapPlot(coastlineWorld)
mapPoints(info$LON, info$LAT, pch=19, col=rgb(1, 0, 0, 0.25), cex=0.1)

Lots of stations over the US and Europe! Another approach is to count the stations in 1 degree by 1 degree boxes:

b <- binCount2D(info$LON, info$LAT, seq(-180, 180, 1), seq(-90, 90, 1))
cm <- colormap(log10(b$number), zlim=c(0, 2), zclip=TRUE, missingColor = NA)
drawPalette(colormap = cm, zlab=expression(log[10]*group('[', stations, ']')))
mapPlot(coastlineWorld)
mapImage(b$xmids, b$ymids, log10(b$number), colormap = cm)