diff --git a/NAMESPACE b/NAMESPACE index 072919dafaa86df497c4a62e8cbabcfedf1d4e83..54a74e7e68d411bb1060488800e9a0ebbbab3a94 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,17 +1,19 @@ # Generated by roxygen2: do not edit by hand -S3method(as.data.frame,mfft_deco) -S3method(lines,mfft_deco) +S3method(as.data.frame,discreteSpectrum) +S3method(develop,discreteSpectrum) +S3method(lines,discreteSpectrum) S3method(plot,SSAObject) +S3method(plot,discreteSpectrum) S3method(plot,memObject) -S3method(plot,mfft_deco) S3method(plot,periodogram) S3method(plot,wavelet) -S3method(print,mfft_deco) +S3method(print,discreteSpectrum) export(approx_ts) export(arspec) export(attributeTones) export(cwt_morlet) +export(develop) export(hilbert_extension) export(mem) export(mfft) @@ -21,7 +23,6 @@ export(mfft_real) export(mmfft) export(periodogram) export(powerspectrum.wavelet) -export(reconstruct_mfft) export(reconstruct_morlet) export(ssa) export(toneCombinations) diff --git a/R/data_description.R b/R/data_description.R index 787a357c98109cf69173e34612d1f2df04886c9e..1255adf364ca887bf4a6d975d27c3810b77a444f 100644 --- a/R/data_description.R +++ b/R/data_description.R @@ -2,7 +2,7 @@ #' #' A simple time series for testing purposes #' cos(t*0.13423167+0.00) + 1.3 * cos(t*0.119432+2.314) + 0.134994 + 0.4*cos(t*0.653167) + 0.11 * cos(t*0.78913498) -#' for t from 0 to 1023. It is coded as a `mfft_deco` objcet and self-describes its spectral decomposition. +#' for t from 0 to 1023. It is coded as a `discreteSpectrum` objcet and self-describes its spectral decomposition. #' #' @docType data #' @keywords datasets @@ -14,7 +14,7 @@ NULL #' #' Same has harmonic sample but with an added random time series #' cos(t*0.13423167+0.00) + 1.3 * cos(t*0.119432+2.314) + 0.134994 + 0.4*cos(t*0.653167) + 0.11 * cos(t*0.78913498) + rnorm(1024)*0.12 -#' for t from 0 to 1023. It is coded as a `mfft_deco` objcet and self-describes its spectral decomposition. +#' for t from 0 to 1023. It is coded as a `discreteSpectrum` objcet and self-describes its spectral decomposition. #' #' @docType data #' @keywords datasets @@ -26,7 +26,7 @@ NULL #' #' A truncated version (squared-shape) version of a harmonic signal for mimicking a sedimentary log #' sign( cos(t*0.13423167+0.00) + 1.3 * cos(t*0.119432+2.314) + 0.134994 + 0.4*cos(t*0.653167) + 0.11 * cos(t*0.78913498) + 0.8)) -#' for t from 0 to 8191. It is coded as a `mfft_deco` objcet and self-describes its spectral decomposition. +#' for t from 0 to 8191. It is coded as a `discreteSpectrum` objcet and self-describes its spectral decomposition. #' @docType data #' @keywords datasets #' @name pseudo_log diff --git a/R/mfft_complex.R b/R/mfft_complex.R index d8fedd769d5dbee36e174ecf0a40613655c07a82..78a3305bf4f39d8631fe8de7a0945aaa0a9e64ca 100644 --- a/R/mfft_complex.R +++ b/R/mfft_complex.R @@ -22,7 +22,7 @@ #' not documented for real time series #' @param nfreq is the number of frequencies returned, must be smaller that the length of xdata. #' @param force_complex : use the complex number implementation even if the time series is real. -#' @return a `mfft_deco` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". +#' @return a `discreteSpectrum` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". #' note that because of a language glitch (to be fixed), "Freq" actually means "Rate" #' @author Michel Crucifix for the R code, and David Nesvorny for most of the supporting C code doing the #' Frequency Modified Fourier transform for Complex Numbers @@ -56,7 +56,7 @@ mfft <- function(xdata, nfreq=15, minfreq=NULL, maxfreq=NULL, correction=1, forc #' the third algorithm should be in general much more precise). #' The computed frequencies are in the range given by minfreq and maxfreq. #' @param nfreq is the number of frequencies returned, must be smaller that the length of xdata. -#' @return a `mfft_deco` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". +#' @return a `discreteSpectrum` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". #' note that because of a language glitch (to be fixed), "Freq" actually means "Rate" #' @author Michel Crucifix for the R code, and David Nesvorny for most of the supporting C code doing the #' actual computations @@ -104,7 +104,7 @@ mfft_complex <- function(xdata, nfreq=30, minfreq=NULL, maxfreq=NULL, correctio # corresponding to the same amplitude OUT <- data.frame(Freq=Freq, Amp=Ampl, Phases=Phase) - class(OUT) <- c("mfft_deco", "data.frame") + class(OUT) <- c("discreteSpectrum", "data.frame") attr(OUT,"nfreq") <- nfreq attr(OUT,"data") <- xdata diff --git a/R/mfft_real.R b/R/mfft_real.R index b8f60c3f6e7f6536c55045f22366f011b647f08c..2f9fc873d5f7bce45e9b5049b263bf54133d408d 100644 --- a/R/mfft_real.R +++ b/R/mfft_real.R @@ -342,15 +342,15 @@ mfft_analyse <- function(xdata, nfreq, fast = TRUE, nu = NULL, minfreq=NULL, max #' @param nfreq is the number of frequencies returned, must be smaller that the length of xdata. #' @param fast (default = TRUE) uses analytical formulations for the crossproducts involving sines and cosines. #' note: this is not really faster because the bottleneck is actually the goden section search. But more elegant. -#' @return a `mfft_deco` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". +#' @return a `discreteSpectrum` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". #' @author Michel Crucifix #' @references #' \insertRef{sidlichovsky97aa}{gtseries} #' @examples #' #' data(harmonic_sample) -#' spectrum <- mfft_real(harmonic_sample$data) -#' print(spectrum) +#' spec <- mfft_real(harmonic_sample$data) +#' print(spec) #' #' @export mfft_real mfft_real <- function(xdata, nfreq=5, minfreq=NULL, maxfreq=NULL, correction = 1 , fast=TRUE){ @@ -430,7 +430,7 @@ mfft_real <- function(xdata, nfreq=5, minfreq=NULL, maxfreq=NULL, correction = # rename for class compatibility names(OUT) <- c("Freq","Amp","Phases") - class(OUT) <- c("mfft_deco", "data.frame") + class(OUT) <- c("discreteSpectrum", "data.frame") attr(OUT, "data") <- xdata attr(OUT, "nfreq") <- nfreq return(OUT) diff --git a/R/mfft_support.R b/R/mfft_support.R index cb8dfc7839cb3f555f7e9863593c11e5cfbf7615..b468684aa561f55a10330b8831e989ded45da334 100644 --- a/R/mfft_support.R +++ b/R/mfft_support.R @@ -1,27 +1,54 @@ #' MFFT reconstruction #' @importFrom stats deltat start -#' @rdname mfft_deco -#' @param M : mfft_deco object -#' @param sum : TRUE if user wants to sum components in the reconstruction -#' @export reconstruct_mfft +#' @param M : discreteSpectrum object +#' @param times: if supplied, times of the decomposition +#' @param start: if supplied, overrides time and will generate a time series with start and deltat, which must then +#' be supplied as well +#' @param deltat : see start. +#' @param sum : TRUE if user wants to sum components %in% the reconstruction +#' @note if none if times, start and deltat are supplied, will reconstruct based on the attribute `xdata` +#' which must then be present. If no `xdata` is availble, return an error. #' @return list of reconstructed components if sum=FALSE, full #' reconstructed time series otherwise -reconstruct_mfft <- function(M, sum=TRUE){ - if (!(attr(M,"class")=="mfft_deco")) stop ("object is not a MFFT decomposition") +#' @method develop discreteSpectrum +#' @export +develop.discreteSpectrum <- function(M, times=NULL, start=NULL, end=NULL, deltat=NULL, sum=TRUE){ + if (!("discreteSpectrum" %in% class(M))) stop ("object is not a discreteSpectrum decomposition") + + timesIsATseries = FALSE + if (!is.null(start)){ + if (is.null(deltat) || is.null(end)) stop ("if you supply start, you must also supply deltat and end"); + times <- start + seq((end - start) %/% deltat) * deltat + timesIsATseries = TRUE + } + + if (is.null(times)){ + if (is.null(attr(M,"data"))) stop ("if you do not supply any time argument (times, or (start, end, deltat)), then object must have a valid data attribute") xdata <- attr(M,"data") + start <- stats::start(xdata) + deltat <- stats::deltat(deltat) + times <- (seq(length(xdata))-1) * stats::deltat(xdata) + stats::start(xdata)[1] + timesIsATseries = TRUE + } + nfreq <- attr(M,"nfreq") - times <- seq(length(xdata))*stats::deltat(xdata) + stats::start(xdata)[1] - reconstructed <- lapply(seq(nfreq), function(i) ts( M$Amp[i]*cos(M$Freq[i]*times + M$Phase[i]), start=stats::start(xdata), deltat=stats::deltat(xdata)) ) + if (is.null(nfreq)) nfreq <- length(M$Amp) + if (timesIsATseries){ + reconstructed <- lapply(seq(nfreq), function(i) ts( M$Amp[i] * cos(M$Freq[i] * times + M$Phase[i]), start=start, deltat=deltat) )} + else { + reconstructed <- lapply(seq(nfreq), function(i) M$Amp[i] * cos(M$Freq[i] * times + M$Phase[i])) + } if ( sum ) reconstructed <- ts(apply(simplify2array(reconstructed), 1 , sum), start=stats::start(xdata), deltat=stats::deltat(xdata)) + return(reconstructed) } #' MFFT ANOVA #' not ready. do not use. -#' @rdname mfft_deco +#' @rdname discreteSpectrum #' @export mfft_anova mfft_anova <- function(M){ - if (!(attr(M,"class")=="mfft_deco")) stop ("object is not a MFFT decomposition") + if (!("discreteSpectrum" %in% class(M))) stop ("object is not a discreteSpectrum decomposition") xdata <- attr(M,"xdata") nfreq <- attr(M,"nfreq") N <- length(xdata) @@ -37,17 +64,17 @@ mfft_anova <- function(M){ } -#' @rdname mfft_deco +#' @rdname discreteSpectrum #' @export -as.data.frame.mfft_deco <- function(x) {data.frame(Freq=x$Freq, Amp=x$Amp, Phases=x$Phases)} +as.data.frame.discreteSpectrum <- function(x) {data.frame(Freq=x$Freq, Amp=x$Amp, Phases=x$Phases)} -#' @rdname mfft_deco -#' @param a `mfft_deco` object, typically the output of a `mfft` call. +#' @rdname discreteSpectrum +#' @param a `discreteSpectrum` object, typically the output of a `mfft` call. #' @param labels to be set above the frequency peaks. Can be the output of `attributeTone` #' @param periods if TRUE will add a lower axis with period labels #' @export -plot.mfft_deco <- function (M,periods=FALSE,labels=NULL,...){ +plot.discreteSpectrum <- function (M,periods=FALSE,labels=NULL,...){ # O <- order(M$Freq) plot(abs(M$Freq), abs(M$Amp),'h',ylab="Amplitudes", xlab="", ...) if (periods) { @@ -60,25 +87,25 @@ plot.mfft_deco <- function (M,periods=FALSE,labels=NULL,...){ } else { mtext("Rate", 1, 3) } - points(abs(M$Freq), abs(M$Amp),'p',...) + # points(abs(M$Freq), abs(M$Amp),'p',...) if (!is.null(labels)) { yshift <- 0.05*diff(range(M$Amp)) text(M$Freq, M$Amp, labels, srt=90, , adj=-0.4) } } -#' @rdname mfft_deco +#' @rdname discreteSpectrum #' @export -lines.mfft_deco <- function (M,...){ +lines.discreteSpectrum <- function (M,...){ # O <- order(M$Freq) lines(abs(M$Freq), abs(M$Amp),'h',...) points(abs(M$Freq), abs(M$Amp),'p',...) } -#' @rdname mfft_deco +#' @rdname discreteSpectrum #' @export -print.mfft_deco <- function (M,...){ +print.discreteSpectrum <- function (M,...){ print.data.frame(cbind(as.data.frame(M), Period=2*pi/M$Freq)) } diff --git a/data-raw/harmonic_sample.R b/data-raw/harmonic_sample.R index bebfb49c5b80d71e55b4716879ad4b31e6a0e9b7..4c2eaff0b84044dbe153fdfceca98435cdb9d51f 100644 --- a/data-raw/harmonic_sample.R +++ b/data-raw/harmonic_sample.R @@ -10,12 +10,10 @@ harmonic_sample_spectrum <- list( Freq = c(0.13423167, 0.119432, 0, 0.653167, 0.78913498), Phases = c(0, 2.314, 0, 0.653167, 0)) -attr(harmonic_sample_spectrum,"class") <- 'mfft_deco' - -harmonic_sample_noisy$data <- harmonic_sample_data -harmonic_sample_noisy$spectrum <- harmonic_sample_spectrum - +class(harmonic_sample_spectrum) <- 'discreteSpectrum' +harmonic_sample_noisy <- harmonic_sample_spectrum +attr(harmonic_sample_noisy,"data") <- harmonic_sample_data usethis::use_data(harmonic_sample_noisy, overwrite = TRUE) @@ -29,10 +27,9 @@ harmonic_sample_spectrum <- list( Freq = c(0.13423167, 0.119432, 0, 0.653167, 0.78913498), Phases = c(0, 2.314, 0, 0.653167, 0)) -attr(harmonic_sample_spectrum,"class") <- 'mfft_deco' - -harmonic_sample$data <- harmonic_sample_data -harmonic_sample$spectrum <- harmonic_sample_spectrum +class(harmonic_sample_spectrum) <- 'discreteSpectrum' +harmonic_sample <- harmonic_sample_spectrum +attr(harmonic_sample,"data") <- harmonic_sample_data usethis::use_data(harmonic_sample, overwrite = TRUE) diff --git a/data-raw/pseudo_log.R b/data-raw/pseudo_log.R index 87bfda7b03504d6b5a84121c240f5de3cd259e7f..b9af3bb12499a32041c37d546a19b1e0472c0cdc 100644 --- a/data-raw/pseudo_log.R +++ b/data-raw/pseudo_log.R @@ -13,10 +13,10 @@ pseudo_log_spectrum <- list( Phases = c(0, 2.314, 0, 0.653167, 0)) -attr(pseudo_log_spectrum,"class") <- 'mfft_deco' +class(pseudo_log_spectrum) <- 'discreteSpectrum' -pseudo_log$data <- pseudo_log_data -pseudo_log$spectrum <- pseudo_log_spectrum +pseudo_log_spectrum <- pseudo_log_spectrum +attr(pseudo_log,"data") <- pseudo_log_data usethis::use_data(pseudo_log, overwrite = TRUE) diff --git a/data/harmonic_sample.rda b/data/harmonic_sample.rda index aef1bb47e236ae0c4327aa3f946e5de953cbf3a7..9ed91e93363079b3f9c90940999733c0fe462705 100644 Binary files a/data/harmonic_sample.rda and b/data/harmonic_sample.rda differ diff --git a/data/harmonic_sample_noisy.rda b/data/harmonic_sample_noisy.rda index 7c8b13fda68bcc9b72f28963720e56a215cc48a7..5e41cdae53d231f3a67cd6bfede311246d95af15 100644 Binary files a/data/harmonic_sample_noisy.rda and b/data/harmonic_sample_noisy.rda differ diff --git a/data/pseudo_log.rda b/data/pseudo_log.rda index 8c9055c8e5cc584113d980757aca23b3183f6f7a..e3c8be88bda306321c36b06b10da5bb57b7027ec 100644 Binary files a/data/pseudo_log.rda and b/data/pseudo_log.rda differ diff --git a/man/harmonic_sample.Rd b/man/harmonic_sample.Rd index beaa588f8f5c64ae7b3da76243c789664baaad94..21e0e93968f7a1b9e558c01218495775402e02c1 100644 --- a/man/harmonic_sample.Rd +++ b/man/harmonic_sample.Rd @@ -10,6 +10,6 @@ data(harmonic_sample) \description{ A simple time series for testing purposes cos(t*0.13423167+0.00) + 1.3 * cos(t*0.119432+2.314) + 0.134994 + 0.4*cos(t*0.653167) + 0.11 * cos(t*0.78913498) -for t from 0 to 1023. It is coded as a `mfft_deco` objcet and self-describes its spectral decomposition. +for t from 0 to 1023. It is coded as a `discreteSpectrum` objcet and self-describes its spectral decomposition. } \keyword{datasets} diff --git a/man/harmonic_sample_noisy.Rd b/man/harmonic_sample_noisy.Rd index cc57b10eefc420891cd9ebac442bc7fa47b7115d..adda2d507109568dc616b2e458711be621510357 100644 --- a/man/harmonic_sample_noisy.Rd +++ b/man/harmonic_sample_noisy.Rd @@ -10,6 +10,6 @@ data(harmonic_sample_noisy) \description{ Same has harmonic sample but with an added random time series cos(t*0.13423167+0.00) + 1.3 * cos(t*0.119432+2.314) + 0.134994 + 0.4*cos(t*0.653167) + 0.11 * cos(t*0.78913498) + rnorm(1024)*0.12 -for t from 0 to 1023. It is coded as a `mfft_deco` objcet and self-describes its spectral decomposition. +for t from 0 to 1023. It is coded as a `discreteSpectrum` objcet and self-describes its spectral decomposition. } \keyword{datasets} diff --git a/man/mfft.Rd b/man/mfft.Rd index 5d229ac6cb71dd630d397f61d14ce365ed94b7a9..6ab12cd6cd083da11096626fb77a6c4a973965a0 100644 --- a/man/mfft.Rd +++ b/man/mfft.Rd @@ -31,7 +31,7 @@ The computed frequencies are in the range given by minfreq and maxfreq.} not documented for real time series} } \value{ -a `mfft_deco` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". +a `discreteSpectrum` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". note that because of a language glitch (to be fixed), "Freq" actually means "Rate" } \description{ diff --git a/man/mfft_complex.Rd b/man/mfft_complex.Rd index 0afdc3b2a1ad597e9b41fc6e7bde0de3d8a2a463..92d92014ca326bc96ddd4d7f1d89f629e075dc8e 100644 --- a/man/mfft_complex.Rd +++ b/man/mfft_complex.Rd @@ -25,7 +25,7 @@ the third algorithm should be in general much more precise). The computed frequencies are in the range given by minfreq and maxfreq.} } \value{ -a `mfft_deco` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". +a `discreteSpectrum` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". note that because of a language glitch (to be fixed), "Freq" actually means "Rate" } \description{ diff --git a/man/mfft_deco.Rd b/man/mfft_deco.Rd deleted file mode 100644 index b85616c3c53aa4e8656de358b086addb15fef44e..0000000000000000000000000000000000000000 --- a/man/mfft_deco.Rd +++ /dev/null @@ -1,44 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/mfft_support.R -\name{reconstruct_mfft} -\alias{reconstruct_mfft} -\alias{mfft_anova} -\alias{as.data.frame.mfft_deco} -\alias{plot.mfft_deco} -\alias{lines.mfft_deco} -\alias{print.mfft_deco} -\title{MFFT reconstruction} -\usage{ -reconstruct_mfft(M, sum = TRUE) - -mfft_anova(M) - -\method{as.data.frame}{mfft_deco}(x) - -\method{plot}{mfft_deco}(M, periods = FALSE, labels = NULL, ...) - -\method{lines}{mfft_deco}(M, ...) - -\method{print}{mfft_deco}(M, ...) -} -\arguments{ -\item{M}{: mfft_deco object} - -\item{sum}{: TRUE if user wants to sum components in the reconstruction} - -\item{periods}{if TRUE will add a lower axis with period labels} - -\item{labels}{to be set above the frequency peaks. Can be the output of `attributeTone`} - -\item{a}{`mfft_deco` object, typically the output of a `mfft` call.} -} -\value{ -list of reconstructed components if sum=FALSE, full - reconstructed time series otherwise -} -\description{ -MFFT reconstruction - -MFFT ANOVA -not ready. do not use. -} diff --git a/man/mfft_real.Rd b/man/mfft_real.Rd index 791d76b3b0be396784f60e261d2a21a9f56d01d4..4a61373325f7a3ee819d7c881f75eb8b3b858e91 100644 --- a/man/mfft_real.Rd +++ b/man/mfft_real.Rd @@ -29,7 +29,7 @@ note: this is not really faster because the bottleneck is actually the goden sec 3: second order-correction using synthetic data (all documented in the Sidlichovsky and Nesvorny reference)} } \value{ -a `mfft_deco` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". +a `discreteSpectrum` object, based on a data.frame with columns "Freq", "Ampl" and "Phases". } \description{ R-coded version of the Modified Fourier Transform @@ -43,8 +43,8 @@ A C-version should be supplied one day. \examples{ data(harmonic_sample) -spectrum <- mfft_real(harmonic_sample$data) -print(spectrum) +spec <- mfft_real(harmonic_sample$data) +print(spec) } \references{ diff --git a/man/pseudo_log.Rd b/man/pseudo_log.Rd index 590aef4459912c92e569115cbc9bff79d6bdb36d..83842a9cbbd91b79c1ccb65409ed466bfd93c9b7 100644 --- a/man/pseudo_log.Rd +++ b/man/pseudo_log.Rd @@ -10,6 +10,6 @@ data(pseudo_log) \description{ A truncated version (squared-shape) version of a harmonic signal for mimicking a sedimentary log sign( cos(t*0.13423167+0.00) + 1.3 * cos(t*0.119432+2.314) + 0.134994 + 0.4*cos(t*0.653167) + 0.11 * cos(t*0.78913498) + 0.8)) -for t from 0 to 8191. It is coded as a `mfft_deco` objcet and self-describes its spectral decomposition. +for t from 0 to 8191. It is coded as a `discreteSpectrum` objcet and self-describes its spectral decomposition. } \keyword{datasets} diff --git a/notebook/mfft_random_test.ipynb b/notebook/mfft_random_test.ipynb index c11aad16ac9b5bd6f3061bb07489a5edfb36e4e0..9b54a3b251ab62a1a92c6509048a7c4a25783c44 100644 --- a/notebook/mfft_random_test.ipynb +++ b/notebook/mfft_random_test.ipynb @@ -259,7 +259,7 @@ "[1] 0.000000 2.314000 0.000000 0.653167 0.000000\n", "\n", "attr(,\"class\")\n", - "[1] \"mfft_deco\"\n" + "[1] \"spectrum\"\n" ] }, "metadata": {}, @@ -3832,7 +3832,7 @@ "[1] -5.4832704 0.4030191 -1.6431082 2.3957822 -1.1782940 0.0000000\n", "\n", "attr(,\"class\")\n", - "[1] \"mfft_deco\"\n", + "[1] \"spectrum\"\n", "attr(,\"row.names\")\n", "[1] 1 2 3 4 5 6\n", "attr(,\"nfreq\")\n", @@ -4269,7 +4269,7 @@ "[1] -0.6804843 -1.4340620 -2.6011862 0.2136524 0.0000000 -5.6976628 0.5740075\n", "\n", "attr(,\"class\")\n", - "[1] \"mfft_deco\"\n", + "[1] \"spectrum\"\n", "attr(,\"row.names\")\n", "[1] 1 2 3 4 5 6 7\n", "attr(,\"nfreq\")\n", @@ -4481,7 +4481,7 @@ "[1] 0.0000000 -1.6515515 -1.5786453 -1.4144484 2.8983504 -0.8902084\n", "\n", "attr(,\"class\")\n", - "[1] \"mfft_deco\"\n", + "[1] \"spectrum\"\n", "attr(,\"row.names\")\n", "[1] 1 2 3 4 5 6\n", "attr(,\"nfreq\")\n",