Commit | Line | Data |
---|---|---|
dc1aa85a BA |
1 | #Entrée : courbes synchrones, soit après étape 1 itérée, soit après chaqure étape 1 |
2 | ||
3 | #(Benjamin) | |
4 | #à partir de là, "conso" == courbes synchrones | |
5 | n <- nrow(conso) | |
6 | delta <- ncol(conso) | |
7 | ||
8 | ||
9 | #17000 colonnes coeff 1, puis 17000 coeff 2... [non : dans chaque tranche du cube] | |
10 | ||
11 | #TODO: une fonction qui fait lignes 59 à 91 | |
12 | ||
13 | #cube: | |
14 | # Xcwt4 <- toCWT(conso, noctave = noctave4, dt = 1, | |
15 | # scalevector = scalevector4, | |
16 | # lt = delta, smooth = FALSE, | |
17 | # nvoice = nvoice) # observations node with CWT | |
18 | # | |
19 | # #matrix: | |
20 | # ############Xcwt2 <- matrix(0.0, nrow= n, ncol= 2 + delta * lscvect) | |
21 | # #Xcwt2 <- matrix(NA_complex_, nrow= n, ncol= 2 + length((c(Xcwt4[,,1])))) | |
22 | # | |
23 | # #NOTE: delta et lscvect pourraient etre gardés à part (communs) | |
24 | # for(i in 1:n) | |
25 | # Xcwt2[i,] <- c(delta, lscvect, Xcwt4[,,i] / max(Mod(Xcwt4[,,i])) ) | |
26 | # | |
27 | # #rm(conso, Xcwt4); gc() | |
28 | # | |
29 | # ## _.b WER^2 distances ######## | |
30 | # Xwer_dist <- matrix(0.0, n, n) | |
31 | # for(i in 1:(n - 1)){ | |
32 | # mat1 <- vect2mat(Xcwt2[i,]) | |
33 | # for(j in (i + 1):n){ | |
34 | # mat2 <- vect2mat(Xcwt2[j,]) | |
35 | # num <- Mod(mat1 * Conj(mat2)) | |
36 | # WX <- Mod(mat1 * Conj(mat1)) | |
37 | # WY <- Mod(mat2 * Conj(mat2)) | |
38 | # smsmnum <- smCWT(num, scalevector = scalevector4) | |
39 | # smsmWX <- smCWT(WX, scalevector = scalevector4) | |
40 | # smsmWY <- smCWT(WY, scalevector = scalevector4) | |
41 | # wer2 <- sum(colSums(smsmnum)^2) / | |
42 | # sum( sum(colSums(smsmWX) * colSums(smsmWY)) ) | |
43 | # Xwer_dist[i, j] <- sqrt(delta * lscvect * (1 - wer2)) | |
44 | # Xwer_dist[j, i] <- Xwer_dist[i, j] | |
45 | # } | |
46 | # } | |
47 | # diag(Xwer_dist) <- numeric(n) | |
48 | # | |
49 | # save(Xwer_dist, file = "../res/2009_synchros200WER.Rdata") | |
50 | # save(Xwer_dist, file = "../res/2009_synchros200-randomWER.Rdata") | |
51 | ||
52 | ||
53 | ||
54 | #lignes 59 à 91 "dépliées" : | |
55 | Xcwt4 <- toCWT(conso, noctave = noctave4, dt = 1, | |
56 | scalevector = scalevector4, | |
57 | lt = delta, smooth = FALSE, | |
58 | nvoice = nvoice) # observations node with CWT | |
59 | ||
60 | #matrix: | |
61 | ############Xcwt2 <- matrix(0.0, nrow= n, ncol= 2 + delta * lscvect) | |
62 | Xcwt2 <- matrix(NA_complex_, nrow= n, ncol= 2 + length((c(Xcwt4[,,1])))) | |
63 | ||
64 | #NOTE: delta et lscvect pourraient etre gardés à part (communs) | |
65 | for(i in 1:n) | |
66 | Xcwt2[i,] <- c(delta, lscvect, Xcwt4[,,i] / max(Mod(Xcwt4[,,i])) ) | |
67 | ||
68 | #rm(conso, Xcwt4); gc() | |
69 | ||
70 | ## _.b WER^2 distances ######## | |
71 | Xwer_dist <- matrix(0.0, n, n) | |
72 | for(i in 1:(n - 1)){ | |
73 | mat1 <- vect2mat(Xcwt2[i,]) | |
74 | ||
75 | #NOTE: vect2mat = as.matrix ?! (dans aux.R) | |
76 | vect2mat <- function(vect){ | |
77 | vect <- as.vector(vect) | |
78 | matrix(vect[-(1:2)], delta, lscvect) | |
79 | } | |
80 | ||
81 | for(j in (i + 1):n){ | |
82 | mat2 <- vect2mat(Xcwt2[j,]) | |
83 | num <- Mod(mat1 * Conj(mat2)) | |
84 | WX <- Mod(mat1 * Conj(mat1)) | |
85 | WY <- Mod(mat2 * Conj(mat2)) | |
86 | smsmnum <- smCWT(num, scalevector = scalevector4) | |
87 | smsmWX <- smCWT(WX, scalevector = scalevector4) | |
88 | smsmWY <- smCWT(WY, scalevector = scalevector4) | |
89 | wer2 <- sum(colSums(smsmnum)^2) / | |
90 | sum( sum(colSums(smsmWX) * colSums(smsmWY)) ) | |
91 | Xwer_dist[i, j] <- sqrt(delta * lscvect * (1 - wer2)) | |
92 | Xwer_dist[j, i] <- Xwer_dist[i, j] | |
93 | } | |
94 | } | |
95 | diag(Xwer_dist) <- numeric(n) | |
96 | ||
97 | #fonction smCWT (dans aux.R) | |
98 | smCWT <- function(CWT, sw= 0, tw= 0, swabs= 0, | |
99 | nvoice= 12, noctave= 2, s0= 2, w0= 2*pi, | |
100 | lt= 24, dt= 0.5, scalevector ) | |
101 | { | |
102 | # noctave <- adjust.noctave(lt, dt, s0, tw, noctave) | |
103 | # scalevector <- 2^(0:(noctave * nvoice) / nvoice) * s0 | |
104 | wsp <- Mod(CWT) | |
105 | smwsp <- smooth.matrix(wsp, swabs) | |
106 | smsmwsp <- smooth.time(smwsp, tw, dt, scalevector) | |
107 | smsmwsp | |
108 | } | |
109 | ||
110 | #dans sowas.R | |
111 | smooth.matrix <- function(wt,swabs){ | |
112 | ||
113 | if (swabs != 0) | |
114 | smwt <- t(filter(t(wt),rep(1,2*swabs+1)/(2*swabs+1))) | |
115 | else | |
116 | smwt <- wt | |
117 | ||
118 | smwt | |
119 | ||
120 | } | |
121 | smooth.time <- function(wt,tw,dt,scalevector){ | |
122 | ||
123 | smwt <- wt | |
124 | ||
125 | if (tw != 0){ | |
126 | for (i in 1:length(scalevector)){ | |
127 | ||
128 | twi <- as.integer(scalevector[i]*tw/dt) | |
129 | smwt[,i] <- filter(wt[,i],rep(1,2*twi+1)/(2*twi+1)) | |
130 | ||
131 | } | |
132 | } | |
133 | smwt | |
134 | } | |
135 | ||
136 | #et filter() est dans stats:: | |
137 | ||
138 | #cf. filters en C dans : https://svn.r-project.org/R/trunk/src/library/stats/src/filter.c | |
139 |