R/smooth.R

Defines functions wav_signal.wthr wav_signal.ebthr_modwt wav_signal.ebthr_dwt wav_signal.default wav_signal

# #' Wavelet-Based Signal Estimation
# #'
# #' This function is a generic function that estimates a signal by using
# #'  several methods of wavelet shirinkage and thresholding.
# #'
# #' @param x An object generated by the function \code{\link{wav_args}}.
# #' @param ... additional parameters
# #'
# #' @return A vector
wav_signal <- function(x, ...) {
  UseMethod("wav_signal")
}

wav_signal.default <- function(x, ...) {
  wmtsa::wavShrink(x, ...)
}

wav_signal.ebthr_dwt <- function(x, ...) {

  x_nobs <- length(x)

  if (!(log2(x_nobs) %% 1 == 0)) {
    x <- WiSEBoot::padVector(x, pad.direction = "rear", replaceLinearTrend = TRUE)[[1]]
  }

  args <- list(...)

  dwtformals <- names(formals(waveslim::dwt))
  ebthrformals <- names(formals(EbayesThresh::ebayesthresh.wavelet))

  dwtlist <- args[names(args) %in% dwtformals]
  dwtargs <- c(list(x = x), dwtlist)
  xdwt <- do.call(waveslim::dwt, dwtargs)

  ebthrlist <- args[!(names(args) %in% dwtformals)]
  ebthrargs <- c(list(xtr = xdwt), ebthrlist)
  xebthr <- do.call(EbayesThresh::ebayesthresh.wavelet, ebthrargs)

  waveslim::idwt(xebthr)[1:x_nobs]

}

wav_signal.ebthr_modwt <- function(x, ...) {
  args <- list(...)

  modwtformals <- names(formals(waveslim::modwt))
  ebthrformals <- names(formals(EbayesThresh::ebayesthresh.wavelet))

  modwtlist <- args[names(args) %in% modwtformals]
  modwtargs <- c(list(x = x), modwtlist)
  xmodwt <- do.call(waveslim::modwt, modwtargs)

  ebthrlist <- args[!(names(args) %in% modwtformals)]
  ebthrargs <- c(list(xtr = xmodwt), ebthrlist)
  xebthr <- do.call(EbayesThresh::ebayesthresh.wavelet, ebthrargs)

  waveslim::imodwt(xebthr)
}

wav_signal.wthr <- function(x, ...) {

  x_nobs <- length(x)

  if (!(log2(x_nobs) %% 1 == 0)) {
    x <- WiSEBoot::padVector(x, pad.direction = "rear", replaceLinearTrend = TRUE)[[1]]
  }

  args <- list(...)

  wdformals <- names(formals(wavethresh::wd))
  thrformals <- names(formals(wavethresh::threshold.wd))

  type_adjusted <- function(list, fun = "wd") {
    type <- list[["type"]]
    wd_list <- list[!(names(list) == "type")]
    if (fun == "wd") {
      type_fun <-  type[type %in% c("wavelet", "station")]
    } else {
      type_fun <-  type[type %in% c("hard", "soft")]
    }
    c(wd_list, type = type_fun)
  }

  wdlist <- args[names(args) %in% wdformals]

  if (!is.null(wdlist$type)) {
    wdlist <- type_adjusted(wdlist, "wd")
  }

  wdargs <- c(list(data = x), wdlist)
  xwd <- do.call(wavethresh::wd, wdargs)

  thrlist <- args[names(args) %in% thrformals]

  if (!is.null(thrlist$type)) {
    thrlist <- type_adjusted(thrlist, "thr")
  }

  if (is.null(wdlist$type) || wdlist$type == "wavelet") {
    thrargs <- c(list(wd = xwd), thrlist)
    xthr <- do.call(wavethresh::threshold, thrargs)
    wavethresh::wr(xthr)[1:x_nobs]
  } else {
    thrargs <- c(list(wst = wavethresh::convert(xwd)), thrlist)
    xthr <- do.call(wavethresh::threshold, thrargs)
    wavethresh::AvBasis(xthr)[1:x_nobs]
  }
}
nelson16silva/wavsigmap documentation built on March 7, 2023, 10:45 a.m.