+
+#toCWT: (aux)
+##NOTE: renvoie une matrice 3D
+ toCWT <- function(X, sw= 0, tw= 0, swabs= 0,
+ nvoice= 12, noctave= 5,
+ s0= 2, w0= 2*pi, lt= 24, dt= 0.5,
+ spectra = FALSE, smooth = TRUE,
+ scaled = FALSE,
+ scalevector)
+ { noctave <- adjust.noctave(lt, dt, s0, tw, noctave)
+ if(missing(scalevector))
+ scalevector <- 2^(0:(noctave * nvoice) / nvoice) * s0
+ res <- lapply(1:nrow(X), function(n)
+ { tsX <- ts( X[n,] )
+ tsCent <- tsX - mean(tsX)
+ if(scaled) tsCent <- ts(scale(tsCent))
+ tsCent.cwt <- cwt.ts(tsCent, s0, noctave, nvoice, w0)
+ tsCent.cwt
+ } )
+ if( spectra ) res <- lapply(res, function(l) Mod(l)^2 )
+ if( smooth ) res <- lapply(res, smCWT, swabs = swabs,
+ tw = tw, dt = dt,
+ scalevector = scalevector)
+ resArray <- array(NA, c(nrow(res[[1]]), ncol(res[[1]]),
+ length(res)))
+ for( l in 1:length(res) ) resArray[ , , l] <- res[[l]]
+ resArray
+ }
+
+#from sowas
+cwt.ts <- function(ts,s0,noctave=5,nvoice=10,w0=2*pi){
+
+ if (class(ts)!="ts"){
+
+ cat("# This function needs a time series object as input. You may construct this by using the function ts(data,start,deltat). Try '?ts' for help.\n")
+
+ }
+ else{
+
+ t=time(ts)
+ dt=t[2]-t[1]
+
+ s0unit=s0/dt*w0/(2*pi)
+ s0log=as.integer((log2(s0unit)-1)*nvoice+1.5)
+
+ if (s0log<1){
+ cat(paste("# s0unit = ",s0unit,"\n",sep=""))
+ cat(paste("# s0log = ",s0log,"\n",sep=""))
+ cat("# s0 too small for w0! \n")
+ }
+ totnoct=noctave+as.integer(s0log/nvoice)+1
+
+ #cwt from package Rwave
+ totts.cwt=cwt(ts,totnoct,nvoice,w0,plot=0)
+
+ ts.cwt=totts.cwt[,s0log:(s0log+noctave*nvoice)]
+
+ #Normalization
+ sqs <- sqrt(2^(0:(noctave*nvoice)/nvoice)*s0)
+ smat <- matrix(rep(sqs,length(t)),nrow=length(t),byrow=TRUE)
+
+ ts.cwt*smat
+
+ }
+
+}
+