Nothing
#' Calculate the spectrum of a time series
#'
#' The power spectral density estimate of the time series is calculated using
#' different approaches.
#'
#' @param data \code{eseis} object, \code{numeric} vector or list of
#' objects, data set to be processed.
#'
#' @param dt \code{Numeric} value, sampling period. If omitted, \code{dt}
#' is set to 1/200. Only needed if \code{data} is no \code{eseis} object.
#'
#' @param method \code{Character} value, calculation method. One out of
#' \code{"periodogram"} and \code{"autoregressive"}.
#' default is \code{"periodogram"}.
#'
#' @param n \code{Numeric} value, optional number of samples in
#' running window used for smoothing the spectrogram. Only applied if a
#' number is provided. Smoothing is performed as running mean.
#'
#' @param \dots Additional arguments passed to the function.
#'
#' @return \code{Data frame} with frequency and power vector
#'
#' @author Michael Dietze
#'
#' @keywords eseis
#'
#' @examples
#'
#' ## load example data set
#' data(rockfall)
#'
#' ## calculate spectrum with standard setup
#' s <- signal_spectrum(data = rockfall_eseis)
#'
#' ## plot spectrum
#' plot_spectrum(data = s)
#'
#' @export signal_spectrum
#'
signal_spectrum <- function(
data,
dt,
method = "periodogram",
n,
...
) {
## check/set dt
if(missing(dt) == TRUE && class(data)[1] != "eseis") {
## try to get object class
class_data <- try(class(data[[1]])[1] == "eseis",
silent = TRUE)
if(class(class_data)[1] == "try-error " | class_data == FALSE) {
warning("Sampling frequency missing! Set to 1/200")
dt <- 1 / 200
}
} else if(missing(dt) == TRUE){
dt <- NULL
}
## check/set dt
if(missing(n) == TRUE) {
n <- NULL
}
## check data structure
if(class(data)[1] == "list") {
## apply function to list
data_out <- lapply(X = data,
FUN = signal_spectrum,
dt = dt,
method = method,
n = n,
...)
## return output
return(data_out)
} else {
## get start time
eseis_t_0 <- Sys.time()
## collect function arguments
eseis_arguments <- list(data = "",
dt = dt,
method = method)
## check if input object is of class eseis
if(class(data)[1] == "eseis") {
## set eseis flag
eseis_class <- TRUE
## store initial object
eseis_data <- data
## extract signal vector
data <- eseis_data$signal
## extract sampling period
dt <- eseis_data$meta$dt
} else {
## set eseis flag
eseis_class <- FALSE
}
if(method == "periodogram") {
## calculate spectrum
s <- spec.pgram(x = data,
plot = FALSE,
...)
## recalculate frequency vector
s$freq <- seq(from = 0,
to = 1 / dt / 2,
length.out = length(s$freq))
## recompose data set
data_out <- data.frame(frequency = s$freq,
power = s$spec)
} else if(method == "autoregressive") {
## calculate spectrum
s <- spec.ar(x = data,
plot = FALSE,
...)
## recalculate frequency vector
s$freq <- seq(from = 0,
to = 1 / dt / 2,
length.out = length(s$freq))
## recompose data set
data_out <- data.frame(frequency = s$freq,
power = s$spec)
} else {
stop("Keyword for method not supported!")
}
## optionally perform smoothing
if(is.null(n) == FALSE) {
data_out$power <- caTools::runmean(x = data_out$power, k = n)
}
## optionally rebuild eseis object
if(eseis_class == TRUE) {
## assign aggregated signal vector
eseis_data$signal <- data_out
## rename output list element
names(eseis_data)[1] <- "spectrum"
## calculate function call duration
eseis_duration <- as.numeric(difftime(time1 = Sys.time(),
time2 = eseis_t_0,
units = "secs"))
## update object history
eseis_data$history[[length(eseis_data$history) + 1]] <-
list(time = Sys.time(),
call = "signal_spectrum()",
arguments = eseis_arguments,
duration = eseis_duration)
names(eseis_data$history)[length(eseis_data$history)] <-
as.character(length(eseis_data$history))
## update data type
eseis_data$meta$type = "spectrum"
## assign eseis object to output data set
data_out <- eseis_data
}
## return output
return(data_out)
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.