Scanfish from OpenEarth OPeNDAP
Willem Stolte
January 8, 2016
require("RNetCDF")
FALSE Loading required package: RNetCDF
library(rworldmap)
FALSE Loading required package: sp
FALSE ### Welcome to rworldmap ###
FALSE For a short introduction type : vignette('rworldmap')
library(ggplot2)
library(scales)
library(plyr)
Download netcdf scanfish data
Scanfish data are available via the OpenEarth OPeNDAP server at Deltares. The NetCDF file is first downloaded to a local directory before it is read into memory. Direct acces to NetCDF files using R is possible but this is not yet a standard procedure under MS-Windows.
setwd("d:/stolte/temp/")
url_scanfish <-"http://opendap.deltares.nl/thredds/fileServer/opendap/rijkswaterstaat/scanfish/scanfish_temperature.nc"
# download the netcdf file
download.file(url_scanfish, "scanfish_temperature.nc" , method = "auto",
quiet = FALSE, mode="wb", cacheOK = TRUE)
# check directory
list.files()
FALSE [1] "1" "2"
FALSE [3] "scanfish_temperature.nc" "scanfishtemperature.csv"
# open connection to nc file
con <- open.nc("scanfish_temperature.nc")
Before plotting, the data are converted into a dataframe. the RNetCDF package function read.nc is used for this in combination with the standard function data.frame.
# make into dataframe
df <- data.frame(read.nc(con))
# turn time array into something useful
df$date <- as.Date(df$TIME, origin = "1970-01-01 00:00:00", tz = "UTC")
df$datetime <- as.POSIXct(df$date)
df$months <- format(df$date, "%m")
The data are now plotted on a static map, using a background map from the rworldmap package.
# load country polygon data
data(countriesLow)
# turn into dataframe for plotting
world <- fortify(countriesLow)
FALSE Regions defined for each Polygons
# make subset of scanfish data
df2 <- subset(df, df$date >= as.Date("2011-05-01"))
# define bbox for plotting
xxlim <- expand_range(range(df2$lon), mul = 0.1)
yylim <- expand_range(range(df2$lat), mul = 0.1)
p <- ggplot(df2, aes(lon, lat))
p +
geom_polygon(data = world,
aes(x=long, y=lat, group=group),
color = "lightgrey", fill = "darkgrey") +
geom_jitter(aes(color = temperature), order = desc(df$z)) +
facet_wrap(~ months) +
coord_map(xlim = xxlim, ylim = yylim)
time plot
p <- ggplot(df2, aes(datetime, -z))
p + geom_point(aes(color = temperature))
xy plot
p <- ggplot(df2, aes(-z, temperature))
p + geom_point(aes(color = months)) +
coord_flip()