set final draft for package
[epclust.git] / epclust / R / stage2.R
diff --git a/epclust/R/stage2.R b/epclust/R/stage2.R
deleted file mode 100644 (file)
index 3ccbbad..0000000
+++ /dev/null
@@ -1,51 +0,0 @@
-library("Rwave")
-
-#Entrée : courbes synchrones, soit après étape 1 itérée, soit après chaqure étape 1
-step2 = function(conso)
-{
-       n <- nrow(conso)
-       delta <- ncol(conso)
-       #TODO: automatic tune of all these parameters ? (for other users)
-       nvoice   <- 4
-       # noctave = 2^13 = 8192 half hours ~ 180 days ; ~log2(ncol(conso))
-       noctave = 13
-       # 4 here represent 2^5 = 32 half-hours ~ 1 day
-       #NOTE: default scalevector == 2^(0:(noctave * nvoice) / nvoice) * s0 (?)
-       scalevector  <- 2^(4:(noctave * nvoice) / nvoice) * 2
-       #condition: ( log2(s0*w0/(2*pi)) - 1 ) * nvoice + 1.5 >= 1
-       s0=2
-       w0=2*pi
-       scaled=FALSE
-       s0log = as.integer( (log2( s0*w0/(2*pi) ) - 1) * nvoice + 1.5 )
-       totnoct = noctave + as.integer(s0log/nvoice) + 1
-
-       # (normalized) observations node with CWT
-       Xcwt4 <- lapply(seq_len(n), function(i) {
-               ts <- scale(ts(conso[i,]), center=TRUE, scale=scaled)
-               totts.cwt = Rwave::cwt(ts,totnoct,nvoice,w0,plot=0)
-               ts.cwt = totts.cwt[,s0log:(s0log+noctave*nvoice)]
-               #Normalization
-               sqs <- sqrt(2^(0:(noctave*nvoice)/nvoice)*s0)
-               sqres <- sweep(ts.cwt,MARGIN=2,sqs,'*')
-               sqres / max(Mod(sqres))
-       })
-
-       Xwer_dist <- matrix(0., n, n)
-       fcoefs = rep(1/3, 3) #moving average on 3 values (TODO: very slow! correct?!)
-       for (i in 1:(n-1))
-       {
-               for (j in (i+1):n)
-               {
-                       #TODO: later, compute CWT here (because not enough storage space for 32M series)
-                       #      'circular=TRUE' is wrong, should just take values on the sides; to rewrite in C
-                       num <- filter(Mod(Xcwt4[[i]] * Conj(Xcwt4[[j]])), fcoefs, circular=TRUE)
-                       WX <- filter(Mod(Xcwt4[[i]] * Conj(Xcwt4[[i]])), fcoefs, circular=TRUE)
-                       WY <- filter(Mod(Xcwt4[[j]] * Conj(Xcwt4[[j]])), fcoefs, circular=TRUE)
-                       wer2    <- sum(colSums(num)^2) / sum( sum(colSums(WX) * colSums(WY)) )
-                       Xwer_dist[i,j] <- sqrt(delta * ncol(Xcwt4[[1]]) * (1 - wer2))
-                       Xwer_dist[j,i] <- Xwer_dist[i,j]
-               }
-       }
-       diag(Xwer_dist) <- numeric(n)
-       Xwer_dist
-}