R/smooth_percentile.R

In robwschlegel/heatwaveR: Detect Heatwaves and Cold-Spells

#' Detect the climatology for a time series.
#'
#' An internal function that helps to create climatologies in
#' accordance with the Hobday et al. (2016) standard.
#'
#' @keywords internal
#'
#' @param data The data given to this function during the calculations
#' performed by \code{\link{ts2clm}}.
#' @param smoothPercentileWidth The width of the smoothing window
#' to be applied. The default is \code{31} days.
#' @param var_calc This is passed from the ts2clm argument \code{var}
#' and tells the function if a var column exists that needs to be smoothed.
#'
#' @return The function returns the data in the same format it was
#' input as, with the climatology values smoothed as desired.
#'
#' @author Smit, A. J., Schlegel, R. W.
#'
smooth_percentile <- function(data, smoothPercentileWidth, var_calc) {

  seas <- thresh <- NULL

  # testing...
  # data <- ts_mat

  if("hoy" %in% colnames(data)) {
    prep <- rbind(utils::tail(data[,-c(1, 2)], smoothPercentileWidth),
                  data[,-c(1, 2)],
                  utils::head(data[,-c(1, 2)], smoothPercentileWidth))
  } else {
    prep <- rbind(utils::tail(data[,-1], smoothPercentileWidth),
                  data[,-1],
                  utils::head(data[,-1], smoothPercentileWidth))
  }

  len_clim_year <- 366
  len_hour_clim_year <- 8784

  seas <- RcppRoll::roll_mean(as.numeric(prep[,1]), n = smoothPercentileWidth, na.rm = FALSE)
  thresh <- RcppRoll::roll_mean(as.numeric(prep[,2]), n = smoothPercentileWidth, na.rm = FALSE)

  if("hoy" %in% colnames(data)) {
    clim <- data.table::data.table(doy = rep(seq_len(len_clim_year), each = 24),
                                   hoy = seq_len(len_hour_clim_year),
                                   seas = seas[(smoothPercentileWidth/2 + 2):((smoothPercentileWidth/2 + 1) + len_hour_clim_year)],
                                   thresh = thresh[(smoothPercentileWidth/2 + 2):((smoothPercentileWidth/2 + 1) + len_hour_clim_year)])
  } else {
    clim <- data.table::data.table(doy = seq_len(len_clim_year),
                                   seas = seas[(smoothPercentileWidth/2 + 2):((smoothPercentileWidth/2 + 1) + len_clim_year)],
                                   thresh = thresh[(smoothPercentileWidth/2 + 2):((smoothPercentileWidth/2 + 1) + len_clim_year)])
  }

  if (var_calc) {
    var <- NULL
    var <- RcppRoll::roll_mean(as.numeric(prep[,3]), n = smoothPercentileWidth, na.rm = FALSE)
    if("hoy" %in% colnames(data)) {
      clim$var <- var[(smoothPercentileWidth/2 + 2):((smoothPercentileWidth/2 + 1) + len_hour_clim_year)]
    } else {
      clim$var <- var[(smoothPercentileWidth/2 + 2):((smoothPercentileWidth/2 + 1) + len_clim_year)]
    }
  }
  rm(data, prep)

  return(clim)
}