Skip to content

Commit

Permalink
Merge pull request #276 from geco-bern/275-repeat-last-year-of-forcing
Browse files Browse the repository at this point in the history 
BiomeE: Repeat last year forcing if required by the specified `nyeartrend`
  • Loading branch information
@fabern
fabern authored Dec 12, 2024
2 parents 1dbfe71 + 762a2d6 commit c50d858
Show file tree
Hide file tree
Showing 23 changed files with 318 additions and 302 deletions.
4 changes: 3 additions & 1 deletion NEWS.md
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,7 @@
# rsofun (development version)
# rsofun (development version v5.1.0)

* new BiomeE behavior to recycle last year of forcing if requested simulation time span (`nyeartrend`) is longer than
available forcing data
* Breaking change: biomee drivers' `init_cohort$init_n_cohorts` column has been phased out and must not be present in
drivers to protect against data corruption.

Expand Down
26 changes: 17 additions & 9 deletions R/data.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,7 @@
#' \describe{
#' \item{spinup}{A logical value indicating whether this simulation does spin-up.}
#' \item{spinupyears}{Number of spin-up years.}
#' \item{recycle}{Length of standard recycling period, in years.}
#' \item{recycle}{Number of first N years of forcing data.frame that are recycled for spin-up.}
#' \item{outdt}{An integer indicating the output periodicity.}
#' \item{ltre}{A logical value, \code{TRUE} if evergreen tree.}
#' \item{ltne}{A logical value, \code{TRUE} if evergreen tree and N-fixing.}
Expand Down Expand Up @@ -166,11 +166,15 @@
#' rsofun P-model output data
#'
#' Example output dataset from a p-model run using \code{\link{p_model_drivers}}
#' See \code{\link{run_pmodel_f_bysite}} for a detailed
#' description of the outputs.
"p_model_output"

#' rsofun P-model output data (using vcmax25 drivers)
#'
#' Example output dataset from a p-model run using \code{\link{p_model_drivers_vcmax25}}
#' See \code{\link{run_pmodel_f_bysite}} for a detailed
#' description of the outputs.
"p_model_output_vcmax25"

#' rsofun BiomeE driver data (Leuning photosynthesis model)
Expand All @@ -188,10 +192,10 @@
#' \describe{
#' \item{spinup}{Flag indicating whether this simulation does spin-up.}
#' \item{spinupyears}{Number of spin-up years.}
#' \item{recycle}{Length of standard recycling period (years).}
#' \item{firstyeartrend}{First transient year.}
#' \item{nyeartrend}{Number of transient years.}
#' \item{steps_per_day}{Time resolution (day-1).}
#' \item{recycle}{Number of first N years of forcing data.frame that are recycled for spin-up.}
#' \item{firstyeartrend}{Year of first transient year (AD) (optional). Is only used to set years in output data frames. Defaults to 0 if not provided.}
#' \item{nyeartrend}{Number of transient years (optional). Determines the length of simulation output after spin-up. Defaults to number of years contained in the forcing data. (If longer than forcing data, last year of forcing is repeated until the end (spin-down).)}
#' \item{steps_per_day}{Time resolution of the forcing (day-1).}
#' \item{do_U_shaped_mortality}{Flag indicating whether U-shaped
#' mortality is used.}
#' \item{update_annualLAImax}{Flag indicating whether updating
Expand All @@ -208,9 +212,9 @@
#' \item{site_info}{Site meta info in a data.frame.
#' This data structure can be freely used for documenting the dataset, but must include at least the following data:
#' \describe{
#' \item{lon}{Longitude of the site location.}
#' \item{lat}{Latitude of the site location.}
#' \item{elv}{Elevation of the site location, in meters.}
#' \item{lon}{Longitude of the site location in degrees east.}
#' \item{lat}{Latitude of the site location in degrees north.}
#' \item{elv}{Elevation of the site location, in meters above sea level.}
#' }}
#' \item{forcing}{Forcing data.frame used as input
#' \describe{
Expand Down Expand Up @@ -256,7 +260,7 @@
#' \item{phenotype}{Integer set to 0 for deciduous and 1 for evergreen.}
#' \item{pt}{Integer indicating the type of plant according to photosynthesis:
#' 0 for C3; 1 for C4}
#' \item{alpha_FR}{Fine root turnonver rate (year\eqn{^{-1}}).}
#' \item{alpha_FR}{Fine root turnover rate (year\eqn{^{-1}}).}
#' \item{rho_FR}{Material density of fine roots (kg C m\eqn{^{-3}}).}
#' \item{root_r}{Radius of the fine roots, in m.}
#' \item{root_zeta}{e-folding parameter of root vertical distribution, in m.}
Expand Down Expand Up @@ -355,9 +359,13 @@
#' rsofun BiomeE (P-model) output data
#'
#' Example output dataset from a BiomeE-model run (p-model)
#' See \code{\link{run_biomee_f_bysite}} for a detailed
#' description of the outputs.
"biomee_p_model_output"

#' rsofun BiomeE (gs_leuning) output data
#'
#' Example output dataset from a BiomeE-model run (gs_leuning)
#' See \code{\link{run_biomee_f_bysite}} for a detailed
#' description of the outputs.
"biomee_gs_leuning_output"
9 changes: 9 additions & 0 deletions R/init_dates_dataframe.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -78,6 +78,15 @@ init_dates_dataframe <- function(
))

# convert to decimal date
numeric_year <- function(x){
y <- as.numeric(format(x, format="%Y"))
doy <- as.numeric(format(x, format="%j")) - 1

ifelse(y %% 4 == 0,
round(y + doy/366, 3),
round(y + doy/365, 3)
)
}
date_range$year_dec <- numeric_year(date_range$date)

# leap year filter
Expand Down
162 changes: 110 additions & 52 deletions R/run_biomee_f_bysite.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,23 +4,17 @@
#'
#' @param sitename Site name.
#' @param params_siml Simulation parameters.
#' See examples \code{\link{biomee_gs_leuning_drivers}} or \code{\link{biomee_p_model_drivers}}
#' @param site_info Site meta info in a data.frame.
#' See examples \code{\link{biomee_gs_leuning_drivers}} or \code{\link{biomee_p_model_drivers}}
#' @param forcing Forcing data.frame used as input.
#' See examples \code{\link{biomee_gs_leuning_drivers}} or \code{\link{biomee_p_model_drivers}}
#' @param forcing A data frame of forcing climate data, used as input.
#' @param params_tile Tile-level model parameters, into a single row data.frame.
#' See examples \code{\link{biomee_gs_leuning_drivers}} or \code{\link{biomee_p_model_drivers}}
#' @param params_species A data.frame containing species-specific model parameters,
#' with one species per row. See examples \code{\link{biomee_gs_leuning_drivers}} or \code{\link{biomee_p_model_drivers}}
#' @param init_cohort A data.frame of initial cohort specifications.
#' See examples \code{\link{biomee_gs_leuning_drivers}} or \code{\link{biomee_p_model_drivers}}
#' @param init_soil A data.frame of initial soil pools.
#' See examples \code{\link{biomee_gs_leuning_drivers}} or \code{\link{biomee_p_model_drivers}}
#' @param makecheck Flag specifying whether checks are performed to verify model inputs and parameters.
#' @param makecheck A logical specifying whether checks are performed
#' to verify forcings and model parameters. \code{TRUE} by default.
#'
#' @export
#' @useDynLib rsofun
#' For further specifications of above inputs and examples see \code{\link{biomee_gs_leuning_drivers}} or \code{\link{biomee_p_model_drivers}}
#'
#' @returns Model output is provided as a list, with elements:
#' \describe{
Expand Down Expand Up @@ -225,7 +219,10 @@
#' \item{deathrate}{Mortality rate of this cohort (yr\eqn{^{-1}}).}
#' }}
#' }
#'
#'
#' @export
#' @useDynLib rsofun
#'
#' @examples
#' \donttest{
#' # Example BiomeE model run
Expand Down Expand Up @@ -260,27 +257,32 @@ run_biomee_f_bysite <- function(

# predefine variables for CRAN check compliance
type <- NULL

forcing_features <- c(
'ppfd',
'temp',
'vpd',
'rain',
'wind',
'patm',
'co2'
)

# select relevant columns of the forcing data
forcing <- forcing %>%
select(
any_of(forcing_features)
)


# base state, always execute the call
continue <- TRUE

# record number of years in forcing data
# frame to use as default values (unless provided othrwise as params_siml$nyeartrend)
ndayyear <- 365
forcing_years <- nrow(forcing)/(ndayyear * params_siml$steps_per_day)

`%nin%` <- Negate(`%in%`)
# Default value for nyeartrend
if ('nyeartrend' %nin% names(params_siml)) {
params_siml$nyeartrend <- forcing_years
}
# Default value for firstyeartrend
# If not provided, we anchor to 0, meaning that spinup years are negative and transient years are positive.
# firstyeartrend is currently not used.
if ('firstyeartrend' %nin% names(params_siml)) {
params_siml$firstyeartrend <- 0
}

runyears <- ifelse(
params_siml$spinup,
(params_siml$spinupyears + params_siml$nyeartrend),
params_siml$nyeartrend)
params_siml$nyeartrend
)

n_daily <- params_siml$nyeartrend * 365

Expand Down Expand Up @@ -325,37 +327,93 @@ run_biomee_f_bysite <- function(
params_siml$method_mortality))
}

# base state, always execute the call
continue <- TRUE
# re-define units and naming of forcing dataframe
# keep the order of columns - it's critical for Fortran (reading by column number)
forcing_features <- c(
'ppfd',
'temp',
'vpd',
'rain',
'wind',
'patm',
'co2'
)

# select relevant columns of the forcing data
forcing <- forcing %>%
select(
any_of(forcing_features)
)

if ('init_n_cohorts' %in% names(init_cohort)) {
warning("Warning: Ignoring column 'init_n_cohorts' under 'init_cohort' in drivers. It has been phased out and should be removed from drivers.")
init_cohort <- select(init_cohort, -'init_n_cohorts')
}


# validate input
if (makecheck){
# Add input and parameter checks here if applicable.
data_integrity <- lapply(
forcing_features,
function(check_var){
if (any(is.nanull(forcing[check_var]))){
warning(
sprintf("Error: Missing forcing %s for site %s",
check_var, sitename))
return(FALSE)
} else {
return(TRUE)
}
})

if ('init_n_cohorts' %in% names(init_cohort)) {
warning("Error: column 'init_n_cohorts' under 'init_cohort' has been phased out and must be removed from the drivers.")
data_integrity <- append(data_integrity, FALSE)
is.nanull <- function(x) ifelse(any(is.null(x), is.na(x)), TRUE, FALSE)

if (params_siml$nyeartrend < forcing_years) {
warning(sprintf(
"Info: provided value of nyeartrend is less than the number of years of forcing data (%i). Only the first %i will be used."
, forcing_years, params_siml$nyeartrend))
}
if (params_siml$nyeartrend > forcing_years) {
warning(sprintf(
"Info: provided value of nyeartrend is greater than the number of years of forcing data (%i). The final year will be repeated as much as needed."
, forcing_years))
}

# create a loop to loop over a list of variables
# to check validity
data_integrity <- lapply(forcing_features, function(check_var){
if (any(is.nanull(forcing[check_var]))){
warning(sprintf("Error: Missing forcing %s for site %s",
check_var, sitename))
return(FALSE)
} else {
return(TRUE)
}
})

# only return true if all checked variables are TRUE
# only run simulation if all checked variables are valid
# suppress warning on coercion of list to single logical
continue <- suppressWarnings(all(as.vector(data_integrity)))
if (suppressWarnings(!all(as.vector(data_integrity)))) {
continue <- FALSE
}

# simulation parameters to check
check_param <- c(
"spinup",
"spinupyears",
"recycle",
"firstyeartrend",
"nyeartrend",
"steps_per_day",
"do_U_shaped_mortality",
"update_annualLAImax",
"do_closedN_run"
)
parameter_integrity <- lapply(check_param, function(check_var){
if (any(is.nanull(params_siml[check_var]))){
warning(sprintf("Error: Missing value in %s for %s",
check_var, sitename))
return(FALSE)
} else {
return(TRUE)
}
})

if (suppressWarnings(!all(parameter_integrity))){
continue <- FALSE
}
}

if (continue) {

## C wrapper call
biomeeout <- .Call(

Expand All @@ -367,6 +425,7 @@ run_biomee_f_bysite <- function(
recycle = as.integer(params_siml$recycle),
firstyeartrend = as.integer(params_siml$firstyeartrend),
nyeartrend = as.integer(params_siml$nyeartrend),
steps_per_day = as.integer(params_siml$steps_per_day),
do_U_shaped_mortality = as.logical(params_siml$do_U_shaped_mortality),
update_annualLAImax = as.logical(params_siml$update_annualLAImax),
do_closedN_run = as.logical(params_siml$do_closedN_run),
Expand Down Expand Up @@ -417,8 +476,7 @@ run_biomee_f_bysite <- function(
n_annual = as.integer(runyears),
n_annual_cohorts = as.integer(params_siml$nyeartrend), # to get cohort outputs after spinup year
#n_annual_cohorts = as.integer(runyears), # to get cohort outputs from year 1
forcing = as.matrix(forcing),
steps_per_day = as.integer(params_siml$steps_per_day)
forcing = as.matrix(forcing)
)

# If simulation is very long, output gets massive.
Expand Down
Loading

0 comments on commit c50d858

Please sign in to comment.