{resourcecode}
The goal of {resourcecode}
is to provide an easy access
to the ResourceCODE hindcast database of sea-states. More information on
the database can be found here. We will explain in this
vignette the way to retrieve data from this database and some the
functionalities offered by this package.
This package depends on data in a data package
{resourcecodedata}
that is available through a
drat
repository on GitHub. To use the
{resourcecode}
package, you will need to install
{resourcecodedata}
on your computer, using the following
install.packages
function (and later update it using the
update.packages
function):
install.packages("resourcecodedata",
repos = "https://resourcecode-project.github.io/drat/",
type = "source"
)
Once the packages are installed, we can proceed as usual:
There are two grids in the ResourceCODE hindcast database: the full, dense, FIELD grid, which is the grid where the numerical WWIII model is run and where sea-state parameters and only 1D spectra are available. The second grid, SPEC, is coarser, but contains the full 2D spectral data at each node.
The variable rscd_field
contains the coordinates of the
FIELD grid, along with the depth of these points and d50
,
the bottom sediment type.
str(rscd_field)
#> 'data.frame': 328030 obs. of 5 variables:
#> $ node : int 1 2 3 4 5 6 7 8 9 10 ...
#> $ longitude: num 6.92 6.91 6.92 6.93 6.95 ...
#> $ latitude : num 62.9 62.9 62.9 62.9 62.9 ...
#> $ depth : num 3 3 3 3 3 3 3 3 3 3 ...
#> $ d50 : num 0 1.32 1.32 1.32 1.32 ...
#> - attr(*, ".internal.selfref")=<externalptr>
head(rscd_field)
#> node longitude latitude depth d50
#> 1 1 6.922498 62.91467 3 0.000
#> 2 2 6.910178 62.91301 3 1.322
#> 3 3 6.920552 62.90438 3 1.322
#> 4 4 6.933673 62.89938 3 1.322
#> 5 5 6.945979 62.90012 3 1.322
#> 6 6 6.958059 62.89899 3 1.322
At each location, a large number of sea-state parameters are
available, described in rscd_variables
:
str(rscd_variables)
#> 'data.frame': 88 obs. of 3 variables:
#> $ name : chr "longitude" "latitude" "tri" "MAPSTA" ...
#> $ longname: chr "longitude" "latitude" "tri" "status map" ...
#> $ unit : chr "degree_east" "degree_north" "" "1" ...
head(rscd_variables)
#> name longname unit
#> 1 longitude longitude degree_east
#> 2 latitude latitude degree_north
#> 3 tri tri
#> 4 MAPSTA status map 1
#> 5 dpt depth m
#> 6 ucur eastward current m s-1
The available location can be plotted on a map using the following
code. It can be noticed here that the variables
rscd_coastline
and rscd_islands
contain the
path of the coast line and the islands, respectively.
lim_lon <- c(-5.25, -4.25)
lim_lat <- c(47.75, 48.75)
field_bzh <- ggplot(rscd_field, aes(x = longitude, y = latitude)) +
geom_point(size = .1, col = "lightblue") +
geom_path(data = rscd_coastline, linewidth = .2) +
geom_path(data = rscd_islands, aes(group = .data$ID), linewidth = .2) +
coord_sf(xlim = lim_lon, ylim = lim_lat, expand = FALSE, crs = sf::st_crs(4326)) +
theme_void()
field_bzh
The spectral grid, or SPEC grid is coarser, with more than 24,000 nodes where the full 2D spectra are available. The grid can be added to the previous plot to see the differences in the spatial coverage.
str(rscd_spectral)
#> 'data.frame': 24276 obs. of 5 variables:
#> $ longitude: num 1.5 -5.5 -6 -6.5 -7 4.5 4 3.5 3 2.5 ...
#> $ latitude : num 52 52 52 52 52 52.5 52.5 52.5 52.5 52.5 ...
#> $ name : chr "E001500N52000" "W005500N52000" "W006000N52000" "W006500N52000" ...
#> $ depth : num 15.4 89.6 95.8 52.7 54.7 ...
#> $ d50 : num 0.737 -0.8265 -1.585 0.737 0.0346 ...
#> - attr(*, ".internal.selfref")=<externalptr>
head(rscd_spectral)
#> longitude latitude name depth d50
#> 1 1.5 52.0 E001500N52000 15.41327 0.73700004
#> 2 -5.5 52.0 W005500N52000 89.64895 -0.82650555
#> 3 -6.0 52.0 W006000N52000 95.83614 -1.58500008
#> 4 -6.5 52.0 W006500N52000 52.67409 0.73700004
#> 5 -7.0 52.0 W007000N52000 54.67464 0.03464983
#> 6 4.5 52.5 E004500N52500 14.18527 0.73700004
When one is interested in some location, we also provide helpers to
find the closest point in each of the grid, via the
closest_point_FIELD()
and closest_point_SPEC()
functions.
We have two different sources of data, the FIELD grid, and the SPEC grid. We will cover each use case here.
get_parameters()
function, which allows to download a time
series easily at a given location. For the moment (Sept. 2023), the 1D
spectra on this grid are not available to direct download.For example, if one is interested in the time series of sea-state parameters at some location, the following code can be adapted:
point_of_interest <- c(longitude = -4.6861533, latitude = 48.3026514)
node <- closest_point_field(point_of_interest)
node
#> $point
#> [1] 134940
#>
#> $distance
#> [1] 288.974
ts <- get_parameters(node = node$point, parameters = c("hs", "tp", "dp", "cge"))
ggplot(tidyr::pivot_longer(ts, -1), aes(x = time, y = value, col = name)) +
geom_line() +
coord_cartesian(expand = FALSE) +
facet_wrap(~name, ncol = 2, scales = "free_y") +
scale_x_datetime(name = NULL, date_breaks = "month") +
scale_y_continuous(name = NULL) +
theme_minimal() +
theme(
legend.position = "none",
axis.text.x = element_text(angle = 60, hjust = 1)
)
get_1Dspectrum()
and
get_2Dspectrum()
. We also provide a plotting function for
the 2D spectrum.Then one can extract the 1D spectrum time-series and see the contents of these downloaded data.
spec_1d <- get_1d_spectrum(node_spectral_grid$point, start = "1994-01-01", end = "1994-02-28")
str(spec_1d)
#> List of 12
#> $ longitude : num -4.68
#> $ latitude : num 48.3
#> $ frequency1: num [1:36(1d)] 0.0339 0.0356 0.0392 0.0431 0.0474 ...
#> $ frequency2: num [1:36(1d)] 0.0356 0.0392 0.0431 0.0474 0.0521 ...
#> $ ef : num [1:36, 1:1416] 0.000526 0.0028 0.009767 0.145907 1.126401 ...
#> ..- attr(*, "dimnames")=List of 2
#> .. ..$ : NULL
#> .. ..$ : NULL
#> $ th1m : num [1:36, 1:1416] 202 203 205 226 230 ...
#> ..- attr(*, "dimnames")=List of 2
#> .. ..$ : NULL
#> .. ..$ : NULL
#> $ th2m : num [1:36, 1:1416] 202 203 205 180 180 ...
#> ..- attr(*, "dimnames")=List of 2
#> .. ..$ : NULL
#> .. ..$ : NULL
#> $ sth1m : num [1:36, 1:1416] 13.5 13.9 14.5 16.6 14.7 ...
#> ..- attr(*, "dimnames")=List of 2
#> .. ..$ : NULL
#> .. ..$ : NULL
#> $ sth2m : num [1:36, 1:1416] 26.1 26.8 28 32.2 28.7 ...
#> ..- attr(*, "dimnames")=List of 2
#> .. ..$ : NULL
#> .. ..$ : NULL
#> $ freq : num [1:36(1d)] 0.0339 0.0373 0.041 0.0451 0.0496 ...
#> $ forcings : tibble [1,416 × 14] (S3: tbl_df/tbl/data.frame)
#> ..$ time : POSIXct[1:1416], format: "1994-01-01 00:00:00" "1994-01-01 01:00:00" ...
#> ..$ dpt : num [1:1416] 27 27.5 29 30.5 31.5 32.5 32.5 32 31 29.5 ...
#> ..$ wnd : num [1:1416] 12.2 11.6 10.8 10.9 10.8 ...
#> ..$ wnddir: num [1:1416] 301 299 298 295 294 ...
#> ..$ cur : num [1:1416] 0.1 0.2 0.3 0.4 0.3 ...
#> ..$ curdir: num [1:1416] 325 27.5 41.2 48.6 55.9 ...
#> ..$ hs : num [1:1416] 3.02 3.11 3.23 3.35 3.43 ...
#> ..$ fp : num [1:1416] 0.0672 0.067 0.0671 0.0672 0.0727 ...
#> ..$ f02 : num [1:1416] 0.132 0.129 0.124 0.121 0.12 ...
#> ..$ f0m1 : num [1:1416] 0.0905 0.0882 0.0864 0.0858 0.0862 ...
#> ..$ th1p : num [1:1416] 234 234 235 237 238 ...
#> ..$ sth1p : num [1:1416] 13.4 13.3 13.1 12.9 12.7 ...
#> ..$ dir : num [1:1416] 239 239 239 240 241 ...
#> ..$ spr : num [1:1416] 19.6 19.5 18.8 18 17.4 ...
#> $ station : chr "W004679N48304"
The same applies to 2D spectra.
spec_2d <- get_2d_spectrum(node_spectral_grid$point, start = "1994-01-01", end = "1994-02-28")
str(spec_2d)
#> List of 9
#> $ longitude : num -4.68
#> $ latitude : num 48.3
#> $ frequency1: num [1:36(1d)] 0.0339 0.0356 0.0392 0.0431 0.0474 ...
#> $ frequency2: num [1:36(1d)] 0.0356 0.0392 0.0431 0.0474 0.0521 ...
#> $ efth : num [1:36, 1:36, 1:1416] 0.00024 0.000629 0.000872 0.000705 0.000423 ...
#> ..- attr(*, "dimnames")=List of 3
#> .. ..$ : NULL
#> .. ..$ : NULL
#> .. ..$ : NULL
#> $ freq : num [1:36(1d)] 0.0339 0.0373 0.041 0.0451 0.0496 ...
#> $ dir : num [1:36(1d)] 0 10 20 30 40 50 60 70 80 90 ...
#> $ forcings : tibble [1,416 × 6] (S3: tbl_df/tbl/data.frame)
#> ..$ time : POSIXct[1:1416], format: "1994-01-01 00:00:00" "1994-01-01 01:00:00" ...
#> ..$ dpt : num [1:1416] 27 27.5 29 30.5 31.5 32.5 32.5 32 31 29.5 ...
#> ..$ wnd : num [1:1416] 12.2 11.6 10.8 10.9 10.8 ...
#> ..$ wnddir: num [1:1416] 301 299 298 295 294 ...
#> ..$ cur : num [1:1416] 0.1 0.2 0.3 0.4 0.3 ...
#> ..$ curdir: num [1:1416] 145 208 221 229 236 ...
#> $ station : chr "W004679N48304"
In addition, it is possible to plot the wave elevation directional spectra for any given time, which can be specified by the time index or directly the date: