2
|
1 beginTreatment <- function(name, Signal_data = NULL, Signal_info = NULL,
|
|
2 force.real = FALSE) {
|
|
3
|
|
4 cat("Begin", name, "\n")
|
|
5
|
|
6
|
|
7 # Formatting the Signal_data and Signal_info -----------------------
|
|
8
|
|
9 vec <- is.vector(Signal_data)
|
|
10 if (vec) {
|
|
11 Signal_data <- vec2mat(Signal_data)
|
|
12 }
|
|
13 if (is.vector(Signal_info)) {
|
|
14 Signal_info <- vec2mat(Signal_info)
|
|
15 }
|
|
16 if (!is.null(Signal_data)) {
|
|
17 if (!is.matrix(Signal_data)) {
|
|
18 stop("Signal_data is not a matrix.")
|
|
19 }
|
|
20 if (!is.complex(Signal_data) && !is.numeric(Signal_data)) {
|
|
21 stop("Signal_data contains non-numerical values.")
|
|
22 }
|
|
23 }
|
|
24 if (!is.null(Signal_info) && !is.matrix(Signal_info)) {
|
|
25 stop("Signal_info is not a matrix.")
|
|
26 }
|
|
27
|
|
28
|
|
29 Original_data <- Signal_data
|
|
30
|
|
31 # Extract the real part of the spectrum ---------------------------
|
|
32
|
|
33 if (force.real) {
|
|
34 if (is.complex(Signal_data)) {
|
|
35 Signal_data <- Re(Signal_data)
|
|
36 } else {
|
|
37 # The signal is numeric Im(Signal_data) is zero anyway so let's avoid
|
|
38 # using complex(real=...,imaginary=0) which would give a complex signal
|
|
39 # in endTreatment()
|
|
40 force.real <- FALSE
|
|
41 }
|
|
42 }
|
|
43
|
|
44
|
|
45 # Return the formatted data and metadata entries --------------------
|
|
46
|
|
47 return(list(start = proc.time(), vec = vec, force.real = force.real,
|
|
48 Original_data = Original_data, Signal_data = Signal_data, Signal_info = Signal_info))
|
|
49 }
|
|
50
|
|
51 binarySearch <- function(a, target, lower = TRUE) {
|
|
52 # search the index i in a such that a[i] == target
|
|
53 # if it doesn't exists and lower, it searches the closer a[i] such that a[i] < target
|
|
54 # if !lower, it seraches the closer a[i] such that a[i] > target
|
|
55 # a should be monotone but can be increasing or decreasing
|
|
56
|
|
57 # if a is increasing INVARIANT: a[amin] < target < a[amax]
|
|
58 N <- length(a)
|
|
59 if ((a[N] - target) * (a[N] - a[1]) <= 0) {
|
|
60 return(N)
|
|
61 }
|
|
62 if ((a[1] - target) * (a[N] - a[1]) >= 0) {
|
|
63 return(1)
|
|
64 }
|
|
65 amin <- 1
|
|
66 amax <- N
|
|
67 while (amin + 1 < amax) {
|
|
68 amid <- floor((amin + amax)/2)
|
|
69 if ((a[amid] - target) * (a[amax] - a[amid]) < 0) {
|
|
70 amin <- amid
|
|
71 } else if ((a[amid] - target) * (a[amax] - a[amid]) > 0) {
|
|
72 amax <- amid
|
|
73 } else {
|
|
74 # a[amid] == a[amax] or a[amid] == target In both cases, a[amid] ==
|
|
75 # target
|
|
76 return(amid)
|
|
77 }
|
|
78 }
|
|
79 if (xor(lower, a[amin] > a[amax])) {
|
|
80 # (lower && a[amin] < a[amax]) || (!lower && a[min] > a[max])
|
|
81 # If increasing and we want the lower, we take amin
|
|
82 # If decreasing and we want the bigger, we take amin too
|
|
83 return(amin)
|
|
84 } else {
|
|
85 return(amax)
|
|
86 }
|
|
87 }
|
|
88
|
|
89 checkArg <- function(arg, checks, can.be.null=FALSE) {
|
|
90 check.list <- list(bool=c(is.logical, "a boolean"),
|
|
91 int =c(function(x){x%%1==0}, "an integer"),
|
|
92 num =c(is.numeric, "a numeric"),
|
|
93 str =c(is.character, "a string"),
|
|
94 pos =c(function(x){x>0}, "positive"),
|
|
95 pos0=c(function(x){x>=0}, "positive or zero"),
|
|
96 l1 =c(function(x){length(x)==1}, "of length 1")
|
|
97 )
|
|
98 if (is.null(arg)) {
|
|
99 if (!can.be.null) {
|
|
100 stop(deparse(substitute(arg)), " is null.")
|
|
101 }
|
|
102 } else {
|
|
103 if (is.matrix(arg)) {
|
|
104 stop(deparse(substitute(arg)), " is not scalar.")
|
|
105 }
|
|
106 for (c in checks) {
|
|
107 if (!check.list[[c]][[1]](arg)) {
|
|
108 stop(deparse(substitute(arg)), " is not ", check.list[[c]][[2]], ".")
|
|
109 }
|
|
110 }
|
|
111 }
|
|
112 }
|
|
113
|
|
114 endTreatment <- function(name, begin_info, Signal_data) {
|
|
115
|
|
116 # begin_info: object outputted from beginTreatment
|
|
117
|
|
118
|
|
119 # Formatting the entries and printing process time -----------------------
|
|
120 end_time <- proc.time() # record it as soon as possible
|
|
121 start_time <- begin_info[["start"]]
|
|
122 delta_time <- end_time - start_time
|
|
123 delta <- delta_time[]
|
|
124 cat("End", name, "\n")
|
|
125 cat("It lasted", round(delta["user.self"], 3), "s user time,", round(delta["sys.self"],3),
|
|
126 "s system time and", round(delta["elapsed"], 3), "s elapsed time.\n")
|
|
127
|
|
128
|
|
129 if (begin_info[["force.real"]]) {
|
|
130 # The imaginary part is left untouched
|
|
131 i <- complex(real = 0, imaginary = 1)
|
|
132 Signal_data <- Signal_data + i * Im(begin_info[["Original_data"]])
|
|
133 }
|
|
134
|
|
135 if (begin_info[["vec"]]) {
|
|
136 Signal_data <- Signal_data[1, ]
|
|
137 }
|
|
138
|
|
139 # Return the formatted data and metadata entries --------------------
|
|
140 return(Signal_data)
|
|
141 }
|
|
142
|
|
143 Interpol <- function(t, y) {
|
|
144 # y: sample
|
|
145 # t : warping function
|
|
146
|
|
147 m <- length(y)
|
|
148 # t <= m-1
|
|
149 # because if t > m-1, y[ti+1] will be NA when we compute g
|
|
150 valid <- 1 <= t & t <= m-1 # FIXME it was '<' in Bubble v2
|
|
151 s <- (1:m)[valid]
|
|
152 ti <- floor(t[s])
|
|
153 tr <- t[s] - ti
|
|
154 g <- y[ti + 1] - y[ti]
|
|
155 f <- y[ti] + tr * g
|
|
156 list(f=f, s=s, g=g)
|
|
157 }
|
|
158
|
|
159 vec2mat <- function(vec) {
|
|
160
|
|
161 return(matrix(vec, nrow = 1, dimnames = list(c(1), names(vec))))
|
|
162
|
|
163 }
|
|
164
|
|
165 getArg <- function(arg, info, argname, can.be.absent=FALSE) {
|
|
166 if (is.null(arg)) {
|
|
167 start <- paste("impossible to get argument", argname, "it was not given directly and");
|
|
168 if (!is.matrix(info)) {
|
|
169 stop(paste(start, "the info matrix was not given"))
|
|
170 }
|
|
171 if (!(argname %in% colnames(info))) {
|
|
172 if (can.be.absent) {
|
|
173 return(NULL)
|
|
174 } else {
|
|
175 stop(paste(start, "is not in the info matrix"))
|
|
176 }
|
|
177 }
|
|
178 if (nrow(info) < 1) {
|
|
179 stop(paste(start, "the info matrix has no row"))
|
|
180 }
|
|
181 arg <- info[1,argname]
|
|
182 if (is.na(arg)) {
|
|
183 stop(paste(start, "it is NA in the info matrix"))
|
|
184 }
|
|
185 }
|
|
186 return(arg)
|
|
187 }
|
|
188
|
|
189 binarySearch <- function(a, target, lower = TRUE) {
|
|
190 # search the index i in a such that a[i] == target
|
|
191 # if it doesn't exists and lower, it searches the closer a[i] such that a[i] < target
|
|
192 # if !lower, it seraches the closer a[i] such that a[i] > target
|
|
193 # a should be monotone but can be increasing or decreasing
|
|
194
|
|
195 # if a is increasing INVARIANT: a[amin] < target < a[amax]
|
|
196 N <- length(a)
|
|
197 if ((a[N] - target) * (a[N] - a[1]) <= 0) {
|
|
198 return(N)
|
|
199 }
|
|
200 if ((a[1] - target) * (a[N] - a[1]) >= 0) {
|
|
201 return(1)
|
|
202 }
|
|
203 amin <- 1
|
|
204 amax <- N
|
|
205 while (amin + 1 < amax) {
|
|
206 amid <- floor((amin + amax)/2)
|
|
207 if ((a[amid] - target) * (a[amax] - a[amid]) < 0) {
|
|
208 amin <- amid
|
|
209 } else if ((a[amid] - target) * (a[amax] - a[amid]) > 0) {
|
|
210 amax <- amid
|
|
211 } else {
|
|
212 # a[amid] == a[amax] or a[amid] == target In both cases, a[amid] ==
|
|
213 # target
|
|
214 return(amid)
|
|
215 }
|
|
216 }
|
|
217 if (xor(lower, a[amin] > a[amax])) {
|
|
218 # (lower && a[amin] < a[amax]) || (!lower && a[min] > a[max])
|
|
219 # If increasing and we want the lower, we take amin
|
|
220 # If decreasing and we want the bigger, we take amin too
|
|
221 return(amin)
|
|
222 } else {
|
|
223 return(amax)
|
|
224 }
|
|
225 }
|
|
226
|
|
227 # returns the discrete borders of the interval for a numeric vector a
|
|
228
|
|
229 indexInterval <- function (a, from, to, inclusive=TRUE) {
|
|
230 # If inclusive and from <= to, we need to take the lower
|
|
231 # If not inclusive and from > to, we need to take the lower too
|
|
232 lowerFrom <- xor(inclusive, from > to)
|
|
233 fromIndex <- binarySearch(a, from, lowerFrom)
|
|
234 toIndex <- binarySearch(a, to, !lowerFrom)
|
|
235 return(fromIndex:toIndex)
|
|
236 }
|