explore_field_2d: Explore 2-dimensional FVCOM fields

View source: R/explore_field_2d.R

explore_field_2dR Documentation

Explore 2-dimensional FVCOM fields

Description

For multiple variables, days and locations, this function extracts model predictions and computes summary statistics and/or creates maps of environmental conditions. The user passes the function a dataframe specifying the variables to be included and their properties alongside various other arguments for data extraction and to customise the summary statistics and plots that are produced. For each variable, the function loads in the FVCOM arrays for each date in turn, computes summary statistics and/or plots a map for the hours specified. Iterative loading in of the files is necessary given the size/memory requirements of FVCOM arrays.

Usage

explore_field_2d(
  field2d,
  dat_ls,
  match_hour = data.frame(hour = 0:23, index = 1:24),
  match_layer = NULL,
  match_mesh_around_nodes = NULL,
  match_mesh_around_elements = NULL,
  corrupt = NULL,
  read_fvcom_ls = list(),
  dir2load_ls = list(),
  mesh_around_nodes = NULL,
  mesh_around_elements = NULL,
  make_plot = TRUE,
  plot_param = list(),
  png_param = list(),
  compute_summary_stats = TRUE,
  summary_stats_param = list(),
  parallelise = "vars",
  cl = NULL,
  pass2varlist = NULL,
  verbose = TRUE,
  ...
)

Arguments

field2d

A dataframe that defines the variables for which model outputs will be explored. This must contain the following columns: 'cov2d', a character vector of variable names; 'mesh_type', a character vector of the mesh type ("element" or "node") across which each variable is structured; 'extension', a character vector which defines the extension of the FVCOM arrays for that variable (see extract); and 'vector_field', a logical vector defining whether or not that variable is a scalar (FALSE) or vector field (TRUE). To save plots directly to file (see make_plot and png_param, below), the dataframe should contain a column 'dir2save' which specifies the folder in which all the plots for each variable will be saved.

dat_ls

A named list of dataframes, with one element for each field2d$cov2d. Each element should comprise a dataframe which defines the FVCOM file date names, hours, and mesh (node or element) IDs for model predictions will be extracted (see extract).

match_hour

A dataframe with two integer columns named 'hour' and 'index' which defines the index in FVCOM arrays (i.e. the row) which corresponds to each hour (see extract).

match_layer

A dataframe with two integer columns named 'layer' and 'index' which defines the index in FVCOM arrays (i.e. the column) which corresponds to each layer (see extract).

match_mesh_around_nodes

A dataframe with two columns named 'mesh' and 'index' which defines the index in FVCOM arrays (columns or sheets for 2- and 3-dimensional arrays respectively) which corresponds to each node cell (see extract). This may be required if field2d contains variables that are resolved at nodes.

match_mesh_around_elements

A dataframe with two columns named 'mesh' and 'index' which defines the index in FVCOM arrays (columns or sheets for 2- and 3-dimensional arrays respectively) which corresponds to each element cell (see extract). This is may be required if field2d contains variables that are resolved at elements.

corrupt

A vector of numbers, representing WeStCOMS date names, which define corrupt files (see extract).

read_fvcom_ls

A named list of functions, with one element for each field2d$cov2d. Each element should contain a function that is used to load files for that variable (see extract).

dir2load_ls

A named list of directories, with one element for each field2d$cov2d. Each element should be a string (for scalar fields) or a vector of strings (for vector fields) that defines the directory from which to load model files (see extract, which is used to load files). For vector fields, the first string should specify the directory of the u component files and the second string should specify the directory of the v component files.

mesh_around_nodes

A mesh, created by build_mesh, that surrounds nodes. This is required for plotting variables that are resolved at nodes.

mesh_around_elements

A mesh, created by build_mesh, that surrounds elements. This is required for plotting variables that are resolved at elements.

make_plot

A logical input defining whether or not to make plots. Plots are either displayed or saved to file (if field2d contains a column named 'dir2save').

plot_param

A list of parameters required to make plots. These parameters either wrap around, or are passed as arguments, to plot_field_2d. The list needs to contain elements with the following names: 'hours4plots', 'par_op', 'coastline', 'zlab', 'zlab_line' and 'vector_scale'. 'hours4plots' is a numeric vector of all the hours for which to create plots on a given day. 'par_op' is an output from par. 'coastline' is an object used to plot coastline (see plot_field_2d); and 'zlab', 'zlab_line' and 'vector_scale' are vectors that define the labels on the z axis, their distance from the z axis and the scale of the arrows used to plot vectors, for each inputted variable. Other graphical parameters are not variable specific and passed as additional arguments outside of this list (see ...).

png_param

A named list of parameters passed to png to customise plots saved as .png files.

compute_summary_stats

A logical input defining whether or not summary statistics should be calculated.

summary_stats_param

A list of parameters necessary to calculate summary statistics. This list should contain elements with the following names: 'hours4stats', 'row_specific' and 'funs' (see summarise_field_2d).

parallelise

A character specifying whether or not to parallelise the algorithm over variables ("vars") or dates ("date_name"). This is only applicable if a cluster is supplied (see below).

cl

(optional) A cluster objected created by the parallel package. If supplied, the algorithm is implemented in parallel. Note that the connection with the cluster is stopped within the function.

pass2varlist

A list containing the names of exported objects. This may be required if cl is supplied. This is passed to the varlist argument of clusterExport. Exported objects must be located in the global environment.

verbose

A logical input which defines whether or not to display messages to the console detailing function progress.

...

Additional graphical parameters that can be passed to plot_field_2d. These affect all plots.

Value

If make_plot = TRUE, the function will produce plots, which are either saved to file or displayed (the latter is only possible if cl = NULL). To save plots, field2d must contain a column named 'dir2save'. Plots are saved in this directory with pre-defined file names pertaining to the date and hour that they represent. The list png_param can be used to customise the saved images. If compute_summary_stats = TRUE, the function will also return a list of dataframes, with one element for each environmental variable. Each dataframe the following columns: a date, hour and a column for each summary statistic specified. This list is returned to the environment.

Author(s)

Edward Lavender

See Also

build_mesh, extract, summarise_field_2d, plot_field_2d

Examples

############################
############################
#### Example (1): Implement explore_field_2d() with a single variable

#### Define field2d, a dataframe which defines variables and details
# ... needed for extraction.
field2d <- data.frame(cov2d = "temp",
                      mesh_type = "element",
                      extension = ".mat",
                      vector_field = FALSE)

#### Define dat_ls, a list of dataframes which defines
# ... the specific FVCOM outputs we'll extract.
# We'll extract all predictions for two days from example files for the top layer:
timestamp <- as.POSIXct(c("2016-03-01", "2016-03-02"), tz = "UTC")
dat <- expand.grid(timestamp = timestamp, mesh_ID = dat_nodexy$node_id)
dat$date_name <- date_name(dat$timestamp)
dat$hour <- lubridate::hour(dat$timestamp)
dat$layer <- 1
dat <- dat[, c("date_name", "hour", "layer", "mesh_ID")]
dat_ls <- list(temp = dat)

#### Define other general parameters for extract()
match_hour <- data.frame(hour = 0:1, index = 1:2)
# match_layer <- NULL # not needed in this example
match_mesh <-
  data.frame(mesh = dat_nodexy$node_id, index = 1:length(dat_nodexy$node_id))
match_mesh_around_nodes <-
  data.frame(mesh = dat_nodexy$node_id, index = 1:length(dat_nodexy$node_id))
match_mesh_around_elements <- NULL # not needed in this example
# corrupt <- NULL # not needed in this example

#### Define a NAMED list of functions to load files
read_fvcom_ls <- list(temp = function(con) R.matlab::readMat(con)$data)

#### Define a NAMED list of directories used to load files
root <- system.file("WeStCOMS_files",
                     package = "fvcom.tbx", mustWork = TRUE)
root <- paste0(root, "/")
dir2load_ls <- list(temp = paste0(root, "temp/"))

#### Define mesh(es)
# We'll use the example data files.
mesh_around_nodes <- dat_mesh_around_nodes
# mesh_around_elements <- NULL

#### Define plot parameters
make_plot <- TRUE
plot_param = list(hours4plots = unique(dat$hour)[1],
                  par_op = par(oma = c(3, 3, 3, 7)),
                  coastline = dat_coast_around_oban,
                  zlab = expression(paste("Temperature (", degree, ")")),
                  zlab_line = 3,
                  vector_scale = NA
)
# png_param <- list() # not needed in this example

#### Define compute stats parameters
compute_summary_stats <- TRUE
summary_stats_param <-
  list(row_specific = TRUE,
       funs = list(mean = function(x) mean(x, na.rm = TRUE),
                   min = function(x) min(x, na.rm = TRUE))
  )

#### Parallelise param
# We'll set inside the function.
# cl <- NULL
# pass2varlist = list()
# parallelise <- "vars"
# verbose <- TRUE

#### Implement explore_field_2d() for a single variable
explore_field_2d(
  field2d = field2d,
  dat_ls = dat_ls,
  match_hour = data.frame(hour = 0:23, index = 1:24),
  match_layer = NULL,
  match_mesh_around_nodes = match_mesh_around_nodes,
  match_mesh_around_elements = NULL,
  corrupt = NULL,
  read_fvcom_ls = read_fvcom_ls,
  dir2load_ls = dir2load_ls,
  mesh_around_nodes = dat_mesh_around_nodes,
  mesh_around_elements = NULL,
  make_plot = TRUE,
  plot_param = plot_param,
  png_param = list(),
  compute_summary_stats = TRUE,
  summary_stats_param = summary_stats_param,
  parallelise = "vars",
  cl = NULL,
  pass2varlist = NULL,
  verbose = TRUE,
  ncols = 50
)


############################
############################
#### Example (2): Implement explore_field_2d() for multiple variables

#### Redefine field2d, a dataframe which defines variables and details needed for extraction.
field2d <- data.frame(cov2d = c("temp", "tidal_elevation", "wind_velocity"),
                      mesh_type = c("element", "element", "node"),
                      extension = c(".mat", ".mat", ".mat"),
                      vector_field = c(FALSE, FALSE, TRUE)
)

#### Redefine dat_ls to include appropriate mesh cells
dat_element <- expand.grid(timestamp = timestamp, mesh_ID = dat_trinodes$element_id)
dat_element$date_name <- date_name(dat_element$timestamp)
dat_element$hour <- lubridate::hour(dat_element$timestamp)
dat_ls <- list(temp = dat,
               tidal_elevation = dat,
               wind_velocity = dat_element)

#### Update dir2load_ls and read_fvcom_ls arguments
dir2load_ls <- list(temp = paste0(root, "temp/"),
                    tidal_elevation = paste0(root, "tidal_elevation/"),
                    wind_velocity = c(paste0(root, "uwind_speed/"), paste0(root, "vwind_speed/")))
read_fvcom_ls <- list(temp = function(con) R.matlab::readMat(con)$data,
                      tidal_elevation = function(con) R.matlab::readMat(con)$data,
                      wind_velocity = function(con) R.matlab::readMat(con)$data)

#### Additional arguments required for vector field
match_mesh_around_elements <-
  data.frame(mesh = dat_trinodes$element_id,
             index = 1:length(dat_trinodes$element_id))

#### Update plot_param
plot_param = list(hours4plots = unique(dat$hour)[1],
                  par_op = par(oma = c(3, 3, 3, 7)),
                  coastline = dat_coast_around_oban,
                  zlab = c(expression(paste("Temperature (", degree, ")")),
                           "Tidal elevation (m)",
                           expression(paste("Wind Velocity (m", s^-1, ")"))
                  ),
                  zlab_line = c(3, 3, 3),
                  vector_scale = c(NA, NA, 20)
)

explore_field_2d(
  field2d = field2d,
  dat_ls = dat_ls,
  match_hour = data.frame(hour = 0:23, index = 1:24),
  match_layer = NULL,
  match_mesh_around_nodes = match_mesh_around_nodes,
  match_mesh_around_elements = match_mesh_around_elements,
  corrupt = NULL,
  read_fvcom_ls = read_fvcom_ls,
  dir2load_ls = dir2load_ls,
  mesh_around_nodes = dat_mesh_around_nodes,
  mesh_around_elements = dat_mesh_around_elements,
  make_plot = TRUE,
  plot_param = plot_param,
  png_param = list(),
  compute_summary_stats = TRUE,
  summary_stats_param = summary_stats_param,
  cl = NULL,
  verbose = TRUE,
  ncols = 50)


############################
############################
#### Example (3): Parallelise Example (2) over variables:

#### These examples take a bit longer to run (a few seconds)
## Not run: 

#### Define 'dir2save' to save figures
# ... which cannot be displayed in parallel.
field2d$dir2save <- paste0(tempdir(), "/")

#### Implement algorithm
explore_field_2d(
  field2d = field2d,
  dat_ls = dat_ls,
  match_hour = data.frame(hour = 0:23, index = 1:24),
  match_layer = NULL,
  match_mesh_around_nodes = match_mesh_around_nodes,
  match_mesh_around_elements = match_mesh_around_elements,
  corrupt = NULL,
  read_fvcom_ls = read_fvcom_ls,
  dir2load_ls = dir2load_ls,
  mesh_around_nodes = dat_mesh_around_nodes,
  mesh_around_elements = dat_mesh_around_elements,
  make_plot = TRUE,
  plot_param = plot_param,
  png_param = list(),
  compute_summary_stats = TRUE,
  summary_stats_param = summary_stats_param,
  verbose = TRUE,
  # Add parallelisation arguments:
  parallelise = "vars",
  cl = parallel::makeCluster(2L),
  pass2varlist = list("dat_mesh_around_nodes",
                      "dat_mesh_around_elements",
                      "match_mesh_around_nodes",
                      "match_mesh_around_elements")
)

#### Examine .png files saved by function
list.files(field2d$dir2save, "*.png")

## End(Not run)


############################
############################
#### Example (4): Parallelise example (3) over dates:

## Not run: 
explore_field_2d(
  field2d = field2d,
  dat_ls = dat_ls,
  match_hour = data.frame(hour = 0:23, index = 1:24),
  match_layer = NULL,
  match_mesh_around_nodes = match_mesh_around_nodes,
  match_mesh_around_elements = match_mesh_around_elements,
  corrupt = NULL,
  read_fvcom_ls = read_fvcom_ls,
  dir2load_ls = dir2load_ls,
  mesh_around_nodes = dat_mesh_around_nodes,
  mesh_around_elements = dat_mesh_around_elements,
  make_plot = TRUE,
  plot_param = plot_param,
  png_param = list(),
  compute_summary_stats = TRUE,
  summary_stats_param = summary_stats_param,
  verbose = TRUE,
  # Add parallelisation arguments:
  parallelise = "date_name",
  cl = parallel::makeCluster(2L),
  pass2varlist = list()
)

## End(Not run)


edwardlavender/fvcom.tbx documentation built on Nov. 26, 2022, 10:28 p.m.