R/dm_add_climate.R

In dendRoAnalyst: A Tool for Processing and Analyzing Dendrometer Data

# =========================================================
# Internal helpers for climate attachment and plotting
# =========================================================

.dm_subset_window <- function(dat, time_col, Year = NULL, DOY = NULL) {
  tt <- dat[[time_col]]

  if (!inherits(tt, "Date") && !inherits(tt, "POSIXct")) {
    tt <- as.POSIXct(tt)
  }

  keep <- rep(TRUE, length(tt))

  if (!is.null(Year)) {
    keep <- keep & lubridate::year(tt) %in% Year
  }

  if (!is.null(DOY)) {
    if (length(DOY) != 2) {
      stop("'DOY' must be a numeric vector of length 2.")
    }
    doy_min <- min(DOY)
    doy_max <- max(DOY)
    keep <- keep & lubridate::yday(tt) >= doy_min & lubridate::yday(tt) <= doy_max
  }

  out <- dat[keep, , drop = FALSE]

  if (nrow(out) == 0) {
    stop("No data available for the selected Year/DOY window.")
  }

  out
}

.dm_label_phase <- function(ph, obj_class = c("ZG", "SC")) {
  obj_class <- match.arg(obj_class)

  if (obj_class == "ZG") {
    factor(ph, levels = c(1, 2), labels = c("TWD", "GRO"))
  } else {
    factor(ph, levels = c(1, 2, 3), labels = c("Shrinkage", "Expansion", "Increment"))
  }
}

.dm_detect_cycle_climate_cols <- function(dat, kind = c("ZG", "SC")) {
  kind <- match.arg(kind)

  core <- if (kind == "ZG") {
    c("Phases", "Start", "End", "Duration_h", "Magnitude", "rate",
      "max.twd", "Max.twd.time", "ABr.value", "Avg.twd", "STD.twd", "DOY")
  } else {
    c("Phases", "Start", "End", "Duration_h", "Duration_m", "Magnitude", "rate", "DOY")
  }

  setdiff(names(dat), core)
}

.dm_detect_point_climate_cols <- function(dat, kind = c("ZG", "SC")) {
  kind <- match.arg(kind)

  core <- if (kind == "ZG") {
    c("TIME", "dm", "Phases", "TWD", "GRO")
  } else {
    c("TIME", "dm", "Phases")
  }

  setdiff(names(dat), core)
}


# =========================================================
# 1) Single helper to add climate
# =========================================================

#' @title Add climate information to dendrometer outputs
#'
#' @description
#' A single helper function that adds climate information to outputs from
#' \code{daily.data()}, \code{phase.zg()}, or \code{phase.sc()}.
#'
#' Depending on the class of \code{x} and the selected \code{scale}, the function:
#' \itemize{
#'   \item joins daily climate to \code{daily.data()} output
#'   \item joins climate summaries to phase windows in \code{ZG_cycle} or \code{SC_cycle}
#'   \item joins subdaily climate features to point-level \code{ZG_phase} or \code{SC_phase}
#' }
#'
#' @param x Object returned by \code{daily.data()}, \code{phase.zg()}, or \code{phase.sc()}.
#' @param clim Climate input. This can be:
#'   \itemize{
#'     \item a standardized object returned by \code{read.climate()}
#'     \item a raw climate data frame
#'     \item a valid climate file path readable by \code{read.climate()}
#'     \item a daily climate table returned by \code{dm_daily_clim()}
#'     \item a subdaily climate-feature table returned by \code{dm_subdaily_clim()}
#'   }
#' @param scale Character string controlling how climate is attached.
#'   One of \code{"auto"}, \code{"daily"}, \code{"phase"}, or \code{"subdaily"}.
#'   \itemize{
#'     \item For \code{daily_output}, only \code{"daily"} is supported.
#'     \item For \code{ZG_output} and \code{SC_output}, \code{"phase"} adds climate
#'           summaries to \code{ZG_cycle}/\code{SC_cycle}, while \code{"subdaily"}
#'           adds timestamp-level climate to \code{ZG_phase}/\code{SC_phase}.
#'   }
#' @param mean_vars,min_vars,max_vars,sum_vars,median_vars Climate variables to
#'   summarize by mean, minimum, maximum, sum, or median.
#'   These are used by \code{dm_daily_clim()} and \code{dm_join_phase_clim()}.
#' @param lag_vars,lagmean_vars,lagsum_vars Variables used by \code{dm_daily_clim()}
#'   to build lagged and antecedent daily climate features.
#' @param lag_days Integer vector of lag/antecedent windows in days for
#'   \code{dm_daily_clim()}.
#' @param sub_mean_vars,sub_sum_vars,sub_lag_vars Variables used by
#'   \code{dm_subdaily_clim()} to build rolling and lagged subdaily features.
#' @param roll_hours,lag_hours Numeric vectors of rolling-window and lag sizes
#'   in hours for \code{dm_subdaily_clim()}.
#' @param suffix Suffix appended to climate summaries added by
#'   \code{dm_join_phase_clim()}.
#'
#' @return
#' The same biological object with climate information added.
#'
#' @examples
#' \donttest{
#' data(nepa17)
#' data(gf_nepa17)
#' data(ktm_clim_hourly)
#'
#' # daily.data() output + daily climate
#' dd <- daily.data(df = nepa17[1:1000, ], TreeNum = 1)
#' dd_clim <- dm_add_climate(
#'   dd,
#'   ktm_clim_hourly,
#'   scale = "daily",
#'   mean_vars = c("temp", "VPD", "RH"),
#'   max_vars  = c("temp", "VPD"),
#'   sum_vars  = c("prec")
#' )
#' head(dd_clim)
#'
#' # phase.zg() output + phase-window climate
#' zg <- phase.zg(df = gf_nepa17[1:800, ], TreeNum = 1)
#' zg_clim <- dm_add_climate(
#'   zg,
#'   ktm_clim_hourly,
#'   scale = "phase",
#'   mean_vars = c("temp", "VPD", "RH"),
#'   max_vars  = c("temp", "VPD"),
#'   sum_vars  = c("prec")
#' )
#' head(zg_clim$ZG_cycle)
#'
#' # phase.sc() output + point-level subdaily climate
#' sc <- phase.sc(df = gf_nepa17[1:800, ], TreeNum = 1, smoothing = 12)
#' sc_point <- dm_add_climate(
#'   sc,
#'   ktm_clim_hourly,
#'   scale = "subdaily",
#'   sub_mean_vars = c("temp", "VPD", "RH"),
#'   sub_sum_vars  = c("prec"),
#'   sub_lag_vars  = c("temp", "VPD", "RH"),
#'   roll_hours = c(1, 3, 6, 24),
#'   lag_hours  = c(1, 3, 6, 24)
#' )
#' head(sc_point$SC_phase)
#' }
#'
#' @export
dm_add_climate <- function(x,
                           clim,
                           scale = c("auto", "daily", "phase", "subdaily"),
                           mean_vars = NULL,
                           min_vars = NULL,
                           max_vars = NULL,
                           sum_vars = NULL,
                           median_vars = NULL,
                           lag_vars = NULL,
                           lagmean_vars = NULL,
                           lagsum_vars = NULL,
                           lag_days = c(1, 3, 7),
                           sub_mean_vars = NULL,
                           sub_sum_vars = NULL,
                           sub_lag_vars = NULL,
                           roll_hours = c(3, 6, 24),
                           lag_hours = c(1, 3, 6, 24),
                           suffix = "_phase") {

  scale <- match.arg(scale)

  # -----------------------------------------------------
  # daily.data() output
  # -----------------------------------------------------
  if (inherits(x, "daily_output")) {
    if (scale %in% c("phase", "subdaily")) {
      stop("For objects from daily.data(), only scale = 'daily' or 'auto' is supported.")
    }

    if (inherits(clim, "daily_clim")) {
      return(dm_join_daily_clim(x, clim))
    }

    clim_day <- dm_daily_clim(
      clim,
      mean_vars = mean_vars,
      min_vars = min_vars,
      max_vars = max_vars,
      sum_vars = sum_vars,
      median_vars = median_vars,
      lag_vars = lag_vars,
      lagmean_vars = lagmean_vars,
      lagsum_vars = lagsum_vars,
      lag_days = lag_days
    )

    return(dm_join_daily_clim(x, clim_day))
  }

  # -----------------------------------------------------
  # phase.zg() / phase.sc() output
  # -----------------------------------------------------
  if (inherits(x, "ZG_output") || inherits(x, "SC_output")) {

    if (scale == "auto") {
      if (inherits(clim, "subdaily_clim")) {
        scale <- "subdaily"
      } else {
        scale <- "phase"
      }
    }

    if (scale == "daily") {
      stop("For phase.zg() and phase.sc() outputs, use scale = 'phase' or scale = 'subdaily'.")
    }

    if (scale == "phase") {
      return(
        dm_join_phase_clim(
          x,
          clim,
          mean_vars = mean_vars,
          min_vars = min_vars,
          max_vars = max_vars,
          sum_vars = sum_vars,
          median_vars = median_vars,
          suffix = suffix
        )
      )
    }

    if (scale == "subdaily") {
      if (inherits(clim, "subdaily_clim")) {
        return(dm_join_subdaily_clim(x, clim))
      }

      clim_sub <- dm_subdaily_clim(
        clim,
        mean_vars = sub_mean_vars,
        sum_vars = sub_sum_vars,
        lag_vars = sub_lag_vars,
        roll_hours = roll_hours,
        lag_hours = lag_hours
      )

      return(dm_join_subdaily_clim(x, clim_sub))
    }
  }

  stop(
    "'x' must be an object returned by daily.data(), phase.zg(), or phase.sc()."
  )
}