diff --git a/.Rbuildignore b/.Rbuildignore index 6bbefad7f80928aa03809bf7fc2cfbe18a748b35..40c7223e81130bf0bdac1b58bab2ee85024a1c15 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -1,6 +1,7 @@ .gitignore test old_man +miscellaneous_notes ^data-raw$ .Rprofile TODO diff --git a/R/haar.R b/R/haar.R new file mode 100644 index 0000000000000000000000000000000000000000..52dcba9796cd87e30e30daf4a29c10cb5964ee7d --- /dev/null +++ b/R/haar.R @@ -0,0 +1,48 @@ +# TO DO: THERE IS A BIAS FOR THE LONG LENGHTSCALES +# I GUESS MUST MULTIPLY BY (N-1)/N where N is +# length / scale + + + +#' Haar fluctuation spectrum +#' +#' @description Haar fluctuation spectrum after Lovejoy +#' still needs to be documented +#' currently requires n=2^i data points + + +#' @param x the input series (typically a numeric) +#' @param q the fluctuation order +#' @examples +#' x = rnorm(2048) +#' xi1 = haar(x, q=1) +#' xi2 = haar(x, q=2) + +haar <- function(x, q=2, discarded_scales = 4) +{ + n = length(x) + nsteps=floor(log(n)/log(2)) + nmax = 2^nsteps + DT = sapply (seq(nsteps), function(i) + { + m = 2^i + nn = nmax / m + M = matrix(x, nmax/m, m, byrow=TRUE) + M1 = M[,1:(m/2), drop=FALSE] + M2 = M[,(m/2+1):m, drop=FALSE] + tM1 = apply(M1, 1, mean) + tM2 = apply(M2, 1, mean) + abs(tM2 - tM1) * ( 2 * nn ) / (2 * nn - 1) + }) + x= 2^seq(nsteps) + y= sapply(DT, function(x) mean(x^q)) + H = data.frame(x = x, y=y, logx=log(x), logy=log(y)) + kept_scales <- seq(max(2,nsteps-discarded_scales)) + out = lm( H$logy[kept_scales] ~ H$logx[kept_scales] )$coefficients[2] + attr(out, "scale") = H + attr(out, "q") = q + attr(out,"class") = "fluctuation spectrum" + out + } + + diff --git a/miscellaneous_notes/fluctuation_spectrum.txt b/miscellaneous_notes/fluctuation_spectrum.txt new file mode 100644 index 0000000000000000000000000000000000000000..8ae1115abd025eee794e3fd6258b301e71d5dd80 --- /dev/null +++ b/miscellaneous_notes/fluctuation_spectrum.txt @@ -0,0 +1,18 @@ +ligne rouge : cumsum de la simulation avec package 'FGN' +ligne noire : simul avec FGN + +x-axis: le Hurst entre dans FGN +y=axis: le H calcule par Haar Wavelet fluct. + +a partir d'ici: on parle du H de Lovejoy + +H s'appelle H is the fluctua- tion (also called “nonconservationâ€) exponent + +H=-0.5 = bruit blanc + +H approx xi(1) +H = xi(2) / 2 + +beta = 1 + xi(2) = 1 + 2H + +donc, bruit blanc : xi(2) = 1 + 2*(-0.5) = 0 diff --git a/src/Makevars b/src/Makevars index c2452b0931ef1dbca5f61c79ad7d5093609b90d3..2871a83d1ccfe66734e662b22634d7c07cbcbb63 100644 --- a/src/Makevars +++ b/src/Makevars @@ -1,7 +1,7 @@ # Kindly supplied by Dirk Eddelbuettel # set by configure -GSL_CFLAGS = -I/usr/local/include -GSL_LIBS = -L/usr/local/lib -lgsl -lgslcblas -lm +GSL_CFLAGS = -I/usr/include +GSL_LIBS = -L/usr/lib/x86_64-linux-gnu -lgsl -lgslcblas -lm # combine to standard arguments for R PKG_CPPFLAGS = $(GSL_CFLAGS) -I.