diff --git a/Makefile b/Makefile index fa508f1c17f..f90819a4fdc 100644 --- a/Makefile +++ b/Makefile @@ -4,7 +4,7 @@ NCPUS ?= 1 BASE := logger utils db settings visualization qaqc remote workflow MODELS := basgra biocro clm45 dalec dvmdostem ed fates gday jules linkages \ - ldndc lpjguess maat maespa sibcasa sipnet stics template + ldndc lpjguess maat maespa rothc sibcasa sipnet stics template MODULES := allometry assim.batch assim.sequential benchmark \ data.atmosphere data.land data.remote \ diff --git a/Makefile.depends b/Makefile.depends index f1a51712e80..9f7518118a3 100644 --- a/Makefile.depends +++ b/Makefile.depends @@ -21,6 +21,7 @@ $(call depends,models/linkages): | .install/base/db .install/base/logger .instal $(call depends,models/lpjguess): | .install/base/logger .install/base/remote .install/base/utils $(call depends,models/maat): | .install/base/logger .install/base/remote .install/base/settings .install/base/utils .install/modules/data.atmosphere $(call depends,models/maespa): | .install/base/logger .install/base/remote .install/base/utils .install/modules/data.atmosphere +$(call depends,models/rothc): | .install/base/logger .install/base/utils $(call depends,models/sibcasa): | .install/base/logger $(call depends,models/sipnet): | .install/base/logger .install/base/remote .install/base/utils .install/modules/data.atmosphere .install/modules/data.land $(call depends,models/stics): | .install/base/logger .install/base/remote .install/base/settings .install/base/utils diff --git a/docker/depends/pecan_package_dependencies.csv b/docker/depends/pecan_package_dependencies.csv index 96dc8b7399e..28b4c59238c 100644 --- a/docker/depends/pecan_package_dependencies.csv +++ b/docker/depends/pecan_package_dependencies.csv @@ -168,6 +168,7 @@ "lubridate",">= 1.6.0","models/lpjguess","Imports",FALSE "lubridate",">= 1.6.0","models/maat","Imports",FALSE "lubridate",">= 1.6.0","models/maespa","Imports",FALSE +"lubridate",">= 1.6.0","models/rothc","Imports",FALSE "lubridate",">= 1.6.0","models/sipnet","Imports",FALSE "lubridate",">= 1.6.0","modules/assim.batch","Imports",FALSE "lubridate",">= 1.6.0","modules/assim.sequential","Imports",FALSE @@ -231,6 +232,7 @@ "ncdf4","*","models/basgra","Imports",FALSE "ncdf4","*","models/dvmdostem","Imports",FALSE "ncdf4","*","models/ldndc","Imports",FALSE +"ncdf4","*","models/rothc","Imports",FALSE "ncdf4","*","models/sibcasa","Imports",FALSE "ncdf4","*","models/stics","Imports",FALSE "ncdf4","*","modules/assim.sequential","Imports",FALSE @@ -338,6 +340,7 @@ "PEcAn.logger","*","models/lpjguess","Imports",TRUE "PEcAn.logger","*","models/maat","Imports",TRUE "PEcAn.logger","*","models/maespa","Imports",TRUE +"PEcAn.logger","*","models/rothc","Imports",TRUE "PEcAn.logger","*","models/sibcasa","Imports",TRUE "PEcAn.logger","*","models/sipnet","Imports",TRUE "PEcAn.logger","*","models/stics","Imports",TRUE @@ -429,6 +432,7 @@ "PEcAn.utils",">= 1.4.8","models/basgra","Imports",TRUE "PEcAn.utils",">= 1.4.8","models/dvmdostem","Imports",TRUE "PEcAn.utils",">= 1.4.8","models/ldndc","Imports",TRUE +"PEcAn.utils",">= 1.4.8","models/rothc","Imports",TRUE "PEcAn.utils",">= 1.4.8","models/stics","Imports",TRUE "PEcAn.utils",">= 1.4.8","models/template","Imports",TRUE "PEcAn.visualization","*","modules/assim.sequential","Suggests",TRUE @@ -536,6 +540,7 @@ "roxygen2","== 7.3.2","models/lpjguess","Roxygen",FALSE "roxygen2","== 7.3.2","models/maat","Roxygen",FALSE "roxygen2","== 7.3.2","models/maespa","Roxygen",FALSE +"roxygen2","== 7.3.2","models/rothc","Roxygen",FALSE "roxygen2","== 7.3.2","models/sibcasa","Roxygen",FALSE "roxygen2","== 7.3.2","models/sipnet","Roxygen",FALSE "roxygen2","== 7.3.2","models/stics","Roxygen",FALSE @@ -571,6 +576,7 @@ "sp","*","modules/data.land","Imports",FALSE "sp","*","modules/data.remote","Imports",FALSE "stats","*","base/qaqc","Imports",FALSE +"stats","*","models/rothc","Imports",FALSE "stats","*","models/sipnet","Imports",FALSE "stats","*","modules/allometry","Imports",FALSE "stats","*","modules/assim.batch","Imports",FALSE @@ -630,6 +636,7 @@ "testthat",">= 2.0.0","base/utils","Suggests",FALSE "testthat",">= 2.0.0","models/biocro","Suggests",FALSE "testthat",">= 2.0.0","modules/benchmark","Suggests",FALSE +"testthat",">= 3.0.0","models/rothc","Suggests",FALSE "testthat",">= 3.0.0","models/sibcasa","Suggests",FALSE "testthat",">= 3.0.4","base/qaqc","Suggests",FALSE "testthat",">= 3.1.7","modules/data.atmosphere","Suggests",FALSE @@ -664,6 +671,7 @@ "utils","*","models/ed","Imports",FALSE "utils","*","models/linkages","Imports",FALSE "utils","*","models/lpjguess","Imports",FALSE +"utils","*","models/rothc","Imports",FALSE "utils","*","modules/allometry","Imports",FALSE "utils","*","modules/assim.batch","Imports",FALSE "utils","*","modules/assim.sequential","Suggests",FALSE @@ -682,6 +690,7 @@ "withr","*","base/workflow","Suggests",FALSE "withr","*","models/basgra","Suggests",FALSE "withr","*","models/ed","Suggests",FALSE +"withr","*","models/rothc","Suggests",FALSE "withr","*","models/sibcasa","Suggests",FALSE "withr","*","models/sipnet","Suggests",FALSE "withr","*","modules/allometry","Suggests",FALSE diff --git a/models/rothc/.Rbuildignore b/models/rothc/.Rbuildignore new file mode 100644 index 00000000000..2d28facae40 --- /dev/null +++ b/models/rothc/.Rbuildignore @@ -0,0 +1,2 @@ +Dockerfile +model_info.json diff --git a/models/rothc/DESCRIPTION b/models/rothc/DESCRIPTION new file mode 100644 index 00000000000..ba38853860a --- /dev/null +++ b/models/rothc/DESCRIPTION @@ -0,0 +1,30 @@ +Package: PEcAn.RothC +Type: Package +Title: PEcAn Package for Integration of the RothC Model +Version: 0.0.0.9000 +Authors@R: c( + person("Chris", "Black", role = c("aut", "cre"), email = "chris@ckblack.org"), + person("Rothamsted Research", role = c("cph")) + ) +Description: This module provides functions to link the Rothamstead soil carbon model "RothC" to PEcAn. + Uses RothC 1.0.0, made available by Rothamstead Research under the Apache License 2.0. TODO: Determine whether we are including a copy of the RothC sources in this package or connecting to externally installed versions, get copyright notice and license requirements met. +Depends: R (>= 4.1) +Imports: + dplyr, + lubridate (>= 1.6.0), + ncdf4, + PEcAn.logger, + PEcAn.utils (>= 1.4.8), + rlang, + stats, + utils +Suggests: + testthat (>= 3.0.0), + withr +SystemRequirements: RothC_Py +OS_type: unix +License: BSD_3_clause + file LICENSE +Copyright: Authors +Encoding: UTF-8 +RoxygenNote: 7.3.2 +Config/testthat/edition: 3 diff --git a/models/rothc/Dockerfile b/models/rothc/Dockerfile new file mode 100644 index 00000000000..a9b5af2dd24 --- /dev/null +++ b/models/rothc/Dockerfile @@ -0,0 +1,40 @@ +# this needs to be at the top, what version are we building +ARG IMAGE_VERSION="latest" + + +# ---------------------------------------------------------------------- +# BUILD PECAN FOR MODEL +# ---------------------------------------------------------------------- +FROM pecan/models:${IMAGE_VERSION} + +# ---------------------------------------------------------------------- +# INSTALL MODEL SPECIFIC PIECES +# ---------------------------------------------------------------------- + +#RUN apt-get update \ +# && apt-get install -y --no-install-recommends \ +# python3.9 \ +# && rm -rf /var/lib/apt/lists/* + +# ---------------------------------------------------------------------- +# SETUP FOR SPECIFIC MODEL +# ---------------------------------------------------------------------- + +# Some variables that can be used to set control the docker build +ARG MODEL_VERSION=2.1.0 + +RUN mkdir -p /tmp/rothc \ + && curl -sSL https://github.com/Rothamsted-Models/RothC_Code/archive/refs/tags/v${MODEL_VERSION}.tar.gz \ + | tar xzf - -C /tmp/rothc \ + && cd /tmp/rothc/RothC_Code-${MODEL_VERSION} \ + && gfortran RothC.for Shell.for -o rothc_bin \ + && mv rothc_bin /usr/local/bin/RothC_v${MODEL_VERSION} \ + && rm -rf /tmp/rothc + +# TODO hard-coding this path is probably wrong +# Setup model_info file +# @VERSION@ is replaced with model version in the model_info.json file +# @BINARY@ is replaced with model binary in the model_info.json file +COPY model_info.json /work/model.json +RUN sed -i -e "s/@VERSION@/${MODEL_VERSION}/g" \ + -e "s#@BINARY@#/usr/local/bin/RothC_v${MODEL_VERSION}#g" /work/model.json diff --git a/models/rothc/LICENSE b/models/rothc/LICENSE new file mode 100644 index 00000000000..09ef35a60b4 --- /dev/null +++ b/models/rothc/LICENSE @@ -0,0 +1,3 @@ +YEAR: 2024 +COPYRIGHT HOLDER: PEcAn Project +ORGANIZATION: PEcAn Project, authors affiliations diff --git a/models/rothc/NAMESPACE b/models/rothc/NAMESPACE new file mode 100644 index 00000000000..bc4c352e18b --- /dev/null +++ b/models/rothc/NAMESPACE @@ -0,0 +1,9 @@ +# Generated by roxygen2: do not edit by hand + +export(met2model.RothC) +export(model2netcdf.RothC) +export(read_restart.ModelName) +export(write.config.RothC) +export(write_restart.ModelName) +importFrom(rlang,.data) +importFrom(rlang,.env) diff --git a/models/rothc/NEWS.md b/models/rothc/NEWS.md new file mode 100644 index 00000000000..cc686419188 --- /dev/null +++ b/models/rothc/NEWS.md @@ -0,0 +1,3 @@ +# PEcAn.RothC 0.0.0.9000 + +Initial development version \ No newline at end of file diff --git a/models/rothc/R/met2model.RothC.R b/models/rothc/R/met2model.RothC.R new file mode 100644 index 00000000000..5c9745456ea --- /dev/null +++ b/models/rothc/R/met2model.RothC.R @@ -0,0 +1,180 @@ +#' Extract monthly weather from CF file for input to RothC +#' +#' Input files need to be named `/.YYYY.nc` +#' +#' Output files are named `/.YY-mm.YY-mm.dat` +#' with one line per month and columns for temperature, rainfall, +#' and evaporation. +#' +#' Note that the created file contains only weather data and not any of the +#' soil or management data needed for RothC's single combined input file. +#' See `write.config.RothC()` for assembly into a model-ready RothC_input.dat`. +#' +#' @param in.path path on disk where CF files live +#' @param in.prefix prefix for each file +#' @param outfolder location where model specific output is written. +#' @param start_date,end_date When to start and end output. +#' Specify as exact dates, but output will be padded to whole months. +#' @param overwrite logical: replace output files if they already exist? +#' @return data frame summarizing file metadata +#' @export +#' @author Chris Black +met2model.RothC <- function(in.path, + in.prefix, + outfolder, + start_date, + end_date, + overwrite = FALSE) { + + PEcAn.logger::logger.info("START met2model.RothC") + + start_date <- as.Date(start_date) + end_date <- as.Date(end_date) + start_year <- strftime(start_date, "%Y") + end_year <- strftime(end_date, "%Y") + year_regex <- paste(start_year:end_year, collapse = "|") + + if (grepl("\\.nc$", in.prefix)) { + # Assume it's the full filename rather than a prefix + # NB also means we assume it contains the whole requested date range + name_pattern <- in.prefix + } else { + name_pattern <- paste0(in.prefix, "\\.(", year_regex, ")\\.nc$") + } + + nc_files <- list.files(in.path, pattern = name_pattern, full.names = TRUE) + + if (length(nc_files) == 0) { + PEcAn.logger::logger.severe( + "No files found matching ", in.prefix, + "for years", start_year, ":", end_year, + "; cannot process data." + ) + } + + # TODO complain (fail?) here if some but not all years found + + out_filename <- paste( + in.prefix, + strftime(start_date, "%Y-%m"), + strftime(end_date, "%Y-%m"), + "dat", + sep = "." + ) + out_path <- file.path(outfolder, out_filename) + results <- data.frame(file = out_path, + host = Sys.getenv( + "FQDN", + unset = Sys.info()[["nodename"]] + ), + mimetype = "text/tab-separated-values", + formatname = "RothC.dat", + startdate = start_date, + enddate = end_date, + dbfile.name = out_filename, + stringsAsFactors = FALSE) + PEcAn.logger::logger.info("internal results") + PEcAn.logger::logger.info(results) + + if (file.exists(out_path) && !overwrite) { + PEcAn.logger::logger.debug( + "File '", out_path, "' already exists, skipping to next file." + ) + return(invisible(results)) + } + + if (!file.exists(outfolder)) { + dir.create(outfolder) + } + + met <- nc_files |> + lapply( + read_nc, + varnames = c("air_temperature", "precipitation_flux", "specific_humidity") + ) |> + do.call(what = "rbind") + + # TODO probably need more care with partial months here: + # we check if data extends to start/end, but not whether that includes enough + # days to treat as a whole month. + # e.g. if start_date = YYYY-01-31 and data starts YYYY-01-30 -> + # current code will aggregate those two days as if they were all of January. + # Consider failing if >n days missing in any output month? + first_month <- lubridate::floor_date(start_date, unit = "month") + last_month <- lubridate::ceiling_date(end_date, unit = "month") + met <- met[(met$timestamp >= first_month) & (met$timestamp < last_month), ] + if (as.Date(min(met$timestamp)) > start_date + || as.Date(max(met$timestamp)) < end_date) { + PEcAn.logger::logger.severe( + "input (", + paste(range(met$timestamp), collapse = " to "), + ") does not cover requested time window (", + start_date, "to", end_date, ")" + ) + } + + met$year <- lubridate::year(met$timestamp) + met$month <- lubridate::month(met$timestamp) + met$Tmp <- met$air_temperature |> + PEcAn.utils::ud_convert("K", "degC") + met$Rain <- 0# TODO... sum up to convert from flux to accumulation, right? + met$Evap <- 0# TODO... how to convert Qair to pan evaporation? + + met_monthly <- merge( + stats::aggregate(met, Tmp ~ year + month, mean), + stats::aggregate(met, cbind(Rain, Evap) ~ year + month, sum), + sort = FALSE # would treat months as strings; sort as numbers below instead + ) + met_monthly <- met_monthly[order(met_monthly$year, met_monthly$month), ] + + utils::write.table( + # as.data.frame to write integer columns as ints not floats + x = format(as.data.frame(met_monthly), digits = 4), + file = out_path, + quote = FALSE, + sep = "\t", + row.names = FALSE, + col.names = TRUE + ) + + results +} + + + + + +# slurp named vars from one PEcAn nc into a dataframe with timestamp +# +# TODO could read other dimensions if present too, but consider if worth it -- +# maybe this function is better for files where only the time dimension varies +# +# if vars = NULL, read all of them +read_nc <- function(ncfile, varnames = NULL) { + + nc <- ncdf4::nc_open(ncfile) + on.exit(ncdf4::nc_close(nc), add = TRUE) + + timestamps <- PEcAn.utils::cf2datetime( + nc$dim$time$vals, + nc$dim$time$units + ) + + if (is.null(varnames)) { + varnames <- names(nc$var) + } + + var_values <- lapply( + varnames, + ncdf4::ncvar_get, + nc = nc + ) + + # todo handle this case (multi-loc files?) + stopifnot(all(sapply(var_values, length) == length(timestamps))) + + var_values |> + stats::setNames(varnames) |> + as.data.frame() |> + transform(timestamp = timestamps) +} diff --git a/models/rothc/R/model2netcdf.RothC.R b/models/rothc/R/model2netcdf.RothC.R new file mode 100644 index 00000000000..9d53f653b36 --- /dev/null +++ b/models/rothc/R/model2netcdf.RothC.R @@ -0,0 +1,117 @@ +#' Convert RothC output to PEcAn-formatted netCDF +#' +#' @param outdir Location of model output +#' @param sitelat Latitude of the site +#' @param sitelon Longitude of the site +#' @param start_date Start time of the simulation +#' @param end_date End time of the simulation +#' @export +#' @importFrom rlang .data .env +#' +#' @author Chris Black +model2netcdf.RothC <- function(outdir, sitelat, sitelon, start_date, end_date) { + start_date <- as.Date(start_date) + end_date <- as.Date(end_date) + + var_map <- rothc_varname_map |> + dplyr::filter(!is.na(.data$pecan_name)) |> + dplyr::left_join( + PEcAn.utils::standard_vars, + by = c("pecan_name" = "Variable.Name") + ) |> + dplyr::rename("pecan_unit" = "Units", "pecan_longname" = "Long.name") + + # used below to rename variables after unit conversion + var_rename_list <- var_map$rothc_name |> + stats::setNames(var_map$pecan_name) + + + # Year, Month, C_Inp_t_C_ha, OA_Inp_t_C_ha, TEMP_C, RM_TMP, RAIN_mm, + # PEVAP_mm, SMD_mm, RM_Moist, PC, RM_PC, DPM_t_C_ha, RPM_t_C_ha, Bio_t_C_ha, + # Hum_t_C_ha, IOM_t_C_ha, SOC_t_C_ha, CO2_t_C_ha + res <- utils::read.csv(file.path(outdir, "month_results.out")) |> + # Remove spinup lines. TODO handle better? + dplyr::filter(.data$Year > 1) |> + # Convert time formats + dplyr::mutate( + midpoint_date = + ISOdate(year = .data$Year, month = .data$Month, day = 15, tz = "UTC"), + month_began = .data$midpoint_date |> + lubridate::rollback(preserve_hms = FALSE, roll_to_first = TRUE), + month_ended = .data$midpoint_date |> + lubridate::rollforward(preserve_hms = FALSE, roll_to_first = FALSE), + time = difftime(.data$midpoint_date, .env$start_date, units = "days") + ) + + # ensure all requested dates are present and then subset, + # padding start/end to full months + stopifnot( + start_date >= min(res$month_began), + end_date <= max(res$month_ended) + ) + res <- res |> + dplyr::filter( + .data$month_began >= lubridate::rollback(start_date, roll_to_first = TRUE), + .data$month_ended <= lubridate::rollforward(end_date, roll_to_first = FALSE) + ) |> + # Convert output columns to PEcAn std units/names + dplyr::mutate( + dplyr::across( + dplyr::all_of(var_map$rothc_name), + \(x) convert_to_pecan_units(x, varmap = var_map) + ) + ) |> + dplyr::rename(dplyr::all_of(var_rename_list)) + + + # Loop over years of output writing netCDFs + # TODO this is consistent with PEcAn practice for other models, + # but feels a bit silly with monthly output. + # Consider writing as one multi-year file instead? + for (yrdat in split(res, ~Year)) { + yr <- yrdat$Year[[1]] + # time_bounds <- array(data = NA, dim = c(nrow(yrdat), 2)) + # time_bounds[, 1] <- yrdat$month_start_date |> + # difftime(start_date, units = "days") + # time_bounds[, 2] <- yrdat$month_end_date |> + # difftime(start_date, units = "days") + nc_dims <- list( + lon = PEcAn.utils::to_ncdim("lon", sitelon), + lat = PEcAn.utils::to_ncdim("lat", sitelat), + time = ncdf4::ncdim_def( + name = "time", + longname = "time", + units = paste0("days since ", start_date, " 00:00:00"), + vals = as.integer(yrdat$time), + calendar = "standard", + unlim = TRUE + ) + # time_bounds = PEcAn.utils::to_ncdim("time_bounds", time_bounds) # TODO should this be "time_bnds"? + ) + nc_vars <- lapply(var_map$pecan_name, PEcAn.utils::to_ncvar, nc_dims) |> + stats::setNames(var_map$pecan_name) + # nc_var$time_bounds <- ??? #TODO + + nc <- ncdf4::nc_create( + file.path(outdir, paste0(yr, ".nc")), + nc_vars + ) + # ncdf4::ncatt_put(nc, "time", "bounds", "time_bounds", prec=NA) # TODO need this line? + for (vn in var_map$pecan_name) { + ncdf4::ncvar_put(nc, nc_vars[[vn]], yrdat[[vn]]) + } + ncdf4::nc_close(nc) + } +} + + + +# Helper to change column units by lookup in var map +# NB `cur_column()` only works inside `dplyr::across()` -- +# it returns the column name `x` has in the enclosing dataframe. +convert_to_pecan_units <- function(x, varmap) { + var_row <- varmap[varmap$rothc_name == dplyr::cur_column(), ] + stopifnot(nrow(var_row) == 1) + + PEcAn.utils::ud_convert(x, var_row$rothc_unit, var_row$pecan_unit) +} diff --git a/models/rothc/R/read_restart.RothC.R b/models/rothc/R/read_restart.RothC.R new file mode 100644 index 00000000000..f8bd824e06c --- /dev/null +++ b/models/rothc/R/read_restart.RothC.R @@ -0,0 +1,21 @@ +#' @title Read restart template for SDA +#' +#' @author Alexey Shiklomanov +#' +#' @param outdir Output directory +#' @param runid Run ID +#' @param stop.time Year that is being read +#' @param settings PEcAn settings object +#' @param var.names Variable names to be extracted +#' @param params Any parameters required for state calculations +#' +#' @description Read restart files from model. +#' +#' @return Forecast numeric matrix +#' @export +read_restart.ModelName <- function(outdir, + runid, + stop.time, + settings, + var.names, + params) {} diff --git a/models/rothc/R/rothc_varname_map.R b/models/rothc/R/rothc_varname_map.R new file mode 100644 index 00000000000..3adf9f66602 --- /dev/null +++ b/models/rothc/R/rothc_varname_map.R @@ -0,0 +1,23 @@ +# Mapping from RothC to PEcAn standard names +# Not complete yet -- missing pecan_name means I haven't looked yet, +# not that we know no equivalent exists +rothc_varname_map <- dplyr::tribble( + ~rothc_name, ~rothc_description, ~rothc_unit, ~pecan_name, + "C_Inp_t_C_ha", "C input", "t C ha-1", NA, + "OA_Inp_t_C_ha", "Farmyard manure", "t C ha-1", NA, + "TEMP_C", "Air temperature", "degC", "Tair", + "RM_TMP", "Rate modifier for temp", "1", NA, + "RAIN_mm", "Rainfall", "mm", NA, + "PEVAP_mm", "Open pan evaporation", "mm", NA, # probably matches to "Evap" -- but careful with units _and_ with conceptual difference between pan evap and total evap + "SWC_mm", "Accumulated soil water deficit", "mm", NA, # TODO may be called "SMD_mm" now? + "RM_Moist", "Rate modifier for soil moist", "1", NA, + "PC", "Soil plant cover (0=bare 1=covered)", "1", NA, + "RM_PC", "rate modifier for crop cover", "1", NA, + "DPM_t_C_ha", "Decomposable plant material", "t C ha-1", NA, + "RPM_t_C_ha", "Resistant plant material", "t C ha-1", NA, + "Bio_t_C_ha", "Microbial biomass", "t C ha-1", NA, + "Hum_t_C_ha", "Humified organic matter", "t C ha-1", "slow_soil_pool_carbon_content", + "IOM_t_C_ha", "Inert organic matter", "t C ha-1", "structural_soil_pool_carbon_content", + "SOC_t_C_ha", "Total soil organic carbon", "t C ha-1", "TotSoilCarb", + "CO2_t_C_ha", "Accumulated CO2", "t C ha-1", NA # needs time dimension to match to "TotalResp" which is kg C m-2 sec-1 +) diff --git a/models/rothc/R/write.config.RothC.R b/models/rothc/R/write.config.RothC.R new file mode 100644 index 00000000000..f4a1a4f2483 --- /dev/null +++ b/models/rothc/R/write.config.RothC.R @@ -0,0 +1,149 @@ +##' Writes a RothC config file. +##' +##' Requires a pft xml object, a list of trait values for a single model run, +##' and the name of the file to create +##' +##' @param defaults list of defaults to process +##' @param trait.values vector of samples for a given trait +##' @param settings list of settings from pecan settings file +##' @param run.id id of run +##' @return configuration file for MODEL for given run +##' @export +##' @author Rob Kooper +write.config.RothC <- function(defaults, trait.values, settings, run.id) { + + # find out where to write run/ouput + rundir <- file.path(settings$host$rundir, run.id) + outdir <- file.path(settings$host$outdir, run.id) + + #----------------------------------------------------------------------- + # create launch script (which will create symlink) + if (!is.null(settings$model$jobtemplate) && file.exists(settings$model$jobtemplate)) { + jobsh <- readLines(con = settings$model$jobtemplate, n = -1) + } else { + jobsh <- readLines(con = system.file("template.job", package = "PEcAn.RothC"), n = -1) + } + + # create host specific setttings + hostsetup <- "" + if (!is.null(settings$model$prerun)) { + hostsetup <- paste(hostsetup, sep = "\n", paste(settings$model$prerun, collapse = "\n")) + } + if (!is.null(settings$host$prerun)) { + hostsetup <- paste(hostsetup, sep = "\n", paste(settings$host$prerun, collapse = "\n")) + } + + hostteardown <- "" + if (!is.null(settings$model$postrun)) { + hostteardown <- paste(hostteardown, sep = "\n", paste(settings$model$postrun, collapse = "\n")) + } + if (!is.null(settings$host$postrun)) { + hostteardown <- paste(hostteardown, sep = "\n", paste(settings$host$postrun, collapse = "\n")) + } + + # create job.sh + jobsh <- gsub("@HOST_SETUP@", hostsetup, jobsh) + jobsh <- gsub("@HOST_TEARDOWN@", hostteardown, jobsh) + + jobsh <- gsub("@SITE_LAT@", settings$run$site$lat, jobsh) + jobsh <- gsub("@SITE_LON@", settings$run$site$lon, jobsh) + jobsh <- gsub("@SITE_MET@", settings$run$site$met, jobsh) + + jobsh <- gsub("@START_DATE@", settings$run$start.date, jobsh) + jobsh <- gsub("@END_DATE@", settings$run$end.date, jobsh) + jobsh <- gsub("@OUTDIR@", outdir, jobsh) + jobsh <- gsub("@RUNDIR@", rundir, jobsh) + + jobsh <- gsub("@BINARY@", settings$model$binary, jobsh) + + writeLines(jobsh, con = file.path(settings$rundir, run.id, "job.sh")) + Sys.chmod(file.path(settings$rundir, run.id, "job.sh")) + + #----------------------------------------------------------------------- + ### Edit a templated config file for runs + if (!is.null(settings$model$config) && file.exists(settings$model$config)) { + config.text <- readLines(con = settings$model$config, n = -1) + } else { + filename <- system.file(settings$model$config, package = "PEcAn.MODEL") + if (filename == "") { + if (!is.null(settings$model$revision)) { + filename <- system.file(paste0("config.", settings$model$revision), package = "PEcAn.MODEL") + } else { + model <- PEcAn.DB::db.query(paste("SELECT * FROM models WHERE id =", settings$model$id), params = settings$database$bety) + filename <- system.file(paste0("config.r", model$revision), package = "PEcAn.MODEL") + } + } + if (filename == "") { + PEcAn.logger::logger.severe("Could not find config template") + } + PEcAn.logger::logger.info("Using", filename, "as template") + config.text <- readLines(con = filename, n = -1) + } + + config.text <- gsub("@SITE_LAT@", settings$run$site$lat, config.text) + config.text <- gsub("@SITE_LON@", settings$run$site$lon, config.text) + config.text <- gsub("@SITE_MET@", settings$run$inputs$met$path, config.text) + config.text <- gsub("@MET_START@", settings$run$site$met.start, config.text) + config.text <- gsub("@MET_END@", settings$run$site$met.end, config.text) + config.text <- gsub("@START_MONTH@", format(settings$run$start.date, "%m"), config.text) + config.text <- gsub("@START_DAY@", format(settings$run$start.date, "%d"), config.text) + config.text <- gsub("@START_YEAR@", format(settings$run$start.date, "%Y"), config.text) + config.text <- gsub("@END_MONTH@", format(settings$run$end.date, "%m"), config.text) + config.text <- gsub("@END_DAY@", format(settings$run$end.date, "%d"), config.text) + config.text <- gsub("@END_YEAR@", format(settings$run$end.date, "%Y"), config.text) + config.text <- gsub("@OUTDIR@", settings$host$outdir, config.text) + config.text <- gsub("@ENSNAME@", run.id, config.text) + config.text <- gsub("@OUTFILE@", paste0("out", run.id), config.text) + +# TODO make these editable -- hard-coding for MVP +# OPT_RMMOIST: soil water parameterization. +# "1: Standard RothC soil water parameters" +# "2: Van Genuchten soil properties and soil is allowed to be drier (ie hygroscopic / capillary water, -1000bar)" +# "3: Van Genuchten soil properties, but uses the Standard RothC soil water function" +config.text <- gsub("@OPT_RMMOIST@", "1", config.text) +# Bare SMD: wilting point configuration +# "1: Standard RothC bareSMD" +# "2: bareSMD is set to wilting point -15bar (could be better for dry soils)" +config.text <- gsub("@OPT_SDDBARE@", "1", config.text) + +# Soil parameters -- TODO read from run$inputs$soil_physics +# clay_pct, depth_cm, iom_tC_ha, nsteps, silt_pct, bulkdens_g_m3, org_C_pct, min_RM_moist +config.text <- gsub("@SOIL_PARAMS@", "23.4 23.0 3.0041 840 58.6 1.27 0.94 0.2", config.text) + +# Climate data + management inputs +# TODO all managements hardcoded for MVP +met_path <- inputs$met$path %||% settings$run$inputs$met$path +met_in <- read.table(met_path, header = FALSE) +zros <- rep(0, nrow(met_in)) +inputs <- data.frame( + C_inp_tC_ha = zros, + FYM_tC_ha = zros, + PC = zros, + PL_DPM_f = zros, + PL_RPM_f = zros, + OA_DPM_f = zros, + OA_RPM_f = zros, + OA_BIO_f = zros, + OA_HUM_f = zros +) + +input_rows <- met_in |> + bind_cols(inputs) |> + dplyr::select(all_of( + "year", "month", + "modern_pct", + "Tmp_C", "Rain_mm", "Evap_mm", + "C_inp_tC_ha", "FYM_tC_ha", "PC", + "PL_DPM_f", "PL_RPM_f", + "OA_DPM_f", "OA_RPM_f", "OA_BIO_f", "OA_HUM_f" + )) |> + # Kinda ugly: Convert to one string to cram it into the template via gsub + format() |> + apply(1, paste, collapse = " ") |> + paste(collapse = "\n") + + config.text <- gsub("@CLIM_DATA@", input_rows, config.text) + + config.file.name <- "RothC_input.dat" + writeLines(config.text, con = file.path(outdir, config.file.name)) +} # write.config.MODEL diff --git a/models/rothc/R/write_restart.RothC.R b/models/rothc/R/write_restart.RothC.R new file mode 100644 index 00000000000..86f4b0cd48a --- /dev/null +++ b/models/rothc/R/write_restart.RothC.R @@ -0,0 +1,28 @@ +#' @title Write restart template for SDA +#' +#' @author Alexey Shiklomanov +#' +#' @param outdir outout directory +#' @param runid run id +#' @param start.time Time of current assimilation step +#' @param stop.time Time of next assimilation step +#' @param settings pecan settings list +#' @param new.state Analysis state matrix returned by \code{sda.enkf} +#' @param RENAME flag to either rename output file or not +#' @param new.params optional, additionals params to pass write.configs that +#' are deterministically related to the parameters updated by the analysis +#' @param inputs new input paths updated by the SDA workflow, +#' will be passed to write.configs +#' +#' @description Write restart files for model +#' +#' @export +write_restart.ModelName <- function(outdir, + runid, + start.time, + stop.time, + settings, + new.state, + RENAME, + new.params, + inputs) {} diff --git a/models/rothc/README.md b/models/rothc/README.md new file mode 100644 index 00000000000..79ecc18223e --- /dev/null +++ b/models/rothc/README.md @@ -0,0 +1,100 @@ +A generic template for adding a new model to PEcAn +========================================================================== + +Adding a new model to PEcAn in a few easy steps: + +1. add modeltype to BETY +2. add a model and PFT to BETY for use with modeltype +3. implement 3 functions as described below +4. Add tests to `tests/testthat` +5. Update README, documentation +6. Update Dockerfile and model_info.json +7. execute pecan with new model + + +### Three Functions + +There are 3 functions that will need to be implemented, each of these +functions will need to have MODEL be replaced with the actual modeltype as +it is defined in the BETY database. + +* `write.config.MODEL.R` + + This will write the configuratin file as well as the job launcher used by + PEcAn. There is an example of the job execution script in the template + folder. The configuration file can also be a template that is found based + on the revision number of the model. This should use the computed results + specified in defaults and trait.values to write a configuration file + based on the PFT and traits found. + +* `met2model.MODEL.R` + + This will convert the standard Met CF file to the model specific file + format. This will allow PEcAn to create metereological files for the + specific site and model. This will only be called if no meterological + data is found for that specific site and model combination. + +* `model2netcdf.MODEL.R` + + This will convert the model specific output to NACP Intercomparison + format. After this function is finished PEcAn will use the generated + output and not use the model specific outputs. The outputs should be + named YYYY.nc + +### Dockerization + +The PEcAn system is leveraging Docker to encapsulate most of the code. +This will make it easier to share new model with others, without them +having to compile the models. The goal is for people to be able to +launch the model in docker, and it will register with PEcAn and is +almost immediatly available to be used. To accomplish this you will need to modify two files. + +* `Dockerfile` + + The [Dockerfile](https://docs.docker.com/engine/reference/builder/) is + like the Makefile for docker. This file is split in two pieces the + part at the top is to actually build the binary. This is where you + specify all the libraries that are needed, as well as all the build + tools to compile your model. The second part, starting at the second + `FROM` line, is where you will install only the libraries needed to + run the binary and copy the binary from the build stage, using the + `COPY --from` line. + +* `model_info.json` + + The model_info.json describes the model and is used to register the + model with PEcAn. In the model_info.json the only fields that are + really required are those at the top: `name`, `type`, `version` and + `binary`. All other fields are optional but are good to be filled + out. You can leave `version` and `binary` with the special values + which will be updated by the Dockerfile. + +Once the image can be build it can be pushed so others can leverage +of the model. For PEcAn we have been using the following naming scheme +for the docker images: `pecan/model--:` +where the `model` and `model_version` are the same as those used to +build the model, and `pecan_version` is the version of PEcAn this +model is compiled for. + +### Additional Changes + +* `README.md` + +This file should contain basic background information about the model. +At a minimum, this should include the scientific motivation and scope, +name(s) of maintainer(s), links to project homepage, and a list of a few +key publications. +relevant publications. + +* `/tests/testthat/` + +Each package should have tests that cover the key functions of the package, +at a minimum, the three functions above. + +* documentation + +Update the `NAMESPACE`, `DESCRIPTION` and `man/*.Rd` files by running + +```r +devtools("models//") +``` diff --git a/models/rothc/inst/RothC_input_template.dat b/models/rothc/inst/RothC_input_template.dat new file mode 100644 index 00000000000..e005b1fbbb1 --- /dev/null +++ b/models/rothc/inst/RothC_input_template.dat @@ -0,0 +1,11 @@ +# RothC 2.0 input for @ENSNAME@ (@SITE_LAT@, @SITE_LON@) +# @START_YEAR@-@START_MONTH@-@START_DAY@ to @END_YEAR@-@END_MONTH@-@END_DAY@ +# +# opt_RMmoist opt_SMDbare + @OPT_RMMOIST@ @OPT_SDDBARE@ +# (%) (cm) (t C/ha) (-) (%) (g/m3) (%) (-) +# clay depth iom nsteps siltper BD OC minRM_Moist +@SOIL_PARAMS@ +# (-) (-) (%) (C) (mm) (mm) (t C/ha) (t C/ha) (-) (-) (-) (-) (-) (-) (-) +# year month modern Tmp Rain Evap C_inp FYM PC PL_DPM_f PL_RPM_f OA_DPM_f OA_RPM_f OA_BIO_f OA_HUM_f +@CLIM_DATA@ \ No newline at end of file diff --git a/models/rothc/inst/template.job b/models/rothc/inst/template.job new file mode 100644 index 00000000000..52add090821 --- /dev/null +++ b/models/rothc/inst/template.job @@ -0,0 +1,42 @@ +#!/bin/bash + +# redirect output +exec 3>&1 +exec &> "@OUTDIR@/logfile.txt" + +# host specific setup +@HOST_SETUP@ +@CDO_SETUP@ + +# create output folder +mkdir -p "@OUTDIR@" + +# see if application needs running +if [ ! -e "@OUTDIR@/results.csv" ]; then + cd "@RUNDIR@" + + "@BINARY@" + STATUS=$? + + # check the status + if [ $STATUS -ne 0 ]; then + echo -e "ERROR IN MODEL RUN\nLogfile is located at '@OUTDIR@/logfile.txt'" >&3 + exit $STATUS + fi + + # convert to MsTMIP + echo "require (PEcAn.MODEL) +model2netcdf.MODEL('@OUTDIR@', @SITE_LAT@, @SITE_LON@, '@START_DATE@', '@END_DATE@') +" | R --vanilla +fi + +# copy readme with specs to output +cp "@RUNDIR@/README.txt" "@OUTDIR@/README.txt" + +# run getdata to extract right variables + +# host specific teardown +@HOST_TEARDOWN@ + +# all done +echo -e "MODEL FINISHED\nLogfile is located at '@OUTDIR@/logfile.txt'" >&3 diff --git a/models/rothc/inst/template_geo.job b/models/rothc/inst/template_geo.job new file mode 100644 index 00000000000..e4fe27e04a7 --- /dev/null +++ b/models/rothc/inst/template_geo.job @@ -0,0 +1,54 @@ +#!/bin/bash -l + +# redirect output +exec 3>&1 +exec &> "@OUTDIR@/logfile.txt" + +# host specific setup +@HOST_SETUP@ +@CDO_SETUP@ + +# create output folder +mkdir -p "@OUTDIR@" + +# see if application needs running +if [ ! -e "@OUTDIR@/sipnet.out" ]; then + cd "@RUNDIR@" + ln -s "@SITE_MET@" sipnet.clim + + "@BINARY@" + STATUS=$? + + # copy output + mv "@RUNDIR@/sipnet.out" "@OUTDIR@" + + # check the status + if [ $STATUS -ne 0 ]; then + echo -e "ERROR IN MODEL RUN\nLogfile is located at '@OUTDIR@/logfile.txt'" >&3 + exit $STATUS + fi + + # convert to MsTMIP + echo "require (PEcAn.SIPNET) + model2netcdf.SIPNET('@OUTDIR@', @SITE_LAT@, @SITE_LON@, '@START_DATE@', '@END_DATE@', @DELETE.RAW@, '@REVISION@', '@PREFIX@', @OVERWRITE@, @CONFLICT@) + " | R --no-save +fi + +# copy readme with specs to output +cp "@RUNDIR@/README.txt" "@OUTDIR@/README.txt" + +# run getdata to extract right variables + +# host specific teardown +@HOST_TEARDOWN@ + +#copy files back. +@CPRUNCMD@ +@CPOUTCMD@ + +#delete files in the run and out folder. +@RMRUNDIRCMD@ +@RMOUTDIRCMD@ + +# all done +echo -e "MODEL FINISHED\nLogfile is located at '@OUTDIR@/logfile.txt'" >&3 \ No newline at end of file diff --git a/models/rothc/man/met2model.RothC.Rd b/models/rothc/man/met2model.RothC.Rd new file mode 100644 index 00000000000..1d45c014888 --- /dev/null +++ b/models/rothc/man/met2model.RothC.Rd @@ -0,0 +1,45 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/met2model.RothC.R +\name{met2model.RothC} +\alias{met2model.RothC} +\title{Extract monthly weather from CF file for input to RothC} +\usage{ +met2model.RothC( + in.path, + in.prefix, + outfolder, + start_date, + end_date, + overwrite = FALSE +) +} +\arguments{ +\item{in.path}{path on disk where CF files live} + +\item{in.prefix}{prefix for each file} + +\item{outfolder}{location where model specific output is written.} + +\item{start_date, end_date}{When to start and end output. +Specify as exact dates, but output will be padded to whole months.} + +\item{overwrite}{logical: replace output files if they already exist?} +} +\value{ +data frame summarizing file metadata +} +\description{ +Input files need to be named `/.YYYY.nc` +} +\details{ +Output files are named `/.YY-mm.YY-mm.dat` +with one line per month and columns for temperature, rainfall, +and evaporation. + +Note that the created file contains only weather data and not any of the +soil or management data needed for RothC's single combined input file. +See `write.config.RothC()` for assembly into a model-ready RothC_input.dat`. +} +\author{ +Chris Black +} diff --git a/models/rothc/man/model2netcdf.RothC.Rd b/models/rothc/man/model2netcdf.RothC.Rd new file mode 100644 index 00000000000..857e5da161c --- /dev/null +++ b/models/rothc/man/model2netcdf.RothC.Rd @@ -0,0 +1,25 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/model2netcdf.RothC.R +\name{model2netcdf.RothC} +\alias{model2netcdf.RothC} +\title{Convert RothC output to PEcAn-formatted netCDF} +\usage{ +model2netcdf.RothC(outdir, sitelat, sitelon, start_date, end_date) +} +\arguments{ +\item{outdir}{Location of model output} + +\item{sitelat}{Latitude of the site} + +\item{sitelon}{Longitude of the site} + +\item{start_date}{Start time of the simulation} + +\item{end_date}{End time of the simulation} +} +\description{ +Convert RothC output to PEcAn-formatted netCDF +} +\author{ +Chris Black +} diff --git a/models/rothc/man/read_restart.ModelName.Rd b/models/rothc/man/read_restart.ModelName.Rd new file mode 100644 index 00000000000..45baa5f051b --- /dev/null +++ b/models/rothc/man/read_restart.ModelName.Rd @@ -0,0 +1,30 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/read_restart.RothC.R +\name{read_restart.ModelName} +\alias{read_restart.ModelName} +\title{Read restart template for SDA} +\usage{ +read_restart.ModelName(outdir, runid, stop.time, settings, var.names, params) +} +\arguments{ +\item{outdir}{Output directory} + +\item{runid}{Run ID} + +\item{stop.time}{Year that is being read} + +\item{settings}{PEcAn settings object} + +\item{var.names}{Variable names to be extracted} + +\item{params}{Any parameters required for state calculations} +} +\value{ +Forecast numeric matrix +} +\description{ +Read restart files from model. +} +\author{ +Alexey Shiklomanov +} diff --git a/models/rothc/man/write.config.RothC.Rd b/models/rothc/man/write.config.RothC.Rd new file mode 100644 index 00000000000..03c47aa2d52 --- /dev/null +++ b/models/rothc/man/write.config.RothC.Rd @@ -0,0 +1,27 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/write.config.RothC.R +\name{write.config.RothC} +\alias{write.config.RothC} +\title{Writes a RothC config file.} +\usage{ +write.config.RothC(defaults, trait.values, settings, run.id) +} +\arguments{ +\item{defaults}{list of defaults to process} + +\item{trait.values}{vector of samples for a given trait} + +\item{settings}{list of settings from pecan settings file} + +\item{run.id}{id of run} +} +\value{ +configuration file for MODEL for given run +} +\description{ +Requires a pft xml object, a list of trait values for a single model run, +and the name of the file to create +} +\author{ +Rob Kooper +} diff --git a/models/rothc/man/write_restart.ModelName.Rd b/models/rothc/man/write_restart.ModelName.Rd new file mode 100644 index 00000000000..c6cdccf73c2 --- /dev/null +++ b/models/rothc/man/write_restart.ModelName.Rd @@ -0,0 +1,45 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/write_restart.RothC.R +\name{write_restart.ModelName} +\alias{write_restart.ModelName} +\title{Write restart template for SDA} +\usage{ +write_restart.ModelName( + outdir, + runid, + start.time, + stop.time, + settings, + new.state, + RENAME, + new.params, + inputs +) +} +\arguments{ +\item{outdir}{outout directory} + +\item{runid}{run id} + +\item{start.time}{Time of current assimilation step} + +\item{stop.time}{Time of next assimilation step} + +\item{settings}{pecan settings list} + +\item{new.state}{Analysis state matrix returned by \code{sda.enkf}} + +\item{RENAME}{flag to either rename output file or not} + +\item{new.params}{optional, additionals params to pass write.configs that +are deterministically related to the parameters updated by the analysis} + +\item{inputs}{new input paths updated by the SDA workflow, +will be passed to write.configs} +} +\description{ +Write restart files for model +} +\author{ +Alexey Shiklomanov +} diff --git a/models/rothc/model_info.json b/models/rothc/model_info.json new file mode 100644 index 00000000000..03de0f14082 --- /dev/null +++ b/models/rothc/model_info.json @@ -0,0 +1,21 @@ +{ + "name": "RothC", + "type": "Soil C", + "version": "@VERSION@", + "binary": "@BINARY@", + "description": "RothC models the turnover of organic carbon in non-waterlogged top-soil. It accounts for the effects of soil texture, temperature, moisture content and plant cover on the turnover process. It uses a monthly time step to calculate total organic carbon (t ha-1), microbial biomass carbon (t ha-1) and d14C (from which the equivalent radiocarbon age of the soil can be calculated).", + "creator": "Model by D.S. Jenkinson and others at Rothamsted Research. PEcAn coupler by Chris Black", + "contributors": ["Kevin Coleman", "Jonah Prout", "Alice Milne"], + "links": [ + { + "type": "git", + "description": "Source code for RothC", + "url": "https://github.com/Rothamsted-Models/RothC_Py" + }, + ], + "inputs": {}, + "bibtex": [ + "How to cite the model, paper 1", + "How to cite the model, paper 2" + ] +} diff --git a/models/rothc/tests/Rcheck_reference.log b/models/rothc/tests/Rcheck_reference.log new file mode 100644 index 00000000000..ef66dc4460f --- /dev/null +++ b/models/rothc/tests/Rcheck_reference.log @@ -0,0 +1,108 @@ +* using log directory ‘/home/tanishq010/pecan/models/PEcAn.ModelName.Rcheck’ +* using R version 4.2.1 (2022-06-23) +* using platform: x86_64-pc-linux-gnu (64-bit) +* using session charset: UTF-8 +* using options ‘--no-manual --as-cran’ +* checking for file ‘PEcAn.ModelName/DESCRIPTION’ ... OK +* checking extension type ... Package +* this is package ‘PEcAn.ModelName’ version ‘1.7.2’ +* package encoding: UTF-8 +* checking CRAN incoming feasibility ... WARNING +Maintainer: ‘Jane Doe ’ + +New submission + +License components with restrictions and base license permitting such: + BSD_3_clause + file LICENSE +File 'LICENSE': + ## This is the master copy of the PEcAn License + + University of Illinois/NCSA Open Source License + + Copyright (c) 2012, University of Illinois, NCSA. All rights reserved. + + PEcAn project + www.pecanproject.org + + Permission is hereby granted, free of charge, to any person obtaining + a copy of this software and associated documentation files (the + "Software"), to deal with the Software without restriction, including + without limitation the rights to use, copy, modify, merge, publish, + distribute, sublicense, and/or sell copies of the Software, and to + permit persons to whom the Software is furnished to do so, subject to + the following conditions: + + - Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimers. + - Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimers in the + documentation and/or other materials provided with the distribution. + - Neither the names of University of Illinois, NCSA, nor the names + of its contributors may be used to endorse or promote products + derived from this Software without specific prior written permission. + + THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. + IN NO EVENT SHALL THE CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR + ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF + CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION + WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS WITH THE SOFTWARE. + +Strong dependencies not in mainstream repositories: + PEcAn.DB, PEcAn.logger, PEcAn.utils + +The Date field is over a month old. +* checking package namespace information ... OK +* checking package dependencies ... OK +* checking if this is a source package ... OK +* checking if there is a namespace ... OK +* checking for executable files ... OK +* checking for hidden files and directories ... OK +* checking for portable file names ... OK +* checking for sufficient/correct file permissions ... OK +* checking serialization versions ... OK +* checking whether package ‘PEcAn.ModelName’ can be installed ... OK +* checking installed package size ... OK +* checking package directory ... OK +* checking for future file timestamps ... OK +* checking DESCRIPTION meta-information ... OK +* checking top-level files ... OK +* checking for left-over files ... OK +* checking index information ... OK +* checking package subdirectories ... OK +* checking R files for non-ASCII characters ... OK +* checking R files for syntax errors ... OK +* checking whether the package can be loaded ... OK +* checking whether the package can be loaded with stated dependencies ... OK +* checking whether the package can be unloaded cleanly ... OK +* checking whether the namespace can be loaded with stated dependencies ... OK +* checking whether the namespace can be unloaded cleanly ... OK +* checking loading without being on the library search path ... OK +* checking use of S3 registration ... OK +* checking dependencies in R code ... NOTE +Namespace in Imports field not imported from: ‘PEcAn.utils’ + All declared Imports should be used. +* checking S3 generic/method consistency ... OK +* checking replacement functions ... OK +* checking foreign function calls ... OK +* checking R code for possible problems ... OK +* checking Rd files ... OK +* checking Rd metadata ... OK +* checking Rd line widths ... OK +* checking Rd cross-references ... OK +* checking for missing documentation entries ... OK +* checking for code/documentation mismatches ... OK +* checking Rd \usage sections ... OK +* checking Rd contents ... OK +* checking for unstated dependencies in examples ... OK +* checking examples ... NONE +* checking for unstated dependencies in ‘tests’ ... OK +* checking tests ... + Running ‘testthat.R’ + OK +* checking for non-standard things in the check directory ... OK +* checking for detritus in the temp directory ... OK +* DONE + +Status: 1 WARNING, 1 NOTE diff --git a/models/rothc/tests/testthat.R b/models/rothc/tests/testthat.R new file mode 100644 index 00000000000..86501762e35 --- /dev/null +++ b/models/rothc/tests/testthat.R @@ -0,0 +1,4 @@ +library(testthat) +library(PEcAn.utils) + +test_check("PEcAn.RothC") diff --git a/models/rothc/tests/testthat/README.txt b/models/rothc/tests/testthat/README.txt new file mode 100644 index 00000000000..1586c8c86c3 --- /dev/null +++ b/models/rothc/tests/testthat/README.txt @@ -0,0 +1,6 @@ +# Tests for the PEcAn.RothC package + +## Contents of the data folder: +* month_results.out: + First 3 years of the sample monthly output file provided with the RothC model code. + The full file (1939-2007) is available at https://github.com/Rothamsted-Models/RothC_Code/blob/main/month_results.out diff --git a/models/rothc/tests/testthat/data/month_results.out b/models/rothc/tests/testthat/data/month_results.out new file mode 100644 index 00000000000..a37fc1e77cf --- /dev/null +++ b/models/rothc/tests/testthat/data/month_results.out @@ -0,0 +1,39 @@ +Year, Month, C_Inp_t_C_ha, OA_Inp_t_C_ha, TEMP_C, RM_TMP, RAIN_mm, PEVAP_mm, SMD_mm, RM_Moist, PC, RM_PC, DPM_t_C_ha, RPM_t_C_ha, Bio_t_C_ha, Hum_t_C_ha, IOM_t_C_ha, SOC_t_C_ha, CO2_t_C_ha + 1, 0, , , , , , , , , , , 0.0000, 0.0000, 0.0000, 0.0000, 3.0041, 3.0041, 0.0000 + 1, 20916, , , , , , , , , , , 0.1606, 5.8208, 0.8717, 32.6202, 3.0041, 42.4774, 0.0000 + 1939, 1, 0.000, 0.000, 4.0, 0.4050, 114.5, 8.2, 0.00, 1.0000, 1, 0.6, 0.1311, 5.7856, 0.8693, 32.6177, 3.0041, 42.4078, 0.0697 + 1939, 2, 0.000, 0.000, 5.4, 0.5335, 30.7, 17.2, 0.00, 1.0000, 1, 0.6, 0.1004, 5.7395, 0.8653, 32.6133, 3.0041, 42.3226, 0.1548 + 1939, 3, 0.000, 0.000, 5.5, 0.5425, 51.6, 43.1, 0.00, 1.0000, 1, 0.6, 0.0766, 5.6929, 0.8605, 32.6080, 3.0041, 42.2421, 0.2353 + 1939, 4, 0.000, 0.000, 8.3, 0.8688, 90.3, 68.2, 0.00, 1.0000, 1, 0.6, 0.0496, 5.6192, 0.8518, 32.5980, 3.0041, 42.1227, 0.3547 + 1939, 5, 0.000, 0.000, 10.5, 1.1757, 40.4, 99.3, -34.07, 0.5480, 1, 0.6, 0.0359, 5.5652, 0.8447, 32.5898, 3.0041, 42.0397, 0.4377 + 1939, 6, 0.000, 0.000, 13.6, 1.6506, 74.9, 112.5, -43.55, 0.2446, 1, 0.6, 0.0294, 5.5316, 0.8401, 32.5844, 3.0041, 41.9895, 0.4879 + 1939, 7, 0.000, 0.000, 15.6, 1.9970, 59.8, 128.9, -44.94, 0.2000, 1, 0.6, 0.0241, 5.4985, 0.8355, 32.5788, 3.0041, 41.9410, 0.5364 + 1939, 8, 1.464, 0.000, 16.1, 2.0863, 96.4, 104.4, -26.84, 0.7795, 1, 0.6, 0.8749, 5.9661, 0.8165, 32.5549, 3.0041, 43.2165, 0.7252 + 1939, 9, 0.000, 0.000, 14.1, 1.7382, 25.2, 64.2, -26.84, 0.7795, 0, 1.0, 0.2829, 5.7674, 0.8519, 32.5918, 3.0041, 42.4982, 1.4435 + 1939, 10, 0.000, 0.000, 7.9, 0.8195, 103.6, 27.7, 0.00, 1.0000, 0, 1.0, 0.1429, 5.6505, 0.8489, 32.5879, 3.0041, 42.2343, 1.7074 + 1939, 11, 0.000, 0.000, 8.5, 0.8907, 105.4, 12.9, 0.00, 1.0000, 1, 0.6, 0.0915, 5.5755, 0.8427, 32.5804, 3.0041, 42.0942, 1.8475 + 1939, 12, 0.000, 0.000, 2.4, 0.2787, 58.0, 3.6, 0.00, 1.0000, 1, 0.6, 0.0796, 5.5522, 0.8402, 32.5776, 3.0041, 42.0538, 1.8880 + 1940, 1, 0.000, 0.000, -2.6, 0.0550, 56.1, 2.4, 0.00, 1.0000, 1, 0.6, 0.0775, 5.5476, 0.8397, 32.5770, 3.0041, 42.0459, 1.8958 + 1940, 2, 0.000, 0.000, 1.6, 0.2298, 52.9, 18.2, 0.00, 1.0000, 1, 0.6, 0.0691, 5.5286, 0.8376, 32.5744, 3.0041, 42.0138, 1.9280 + 1940, 3, 0.000, 0.000, 6.0, 0.5954, 79.1, 47.1, 0.00, 1.0000, 1, 0.6, 0.0513, 5.4794, 0.8318, 32.5673, 3.0041, 41.9339, 2.0078 + 1940, 4, 0.057, 0.000, 8.5, 0.8972, 69.7, 78.9, 0.00, 1.0000, 1, 0.6, 0.0664, 5.4296, 0.8223, 32.5555, 3.0041, 41.8779, 2.1209 + 1940, 5, 0.057, 0.000, 12.1, 1.4213, 45.8, 105.3, -33.17, 0.5768, 1, 0.6, 0.0778, 5.3866, 0.8143, 32.5452, 3.0041, 41.8281, 2.2278 + 1940, 6, 0.057, 0.000, 15.7, 2.0148, 35.4, 120.3, -44.94, 0.2000, 1, 0.6, 0.0973, 5.3775, 0.8108, 32.5406, 3.0041, 41.8303, 2.2827 + 1940, 7, 0.057, 0.000, 14.8, 1.8584, 82.2, 124.7, -44.94, 0.2000, 1, 0.6, 0.1145, 5.3710, 0.8079, 32.5366, 3.0041, 41.8342, 2.3360 + 1940, 8, 0.057, 0.000, 15.5, 1.9811, 6.0, 111.3, -44.94, 0.2000, 1, 0.6, 0.1276, 5.3626, 0.8051, 32.5328, 3.0041, 41.8323, 2.3950 + 1940, 9, 0.285, 0.000, 12.6, 1.4979, 29.7, 62.9, -44.94, 0.2000, 1, 0.6, 0.2783, 5.4555, 0.8032, 32.5302, 3.0041, 42.0714, 2.4413 + 1940, 10, 0.000, 0.000, 8.8, 0.9339, 81.4, 28.1, 0.00, 1.0000, 1, 0.6, 0.1745, 5.3796, 0.8027, 32.5279, 3.0041, 41.8888, 2.6238 + 1940, 11, 0.000, 0.000, 6.4, 0.6408, 189.3, 6.0, 0.00, 1.0000, 1, 0.6, 0.1266, 5.3282, 0.7999, 32.5234, 3.0041, 41.7822, 2.7304 + 1940, 12, 0.000, 0.000, 3.0, 0.3273, 46.2, 3.0, 0.00, 1.0000, 1, 0.6, 0.1075, 5.3021, 0.7979, 32.5205, 3.0041, 41.7321, 2.7806 + 1941, 1, 0.000, 0.000, -0.0, 0.1425, 74.8, 4.0, 0.00, 1.0000, 0, 1.0, 0.0955, 5.2832, 0.7962, 32.5181, 3.0041, 41.6971, 2.8155 + 1941, 2, 0.000, 0.000, 3.3, 0.3466, 55.9, 16.1, 0.00, 1.0000, 0, 1.0, 0.0715, 5.2376, 0.7917, 32.5117, 3.0041, 41.6167, 2.8959 + 1941, 3, 0.000, 0.000, 4.5, 0.4502, 88.3, 43.0, 0.00, 1.0000, 0, 1.0, 0.0492, 5.1790, 0.7851, 32.5022, 3.0041, 41.5196, 2.9930 + 1941, 4, 0.000, 0.000, 6.2, 0.6179, 45.6, 73.9, -9.83, 1.0000, 0, 1.0, 0.0294, 5.0996, 0.7750, 32.4878, 3.0041, 41.3959, 3.1167 + 1941, 5, 0.000, 0.000, 8.4, 0.8842, 58.6, 104.3, -29.45, 0.6960, 1, 0.6, 0.0216, 5.0528, 0.7686, 32.4786, 3.0041, 41.3257, 3.1869 + 1941, 6, 0.000, 0.000, 15.0, 1.9002, 40.9, 120.8, -44.94, 0.2000, 1, 0.6, 0.0179, 5.0240, 0.7646, 32.4728, 3.0041, 41.2834, 3.2293 + 1941, 7, 0.000, 0.000, 17.5, 2.3546, 67.8, 134.9, -44.94, 0.2000, 1, 0.6, 0.0141, 4.9887, 0.7595, 32.4654, 3.0041, 41.2318, 3.2808 + 1941, 8, 0.000, 0.000, 14.5, 1.8117, 105.0, 102.9, -17.12, 1.0000, 1, 0.6, 0.0057, 4.8549, 0.7404, 32.4359, 3.0041, 41.0411, 3.4715 + 1941, 9, 1.327, 0.000, 13.9, 1.7178, 25.5, 65.7, -40.89, 0.3297, 1, 0.6, 0.7873, 5.3575, 0.7343, 32.4265, 3.0041, 42.3097, 3.5295 + 1941, 10, 0.000, 0.000, 9.8, 1.0763, 52.5, 29.5, -10.52, 1.0000, 0, 1.0, 0.3211, 5.2153, 0.7643, 32.4532, 3.0041, 41.7579, 4.0813 + 1941, 11, 0.000, 0.000, 6.2, 0.6255, 69.6, 10.2, 0.00, 1.0000, 0, 1.0, 0.1906, 5.1344, 0.7661, 32.4518, 3.0041, 41.5470, 4.2923 + 1941, 12, 0.000, 0.000, 4.9, 0.4914, 42.8, 4.4, 0.00, 1.0000, 0, 1.0, 0.1266, 5.0717, 0.7634, 32.4460, 3.0041, 41.4117, 4.4275 diff --git a/models/rothc/tests/testthat/setup.R b/models/rothc/tests/testthat/setup.R new file mode 100644 index 00000000000..d5829d4d6ad --- /dev/null +++ b/models/rothc/tests/testthat/setup.R @@ -0,0 +1,11 @@ + +# avoid exiting on error +oq <- PEcAn.logger::logger.setQuitOnSevere(FALSE) +withr::defer(PEcAn.logger::logger.setQuitOnSevere(oq), teardown_env()) + +# write to stdout so expect_output() works on logger messages +oc <- PEcAn.logger::logger.setUseConsole(console = TRUE, stderr = FALSE) +withr::defer( + PEcAn.logger::logger.setUseConsole(console = oc$console, stderr = oc$stderr), + teardown_env() +) diff --git a/models/rothc/tests/testthat/test-model2netcdf.RothC.R b/models/rothc/tests/testthat/test-model2netcdf.RothC.R new file mode 100644 index 00000000000..40879195957 --- /dev/null +++ b/models/rothc/tests/testthat/test-model2netcdf.RothC.R @@ -0,0 +1,51 @@ +test_that("model2netcdf", { + outdir <- withr::local_tempdir() + file.copy("data/month_results.out", outdir) + ref_output <- read.csv("data/month_results.out") + + + # happy path + expect_silent( + model2netcdf.RothC( + outdir = outdir, + sitelat = 0, + sitelon = 0, + start_date = "1939-01-01", + end_date = "1941-12-31" + ) + ) + expected_files <- file.path(outdir, c("1939.nc", "1940.nc", "1941.nc")) + expect_true(all(file.exists(expected_files))) + + pred_soc <- PEcAn.utils::read.output( + ncfiles = expected_files, + variables = "TotSoilCarb", + dataframe = TRUE, + verbose = FALSE, + print_summary = FALSE + ) + expect_equal(nrow(pred_soc), 3 * 12) + expect_true(min(pred_soc$posix) == as.Date("1939-01-15")) + expect_true(max(pred_soc$posix) == as.Date("1941-12-15")) + expect_equal( + pred_soc$TotSoilCarb, + ref_output |> + dplyr::filter(.data$Year > 1) |> + dplyr::pull("SOC_t_C_ha") |> + PEcAn.utils::ud_convert("t C ha-1", "kg C m-2") + ) + + + # date ranges enforced + expect_error( + model2netcdf.RothC( + outdir = outdir, + sitelat = 0, + sitelon = 0, + start_date = "1939-01-01", + end_date = "1945-12-31" + ), + "end_date <= max(res$month_ended) is not TRUE", + fixed = TRUE + ) +}) diff --git a/models/rothc/tests/testthat/test.met2model.R b/models/rothc/tests/testthat/test.met2model.R new file mode 100644 index 00000000000..37195568c0e --- /dev/null +++ b/models/rothc/tests/testthat/test.met2model.R @@ -0,0 +1,85 @@ +test_that("met2model", { + + test_dir <- withr::local_tempdir("rothc_met_tmp") + nc_path <- system.file( + "test-data", "CRUNCEP.2000.nc", + package = "PEcAn.utils" + ) + indir <- dirname(nc_path) + + + # full year + res_1 <- met2model.RothC( + in.path = indir, + in.prefix = "CRUNCEP", + outfolder = test_dir, + start_date = "2000-01-01", + end_date = "2000-12-31" + ) + expect_s3_class(res_1, "data.frame") + datfile <- res_1[["file"]][[1]] + expect_true(file.exists(datfile)) + means_full <- read.table(datfile, header = TRUE) + expect_identical(dim(means_full), c(12L, 5L)) + expect_identical(means_full$month, 1:12) + expect_false(anyNA(means_full)) + + + # skips rewriting if overwrite not passed + mtime_1 <- file.mtime(datfile) + expect_output( + met2model.RothC( + indir, "CRUNCEP", test_dir, + "2000-01-01", "2000-12-31" + ), + "already exists, skipping to next file" + ) + expect_equal(file.mtime(datfile), mtime_1) + met2model.RothC( + indir, "CRUNCEP", test_dir, + "2000-01-01", "2000-12-31", + overwrite = TRUE + ) + expect_gt(file.mtime(datfile), mtime_1) + + + # partial year + res_2 <- met2model.RothC( + in.path = indir, + in.prefix = "CRUNCEP", + outfolder = test_dir, + start_date = "2000-03-15", + end_date = "2000-06-01" + ) + datfile <- res_2[["file"]][[1]] + expect_match(datfile, "CRUNCEP.2000-03.2000-06") + means_summer <- read.table(datfile, header = TRUE) + expect_identical(dim(means_summer), c(4L, 5L)) + expect_equal(means_summer, means_full[3:6, ], ignore_attr = "row.names") + + + # all years unavailable + expect_error( + met2model.RothC( + in.path = indir, + in.prefix = "CRUNCEP", + outfolder = test_dir, + start_date = "2002-01-01", + end_date = "2003-12-31" + ), + "No files found" + ) + + # some years unavailable + expect_error( + met2model.RothC( + in.path = indir, + in.prefix = "CRUNCEP", + outfolder = test_dir, + start_date = "2000-01-01", + end_date = "2002-12-31" + ), + "does not cover" + ) + +})