diff --git a/DESCRIPTION b/DESCRIPTION index abe9adb3a..2cf76a6a7 100755 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -30,6 +30,6 @@ URL: https://github.com/wadpac/GGIR/, https://groups.google.com/forum/#!forum/Rp BugReports: https://github.com/wadpac/GGIR/issues License: Apache License (== 2.0) | file LICENSE Suggests: testthat, covr, knitr, rmarkdown, actilifecounts, ActCR, GGIRread, read.gt3x, readxl -Imports: data.table, foreach, doParallel, signal, zoo, unisensR, ineq, methods, psych, irr +Imports: data.table, foreach, doParallel, signal, zoo, unisensR, ineq, methods, psych, irr, lubridate Depends: stats, utils, R (>= 3.5) VignetteBuilder: knitr diff --git a/NAMESPACE b/NAMESPACE index 6355ea17b..e99869c32 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -7,7 +7,7 @@ export(g.analyse, g.calibrate, GGIR, g.shell.GGIR, g.getbout, g.abr.day.names, g.applymetrics, g.create.sp.mat, g.detecmidnight, - g.dotorcomma, g.downsample, g.extractheadervars, + g.dotorcomma, g.extractheadervars, g.getM5L5, g.getstarttime, g.loadlog, g.plot5, g.readaccfile, @@ -24,7 +24,7 @@ export(g.analyse, g.calibrate, g.part5.addsib, g.part5.definedays, g.part5.fixmissingnight, g.part5.onsetwaketiming, g.part5.wakesleepwindows, g.part5.savetimeseries, - get_nw_clip_block_params, get_starttime_weekday_meantemp_truncdata, + get_nw_clip_block_params, get_starttime_weekday_truncdata, ismovisens, g.readtemp_movisens, g.conv.actlog, g.fragmentation, g.intensitygradient, g.part5.handle_lux_extremes, diff --git a/NEWS.md b/NEWS.md index 75123452a..ed7f0f48e 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,15 @@ +# CHANGES IN GGIR VERSION 3.0-6 + +Part 1: + +- Improved readability and maintainability of the code #1027 + +- Improved processing speed for Axivity .cwa, GENEActiv .bin, and Movisens files + +- Made sure that g.readaccfile() reads timestamps in the correct timezone, configtz, for all monitor types + +- Note: there will be small differences in both metalong and metashort metrics calculated by this GGIR version, compared to prior versions. This is due to small improvements in the management of timestamps, calibration coefficients, and input data block boundaries. + # CHANGES IN GGIR VERSION 3.0-5 - Part 5: Fix bug in functionality for Sensewear data (externally derived epoch data) #1030 diff --git a/R/check_params.R b/R/check_params.R index 8a0148cad..58d16c4c1 100644 --- a/R/check_params.R +++ b/R/check_params.R @@ -72,6 +72,9 @@ check_params = function(params_sleep = c(), params_metrics = c(), check_class("Raw data", params = params_rawdata, parnames = numeric_params, parclass = "numeric") check_class("Raw data", params = params_rawdata, parnames = boolean_params, parclass = "boolean") check_class("Raw data", params = params_rawdata, parnames = character_params, parclass = "character") + + if (params_rawdata[["chunksize"]] > 1.5) params_rawdata[["chunksize"]] = 1.5 + if (params_rawdata[["chunksize"]] < 0.2) params_rawdata[["chunksize"]] = 0.2 } if (length(params_247) > 0) { # iglevels and qwindow can be numeric or character, so not tested @@ -127,6 +130,29 @@ check_params = function(params_sleep = c(), params_metrics = c(), check_class("general", params = params_general, parnames = numeric_params, parclass = "numeric") check_class("general", params = params_general, parnames = boolean_params, parclass = "boolean") check_class("general", params = params_general, parnames = character_params, parclass = "character") + + ws3 = params_general[["windowsizes"]][1]; ws2 = params_general[["windowsizes"]][2]; ws = params_general[["windowsizes"]][3] + if (ws2/60 != round(ws2/60)) { + ws2 = as.numeric(60 * ceiling(ws2/60)) + warning(paste0("The long windowsize needs to be a multitude of 1 minute periods.\n", + "Long windowsize has now been automatically adjusted to ", + ws2, " seconds in order to meet this criteria."), call. = FALSE) + } + if (ws2/ws3 != round(ws2/ws3)) { + def = c(1,5,10,15,20,30,60) + def2 = abs(def - ws3) + ws3 = as.numeric(def[which(def2 == min(def2))]) + warning(paste0("The long windowsize needs to be a multitude of short windowsize.\n", + "The short windowsize has now been automatically adjusted to ", + ws3, " seconds in order to meet this criteria.\n"), call. = FALSE) + } + if (ws/ws2 != round(ws/ws2)) { + ws = ws2 * ceiling(ws/ws2) + warning(paste0("The third value of parameter windowsizes needs to be a multitude of the second value.\n", + "The third value has been automatically adjusted to ", + ws, " seconds in order to meet this criteria.\n"), call. = FALSE) + } + params_general[["windowsizes"]] = c(ws3, ws2, ws) } #----------------------------------------------------------------------------------- # Check value combinations and apply corrections if not logical diff --git a/R/datadir2fnames.R b/R/datadir2fnames.R index 27ba9d119..f225cfcbf 100644 --- a/R/datadir2fnames.R +++ b/R/datadir2fnames.R @@ -13,14 +13,24 @@ datadir2fnames = function(datadir,filelist) { } } else { if (ismovisens(datadir)) { - fnamesfull = dir(datadir, recursive = TRUE, pattern = "acc.bin", full.names = TRUE) - fnames = basename(dirname(fnamesfull)) - nfolders = length(dir(datadir)) - if (nfolders > length(fnamesfull)) { # meaning there are movisens data folder without acc.bin - # folders without acc.bin - allfolders = dir(datadir, full.names = TRUE) - foldersWithAccBin = dirname(fnamesfull) - noAccBin = allfolders[which(!allfolders %in% foldersWithAccBin)] + fnamesfull = dir(datadir, recursive = TRUE, pattern = "acc.bin$", full.names = TRUE) + foldersWithAccBin = dirname(fnamesfull) + fnames = basename(foldersWithAccBin) + + allfolders = list.dirs(datadir, recursive=FALSE, full.names = TRUE) + + # are there any folders that don't contain acc.bin (directly, or in any subfolders)? + noAccBin = c() + for (fld in allfolders) { + if (fld %in% foldersWithAccBin) next # folder contains acc.bin + + # do any subfolders contain acc.bin? + if (!any(grepl(paste(fld, '/', sep = ''), foldersWithAccBin))) { + noAccBin = c(noAccBin, fld) + } + } + + if (length(noAccBin) > 0) { warning(paste0("The following movisens data folders do not contain the ", "acc.bin file with the accelerometer recording, and ", "therefore cannot be processed in GGIR: ", diff --git a/R/detect_nonwear_clipping.R b/R/detect_nonwear_clipping.R index 6e5db97bf..80c16a026 100644 --- a/R/detect_nonwear_clipping.R +++ b/R/detect_nonwear_clipping.R @@ -11,11 +11,6 @@ detect_nonwear_clipping = function(data = c(), windowsizes = c(5, 900, 3600), sf CWav = NWav = rep(0, nmin) crit = ((window/window2)/2) + 1 - if (is.numeric(data[,1]) == FALSE) { - data[,1] = as.numeric(data[,1]) - data[,2] = as.numeric(data[,2]) - data[,3] = as.numeric(data[,3]) - } if (nonwear_approach %in% c("2013", "2023")) { # define windows to check: for (h in 1:nmin) { #number of windows @@ -43,25 +38,24 @@ detect_nonwear_clipping = function(data = c(), windowsizes = c(5, 900, 3600), sf NWflag = h:(h + window/window2 - 1) if (NWflag[length(NWflag)] > nmin) NWflag = NWflag[-which(NWflag > nmin)] # window to check (not aggregated values) - hoc1 = h * window2 - window2 + 1 - hoc2 = hoc1 + window - 1 + hoc1 = h * window2 - window2 + hoc2 = hoc1 + window if (hoc2 > nrow(data)) { hoc2 = nrow(data) } } # --- - if (length(params_rawdata[["rmc.col.wear"]]) > 0) { - wearcol = as.character(data[, which(colnames(data) == "wear")]) - suppressWarnings(storage.mode(wearcol) <- "logical") - wearTable = table(wearcol[(1 + hoc1):hoc2], useNA = FALSE) - NWav[h] = as.logical(tail(names(sort(wearTable)), 1)) * 3 # times 3 to simulate heuristic approach + if ("wear" %in% colnames(data)) { + wearTable = table(data[(1 + hoc1):hoc2, "wear"], useNA = "no") + NWav[h] = as.numeric(tail(names(sort(wearTable)), 1)) * 3 # times 3 to simulate heuristic approach } - for (jj in 1:3) { + xyzCol = which(colnames(data) %in% c("x", "y", "z")) + for (jj in seq(3)) { # Clipping - aboveThreshold = which(abs(data[(1 + cliphoc1):cliphoc2, jj]) > clipthres) + aboveThreshold = which(abs(data[(1 + cliphoc1):cliphoc2, xyzCol[jj]]) > clipthres) CW[h, jj] = length(aboveThreshold) if (length(aboveThreshold) > 0) { - if (length(which(abs(data[c((1 + cliphoc1):cliphoc2)[aboveThreshold],jj]) > clipthres * 1.5)) > 0) { + if (length(which(abs(data[c((1 + cliphoc1):cliphoc2)[aboveThreshold], xyzCol[jj]]) > clipthres * 1.5)) > 0) { CW[h, jj] = window2 # If there is a a value that is more than 150% the dynamic range then ignore entire block. } } @@ -73,18 +67,18 @@ detect_nonwear_clipping = function(data = c(), windowsizes = c(5, 900, 3600), sf } else if (nonwear_approach == "2023") { indices = seq((1 + hoc1), hoc2, by = ceiling(sf / 5)) } - maxwacc = max(data[indices, jj], na.rm = TRUE) - minwacc = min(data[indices, jj], na.rm = TRUE) + maxwacc = max(data[indices, xyzCol[jj]], na.rm = TRUE) + minwacc = min(data[indices, xyzCol[jj]], na.rm = TRUE) absrange = abs(maxwacc - minwacc) if (absrange < racriter) { - sdwacc = sd(data[indices,jj], na.rm = TRUE) + sdwacc = sd(data[indices, xyzCol[jj]], na.rm = TRUE) if (sdwacc < sdcriter) { NW[NWflag,jj] = 1 } } } CW = CW / (window2) - if (length(params_rawdata[["rmc.col.wear"]]) == 0) { + if (!("wear" %in% colnames(data))) { NWav[h] = (NW[h,1] + NW[h,2] + NW[h,3]) #indicator of non-wear } CWav[h] = max(c(CW[h, 1], CW[h, 2], CW[h, 3])) #indicator of clipping diff --git a/R/g.applymetrics.R b/R/g.applymetrics.R index 14db22917..9af3a17f3 100644 --- a/R/g.applymetrics.R +++ b/R/g.applymetrics.R @@ -3,8 +3,10 @@ g.applymetrics = function(data, sf, ws3, metrics2do, zc.lb = 0.25, zc.hb = 3, zc.sb = 0.01, zc.order = 2, actilife_LFE = FALSE){ + # data should be a 3 column matrix with the x, y, and z acceleration + data = data[, c("x", "y", "z")] + epochsize = ws3 #epochsize in seconds - # data is a 3 column matrix with the x, y, and z acceleration do.bfen = metrics2do$do.bfen do.enmo = metrics2do$do.enmo do.lfenmo = metrics2do$do.lfenmo diff --git a/R/g.calibrate.R b/R/g.calibrate.R index f6888696c..22f00bc17 100644 --- a/R/g.calibrate.R +++ b/R/g.calibrate.R @@ -22,16 +22,26 @@ g.calibrate = function(datafile, params_rawdata = c(), params_cleaning = params$params_cleaning params_general = params$params_general } - - use.temp = TRUE - filename = unlist(strsplit(as.character(datafile),"/")) - filename = filename[length(filename)] - # set parameters + #----------------- + # Define local functions: + rollMean = function(x, fs, epochSize) { + x = cumsum(c(0, x)) + select = seq(1, length(x), by = fs * epochSize) + y = diff(x[round(select)]) / abs(diff(round(select))) + return(y) + } + rollSD = function(x, sf, epochSize) { + dim(x) = c(sf * epochSize, ceiling(length(x) / (sf * epochSize))) + y = apply(x, 2, sd) + return(y) + } + #----------------- + use.temp = temp.available = TRUE + filequality = data.frame(filetooshort = FALSE, filecorrupt = FALSE, - filedoesnotholdday = FALSE, stringsAsFactors = TRUE) - ws4 = 10 #epoch for recalibration, don't change - ws2 = params_general[["windowsizes"]][2] #dummy variable - ws = params_general[["windowsizes"]][3] # window size for assessing non-wear time (seconds) + filedoesnotholdday = FALSE, stringsAsFactors = FALSE) + calibEpochSize = 10 # epoch for recalibration as used in the 2014 paper + blockResolution = 3600 # resolution at which raw data is loaded and processed (seconds) cal.error.start = cal.error.end = c() spheredata = c() tempoffset = c() @@ -60,27 +70,10 @@ g.calibrate = function(datafile, params_rawdata = c(), return() } - # if GENEActiv csv, deprecated function - if (mon == MONITOR$GENEACTIV && dformat == FORMAT$CSV && length(params_rawdata[["rmc.firstrow.acc"]]) == 0) { - stop(paste0("The GENEActiv csv reading functionality is deprecated in", - " GGIR from the version 2.6-4 onwards. Please, use either", - " the GENEActiv bin files or try to read the csv files with", - " GGIR::read.myacc.csv"), call. = FALSE) - } - if (sf == 0) { - stop(paste0("\nSample frequency not recognised in ", datafile, - " calibration procedure stopped."), call. = FALSE) - } #creating matrices for storing output S = matrix(0,0,4) #dummy variable needed to cope with head-tailing succeeding blocks of data - NR = ceiling((90*10^6) / (sf*ws4)) + 1000 #NR = number of 'ws4' second rows (this is for 10 days at 80 Hz) - if (mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) || - mon == MONITOR$MOVISENS || (mon == MONITOR$AD_HOC && length(params_rawdata[["rmc.col.temp"]]) > 0)) { - meta = matrix(99999,NR,8) #for meta data - } else if (mon == MONITOR$ACTIGRAPH || (mon == MONITOR$AXIVITY && dformat == FORMAT$CSV) || - (mon == MONITOR$AD_HOC && length(params_rawdata[["rmc.col.temp"]]) == 0)) { - meta = matrix(99999,NR,7) - } + NR = ceiling((90*10^6) / (sf*calibEpochSize)) + 1000 #NR = number of rows to initialise features matrix with + # setting size of blocks that are loaded (too low slows down the process) # the setting below loads blocks size of 12 hours (modify if causing memory problems) blocksize = round((14512 * (sf/50)) * (params_rawdata[["chunksize"]]*0.5)) @@ -97,9 +90,9 @@ g.calibrate = function(datafile, params_rawdata = c(), i = 1 #counter to keep track of which binary block is being read count = 1 #counter to keep track of the number of seconds that have been read LD = 2 #dummy variable used to identify end of file and to make the process stop - switchoffLD = 0 #dummy variable part of "end of loop mechanism" + isLastBlock = FALSE # dummy variable part of "end of loop mechanism" + header = NULL while (LD > 1) { - P = c() if (verbose == TRUE) { if (i == 1) { cat(paste("\nLoading chunk: ",i,sep="")) @@ -110,180 +103,95 @@ g.calibrate = function(datafile, params_rawdata = c(), #try to read data blocks based on monitor type and data format options(warn=-1) #turn off warnings (code complains about unequal rowlengths accread = g.readaccfile(filename = datafile, blocksize = blocksize, blocknumber = i, - filequality = filequality, ws = ws, + filequality = filequality, ws = blockResolution, PreviousEndPage = PreviousEndPage, inspectfileobject = INFI, - params_rawdata = params_rawdata, params_general = params_general) - P = accread$P - switchoffLD = accread$switchoffLD + params_rawdata = params_rawdata, params_general = params_general, + header = header) + header = accread$header + isLastBlock = accread$isLastBlock PreviousEndPage = accread$endpage - PreviousStartPage = accread$startpage - rm(accread); + + if (i == 1) { + use.temp = temp.available = ("temperature" %in% colnames(accread$P$data)) + if (use.temp) { + features = matrix(99999, NR, 8) + } else { + features = matrix(99999, NR, 7) + } + } + options(warn = 0) #turn on warnings #process data as read from binary file - if (length(P) > 0) { #would have been set to zero if file was corrupt or empty - if (mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) { - data = P$data.out - } else if (dformat == FORMAT$CSV || mon == MONITOR$MOVISENS) { - data = as.matrix(P) - } else if (dformat == FORMAT$CWA) { - if (P$header$hardwareType == "AX6") { # cwa AX6 - # Note 18-Feb-2020: For the moment GGIR ignores the AX6 gyroscope signals until robust sensor - # fusion algorithms and gyroscope metrics have been prepared. - # Note however that while AX6 is able to collect gyroscope data, it can also be configured - # to only collect accelerometer data, so only remove gyro data if it's present. - if (ncol(P$data) == 10) { - data = P$data[,-c(2:4)] - } else { - data = P$data - } - } else { - # cwa AX3 - data = P$data - } - } else if (dformat == FORMAT$AD_HOC_CSV) { - if ("timestamp" %in% colnames(P$data)) { - columns_to_use = 2:4 - } else { - columns_to_use = 1:3 - } - if ("temperature" %in% colnames(P$data)) { - columns_to_use = c(columns_to_use, which(colnames(P$data) == "temperature")) - use.temp = TRUE - } else { - use.temp = FALSE - } - data = P$data[, columns_to_use] - } else if (dformat == FORMAT$GT3X) { - data = as.matrix(P[,2:4]) + if (length(accread$P) > 0) { # would have been set to zero if file was corrupt or empty + data = as.matrix(accread$P$data[,which(colnames(accread$P$data) %in% c("x", "y", "z", "time", "temperature"))]) # convert to matrix b/c rbind is much faster on matrices + if (exists("accread")) { + rm(accread) } - rm(P) - #add left over data from last time + # add leftover data from last time if (min(dim(S)) > 1) { data = rbind(S,data) } # remove 0s if ActiGraph csv (idle sleep mode) OR if similar imputation done in ad-hoc csv # current ActiGraph csv's are not with zeros but with last observation carried forward - zeros = which(data[,1] == 0 & data[,2] == 0 & data[,3] == 0) + zeros = c() + if (!(mon == MONITOR$ACTIGRAPH && dformat == FORMAT$CSV)) { + zeros = which(data[, "x"] == 0 & data[, "y"] == 0 & data[, "z"] == 0) + } if ((mon == MONITOR$ACTIGRAPH && dformat == FORMAT$CSV) || length(zeros) > 0) { - data = g.imputeTimegaps(x = as.data.frame(data), xyzCol = 1:3, timeCol = c(), sf = sf, impute = FALSE) + data = g.imputeTimegaps(x = as.data.frame(data), sf = sf, impute = FALSE) data = as.matrix(data$x) } LD = nrow(data) #store data that could not be used for this block, but will be added to next block - use = (floor(LD / (ws*sf))) * (ws*sf) #number of datapoint to use + use = (floor(LD / (blockResolution*sf))) * (blockResolution*sf) #number of datapoint to use if (length(use) > 0) { if (use > 0) { if (use != LD) { - S = as.matrix(data[(use + 1):LD,]) #store left over # as.matrix removed on 22May2019 because redundant - #S = data[(use+1):LD,] #store left over - } - data = as.matrix(data[1:use,]) - LD = nrow(data) #redefine LD because there is less data - ##================================================== - # Initialization of variables - if (dformat != FORMAT$AD_HOC_CSV) { - suppressWarnings(storage.mode(data) <- "numeric") - } - if (dformat == FORMAT$CWA && mon == MONITOR$AXIVITY) { - Gx = data[,2]; Gy = data[,3]; Gz = data[,4] - use.temp = TRUE - } else if (dformat == FORMAT$CSV && mon == MONITOR$AXIVITY) { - Gx = data[,2]; Gy = data[,3]; Gz = data[,4] - use.temp = FALSE - } else if (dformat == FORMAT$GT3X && mon == MONITOR$ACTIGRAPH) { - Gx = data[,1]; Gy = data[,2]; Gz = data[,3] - use.temp = FALSE - } else if (mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) { - Gx = data[,2]; Gy = data[,3]; Gz = data[,4] - use.temp = TRUE - } else if (mon == MONITOR$MOVISENS) { - Gx = data[,1]; Gy = data[,2]; Gz = data[,3] - use.temp = TRUE - } else if (dformat == FORMAT$AD_HOC_CSV) { - Gx = as.numeric(data[,1]); Gy = as.numeric(data[,2]); Gz = as.numeric(data[,3]) - } else if (dformat == FORMAT$CSV && mon != MONITOR$AXIVITY) { # csv and not AX (so, GENEAcitv) - data2 = matrix(NA,nrow(data),3) - if (ncol(data) == 3) extra = 0 - if (ncol(data) >= 4) extra = 1 - for (jij in 1:3) { - data2[,jij] = data[,(jij+extra)] - } - Gx = data[,1]; Gy = data[,2]; Gz = data[,3] - } + S = data[(use + 1):LD,] # store leftover data - if (mon == MONITOR$GENEACTIV) { - if ("temperature" %in% colnames(data)) { - temperaturecolumn = which(colnames(data) == "temperature") #GGIRread - } else { - temperaturecolumn = 7 + # if S is only one row, it gets coersed to a vector, and what we need is a 1-row matrix + if (is.vector(S)) { + S = t(S) } - temperature = as.numeric(data[,temperaturecolumn]) - } else if (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) { - temperaturecolumn = 5 - temperature = as.numeric(data[,temperaturecolumn]) - } else if (dformat == FORMAT$AD_HOC_CSV && use.temp == TRUE) { - temperaturecolumn = 4 - temperature = as.numeric(data[,temperaturecolumn]) - } else if (mon == MONITOR$ACTIGRAPH) { - use.temp = FALSE - } else if (mon == MONITOR$MOVISENS) { - temperature = g.readtemp_movisens(datafile, params_general[["desiredtz"]], PreviousStartPage, PreviousEndPage, - interpolationType = params_rawdata[["interpolationType"]]) - data = cbind(data, temperature[1:nrow(data)]) - colnames(data)[4] = "temp" - temperaturecolumn = 4 } - if ((mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) || mon == MONITOR$MOVISENS || mon == MONITOR$AD_HOC) && - use.temp == TRUE) { - # ignore temperature if the values are unrealisticly high or NA - if (length(which(is.na(mean(as.numeric(data[1:10,temperaturecolumn]))) == T)) > 0) { - warning("\ntemperature ignored for auto-calibration because values are NA\n") - use.temp = FALSE - } else if (length(which(mean(as.numeric(data[1:10,temperaturecolumn])) > 120)) > 0) { + data = data[1:use,] + LD = nrow(data) #redefine LD because there is less data + + Gx = data[, "x"] + Gy = data[, "y"] + Gz = data[, "z"] + + if (use.temp) { + if (mean(data[1:10, "temperature"], na.rm = TRUE) > 120) { warning("\ntemperature ignored for auto-calibration because values are too high\n") use.temp = FALSE - } else if (sd(data[,temperaturecolumn], na.rm = TRUE) == 0) { + } else if (sd(data[, "temperature"], na.rm = TRUE) < 0.01) { warning("\ntemperature ignored for auto-calibration because no variance in values\n") use.temp = FALSE } } #============================================= - # non-integer sample frequency is a pain for deriving epoch based sd - # however, with an epoch of 10 seconds it is an integer number of samples per epoch - EN = sqrt(Gx^2 + Gy^2 + Gz^2) - D1 = g.downsample(EN,sf,ws4,ws2) - EN2 = D1$var2 - #mean acceleration - D1 = g.downsample(Gx,sf,ws4,ws2); GxM2 = D1$var2 - D1 = g.downsample(Gy,sf,ws4,ws2); GyM2 = D1$var2 - D1 = g.downsample(Gz,sf,ws4,ws2); GzM2 = D1$var2 - if (use.temp == TRUE) { - D1 = g.downsample(temperature,sf,ws4,ws2); - TemperatureM2 = D1$var2 - } - #sd acceleration - dim(Gx) = c(sf*ws4,ceiling(length(Gx)/(sf*ws4))); GxSD2 = apply(Gx,2,sd) - dim(Gy) = c(sf*ws4,ceiling(length(Gy)/(sf*ws4))); GySD2 = apply(Gy,2,sd) - dim(Gz) = c(sf*ws4,ceiling(length(Gz)/(sf*ws4))); GzSD2 = apply(Gz,2,sd) - #----------------------------------------------------- #expand 'out' if it is expected to be too short - if (count > (nrow(meta) - (2.5 * (3600/ws4) * 24))) { - extension = matrix(99999, ((3600/ws4) * 24), ncol(meta)) - meta = rbind(meta,extension) - if (verbose == TRUE) cat("\nvariabel meta extended\n") + if (count > (nrow(features) - (2.5 * (3600/calibEpochSize) * 24))) { + extension = matrix(99999, ((3600/calibEpochSize) * 24), ncol(features)) + features = rbind(features, extension) } - #storing in output matrix - meta[count:(count - 1 + length(EN2)), 1] = EN2 - meta[count:(count - 1 + length(EN2)), 2] = GxM2 - meta[count:(count - 1 + length(EN2)), 3] = GyM2 - meta[count:(count - 1 + length(EN2)), 4] = GzM2 - meta[count:(count - 1 + length(EN2)), 5] = GxSD2 - meta[count:(count - 1 + length(EN2)), 6] = GySD2 - meta[count:(count - 1 + length(EN2)), 7] = GzSD2 + # mean acceleration for EN, x, y, and z + EN = sqrt(Gx^2 + Gy^2 + Gz^2) + EN2 = rollMean(EN, sf, calibEpochSize) + endCount = count - 1 + length(EN2) + features[count:endCount, 1] = EN2 + features[count:endCount, 2] = rollMean(Gx, sf, calibEpochSize) + features[count:endCount, 3] = rollMean(Gy, sf, calibEpochSize) + features[count:endCount, 4] = rollMean(Gz, sf, calibEpochSize) + # sd acceleration + features[count:endCount, 5] = rollSD(Gx, sf, calibEpochSize) + features[count:endCount, 6] = rollSD(Gy, sf, calibEpochSize) + features[count:endCount, 7] = rollSD(Gz, sf, calibEpochSize) if (use.temp == TRUE) { - meta[count:(count - 1 + length(EN2)), 8] = TemperatureM2 + features[count:endCount, 8] = rollMean(data[, "temperature"], sf, calibEpochSize) } - count = count + length(EN2) #increasing "count": the indicator of how many seconds have been read + count = endCount + 1 # increasing count, the indicator of how many timesteps (calibEpochSize) have been read rm(Gx); rm(Gy); rm(Gz) # Update blocksize depending on available memory: BlocksizeNew = updateBlocksize(blocksize = blocksize, bsc_qc = bsc_qc) @@ -294,20 +202,19 @@ g.calibrate = function(datafile, params_rawdata = c(), } } else { LD = 0 #once LD < 1 the analysis stops, so this is a trick to stop it - # cat("\nstop reading because there is not enough data in this block\n") } spherepopulated = 0 - if (switchoffLD == 1) { + if (isLastBlock) { LD = 0 } - meta_temp = data.frame(V = meta, stringsAsFactors = FALSE) - cut = which(meta_temp[,1] == 99999) + features_temp = data.frame(V = features, stringsAsFactors = FALSE) + cut = which(features_temp[,1] == 99999) if (length(cut) > 0) { - meta_temp = meta_temp[-cut,] + features_temp = features_temp[-cut,] } - nhoursused = (nrow(meta_temp) * 10) / 3600 - if (nrow(meta_temp) > (params_rawdata[["minloadcrit"]] - 21)) { # enough data for the sphere? - meta_temp = meta_temp[-1,] + nhoursused = (nrow(features_temp) * 10) / 3600 + if (nrow(features_temp) > (params_rawdata[["minloadcrit"]] - 21)) { # enough data for the sphere? + features_temp = features_temp[-1,] #select parts with no movement if (mon == MONITOR$AD_HOC) { if (length(params_rawdata[["rmc.noise"]]) == 0) { @@ -322,28 +229,28 @@ g.calibrate = function(datafile, params_rawdata = c(), } else { sdcriter = 0.013 } - nomovement = which(meta_temp[,5] < sdcriter & meta_temp[,6] < sdcriter & meta_temp[,7] < sdcriter & - abs(as.numeric(meta_temp[,2])) < 2 & abs(as.numeric(meta_temp[,3])) < 2 & - abs(as.numeric(meta_temp[,4])) < 2) #the latter three are to reduce chance of including clipping periods + nomovement = which(features_temp[,5] < sdcriter & features_temp[,6] < sdcriter & features_temp[,7] < sdcriter & + abs(as.numeric(features_temp[,2])) < 2 & abs(as.numeric(features_temp[,3])) < 2 & + abs(as.numeric(features_temp[,4])) < 2) #the latter three are to reduce chance of including clipping periods if (length(nomovement) < 10) { # take only one row to trigger that autocalibration is skipped # with the QCmessage that there is not enough data - meta_temp = meta_temp[1, ] + features_temp = features_temp[1, ] } else { - meta_temp = meta_temp[nomovement,] + features_temp = features_temp[nomovement,] } - dup = which(rowSums(meta_temp[1:(nrow(meta_temp) - 1), 2:7] == meta_temp[2:nrow(meta_temp), 2:7]) == 3) # remove duplicated values - if (length(dup) > 0) meta_temp = meta_temp[-dup,] + dup = which(rowSums(features_temp[1:(nrow(features_temp) - 1), 2:7] == features_temp[2:nrow(features_temp), 2:7]) == 3) # remove duplicated values + if (length(dup) > 0) features_temp = features_temp[-dup,] rm(nomovement, dup) - if (min(dim(meta_temp)) > 1) { - meta_temp = meta_temp[(is.na(meta_temp[,4]) == F & is.na(meta_temp[,1]) == F),] - npoints = nrow(meta_temp) - cal.error.start = sqrt(as.numeric(meta_temp[,2])^2 + as.numeric(meta_temp[,3])^2 + as.numeric(meta_temp[,4])^2) + if (min(dim(features_temp)) > 1) { + features_temp = features_temp[(is.na(features_temp[,4]) == F & is.na(features_temp[,1]) == F),] + npoints = nrow(features_temp) + cal.error.start = sqrt(as.numeric(features_temp[,2])^2 + as.numeric(features_temp[,3])^2 + as.numeric(features_temp[,4])^2) cal.error.start = round(mean(abs(cal.error.start - 1)), digits = 5) #check whether sphere is well populated tel = 0 for (axis in 2:4) { - if ( min(meta_temp[,axis]) < -params_rawdata[["spherecrit"]] & max(meta_temp[,axis]) > params_rawdata[["spherecrit"]]) { + if ( min(features_temp[,axis]) < -params_rawdata[["spherecrit"]] & max(features_temp[,axis]) > params_rawdata[["spherecrit"]]) { tel = tel + 1 } } @@ -354,22 +261,20 @@ g.calibrate = function(datafile, params_rawdata = c(), QC = "recalibration not done because not enough points on all sides of the sphere" } } else { - if (verbose == TRUE) cat(" No non-movement found\n") QC = "recalibration not done because no non-movement data available" - meta_temp = c() + features_temp = c() } } else { QC = "recalibration not done because not enough data in the file or because file is corrupt" } if (spherepopulated == 1) { #only try to improve calibration if there are enough datapoints around the sphere #--------------------------------------------------------------------------- - # START of Zhou Fang's code (slightly edited by vtv21 to use matrix meta_temp from above + # START of Zhou Fang's code (slightly edited by vtv21 to use matrix features_temp from above # instead the similar matrix generated by Zhou Fang's original code. This to allow for - # more data to be used as meta_temp can now be based on 10 or more days of raw data - input = meta_temp[,2:4] #as.matrix() - if (use.temp == TRUE) { #at the moment i am always using temperature if mon == MONITOR$GENEACTIV - # mon == MONITOR$GENEACTIV & removed 19-11-2014 because it is redundant and would not allow for newer monitors to use it - inputtemp = cbind(as.numeric(meta_temp[,8]),as.numeric(meta_temp[,8]),as.numeric(meta_temp[,8])) #temperature + # more data to be used as features_temp can now be based on 10 or more days of raw data + input = features_temp[,2:4] + if (use.temp == TRUE) { + inputtemp = cbind(as.numeric(features_temp[,8]),as.numeric(features_temp[,8]),as.numeric(features_temp[,8])) #temperature } else { inputtemp = matrix(0, nrow(input), ncol(input)) #temperature, here used as a dummy variable } @@ -388,7 +293,6 @@ g.calibrate = function(datafile, params_rawdata = c(), scale(inputtemp, center = F, scale = 1/tempoffset)}, silent = TRUE) if (length(curr) == 0) { # set coefficients to default, because it did not work. - if (verbose == TRUE) cat("\nObject curr has length zero.") break } closestpoint = curr / sqrt(rowSums(curr^2)) @@ -429,28 +333,23 @@ g.calibrate = function(datafile, params_rawdata = c(), if (abs(res[iter + 1] - res[iter]) < tol) break } if (use.temp == FALSE) { - meta_temp2 = scale(as.matrix(meta_temp[,2:4]),center = -offset, scale = 1/scale) + features_temp2 = scale(as.matrix(features_temp[,2:4]),center = -offset, scale = 1/scale) } else { - yy = as.matrix(cbind(as.numeric(meta_temp[,8]),as.numeric(meta_temp[,8]),as.numeric(meta_temp[,8]))) - meta_temp2 = scale(as.matrix(meta_temp[,2:4]),center = -offset, scale = 1/scale) + + yy = as.matrix(cbind(as.numeric(features_temp[,8]),as.numeric(features_temp[,8]),as.numeric(features_temp[,8]))) + features_temp2 = scale(as.matrix(features_temp[,2:4]),center = -offset, scale = 1/scale) + scale(yy, center = rep(meantemp,3), scale = 1/tempoffset) } #equals: D2[,axis] = (D[,axis] + offset[axis]) / (1/scale[axis]) # END of Zhou Fang's code #------------------------------------------- - cal.error.end = sqrt(meta_temp2[,1]^2 + meta_temp2[,2]^2 + meta_temp2[,3]^2) - rm(meta_temp2) + cal.error.end = sqrt(features_temp2[,1]^2 + features_temp2[,2]^2 + features_temp2[,3]^2) + rm(features_temp2) cal.error.end = round(mean(abs(cal.error.end - 1)), digits = 5) # assess whether calibration error has sufficiently been improved if (cal.error.end < cal.error.start & cal.error.end < 0.01 & nhoursused > params_rawdata[["minloadcrit"]]) { #do not change scaling if there is no evidence that calibration improves - if (use.temp == TRUE && (mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) || - mon == MONITOR$MOVISENS || (mon == MONITOR$AD_HOC && use.temp == FALSE))) { + if (use.temp == temp.available) { QC = "recalibration done, no problems detected" - } else if (use.temp == FALSE && (mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) || - (mon == MONITOR$AXIVITY && dformat == FORMAT$CSV) || mon == MONITOR$MOVISENS || - (mon == MONITOR$AD_HOC && use.temp == TRUE))) { + } else { QC = "recalibration done, but temperature values not used" - } else if (mon != MONITOR$GENEACTIV) { - QC = "recalibration done, no problems detected" } LD = 0 #stop loading } else { @@ -471,8 +370,8 @@ g.calibrate = function(datafile, params_rawdata = c(), QC = "recalibration not done because recalibration does not decrease error" } } - if (all(dim(meta_temp)) != 0) { # change 2022-08-18 to handle when filetooshort = TRUE (7 columns, empty rows) - spheredata = data.frame(A = meta_temp, stringsAsFactors = TRUE) + if (all(dim(features_temp)) != 0) { # change 2022-08-18 to handle when filetooshort = TRUE (7 columns, empty rows) + spheredata = features_temp if (use.temp == TRUE) { names(spheredata) = c("Euclidean Norm","meanx","meany","meanz","sdx","sdy","sdz","temperature") } else { @@ -481,7 +380,7 @@ g.calibrate = function(datafile, params_rawdata = c(), } else { spheredata = c() } - rm(meta_temp) + rm(features_temp) QCmessage = QC if (params_rawdata[["printsummary"]] == TRUE & verbose == TRUE) { cat("\nSummary of autocalibration procedure:") diff --git a/R/g.dotorcomma.R b/R/g.dotorcomma.R index 7894c2bfb..011991994 100644 --- a/R/g.dotorcomma.R +++ b/R/g.dotorcomma.R @@ -1,4 +1,4 @@ -g.dotorcomma = function(inputfile, dformat, mon, desiredtz = "", loadGENEActiv = "GGIRread", ...) { +g.dotorcomma = function(inputfile, dformat, mon, ...) { #get input variables (relevant when read.myacc.csv is used) input = list(...) decn = getOption("OutDec") # extract system decimal separator @@ -17,7 +17,7 @@ g.dotorcomma = function(inputfile, dformat, mon, desiredtz = "", loadGENEActiv = if (exists("rmc.dec") == FALSE) rmc.dec = decn if (length(rmc.firstrow.acc) == 1) { dformat = FORMAT$AD_HOC_CSV - mon = MONITOR$MOVISENS + mon = MONITOR$AD_HOC decn = rmc.dec } if (dformat == FORMAT$CSV) { @@ -26,16 +26,17 @@ g.dotorcomma = function(inputfile, dformat, mon, desiredtz = "", loadGENEActiv = # lot of zeros, which makes it impossible to detect decimal separator # "." will then be the default, which is not correct for "," systems. while (skiprows < 1000000) { #foundnonzero == FALSE & - deci = as.matrix(read.csv(inputfile, skip = skiprows, nrow = 10)) + tmp = try(expr = {as.matrix(read.csv(inputfile, skip = skiprows, nrow = 10))}, silent = TRUE) + if (inherits(tmp, "try-error")) break # nothing left in the file to read + deci = tmp + skiprows = skiprows + 10000 if (length(unlist(strsplit(as.character(deci[2,2]), ","))) > 1) { decn = "," - break() + break } numtemp = as.numeric(deci[2,2]) - if (is.na(numtemp) == FALSE) { - if (numtemp != 0) break() - } + if (is.na(numtemp) == FALSE && numtemp != 0) break } if (!exists("deci")) stop("Problem with reading .csv file in GGIR function dotorcomma") if (is.na(suppressWarnings(as.numeric(deci[2,2]))) == T & decn == ".") decn = "," @@ -47,8 +48,8 @@ g.dotorcomma = function(inputfile, dformat, mon, desiredtz = "", loadGENEActiv = if (is.na(as.numeric(deci$data.out[2, 2])) == T & decn == ".") decn = "," } } else if (dformat == FORMAT$CWA) { - try(expr = {deci = GGIRread::readAxivity(filename = inputfile,start = 1, end = 10, desiredtz = desiredtz, - interpolationType = 1)$data},silent = TRUE) + try(expr = {deci = GGIRread::readAxivity(filename = inputfile, start = 1, end = 10, + interpolationType = 1)$data}, silent = TRUE) if (!exists("deci")) stop("Problem with reading .cwa file in GGIR function dotorcomma") if (is.na(suppressWarnings(as.numeric(deci[2,2]))) == T & decn == ".") decn = "," } else if (dformat == FORMAT$GT3X) { diff --git a/R/g.downsample.R b/R/g.downsample.R deleted file mode 100644 index ebd217c3a..000000000 --- a/R/g.downsample.R +++ /dev/null @@ -1,13 +0,0 @@ -g.downsample = function(sig,fs,ws3,ws2) { - #averaging per second => var1 - sig2 =cumsum(c(0,sig)) - select = seq(1,length(sig2),by=fs) - var1 = diff(sig2[round(select)]) / abs(diff(round(select[1:(length(select))]))) - #averaging per ws3 => var2 (e.g. 5 seconds) - select = seq(1,length(sig2),by=fs*ws3) - var2 = diff(sig2[round(select)]) / abs(diff(round(select[1:(length(select))]))) - #averaging per ws2 => var3 (e.g. 15 minutes) - select = seq(1,length(sig2),by=fs*ws2) - var3 = diff(sig2[round(select)]) / abs(diff(round(select[1:(length(select))]))) - invisible(list(var1=var1,var2=var2,var3=var3)) -} diff --git a/R/g.getmeta.R b/R/g.getmeta.R index 02da72e2d..68408d02f 100644 --- a/R/g.getmeta.R +++ b/R/g.getmeta.R @@ -68,11 +68,8 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), do.zcz = params_metrics[["do.zcz"]], do.brondcounts = params_metrics[["do.brondcounts"]], do.neishabouricounts = params_metrics[["do.neishabouricounts"]], - stringsAsFactors = TRUE) - if (length(params_rawdata[["chunksize"]]) == 0) params_rawdata[["chunksize"]] = 1 - if (params_rawdata[["chunksize"]] > 1.5) params_rawdata[["chunksize"]] = 1.5 - if (params_rawdata[["chunksize"]] < 0.2) params_rawdata[["chunksize"]] = 0.2 - gyro_available = FALSE + stringsAsFactors = FALSE) + nmetrics = sum(c(params_metrics[["do.bfen"]], params_metrics[["do.enmo"]], params_metrics[["do.lfenmo"]], params_metrics[["do.en"]], params_metrics[["do.hfen"]], params_metrics[["do.hfenplus"]], @@ -94,44 +91,30 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), if (length(myfun) != 0) { nmetrics = nmetrics + length(myfun$colnames) # check myfun object already, because we do not want to discover - # bugs after waiting for the data to be load + # bugs after waiting for the data to load check_myfun(myfun, params_general[["windowsizes"]]) } if (length(nmetrics) == 0) { - if (verbose == TRUE) cat("\nWARNING: No metrics selected\n") + warning("No metrics selected.", call. = FALSE) } - filename = unlist(strsplit(as.character(datafile),"/")) - filename = filename[length(filename)] - # parameters - ws3 = params_general[["windowsizes"]][1]; ws2 = params_general[["windowsizes"]][2]; ws = params_general[["windowsizes"]][3] #window sizes - if ((ws2/60) != round(ws2/60)) { - ws2 = as.numeric(60 * ceiling(ws2/60)) - if (verbose == TRUE) { - cat("\nWARNING: The long windowsize needs to be a multitude of 1 minute periods. The\n") - cat(paste0("\nlong windowsize has now been automatically adjusted to: ", ws2, " seconds in order to meet this criteria.\n")) - } - } - if ((ws2/ws3) != round(ws2/ws3)) { - def = c(1,5,10,15,20,30,60) - def2 = abs(def - ws3) - ws3 = as.numeric(def[which(def2 == min(def2))]) - if (verbose == TRUE) { - cat("\nWARNING: The long windowsize needs to be a multitude of short windowsize. The \n") - cat(paste0("\nshort windowsize has now been automatically adjusted to: ", ws3, " seconds in order to meet this criteria.\n")) - } - } - params_general[["windowsizes"]] = c(ws3,ws2,ws) - data = PreviousEndPage = PreviousStartPage = starttime = wday = wdayname = c() + + ws3 = params_general[["windowsizes"]][1]; ws2 = params_general[["windowsizes"]][2]; ws = params_general[["windowsizes"]][3] + + PreviousEndPage = c() filequality = data.frame(filetooshort = FALSE, filecorrupt = FALSE, - filedoesnotholdday = FALSE, NFilePagesSkipped = 0, stringsAsFactors = TRUE) + filedoesnotholdday = FALSE, NFilePagesSkipped = 0) + filetooshort = FALSE + filecorrupt = FALSE + filedoesnotholdday = FALSE + NFilePagesSkipped = 0 + i = 1 #counter to keep track of which binary block is being read count = 1 #counter to keep track of the number of seconds that have been read count2 = 1 #count number of blocks read with length "ws2" (long epoch, 15 minutes by default) LD = 2 #dummy variable used to identify end of file and to make the process stop bsc_qc = data.frame(time = c(), size = c(), stringsAsFactors = FALSE) - # inspect file if (length(inspectfileobject) > 0) { INFI = inspectfileobject @@ -142,81 +125,45 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), mon = INFI$monc dformat = INFI$dformc sf = INFI$sf - if (is.null(sf)) { - filequality$filecorrupt = TRUE - filecorrupt = TRUE - filetooshort = FALSE - filedoesnotholdday = FALSE - NFilePagesSkipped = 0 - QClog = NULL - LD = 0 # to prevent while loop and skip reading of file - deviceSerialNumber = "notExtracted" - } - if (LD > 1) { - hvars = g.extractheadervars(INFI) - deviceSerialNumber = hvars$deviceSerialNumber - } - # if GENEActiv csv, deprecated function - if (mon == MONITOR$GENEACTIV && dformat == FORMAT$CSV & length(params_rawdata[["rmc.firstrow.acc"]]) == 0) { - stop("The GENEActiv csv reading functionality is deprecated in GGIR from the version 2.6-4 onwards. Please, use either the GENEActiv bin files or try to read the csv files with GGIR::read.myacc.csv") - } - if (mon == MONITOR$ACTIGRAPH) { - # If Actigraph then try to specify dynamic range based on Actigraph model - if (length(grep(pattern = "CLE", x = deviceSerialNumber)) == 1) { - params_rawdata[["dynrange"]] = 6 - } else if (length(grep(pattern = "MOS", x = deviceSerialNumber)) == 1) { - params_rawdata[["dynrange"]] = 8 - } else if (length(grep(pattern = "NEO", x = deviceSerialNumber)) == 1) { - params_rawdata[["dynrange"]] = 6 - } - } - if (LD > 1) { - if (sf == 0) stop("Sample frequency not recognised") #assume 80Hertz in the absense of any other info - header = INFI$header - ID = hvars$ID - - # get now-wear, clip, and blocksize parameters (thresholds) - ncb_params = get_nw_clip_block_params(chunksize = params_rawdata[["chunksize"]], - dynrange = params_rawdata[["dynrange"]], - mon, rmc.noise = params_rawdata[["rmc.noise"]], - sf, dformat, - rmc.dynamic_range = params_rawdata[["rmc.dynamic_range"]]) - clipthres = ncb_params$clipthres - blocksize = ncb_params$blocksize - sdcriter = ncb_params$sdcriter - racriter = ncb_params$racriter - n_decimal_places = 4 # number of decimal places to which features should be rounded - #creating matrixes for storing output - S = matrix(0,0,4) #dummy variable needed to cope with head-tailing succeeding blocks of data - nev = 80*10^7 # number expected values - # NR = ceiling((90*10^6) / (sf*ws3)) + 1000 #NR = number of 'ws3' second rows (this is for 10 days at 80 Hz) - NR = ceiling(nev / (sf*ws3)) + 1000 #NR = number of 'ws3' second rows (this is for 10 days at 80 Hz) - metashort = matrix(" ",NR,(1 + nmetrics)) #generating output matrix for acceleration signal - if (mon == MONITOR$ACTIGRAPH || mon == MONITOR$VERISENSE || (mon == MONITOR$AXIVITY && dformat == FORMAT$CSV) || - (mon == MONITOR$AD_HOC && length(params_rawdata[["rmc.col.temp"]]) == 0)) { - temp.available = FALSE - } else if (mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) || - mon == MONITOR$MOVISENS || (mon == MONITOR$AD_HOC && length(params_rawdata[["rmc.col.temp"]]) > 0)) { - temp.available = TRUE - } - QClog = NULL - if (temp.available == FALSE) { - metalong = matrix(" ", ((nev/(sf*ws2)) + 100), 4) #generating output matrix for 15 minutes summaries - } else if (temp.available == TRUE && mon != MONITOR$MOVISENS && mon != MONITOR$AD_HOC) { - metalong = matrix(" ", ((nev/(sf*ws2)) + 100), 7) #generating output matrix for 15 minutes summaries - } else if (temp.available == TRUE && (mon == MONITOR$MOVISENS || mon == MONITOR$AD_HOC)) { - metalong = matrix(" ", ((nev/(sf*ws2)) + 100), 5) #generating output matrix for 15 minutes summaries - } - #=============================================== - # Read file - switchoffLD = 0 #dummy variable part "end of loop mechanism" - sforiginal = sf - } else { - metalong = NULL - metashort = NULL + + if (is.null(sf)) { # sf is NULL for corrupt files + return(invisible(list(filecorrupt = TRUE, filetooshort = FALSE, NFilePagesSkipped = 0, + metalong = c(), metashort = c(), wday = c(), wdayname = c(), + windowsizes = c(), bsc_qc = bsc_qc, QClog = NULL))) } + + hvars = g.extractheadervars(INFI) + deviceSerialNumber = hvars$deviceSerialNumber + + # get now-wear, clip, and blocksize parameters (thresholds) + ncb_params = get_nw_clip_block_params(chunksize = params_rawdata[["chunksize"]], + dynrange = params_rawdata[["dynrange"]], + monc = mon, dformat = dformat, + deviceSerialNumber = deviceSerialNumber, + rmc.noise = params_rawdata[["rmc.noise"]], + sf = sf, + rmc.dynamic_range = params_rawdata[["rmc.dynamic_range"]]) + clipthres = ncb_params$clipthres + blocksize = ncb_params$blocksize + sdcriter = ncb_params$sdcriter + racriter = ncb_params$racriter + n_decimal_places = 4 # number of decimal places to which features should be rounded + # creating matrices for storing output + S = matrix(0,0,4) #dummy variable needed to cope with head-tailing succeeding blocks of data + nev = 80*10^7 # number expected values + NR = ceiling(nev / (sf*ws3)) + 1000 #NR = number of 'ws3' second rows (this is for 10 days at 80 Hz) + metashort = matrix(" ",NR,(1 + nmetrics)) #generating output matrix for acceleration signal + QClog = NULL + + #=============================================== + # Read file + isLastBlock = FALSE # dummy variable part "end of loop mechanism" + + PreviousLastValue = c(0, 0, 1) + PreviousLastTime = NULL + header = NULL + while (LD > 1) { - P = c() if (verbose == TRUE) { if (i == 1) { cat(paste0("\nLoading chunk: ", i)) @@ -225,121 +172,92 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), } } - options(warn = -1) #turn off warnings (code complains about unequal rowlengths - - if (!exists("PreviousLastValue")) PreviousLastValue = c(0, 0, 1) - if (!exists("PreviousLastTime")) PreviousLastTime = NULL accread = g.readaccfile(filename = datafile, blocksize = blocksize, blocknumber = i, filequality = filequality, ws = ws, PreviousEndPage = PreviousEndPage, inspectfileobject = INFI, PreviousLastValue = PreviousLastValue, PreviousLastTime = PreviousLastTime, - params_rawdata = params_rawdata, params_general = params_general) + params_rawdata = params_rawdata, params_general = params_general, + header = header) + header = accread$header + if ("PreviousLastValue" %in% names(accread$P)) { # output when reading ad-hoc csv - P = accread$P[1:2] PreviousLastValue = accread$P$PreviousLastValue PreviousLastTime = accread$P$PreviousLastTime - } else { - P = accread$P - } + } + filequality = accread$filequality filetooshort = filequality$filetooshort filecorrupt = filequality$filecorrupt filedoesnotholdday = filequality$filedoesnotholdday NFilePagesSkipped = filequality$NFilePagesSkipped - switchoffLD = accread$switchoffLD + isLastBlock = accread$isLastBlock PreviousEndPage = accread$endpage - startpage = accread$startpage - options(warn = -1) # to ignore warnings relating to failed mmap.load attempt - rm(accread); gc() - options(warn = 0) # to ignore warnings relating to failed mmap.load attempt - if (mon == MONITOR$MOVISENS) { # if movisens, then read temperature - PreviousStartPage = startpage - temperature = g.readtemp_movisens(datafile, desiredtz = params_general[["desiredtz"]], PreviousStartPage, - PreviousEndPage, interpolationType = params_rawdata[["interpolationType"]]) - P = cbind(P, temperature[1:nrow(P)]) - colnames(P)[4] = "temp" - } - options(warn = 0) #turn on warnings + #============ #process data as read from binary file - if (length(P) > 0) { #would have been set to zero if file was corrupt or empty - - if (mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) { - data = P$data.out - } else if (dformat == FORMAT$CSV) { #csv (Actigraph or ad-hoc csv format) - if (params_rawdata[["imputeTimegaps"]] == TRUE) { - P = as.data.frame(P) - if (ncol(P) == 3) { - timeCol = c() - xyzCol = names(P)[1:3] - } else if (ncol(P) == 4) { - timeCol = names(P)[1] - xyzCol = names(P)[2:4] + if (length(accread$P) > 0) { # would have been set to zero if file was corrupt or empty + data = accread$P$data + QClog = rbind(QClog, accread$P$QClog) + rm(accread) + + if (i == 1) { + light.available = ("light" %in% colnames(data)) + + use.temp = ("temperature" %in% colnames(data)) + if (use.temp) { + if (mean(data$temperature[1:10], na.rm = TRUE) > 50) { + warning("temperature value is unreaslistically high (> 50 Celcius) and will not be used.", call. = FALSE) + use.temp = FALSE } - if (!exists("PreviousLastValue")) PreviousLastValue = c(0, 0, 1) - if (!exists("PreviousLastTime")) PreviousLastTime = NULL - P = g.imputeTimegaps(P, xyzCol = xyzCol, timeCol = timeCol, sf = sf, k = 0.25, - PreviousLastValue = PreviousLastValue, - PreviousLastTime = PreviousLastTime, - epochsize = c(ws3, ws2)) - QClog = rbind(QClog, P$QClog) - P = P$x - PreviousLastValue = as.numeric(P[nrow(P), xyzCol]) - if (is.null(timeCol)) PreviousLastTime = NULL else PreviousLastTime = as.POSIXct(P[nrow(P), timeCol]) } - data = P - } else if (dformat == FORMAT$CWA) { - if (P$header$hardwareType == "AX6") { - # GGIR now ignores the AX6 gyroscope signals until added value has robustly been demonstrated. - # Note however that while AX6 is able to collect gyroscope data, it can also be configured - # to only collect accelerometer data, so only remove gyro data if it's present. - if (ncol(P$data) == 10) { - data = P$data[,-c(2:4)] - P$data = P$data[1:min(100,nrow(P$data)),-c(2:4)] # trim object, because rest of data is not needed anymore - } else { - data = P$data - P$data = P$data[1:min(100,nrow(P$data)),] # trim object, because rest of data is not needed anymore - } - gyro_available = FALSE - # If we ever want to use gyroscope data then - # comment out this if statement and set gyro_available = TRUE + + # output matrix for 15 minutes summaries + if (!use.temp && !light.available) { + metalong = matrix(" ", ((nev/(sf*ws2)) + 100), 4) + metricnames_long = c("timestamp","nonwearscore","clippingscore","EN") + } else if (use.temp && !light.available) { + metalong = matrix(" ", ((nev/(sf*ws2)) + 100), 5) + metricnames_long = c("timestamp","nonwearscore","clippingscore","temperaturemean","EN") + } else if (!use.temp && light.available) { + metalong = matrix(" ", ((nev/(sf*ws2)) + 100), 6) + metricnames_long = c("timestamp","nonwearscore","clippingscore","lightmean","lightpeak","EN") + } else if (use.temp && light.available) { + metalong = matrix(" ", ((nev/(sf*ws2)) + 100), 7) + metricnames_long = c("timestamp","nonwearscore","clippingscore","lightmean","lightpeak","temperaturemean","EN") + } + } + + if (params_rawdata[["imputeTimegaps"]] && (dformat == FORMAT$CSV || dformat == FORMAT$GT3X)) { + P = g.imputeTimegaps(data, sf = sf, k = 0.25, + PreviousLastValue = PreviousLastValue, + PreviousLastTime = PreviousLastTime, + epochsize = c(ws3, ws2)) + data = P$x + PreviousLastValue = data[nrow(data), c("x", "y", "z")] + if ("time" %in% colnames(data)) { + PreviousLastTime = as.POSIXct(data$time[nrow(data)]) } else { - data = P$data - P$data = P$data[1:min(100,nrow(P$data)),] # trim object, because rest of data is not needed anymore + PreviousLastTime = NULL } QClog = rbind(QClog, P$QClog) - } else if (dformat == FORMAT$AD_HOC_CSV) { - data = P$data - } else if (mon == MONITOR$MOVISENS) { - data = P - } else if (dformat == FORMAT$GT3X) { - if (params_rawdata[["imputeTimegaps"]] == TRUE) { - if (!exists("PreviousLastValue")) PreviousLastValue = c(0, 0, 1) - if (!exists("PreviousLastTime")) PreviousLastTime = NULL - P = g.imputeTimegaps(P, xyzCol = c("X", "Y", "Z"), timeCol = "time", sf = sf, k = 0.25, - PreviousLastValue = PreviousLastValue, - PreviousLastTime = PreviousLastTime, - epochsize = c(ws3, ws2)) - QClog = rbind(QClog, P$QClog) - P = P$x - PreviousLastValue = as.numeric(P[nrow(P), c("X", "Y", "Z")]) - PreviousLastTime = as.POSIXct(P[nrow(P), "time"]) - } - data = P[,2:ncol(P)] + rm(P) } + + gc() + data = as.matrix(data, rownames.force = FALSE) - #add left over data from last time + #add leftover data from last time if (nrow(S) > 0) { - if (params_rawdata[["imputeTimegaps"]] == TRUE) { + if (params_rawdata[["imputeTimegaps"]]) { if ("remaining_epochs" %in% colnames(data)) { if (ncol(S) == (ncol(data) - 1)) { # this block has time gaps while the previous block did not S = cbind(S, 1) colnames(S)[4] = "remaining_epochs" } - } else if ("remaining_epochs" %in% colnames(S) == TRUE) { + } else if ("remaining_epochs" %in% colnames(S)) { if ((ncol(S) - 1) == ncol(data)) { # this block does not have time gaps while the previous blog did data = cbind(data, 1) @@ -347,40 +265,27 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), } } } - data = suppressWarnings(rbind(S,data)) # suppress warnings about string as factor + data = rbind(S,data) } - SWMT = get_starttime_weekday_meantemp_truncdata(temp.available, mon, dformat, - data, - P, header, desiredtz = params_general[["desiredtz"]], - sf, i, datafile, ws2, - starttime, wday, wdayname, configtz = params_general[["configtz"]]) - starttime = SWMT$starttime - meantemp = SWMT$meantemp - use.temp = SWMT$use.temp - wday = SWMT$wday; wdayname = SWMT$wdayname - params_general[["desiredtz"]] = SWMT$desiredtz; data = SWMT$data - - if (mon == MONITOR$ACTIGRAPH || mon == MONITOR$VERISENSE || - (mon == MONITOR$AXIVITY && dformat == FORMAT$CSV) || - (mon == MONITOR$AD_HOC && use.temp == FALSE)) { - metricnames_long = c("timestamp","nonwearscore","clippingscore","en") - } else if (mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA)) { - metricnames_long = c("timestamp","nonwearscore","clippingscore","lightmean","lightpeak","temperaturemean","EN") - } else if (mon == MONITOR$MOVISENS || (mon == MONITOR$AD_HOC & use.temp == TRUE)) { - # at the moment read.myacc.csv does not facilitate extracting light data, so only temperature is used - metricnames_long = c("timestamp","nonwearscore","clippingscore","temperaturemean","EN") + if (i == 1) { + SWMT = get_starttime_weekday_truncdata(mon, dformat, + data, header, desiredtz = params_general[["desiredtz"]], + sf, datafile, ws2, + configtz = params_general[["configtz"]]) + starttime = SWMT$starttime + wday = SWMT$wday; wdayname = SWMT$wdayname + data = SWMT$data + + rm(SWMT) } - rm(SWMT) - if (exists("P")) rm(P); gc() - if (i != 0 & exists("P")) rm(P); gc() + LD = nrow(data) - if (LD < (ws*sf) & i == 1) { + if (LD < (ws*sf) && i == 1) { warning('\nWarning data too short for doing non-wear detection 3\n') - switchoffLD = 1 + isLastBlock = TRUE LD = 0 #ignore rest of the data and store what has been loaded so far. } - #store data that could not be used for this block, but will be added to next block if (LD >= (ws*sf)) { @@ -392,7 +297,7 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), } } if ((LD - use) > 1) { - S = data[(use + 1):LD,] #store left over + S = data[(use + 1):LD,] #store leftover data if (ncol(S) == 1) { S = t(S) } @@ -411,101 +316,24 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), dur = nrow(data) #duration of experiment in data points durexp = nrow(data) / (sf*ws) #duration of experiment in hrs #-------------------------------------------- - if (mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) || mon == MONITOR$MOVISENS || - (mon == MONITOR$AD_HOC && length(params_rawdata[["rmc.col.temp"]]) > 0)) { - if (mon == MONITOR$GENEACTIV) { - temperaturecolumn = 6; lightcolumn = 5 - } else if (mon == MONITOR$AXIVITY) { - temperaturecolumn = 5; lightcolumn = 7 - if (gyro_available == TRUE) { - temperaturecolumn = temperaturecolumn + 3 - lightcolumn = lightcolumn + 3 - } - } else if (mon == MONITOR$MOVISENS) { - temperaturecolumn = 4 - } else if (mon == MONITOR$AD_HOC) { - temperaturecolumn = which(colnames(data) == "temperature") - } - if (mon != MONITOR$AD_HOC && mon != MONITOR$MOVISENS) { - light = as.numeric(data[, lightcolumn]) - } - if (mon == MONITOR$AD_HOC && length(params_rawdata[["rmc.col.wear"]]) > 0) { - wearcol = as.character(data[, which(colnames(data) == "wear")]) - suppressWarnings(storage.mode(wearcol) <- "logical") - } - temperature = as.numeric(data[, temperaturecolumn]) + temperature = light = c() + if (light.available) { + light = data[, "light"] } - # Initialization of variables - data_scaled = FALSE - if (mon == MONITOR$ACTIGRAPH && dformat == FORMAT$GT3X) { - data = data[, 1:3] - data[, 1:3] = scale(data[, 1:3], center = -offset, scale = 1/scale) #rescale data - data_scaled = TRUE - } else if (mon == MONITOR$AXIVITY && (dformat == FORMAT$CWA || dformat == FORMAT$CSV)) { - extraction_succeeded = FALSE - if (gyro_available == TRUE) { - data[,5:7] = scale(data[,5:7],center = -offset, scale = 1/scale) #rescale data - extraction_succeeded = TRUE - data = data[, 2:7] - } - if (extraction_succeeded == FALSE) { - data[, 2:4] = scale(data[, 2:4],center = -offset, scale = 1/scale) #rescale data - data = data[,2:4] - } - data_scaled = TRUE - } else if (mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) { - yy = as.matrix(cbind(as.numeric(data[,temperaturecolumn]), - as.numeric(data[,temperaturecolumn]), - as.numeric(data[,temperaturecolumn]))) - data = data[,2:4] - data[,1:3] = scale(as.matrix(data[,1:3]),center = -offset, scale = 1/scale) + - scale(yy, center = rep(meantemp,3), scale = 1/tempoffset) #rescale data - rm(yy); gc() - data_scaled = TRUE - } else if (mon == MONITOR$MOVISENS) { - yy = as.matrix(cbind(as.numeric(data[,4]),as.numeric(data[,4]),as.numeric(data[,4]))) - data = data[,1:3] - data[,1:3] = scale(as.matrix(data[,1:3]),center = -offset, scale = 1/scale) + - scale(yy, center = rep(meantemp,3), scale = 1/tempoffset) #rescale data - rm(yy); gc() - data_scaled = TRUE - } else if ((dformat == FORMAT$CSV || dformat == FORMAT$AD_HOC_CSV) && (mon != MONITOR$AXIVITY)) { - # Any brand that is not Axivity with csv or Movisense format data - if (mon == MONITOR$GENEACTIV || (mon == MONITOR$AD_HOC && use.temp == TRUE)) { - tempcolumnvalues = as.numeric(as.character(data[,temperaturecolumn])) - yy = as.matrix(cbind(tempcolumnvalues, tempcolumnvalues, tempcolumnvalues)) - meantemp = mean(as.numeric(data[,temperaturecolumn])) - if (length(meantempcal) == 0) meantempcal = meantemp - } - if (ncol(data) == 3) data = data[,1:3] - if (ncol(data) >= 4) { - data = data[,2:4] - if (is(data[,1], "character")) { - data = apply(data, 2,as.numeric) - } - } - suppressWarnings(storage.mode(data) <- "numeric") - if ((mon == MONITOR$ACTIGRAPH || mon == MONITOR$AD_HOC || mon == MONITOR$VERISENSE) && use.temp == FALSE) { - data = scale(data,center = -offset, scale = 1/scale) #rescale data - } else if ((mon == MONITOR$GENEACTIV || mon == MONITOR$AD_HOC) && use.temp == TRUE) { - # meantemp replaced by meantempcal # 19-12-2013 - data = scale(data,center = -offset, scale = 1/scale) + - scale(yy, center = rep(meantempcal,3), scale = 1/tempoffset) #rescale data - rm(yy); gc() - } - data_scaled = TRUE + if (use.temp) { + temperature = data[, "temperature"] } - if (data_scaled == FALSE) { - warning(paste0("\nAutocalibration was not applied, this should", - "not happen please contact GGIR maintainers")) + + # rescale data + data[, c("x", "y", "z")] = scale(data[, c("x", "y", "z")], center = -offset, scale = 1/scale) + if (use.temp && length(meantempcal) > 0) { + yy = cbind(temperature, + temperature, + temperature) + data[, c("x", "y", "z")] = data[, c("x", "y", "z")] + scale(yy, center = rep(meantempcal,3), scale = 1/tempoffset) } - suppressWarnings(storage.mode(data) <- "numeric") - ## resample experiment to see whehter processing time can be much improved if data is resampled - sfold = sforiginal # keep sf, because light, temperature are not resampled at the moment - # STORE THE RAW DATA - # data[,1], data[,2], data[,3], starttime, (temperature, light) - - EN = sqrt(data[,1]^2 + data[,2]^2 + data[,3]^2) # Do not delete Used for long epoch calculation + + EN = sqrt(data[, "x"]^2 + data[, "y"]^2 + data[, "z"]^2) # Do not delete Used for long epoch calculation accmetrics = g.applymetrics(data = data, sf = sf, ws3 = ws3, metrics2do = metrics2do, @@ -551,7 +379,6 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), metashort = rbind(metashort,extension) extension2 = matrix(" ", ((3600/ws2) * 24) + (totalgap * (ws2/ws3)), ncol(metalong)) #add another day to metashort once you reach the end of it metalong = rbind(metalong, extension2) - if (verbose == TRUE) cat("\nvariable metashort extended\n") } col_msi = 2 # Add metric time series to metashort object @@ -593,7 +420,7 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), blocksize = BlocksizeNew$blocksize ##================================================== # MODULE 2 - non-wear time & clipping - NWCW = detect_nonwear_clipping(data = data, windowsizes = c(ws3, ws2, ws), sf = sfold, + NWCW = detect_nonwear_clipping(data = data, windowsizes = c(ws3, ws2, ws), sf = sf, clipthres = clipthres, sdcriter = sdcriter, racriter = racriter, nonwear_approach = params_cleaning[["nonwear_approach"]], params_rawdata = params_rawdata) @@ -601,49 +428,49 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), # metalong col_mli = 2 metalong[count2:((count2 - 1) + length(NWav)),col_mli] = NWav; col_mli = col_mli + 1 - metalong[(count2):((count2 - 1) + length(NWav)),col_mli] = CWav; col_mli = col_mli + 1 - if (mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) || - mon == MONITOR$MOVISENS || (mon == MONITOR$AD_HOC && length(params_rawdata[["rmc.col.temp"]]) != 0)) { # going from sample to ws2 - if (mon == MONITOR$GENEACTIV || mon == MONITOR$AXIVITY) { - #light (running mean) - lightc = cumsum(c(0,light)) - select = seq(1, length(lightc), by = (ws2 * sfold)) - lightmean = diff(lightc[round(select)]) / abs(diff(round(select))) - rm(lightc); gc() - #light (running max) - lightmax = matrix(0, length(lightmean), 1) - for (li in 1:(length(light)/(ws2*sfold))) { - tempm = max(light[((li - 1) * (ws2 * sfold)):(li * (ws2 * sfold))]) - if (length(tempm) > 0) { - lightmax[li] = tempm[1] - } else { - lightmax[li] = max(light[((li - 1) * (ws2 * sfold)):(li * (ws2 * sfold))]) - } + metalong[count2:((count2 - 1) + length(NWav)),col_mli] = CWav; col_mli = col_mli + 1 + + if(light.available) { + #light (running mean) + lightc = cumsum(c(0,light)) + select = seq(1, length(lightc), by = (ws2 * sf)) + lightmean = diff(lightc[round(select)]) / abs(diff(round(select))) + rm(lightc) + #light (running max) + lightmax = matrix(0, length(lightmean), 1) + for (li in 1:(length(light)/(ws2*sf))) { + tempm = max(light[((li - 1) * (ws2 * sf)):(li * (ws2 * sf))]) + if (length(tempm) > 0) { + lightmax[li] = tempm[1] + } else { + lightmax[li] = max(light[((li - 1) * (ws2 * sf)):(li * (ws2 * sf))]) } } + + metalong[(count2):((count2 - 1) + length(NWav)), col_mli] = round(lightmean, digits = n_decimal_places) + col_mli = col_mli + 1 + metalong[(count2):((count2 - 1) + length(NWav)), col_mli] = round(lightmax, digits = n_decimal_places) + col_mli = col_mli + 1 + } + + if(use.temp) { #temperature (running mean) temperaturec = cumsum(c(0, temperature)) - select = seq(1, length(temperaturec), by = (ws2 * sfold)) + select = seq(1, length(temperaturec), by = (ws2 * sf)) temperatureb = diff(temperaturec[round(select)]) / abs(diff(round(select))) - rm(temperaturec); gc() + rm(temperaturec) + + metalong[(count2):((count2 - 1) + length(NWav)), col_mli] = round(temperatureb, digits = n_decimal_places) + col_mli = col_mli + 1 } + #EN going from sample to ws2 ENc = cumsum(c(0, EN)) select = seq(1, length(ENc), by = (ws2 * sf)) #<= EN is derived from data, so it needs the new sf ENb = diff(ENc[round(select)]) / abs(diff(round(select))) - rm(ENc, EN); gc() - if (mon == MONITOR$GENEACTIV || (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA)) { - metalong[(count2):((count2 - 1) + length(NWav)), col_mli] = round(lightmean, digits = n_decimal_places) - col_mli = col_mli + 1 - metalong[(count2):((count2 - 1) + length(NWav)), col_mli] = round(lightmax, digits = n_decimal_places) - col_mli = col_mli + 1 - metalong[(count2):((count2 - 1) + length(NWav)), col_mli] = round(temperatureb, digits = n_decimal_places) - col_mli = col_mli + 1 - } else if (mon == MONITOR$MOVISENS || (mon == MONITOR$AD_HOC && length(params_rawdata[["rmc.col.temp"]]) != 0)) { - metalong[(count2):((count2 - 1) + length(NWav)), col_mli] = round(temperatureb, digits = n_decimal_places) - col_mli = col_mli + 1 - } + rm(ENc, EN) metalong[(count2):((count2 - 1) + length(NWav)), col_mli] = round(ENb, digits = n_decimal_places) + if (exists("remaining_epochs")) { # Impute long gaps at epoch levels, because imputing them at raw level would # be too memory hungry @@ -686,13 +513,13 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), # metalong metalong = impute_at_epoch_level(gapsize = floor(remaining_epochs[gaps_to_fill] * (ws3/ws2)) + 1, # plus 1 needed to count for current epoch timeseries = metalong, - gap_index = floor(gaps_to_fill / (ws2 * sfold)) + count2, # Using floor so that the gap is filled in the epoch in which it is occurring + gap_index = floor(gaps_to_fill / (ws2 * sf)) + count2, # Using floor so that the gap is filled in the epoch in which it is occurring metnames = metricnames_long) # metashort # added epoch-level nonwear to metashort to get it imputed, then remove it metashort = impute_at_epoch_level(gapsize = remaining_epochs[gaps_to_fill], # gapsize in epochs timeseries = metashort, - gap_index = floor(gaps_to_fill / (ws3 * sfold)) + count, # Using floor so that the gap is filled in the epoch in which it is occurring + gap_index = floor(gaps_to_fill / (ws3 * sf)) + count, # Using floor so that the gap is filled in the epoch in which it is occurring metnames = c("timestamp", metnames)) # epoch level index of gap nr_after = c(nrow(metalong), nrow(metashort)) count2 = count2 + (nr_after[1] - nr_before[1]) @@ -711,7 +538,7 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), LD = 0 #once LD < 1 the analysis stops, so this is a trick to stop it # stop reading because there is not enough data in this block } - if (switchoffLD == 1) LD = 0 + if (isLastBlock) LD = 0 if (ceiling(daylimit) != FALSE) { if (i == ceiling(daylimit)) { #to speed up testing only read first 'i' blocks of data LD = 0 #once LD < 1 the analysis stops, so this is a trick to stop it @@ -721,38 +548,30 @@ g.getmeta = function(datafile, params_metrics = c(), params_rawdata = c(), i = i + 1 #go to next block } # deriving timestamps - if (filecorrupt == FALSE & filetooshort == FALSE & filedoesnotholdday == FALSE) { + if (!filecorrupt && !filetooshort && !filedoesnotholdday) { cut = count:nrow(metashort) if (length(cut) > 1) { - tmp = metashort[-cut,] + metashort = metashort[-cut,] # for a very small file, there could be just one row in metashort[-cut,], so it gets coerced to a vector. # But what we actually need is a 1-row matrix. So we need to transpose it. - if(is.vector(tmp)) { - metashort = as.matrix(t(tmp)) - } else { - metashort = as.matrix(tmp) + if(is.vector(metashort)) { + metashort = as.matrix(t(metashort)) } } if (nrow(metashort) > 1) { starttime3 = round(as.numeric(starttime)) #numeric time but relative to the desiredtz time5 = seq(starttime3, (starttime3 + ((nrow(metashort) - 1) * ws3)), by = ws3) - if (length(params_general[["desiredtz"]]) == 0) { - warning("\ndesiredtz not specified, system timezone used as default") - params_general[["desiredtz"]] = "" - } time6 = as.POSIXlt(time5,origin = "1970-01-01", tz = params_general[["desiredtz"]]) time6 = strftime(time6, format = "%Y-%m-%dT%H:%M:%S%z") metashort[,1] = as.character(time6) } cut2 = count2:nrow(metalong) if (length(cut2) > 1) { - tmp = metalong[-cut2,] + metalong = metalong[-cut2,] # for a very small file, there could be just one row in metalong[-cut2,], so it gets coerced to a vector. # But what we actually need is a 1-row matrix. So we need to transpose it. - if(is.vector(tmp)) { - metalong = as.matrix(t(tmp)) - } else { - metalong = as.matrix(tmp) + if(is.vector(metalong)) { + metalong = as.matrix(t(metalong)) } } if (nrow(metalong) > 2) { diff --git a/R/g.getstarttime.R b/R/g.getstarttime.R index f3149ab67..0a143201b 100644 --- a/R/g.getstarttime.R +++ b/R/g.getstarttime.R @@ -1,146 +1,140 @@ -g.getstarttime = function(datafile, P, header, mon, dformat, desiredtz, configtz = NULL) { - #get input variables (relevant when read.myacc.csv is used) - #------------------------------------------------------------ - if (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) { - starttime = P$data[1,1] - starttime = as.POSIXlt(starttime, tz = desiredtz, origin = "1970-01-01") - } else if (mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) { - starttime = as.POSIXlt(P$data.out$time[1], tz = desiredtz, origin = "1970-01-01") - } else if (dformat == FORMAT$CSV && (mon == MONITOR$ACTIGRAPH || mon == MONITOR$AXIVITY || mon == MONITOR$VERISENSE)) { - if (mon == MONITOR$ACTIGRAPH || mon == MONITOR$VERISENSE) { - tmph = read.csv(datafile, nrow = 8, skip = 1) - tmphi = 1 - while (tmphi < 10) { - if (length(unlist(strsplit(format(tmph[tmphi,1]),"Start Time"))) > 1) { - break - } - tmphi = tmphi + 1 +g.getstarttime = function(datafile, data, mon, dformat, desiredtz, configtz = NULL) { + starttime = NULL + if (is.null(configtz)) { + configtz = desiredtz + } + if ("time" %in% colnames(data)) { + starttime = as.POSIXlt(data[1, "time"], tz = desiredtz, origin = "1970-01-01") + } else if (mon == MONITOR$MOVISENS) { + starttime = unisensR::readUnisensStartTime(dirname(datafile)) + + # movisens timestamp is stored in unisens.xml file as a string "%Y-%m-%dT%H:%M:%S", + # without a timezone. + # unisensR::readUnisensStartTime() converts this string into a POSIXct object, + # but since no timezone information was stored, this conversion happens using + # the system timezone. + # If configtz is different from the system timezone (""), we need to force the + # correct timezone (force, while keeping the same hh:mm:ss time, not just convert to congigtz) + if (configtz != "") { + starttime = lubridate::force_tz(starttime, configtz) } - starttime = unlist(strsplit(format(tmph[tmphi,1]),"Start Time"))[2] - #------------------------------- - tmphi = 1 - while (tmphi < 10) { - if (length(unlist(strsplit(format(tmph[tmphi,1]),"Start Date"))) > 1) { - break - } - tmphi = tmphi + 1 + starttime = as.POSIXlt(starttime, tz = desiredtz) + + } else if (dformat == FORMAT$CSV && (mon == MONITOR$ACTIGRAPH || mon == MONITOR$VERISENSE)) { + starttime = startdate = NULL + header_rows = 8 + tmph = read.csv(datafile, nrow = header_rows, skip = 1) + + tmphi = 1 + while (tmphi < header_rows) { + tmp = unlist(strsplit(format(tmph[tmphi,1]),"Start Time")) + if (length(tmp) > 1) { + starttime = tmp[2] + break } - startdate = unlist(strsplit(format(tmph[tmphi,1]), "Start Date"))[2] - startdate = as.character(unlist(strsplit(format(startdate)," "))) - starttime = as.character(unlist(strsplit(format(starttime)," "))) + tmphi = tmphi + 1 } - if (mon == MONITOR$AXIVITY) { - starttime = P[1,1] - starttime = as.POSIXlt(starttime,tz = desiredtz,origin = "1970-01-01") - startdate = unlist(strsplit(format(starttime), " "))[1] - } else { - #----------------------------------------- - #remove possible spaces in date or time - newstarttime = starttime #20-11-2014 - newstartdate = startdate #20-11-2014 - if (length(startdate) > 1) { - for (rpsi in 1:length(startdate)) { - if (length(unlist(strsplit(startdate[rpsi], ""))) > 1) { - newstartdate = startdate[rpsi] - } - } - } - if (length(starttime) > 1) { - for (rpsi in 1:length(starttime)) { - if (length(unlist(strsplit(starttime[rpsi], ""))) > 1) { - newstarttime = starttime[rpsi] - } - } + if (is.null(starttime)) { + stop(paste0("Start Time not found in the header of ", datafile)) + } + #------------------------------- + tmphi = 1 + while (tmphi < header_rows) { + tmp = unlist(strsplit(format(tmph[tmphi,1]),"Start Date")) + if (length(tmp) > 1) { + startdate = tmp[2] + break } - starttime = newstarttime - startdate = newstartdate + tmphi = tmphi + 1 } + if (is.null(startdate)) { + stop(paste0("Start Date not found in the header of ", datafile)) + } + + #----------------------------------------- + # trim any spaces + startdate = gsub(" ", "", startdate, fixed = TRUE) + starttime = gsub(" ", "", starttime, fixed = TRUE) + #----------------------------------------- # flexible four date/time formats starttime = paste0(startdate," ",starttime) getOption("digits.secs") options(digits.secs = 3) - if ((mon == MONITOR$ACTIGRAPH && dformat != FORMAT$GT3X) || mon == MONITOR$VERISENSE) { - options(warn = -1) - topline = as.matrix(colnames(as.matrix(read.csv(datafile, nrow = 1, skip = 0)))) - topline = topline[1] #To avoid dots - options(warn = 0) - # Extraction of date format. - # all formats to consider following R date formatting symbols: - # Y year including century, y year excluding century, b month as character, m month as a number - # Note: year in the middle of a date is assumed to not occur. - allformats = c("mdY", "mdy","bdY", "bdy", - "dmY", "dmy", "dbY", "dby", - "Ymd", "ymd", "Ybd", "ybd", - "Ydm", "ydm", "Ydb", "ydb") - fc = data.frame(matrix(0,1,length(allformats)), stringsAsFactors = TRUE) # data.frame to keep track of formats that have been checked - names(fc) = allformats - # Now systematically go through all possible formats to identify which one it is - fc$mdY = length(grep("MM[.]dd[.]yyyy|M[.]d[.]yyyy|M[.]dd[.]yyyy|MM[.]d[.]yyyy",topline)) - fc$mdy = length(grep("MM[.]dd[.]yy|M[.]d[.]yy|M[.]dd[.]yy|MM[.]d[.]yy",topline)) - fc$bdY = length(grep("MMM[.]dd[.]yyyy|MMM[.]d[.]yyyy",topline)) - fc$bdy = length(grep("MMM[.]dd[.]yy|MMM[.]d[.]yy",topline)) - if (fc$mdY == 1 & fc$mdy == 1) fc$mdy = 0 # not yy when yyyy is also detected - if (fc$bdY == 1 & fc$bdy == 1) fc$bdy = 0 # not yy when yyyy is also detected - if (fc$mdY == 1 & fc$bdY == 1) fc$mdY = 0 # not M(M) when MMM is also detected for yyyy - if (fc$mdy == 1 & fc$bdy == 1) fc$mdy = 0 # not M(M) when MMM is also detected for yy - fc$dmY = length(grep("d[.]M[.]yyyy|d[.]MM[.]yyyy",topline)) - fc$dmy = length(grep("d[.]M[.]yy|d[.]MM[.]yy",topline)) - fc$dbY = length(grep("d[.]MMM[.]yyyy",topline)) - fc$dby = length(grep("d[.]MMM[.]yy",topline)) - if (fc$dmY == 1 & fc$dmy == 1) fc$dmy = 0 # not yy when yyyy is also detected - if (fc$dbY == 1 & fc$dby == 1) fc$dby = 0 # not yy when yyyy is also detected - if (fc$dmY == 1 & fc$dbY == 1) fc$dmY = 0 # not M(M) when MMM is also detected for yyyy - if (fc$dmy == 1 & fc$dby == 1) fc$dmy = 0 # not M(M) when MMM is also detected for yy - fc$Ymd = length(grep("yyyy[.]M[.]d|yyyy[.]MM[.]d",topline)) - fc$ymd = length(grep("yy[.]M[.]d|yy[.]MM[.]d",topline)) - fc$Ybd = length(grep("yyyy[.]MMM[.]d",topline)) - fc$ybd = length(grep("yy[.]MMM[.]d",topline)) - if (fc$Ymd == 1 & fc$ymd == 1) fc$ymd = 0 # not yy when yyyy is also detected - if (fc$Ybd == 1 & fc$ybd == 1) fc$ybd = 0 # not yy when yyyy is also detected - if (fc$Ymd == 1 & fc$Ybd == 1) fc$Ymd = 0 # not M(M) when MMM is also detected for yyyy - if (fc$ymd == 1 & fc$ybd == 1) fc$ymd = 0 # not M(M) when MMM is also detected for yy - fc$Ydm = length(grep("yyyy[.]dd[.]MM|yyyy[.]d[.]M|yyyy[.]d[.]MM|yyyy[.]dd[.]M",topline)) - fc$ydm = length(grep("yy[.]dd[.]MM|yy[.]d[.]M|yy[.]d[.]MM|yy[.]dd[.]M",topline)) - fc$Ydb = length(grep("yyyy[.]dd[.]MMM|yyyy[.]d[.]MMM",topline)) - fc$ydb = length(grep("yy[.]dd[.]MMM|yy[.]d[.]MMM",topline)) - if (fc$Ydm == 1 & fc$ydm == 1) fc$ydm = 0 # not yy when yyyy is also detected - if (fc$Ydb == 1 & fc$ydb == 1) fc$ydb = 0 # not yy when yyyy is also detected - if (fc$Ydm == 1 & fc$Ydb == 1) fc$Ydm = 0 # not M(M) when MMM is also detected for yyyy - if (fc$ydm == 1 & fc$ydb == 1) fc$ydm = 0 # not M(M) when MMM is also detected for yy - # Now we can identify which format it is: - theformat = names(fc)[which(fc == 1)[1]] - if (length(theformat) == 0) warning("date format not recognised") - splitformat = unlist(strsplit(theformat,"")) - # identify separater: - if (length(grep("/",starttime)) > 0) { - sepa = "/" + + options(warn = -1) + topline = as.matrix(colnames(as.matrix(read.csv(datafile, nrow = 1, skip = 0)))) + topline = topline[1] #To avoid dots + options(warn = 0) + # Extraction of date format. + # all formats to consider following R date formatting symbols: + # Y year including century, y year excluding century, b month as character, m month as a number + # Note: year in the middle of a date is assumed to not occur. + allformats = c("mdY", "mdy","bdY", "bdy", + "dmY", "dmy", "dbY", "dby", + "Ymd", "ymd", "Ybd", "ybd", + "Ydm", "ydm", "Ydb", "ydb") + fc = data.frame(matrix(0,1,length(allformats))) # data.frame to keep track of formats that have been checked + names(fc) = allformats + # Now systematically go through all possible formats to identify which one it is + fc$mdY = length(grep("MM[.]dd[.]yyyy|M[.]d[.]yyyy|M[.]dd[.]yyyy|MM[.]d[.]yyyy",topline)) + fc$mdy = length(grep("MM[.]dd[.]yy|M[.]d[.]yy|M[.]dd[.]yy|MM[.]d[.]yy",topline)) + fc$bdY = length(grep("MMM[.]dd[.]yyyy|MMM[.]d[.]yyyy",topline)) + fc$bdy = length(grep("MMM[.]dd[.]yy|MMM[.]d[.]yy",topline)) + if (fc$mdY == 1 & fc$mdy == 1) fc$mdy = 0 # not yy when yyyy is also detected + if (fc$bdY == 1 & fc$bdy == 1) fc$bdy = 0 # not yy when yyyy is also detected + if (fc$mdY == 1 & fc$bdY == 1) fc$mdY = 0 # not M(M) when MMM is also detected for yyyy + if (fc$mdy == 1 & fc$bdy == 1) fc$mdy = 0 # not M(M) when MMM is also detected for yy + fc$dmY = length(grep("d[.]M[.]yyyy|d[.]MM[.]yyyy",topline)) + fc$dmy = length(grep("d[.]M[.]yy|d[.]MM[.]yy",topline)) + fc$dbY = length(grep("d[.]MMM[.]yyyy",topline)) + fc$dby = length(grep("d[.]MMM[.]yy",topline)) + if (fc$dmY == 1 & fc$dmy == 1) fc$dmy = 0 # not yy when yyyy is also detected + if (fc$dbY == 1 & fc$dby == 1) fc$dby = 0 # not yy when yyyy is also detected + if (fc$dmY == 1 & fc$dbY == 1) fc$dmY = 0 # not M(M) when MMM is also detected for yyyy + if (fc$dmy == 1 & fc$dby == 1) fc$dmy = 0 # not M(M) when MMM is also detected for yy + fc$Ymd = length(grep("yyyy[.]M[.]d|yyyy[.]MM[.]d",topline)) + fc$ymd = length(grep("yy[.]M[.]d|yy[.]MM[.]d",topline)) + fc$Ybd = length(grep("yyyy[.]MMM[.]d",topline)) + fc$ybd = length(grep("yy[.]MMM[.]d",topline)) + if (fc$Ymd == 1 & fc$ymd == 1) fc$ymd = 0 # not yy when yyyy is also detected + if (fc$Ybd == 1 & fc$ybd == 1) fc$ybd = 0 # not yy when yyyy is also detected + if (fc$Ymd == 1 & fc$Ybd == 1) fc$Ymd = 0 # not M(M) when MMM is also detected for yyyy + if (fc$ymd == 1 & fc$ybd == 1) fc$ymd = 0 # not M(M) when MMM is also detected for yy + fc$Ydm = length(grep("yyyy[.]dd[.]MM|yyyy[.]d[.]M|yyyy[.]d[.]MM|yyyy[.]dd[.]M",topline)) + fc$ydm = length(grep("yy[.]dd[.]MM|yy[.]d[.]M|yy[.]d[.]MM|yy[.]dd[.]M",topline)) + fc$Ydb = length(grep("yyyy[.]dd[.]MMM|yyyy[.]d[.]MMM",topline)) + fc$ydb = length(grep("yy[.]dd[.]MMM|yy[.]d[.]MMM",topline)) + if (fc$Ydm == 1 & fc$ydm == 1) fc$ydm = 0 # not yy when yyyy is also detected + if (fc$Ydb == 1 & fc$ydb == 1) fc$ydb = 0 # not yy when yyyy is also detected + if (fc$Ydm == 1 & fc$Ydb == 1) fc$Ydm = 0 # not M(M) when MMM is also detected for yyyy + if (fc$ydm == 1 & fc$ydb == 1) fc$ydm = 0 # not M(M) when MMM is also detected for yy + # Now we can identify which format it is: + theformat = names(fc)[which(fc == 1)[1]] + if (is.na(theformat)) warning("date format not recognised") + splitformat = unlist(strsplit(theformat,"")) + # identify separater: + if (length(grep("/",starttime)) > 0) { + sepa = "/" + } else { + if (length(grep("-",starttime)) > 0) { + sepa = "-" } else { - if (length(grep("-",starttime)) > 0) { - sepa = "-" + if (length(grep(".",starttime)) > 0) { + sepa = "." } else { - if (length(grep(".",starttime)) > 0) { - sepa = "." - } else { - warning("separater character for dates not identified") - } + warning("separator character for dates not identified") } } - expectedformat = paste0('%',splitformat[1],sepa,'%',splitformat[2],sepa,'%',splitformat[3],' %H:%M:%S') - Sys.setlocale("LC_TIME", "C") # set language to English because that is what we use elsewhere in GGIR - starttime = as.POSIXlt(starttime, format = expectedformat) - } - } else if (dformat == FORMAT$AD_HOC_CSV && mon == MONITOR$AD_HOC) { - starttime = P$data$timestamp[1] - } else if (mon == MONITOR$MOVISENS) { - starttime = unisensR::readUnisensStartTime(dirname(datafile)) - } else if (dformat == FORMAT$GT3X && mon == MONITOR$ACTIGRAPH) { - if (is.null(configtz)) configtz = desiredtz - starttime = as.POSIXct(format(P[1, 1]), tz = configtz) - if (configtz != desiredtz) { - attr(starttime, "tzone") <- desiredtz } + expectedformat = paste0('%',splitformat[1],sepa,'%',splitformat[2],sepa,'%',splitformat[3],' %H:%M:%S') + Sys.setlocale("LC_TIME", "C") # set language to English because that is what we use elsewhere in GGIR + starttime = as.POSIXct(starttime, format = expectedformat, tz = configtz) starttime = as.POSIXlt(starttime, tz = desiredtz) + } else { + stop(paste0("Timestamps not found for monitor type ", mon, " and file format type ", dformat, + "\nThis should not happen.")) } + return(starttime) } diff --git a/R/g.imputeTimegaps.R b/R/g.imputeTimegaps.R index 17627b25c..11e8d6568 100644 --- a/R/g.imputeTimegaps.R +++ b/R/g.imputeTimegaps.R @@ -1,4 +1,4 @@ -g.imputeTimegaps = function(x, xyzCol, timeCol = c(), sf, k=0.25, impute = TRUE, +g.imputeTimegaps = function(x, sf, k=0.25, impute = TRUE, PreviousLastValue = c(0, 0, 1), PreviousLastTime = NULL, epochsize = NULL) { if (!is.null(epochsize)) { @@ -6,55 +6,48 @@ g.imputeTimegaps = function(x, xyzCol, timeCol = c(), sf, k=0.25, impute = TRUE, longEpochSize = epochsize[2] } # dummy variables to control the process - remove_time_at_end = dummyTime = FirstRowZeros = imputelast = FALSE + remove_time_at_end = FirstRowZeros = imputelast = FALSE # initialize numberofgaps and GapsLength - NumberOfGaps = GapsLength = NULL + NumberOfGaps = GapsLength = 0 + # add temporary timecolumn to enable timegap imputation where there are zeros - if (length(timeCol) == 1) { - if (!(timeCol %in% colnames(x))) dummyTime = TRUE - } - if (length(timeCol) == 0 | dummyTime == TRUE) { - dummytime = Sys.time() - adhoc_time = seq(dummytime, dummytime + (nrow(x) - 1) * (1/sf), by = 1/sf) - if (length(adhoc_time) < nrow(x)) { - NotEnough = nrow(x) - length(adhoc_time) - adhoc_time = seq(dummytime, dummytime + (nrow(x) + NotEnough) * (1/sf), by = 1/sf) - } - x$time = adhoc_time[1:nrow(x)] - timeCol = "time" + if (!("time" %in% colnames(x))) { + x$time = seq(from = Sys.time(), by = 1/sf, length.out = nrow(x)) remove_time_at_end = TRUE } - # define function to imputation at raw level + xyzCol = which(colnames(x) %in% c("x", "y", "z")) + + # define function for imputation at raw level imputeRaw = function(x, sf) { - # impute raw timestamps because timestamps still need to be meaningul when + # impute raw timestamps because timestamps still need to be meaningful when # resampling or when plotting gapp = which(x$gap != 1) if (length(gapp) > 0) { if (gapp[1] > 1) { - newTime = x$timestamp[1:(gapp[1] - 1)] + newTime = x$time[1:(gapp[1] - 1)] } else { newTime = NULL } for (g in 1:length(gapp)) { - newTime = c(newTime, x$timestamp[gapp[g]] + seq(0, (x$gap[gapp[g]] - 1) / sf, by = 1/sf)) + newTime = c(newTime, x$time[gapp[g]] + seq(0, by = 1/sf, length.out = x$gap[gapp[g]])) if (g < length(gapp)) { - newTime = c(newTime, x$timestamp[(gapp[g] + 1):(gapp[g + 1] - 1)]) + newTime = c(newTime, x$time[(gapp[g] + 1):(gapp[g + 1] - 1)]) } } - newTime = c(newTime, x$time[(gapp[g] + 1):length(x$timestamp)]) + newTime = c(newTime, x$time[(gapp[g] + 1):length(x$time)]) } x <- as.data.frame(lapply(x, rep, x$gap)) if (length(gapp) > 0) { - x$timestamp = newTime[1:nrow(x)] + x$time = newTime[1:nrow(x)] } x = x[, which(colnames(x) != "gap")] return(x) } # find zeros and remove them from dataset - zeros = which(x[,xyzCol[1]] == 0 & x[,xyzCol[2]] == 0 & x[,xyzCol[3]] == 0) + zeros = which(x$x == 0 & x$y == 0 & x$z == 0) if (length(zeros) > 0) { # if first value is a zero, remember value from previous chunk to replicate # if chunk = 1, then it will use c(0, 0, 1) @@ -79,32 +72,32 @@ g.imputeTimegaps = function(x, xyzCol, timeCol = c(), sf, k=0.25, impute = TRUE, if (k < 2/sf) { # prevent trying to impute timegaps shorter than 2 samples k = 2/sf } - deltatime = diff(x[, timeCol]) + deltatime = diff(x$time) if (!is.numeric(deltatime)) { # in csv axivity, the time is directly read as numeric (seconds) units(deltatime) = "secs" deltatime = as.numeric(deltatime) } # refill if first value is not consecutive from last value in previous chunk if (!is.null(PreviousLastTime)) { - first_deltatime = diff(c(PreviousLastTime, x[1, timeCol])) + first_deltatime = diff(c(PreviousLastTime, x$time[1])) if (!is.numeric(first_deltatime)) { # in csv axivity, the time is directly read as numeric (seconds) units(first_deltatime) = "secs" first_deltatime = as.numeric(first_deltatime) } - if (first_deltatime >= k) { # prevent trying to impute timegaps shorter than 2 samples + if (first_deltatime >= k) { # don't impute a timegap shorter than the minimum requested x = rbind(x[1,], x) - x[1, timeCol] = PreviousLastTime + x$time[1] = PreviousLastTime x[1, xyzCol] = PreviousLastValue deltatime = c(first_deltatime, deltatime) } } # impute time gaps - gapsi = which(deltatime >= k) # limit imputation to gaps larger than 0.25 seconds + gapsi = which(deltatime >= k) # limit imputation to gaps longer than the minimum requested NumberOfGaps = length(gapsi) if (NumberOfGaps > 0) { x$gap = 1 x$gap[gapsi] = round(deltatime[gapsi] * sf) # as.integer was problematic many decimals close to wholenumbers (but not whole numbers) resulting in 1 row less than expected - GapsLength = sum(x$gap[gapsi]) - NumberOfGaps # - numberOfGaps because x$gap == 1 means no gap + GapsLength = sum(x$gap[gapsi]) # normalisation to 1 G normalise = which(x$gap > 1) for (i_normalise in normalise) { @@ -178,28 +171,25 @@ g.imputeTimegaps = function(x, xyzCol, timeCol = c(), sf, k=0.25, impute = TRUE, } } } else if (impute == FALSE) { - if (FirstRowZeros == TRUE) x = x[-1,] # since zeros[1] was removed in line 21 - if (imputelast == TRUE) x = x[-nrow(x),] # since zeros[length(zeros)] was removed in line 27 + if (FirstRowZeros == TRUE) x = x[-1,] # since zeros[1] was removed + if (imputelast == TRUE) x = x[-nrow(x),] # since zeros[length(zeros)] was removed } # impute last value? if (imputelast) x[nrow(x), xyzCol] = x[nrow(x) - 1, xyzCol] - if (remove_time_at_end == TRUE) { - x = x[, grep(pattern = "time", x = colnames(x), invert = TRUE)] - } - # keep only timestamp column - if (all(c("time", "timestamp") %in% colnames(x))) { - x = x[, grep(pattern = "timestamp", x = colnames(x), invert = TRUE)] - } + # QClog - start = as.numeric(as.POSIXct(x[1,1], origin = "1970-1-1")) + start = as.numeric(as.POSIXct(x$time[1], origin = "1970-1-1")) end = start + nrow(x) - if (is.null(GapsLength)) GapsLength = 0 - if (is.null(NumberOfGaps)) NumberOfGaps = 0 imputed = NumberOfGaps > 0 QClog = data.frame(imputed = imputed, start = start, end = end, blockLengthSeconds = (end - start) / sf, timegaps_n = NumberOfGaps, timegaps_min = GapsLength/sf/60) + + if (remove_time_at_end == TRUE) { + x = x[, grep(pattern = "time", x = colnames(x), invert = TRUE)] + } + # return data and QClog return(list(x = x, QClog = QClog)) } \ No newline at end of file diff --git a/R/g.inspectfile.R b/R/g.inspectfile.R index aa281ee13..17dcb9b4b 100644 --- a/R/g.inspectfile.R +++ b/R/g.inspectfile.R @@ -140,7 +140,7 @@ g.inspectfile = function(datafile, desiredtz = "", params_rawdata = c(), } else if (mon == MONITOR$AXIVITY) { # sample frequency is not stored tmp = read.csv(datafile, nrow = 100000, skip = 0) - tmp = as.numeric(as.POSIXlt(tmp[, 1])) + tmp = as.numeric(as.POSIXct(tmp[, 1], origin = "1970-01-01")) sf = length(tmp) / (tmp[length(tmp)] - tmp[1]) sf = floor((sf) / 5 ) * 5 # round down to nearest integer of 5, we never want to assume that there is more frequency content in a signal than there truly is } @@ -213,11 +213,21 @@ g.inspectfile = function(datafile, desiredtz = "", params_rawdata = c(), sf = params_rawdata[["rmc.sf"]] if (is.null(sf)) { stop("\nFile header doesn't specify sample rate. Please provide rmc.sf value to process ", datafile) + } else if (sf == 0) { + stop("\nFile header doesn't specify sample rate. Please provide a non-zero rmc.sf value to process ", datafile) } } else { sf = Pusercsvformat$header$sample_rate } } + + if (mon == MONITOR$GENEACTIV && dformat == FORMAT$CSV) { + stop(paste0("The GENEActiv csv reading functionality is deprecated in", + " GGIR from version 2.6-4 onwards. Please, use either", + " the GENEActiv bin files or try to read the csv files with", + " GGIR::read.myacc.csv"), call. = FALSE) + } + if (dformat == FORMAT$BIN) { if (mon == MONITOR$GENEACTIV) { H = GGIRread::readGENEActiv(filename = datafile, start = 0, end = 1)$header @@ -244,10 +254,7 @@ g.inspectfile = function(datafile, desiredtz = "", params_rawdata = c(), } } } else if (dformat == FORMAT$CSV) { - if (mon == MONITOR$GENEACTIV) { - H = read.csv(datafile,nrow = 20, skip = 0) #note that not the entire header is copied - # cat("\nGENEACTIV csv files support is deprecated in GGIR v2.6-2 onwards. Please, either use the GENEACTIV bin files or the read.myacc.csv function on the csv files") - } else if (mon == MONITOR$ACTIGRAPH) { + if (mon == MONITOR$ACTIGRAPH) { H = read.csv(datafile, nrow = 9, skip = 0) } else if (mon == MONITOR$AXIVITY) { H = "file does not have header" # these files have no header @@ -263,6 +270,9 @@ g.inspectfile = function(datafile, desiredtz = "", params_rawdata = c(), H = data.frame(name = row.names(header), value = header, stringsAsFactors = TRUE) } sf = params_rawdata[["rmc.sf"]] + if (sf == 0) { + stop("\nPlease provide a non-zero rmc.sf value to process ", datafile) + } } else if (dformat == FORMAT$GT3X) { # gt3x info = try(expr = {read.gt3x::parse_gt3x_info(datafile, tz = desiredtz)},silent = TRUE) if (inherits(info, "try-error") == TRUE || is.null(info)) { @@ -285,6 +295,11 @@ g.inspectfile = function(datafile, desiredtz = "", params_rawdata = c(), sf = as.numeric(H[which(H[,1] == "Sample Rate"), 2]) } } + + if (sf == 0) { + stop(paste0("\nSample frequency not recognised in ", datafile), call. = FALSE) + } + if (is.null(sf) == FALSE) { H = as.matrix(H) if (ncol(H) == 3 && dformat == FORMAT$CSV && mon == MONITOR$ACTIGRAPH) { @@ -336,10 +351,18 @@ g.inspectfile = function(datafile, desiredtz = "", params_rawdata = c(), } } } + + # detect dot or comma separator in the data file + op <- options(warn = -1) # turn off warnings + on.exit(options(op)) + suppressWarnings(expr = {decn = g.dotorcomma(datafile, dformat, mon, + rmc.dec = params_rawdata[["rmc.dec"]])}) + options(warn = 0) # turn on warnings + monc = mon monn = ifelse(mon > 0, monnames[mon], "unknown") dformc = dformat dformn = fornames[dformat] invisible(list(header = header, monc = monc, monn = monn, - dformc = dformc, dformn = dformn, sf = sf, filename = filename)) + dformc = dformc, dformn = dformn, sf = sf, decn = decn, filename = filename)) } diff --git a/R/g.part1.R b/R/g.part1.R index afcf661b5..7f662e611 100644 --- a/R/g.part1.R +++ b/R/g.part1.R @@ -27,11 +27,23 @@ g.part1 = function(datadir = c(), metadatadir = c(), f0 = 1, f1 = c(), myfun = c filesizes = file.size(fnamesfull) # in bytes bigEnough = which(filesizes/1e6 > params_rawdata[["minimumFileSizeMB"]]) fnamesfull = fnamesfull[bigEnough] + + if (length(bigEnough) > 0) { + fnames_toosmall = fnames[-bigEnough] + } else { + fnames_toosmall = fnames + } + fnames = fnames[bigEnough] + if (verbose && length(fnames_toosmall) > 0) { + warning(paste0("\nSkipping files that are too small for analysis: ", toString(fnames_toosmall), + " (configurable with parameter minimumFileSizeMB)."), call. = FALSE) + } + if (length(fnamesfull) == 0) { stop(paste0("\nNo files to analyse. Check that there are accelerometer files ", - "in the directory specified with argument datadir")) + "in the directory specified with argument datadir"), call. = FALSE) } # create output directory if it does not exist @@ -208,6 +220,7 @@ g.part1 = function(datadir = c(), metadatadir = c(), f0 = 1, f1 = c(), myfun = c bcc.scalei = which(colnames(bcc.data) == "scale.x" | colnames(bcc.data) == "scale.y" | colnames(bcc.data) == "scale.z") bcc.offseti = which(colnames(bcc.data) == "offset.x" | colnames(bcc.data) == "offset.y" | colnames(bcc.data) == "offset.z") bcc.temp.offseti = which(colnames(bcc.data) == "temperature.offset.x" | colnames(bcc.data) == "temperature.offset.y" | colnames(bcc.data) == "temperature.offset.z") + bcc.meantempcali = which(colnames(bcc.data) == "meantempcal") bcc.QCmessagei = which(colnames(bcc.data) == "QCmessage") bcc.npointsi = which(colnames(bcc.data) == "n.10sec.windows") bcc.nhoursusedi = which(colnames(bcc.data) == "n.hours.considered") @@ -215,12 +228,13 @@ g.part1 = function(datadir = c(), metadatadir = c(), f0 = 1, f1 = c(), myfun = c C$scale = as.numeric(bcc.data[bcc.i[1],bcc.scalei]) C$offset = as.numeric(bcc.data[bcc.i[1],bcc.offseti]) C$tempoffset = as.numeric(bcc.data[bcc.i[1],bcc.temp.offseti]) + C$meantempcal = bcc.data[bcc.i[1], bcc.meantempcali] C$cal.error.start = as.numeric(bcc.data[bcc.i[1],bcc.cal.error.start]) C$cal.error.end = as.numeric(bcc.data[bcc.i[1],bcc.cal.error.end]) C$QCmessage = bcc.data[bcc.i[1],bcc.QCmessagei] C$npoints = bcc.data[bcc.i[1], bcc.npointsi] C$nhoursused = bcc.data[bcc.i[1], bcc.nhoursusedi] - C$use.temp = bcc.data[bcc.i[1], bcc.nhoursusedi] + C$use.temp = bcc.data[bcc.i[1], bcc.use.tempi] if (verbose == TRUE) { cat(paste0("\nRetrieved Calibration error (g) before: ",as.numeric(bcc.data[bcc.i[1],bcc.cal.error.start]))) cat(paste0("\nRetrieved Callibration error (g) after: ",as.numeric(bcc.data[bcc.i[1],bcc.cal.error.end]))) @@ -384,7 +398,7 @@ g.part1 = function(datadir = c(), metadatadir = c(), f0 = 1, f1 = c(), myfun = c "g.getstarttime", "POSIXtime2iso8601", "iso8601chartime2POSIX", "datadir2fnames", "read.myacc.csv", "get_nw_clip_block_params", - "get_starttime_weekday_meantemp_truncdata", "ismovisens", + "get_starttime_weekday_truncdata", "ismovisens", "g.extractheadervars", "g.imputeTimegaps", "extract_params", "load_params", "check_params", "detect_nonwear_clipping") diff --git a/R/g.readaccfile.R b/R/g.readaccfile.R index 5e642df13..214ae2c8d 100644 --- a/R/g.readaccfile.R +++ b/R/g.readaccfile.R @@ -1,7 +1,7 @@ g.readaccfile = function(filename, blocksize, blocknumber, filequality, ws, PreviousEndPage = 1, inspectfileobject = c(), PreviousLastValue = c(0, 0, 1), PreviousLastTime = NULL, - params_rawdata = c(), params_general = c(), ...) { + params_rawdata = c(), params_general = c(), header = NULL, ...) { #get input variables input = list(...) if (length(input) > 0 || @@ -16,90 +16,73 @@ g.readaccfile = function(filename, blocksize, blocknumber, filequality, params_rawdata = params$params_rawdata params_general = params$params_general } - # function wrapper to read blocks of accelerationd data from various brands - # the code identifies which accelerometer brand and data format it is - # blocksize = number of pages to read at once - # blocknumber = block count relative to beginning of measurement - # sf = sample frequency (Hertz) - # ws = large window size (default 3600 seconds) - switchoffLD = 0 + desiredtz = params_general[["desiredtz"]] + configtz = params_general[["configtz"]] + if (length(configtz) == 0) configtz = desiredtz + I = inspectfileobject mon = I$monc if (mon == MONITOR$VERISENSE) mon = MONITOR$ACTIGRAPH dformat = I$dformc sf = I$sf - - P = c() - updatepageindexing = function(startpage = c(), deltapage = c(), blocknumber = c(), PreviousEndPage = c(), - mon = c(), dformat = c()) { - # This function ensures that startpage is only specified for blocknumber 1. - # The next time (blocknumber > 1) the startpage will be derived from the previous - # endpage and the blocksize. - if (blocknumber != 1 & length(PreviousEndPage) != 0) { - # if ((mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) || dformat == FORMAT$CSV) { # change this line as the csv data do not need to skip one more row (the skip argument in read.csv does not include this row of the dataset) - if ((mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) | dformat == FORMAT$GT3X) { - # only in GENEActiv binary data and for gt3x format data - # page selection is defined from start to end (including end) - startpage = PreviousEndPage + 1 - } else { - # for other monitor brands and data formats - # page selection is defined from start to end (excluding end itself) - # so start page of one block equals the end page of previous block - startpage = PreviousEndPage - } - } - endpage = startpage + deltapage - return(list(startpage = startpage, endpage = endpage)) + decn = I$decn + + if ((mon == MONITOR$ACTIGRAPH && dformat == FORMAT$CSV) || + (mon == MONITOR$AXIVITY && dformat == FORMAT$CSV) || + dformat == FORMAT$AD_HOC_CSV) { + blocksize = blocksize * 300 } - if (mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) { - startpage = blocksize * (blocknumber - 1) + 1 # GENEActiv starts with page 1 - deltapage = blocksize - UPI = updatepageindexing(startpage = startpage, deltapage = deltapage, - blocknumber = blocknumber, PreviousEndPage = PreviousEndPage, mon = mon, dformat = dformat) - startpage = UPI$startpage; endpage = UPI$endpage - - try(expr = {P = GGIRread::readGENEActiv(filename = filename, start = startpage, - end = endpage, desiredtz = params_general[["desiredtz"]], - configtz = params_general[["configtz"]])}, silent = TRUE) - - if (length(P) > 0) { - if (nrow(P$data.out) < (blocksize*300)) { - switchoffLD = 1 #last block - } + + if (blocknumber < 1) blocknumber = 1 + + # startpage should only be specified for blocknumber 1. + # The next time (blocknumber > 1) the startpage will be derived from the previous + # endpage and the blocksize. + + if ((mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) || dformat == FORMAT$GT3X || + (mon == MONITOR$MOVISENS && dformat == FORMAT$BIN)) { + # for GENEActiv binary data, gt3x format data, and Movisens data, + # page selection is defined from start to end **including end** + if (blocknumber > 1 && length(PreviousEndPage) != 0) { + startpage = PreviousEndPage + 1 + } else { + startpage = blocksize * (blocknumber - 1) + 1 # pages are numbered starting with page 1 } - if (length(P) == 0) { #if first block doens't read then probably corrupt - if (blocknumber == 1) { - #try to read without specifying blocks (file too short) - try(expr = { - P = GGIRread::readGENEActiv(filename = filename, desiredtz = params_general[["desiredtz"]], - configtz = params_general[["configtz"]]) - }, silent = TRUE) - if (length(P) == 0) { - warning('\nFile possibly corrupt\n') - P = c(); switchoffLD = 1 - filequality$filecorrupt = TRUE - } #if not then P is now filled with data - } else { - P = c() #just no data in this last block + endpage = startpage + blocksize - 1 # -1 because both startpage and endpage will be read, + # and we want to read blocksize # of samples + } else { + # for other monitor brands and data formats + # page selection is defined from start to end **excluding end itself**, + # so start page of one block equals the end page of previous block + if (blocknumber > 1 && length(PreviousEndPage) != 0) { + startpage = PreviousEndPage + } else { + startpage = blocksize * (blocknumber - 1) # pages are numbered starting with page 0 + + if (mon == MONITOR$ACTIGRAPH && dformat == FORMAT$CSV) { + headerlength = 10 + startpage = startpage + headerlength } } - if (length(P) > 0) { #check whether there is enough data - if (nrow(P$data.out) < ((sf * ws * 2) + 1) & blocknumber == 1) { - P = c(); switchoffLD = 1 - filequality$filetooshort = TRUE + endpage = startpage + blocksize + } + + P = c() + isLastBlock = FALSE + + if (mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) { + try(expr = {P = GGIRread::readGENEActiv(filename = filename, start = startpage, + end = endpage, desiredtz = desiredtz, + configtz = configtz)}, silent = TRUE) + if (length(P) > 0 && ("data.out" %in% names(P))) { + names(P)[names(P) == "data.out"] = "data" + + if (nrow(P$data) < (blocksize*300)) { + isLastBlock = TRUE } } - #=============== - } else if (mon == MONITOR$ACTIGRAPH && dformat == FORMAT$CSV) { - headerlength = 10 - #-------------- - startpage = (headerlength + (blocksize * 300 * (blocknumber - 1))) - deltapage = blocksize * 300 - UPI = updatepageindexing(startpage = startpage, deltapage = deltapage, - blocknumber = blocknumber, PreviousEndPage = PreviousEndPage, mon = mon, dformat = dformat) - startpage = UPI$startpage; endpage = UPI$endpage - + } else if (mon == MONITOR$ACTIGRAPH && dformat == FORMAT$CSV) { # load rows 11:13 to investigate whether the file has a header # invisible because R complains about poor Actigraph file format, # this is an an ActiGraph problem not a GGIR problem, so we ignore it @@ -110,102 +93,63 @@ g.readaccfile = function(filename, blocksize, blocknumber, filequality, invisible(force(x)) } - op <- options(stringsAsFactors = FALSE) - on.exit(options(op)) - options(warn = -1) # turn off warnings - suppressWarnings(expr = {decn = g.dotorcomma(filename, dformat, mon, - desiredtz = params_general[["desiredtz"]], - rmc.dec = params_rawdata[["rmc.dec"]], - loadGENEActiv = params_rawdata[["loadGENEActiv"]])}) # detect dot or comma dataformat - options(warn = 0) #turn on warnings - - testheader = quiet(as.data.frame(data.table::fread(filename, nrows = 2, skip = 10, - dec = decn, showProgress = FALSE, - header = TRUE), - stringsAsFactors = FALSE)) - if (suppressWarnings(is.na(as.numeric(colnames(testheader)[1]))) == FALSE) { # it has no header, first value is a number - freadheader = FALSE - } else { # it has a header, first value is a character - freadheader = TRUE - headerlength = 11 - # skip 1 more row only in the case the file has a header, only in the first chunk of data (when the header needs to be skipped) - if (blocknumber == 1) { + # skip 1 more row only if the file has a header. Only the first chunk of data can have a header. + if (blocknumber == 1) { + testheader = quiet(data.table::fread(filename, nrows = 2, skip = 10, + dec = decn, showProgress = FALSE, + header = TRUE, data.table=FALSE, stringsAsFactors=FALSE)) + if (suppressWarnings(is.na(as.numeric(colnames(testheader)[1])))) { # first value is *not* a number, so file starts with a header startpage = startpage + 1 endpage = endpage + 1 } } - #-------------- try(expr = { - P = quiet(as.data.frame( - data.table::fread(filename, nrows = deltapage, - skip = startpage, - dec = decn, showProgress = FALSE, header = FALSE), # header should always be FALSE to prevent that acceleration values are taken as header when reading chunks 2 onwards - stringsAsFactors = TRUE)) + P$data = quiet(data.table::fread(filename, nrows = blocksize, skip = startpage, + dec = decn, showProgress = FALSE, + header = FALSE, # header should always be FALSE to prevent acceleration values from being mistaken for header when reading chunks 2 and on + data.table=FALSE, stringsAsFactors=FALSE)) }, silent = TRUE) - if (length(P) > 1) { - # data.matrix turnes num to char if there are missing values. - if (ncol(P) == 3) { - P = data.matrix(P) + if (length(P$data) > 0) { + if (ncol(P$data) < 3) { + P$data = c() } else { - P = data.matrix(P[, 2:ncol(P)]) # avoid timestamp column - } - if (nrow(P) < ((sf * ws * 2) + 1) & blocknumber == 1) { - P = c() ; switchoffLD = 1 - filequality$filetooshort = TRUE + if (ncol(P$data) > 3) { + P$data = P$data[, 2:4] # remove timestamp column, keep only XYZ columns + } + colnames(P$data) = c("x", "y", "z") } - } else { - P = c() } } else if (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) { if (utils::packageVersion("GGIRread") < "0.3.0") { # ignore frequency_tol parameter - apply_readAxivity = function(fname = filename, - bstart, bend, - desiredtz = params_general[["desiredtz"]], - configtz = params_general[["configtz"]], - interpolationType = params_rawdata[["interpolationType"]], - frequency_tol = params_rawdata[["frequency_tol"]]) { - try(expr = {P = GGIRread::readAxivity(filename = fname, start = bstart, - end = bend, progressBar = FALSE, desiredtz = desiredtz, + apply_readAxivity = function(bstart, bend) { + try(expr = {P = GGIRread::readAxivity(filename = filename, start = bstart, end = bend, + progressBar = FALSE, + desiredtz = desiredtz, configtz = configtz, - interpolationType = interpolationType)}, silent = TRUE) + interpolationType = params_rawdata[["interpolationType"]], + header = header) + }, silent = TRUE) return(P) } } else { # pass on frequency_tol parameter to GGIRread::readAxivity function - apply_readAxivity = function(fname = filename, - bstart, bend, - desiredtz = params_general[["desiredtz"]], - configtz = params_general[["configtz"]], - interpolationType = params_rawdata[["interpolationType"]], - frequency_tol = params_rawdata[["frequency_tol"]]) { - try(expr = {P = GGIRread::readAxivity(filename = filename, start = bstart, - end = bend, progressBar = FALSE, desiredtz = desiredtz, + apply_readAxivity = function(bstart, bend) { + try(expr = {P = GGIRread::readAxivity(filename = filename, start = bstart, end = bend, + progressBar = FALSE, + desiredtz = desiredtz, configtz = configtz, - interpolationType = interpolationType, - frequency_tol = frequency_tol)}, silent = TRUE) + interpolationType = params_rawdata[["interpolationType"]], + frequency_tol = params_rawdata[["frequency_tol"]], + header = header) + }, silent = TRUE) return(P) } } - - startpage = blocksize * (blocknumber - 1) - deltapage = blocksize - UPI = updatepageindexing(startpage = startpage, deltapage = deltapage, - blocknumber = blocknumber, PreviousEndPage = PreviousEndPage, mon = mon, dformat = dformat) - startpage = UPI$startpage; endpage = UPI$endpage + P = apply_readAxivity(bstart = startpage, bend = endpage) - if (length(P) > 1) { # data reading succesful - if (length(P$data) == 0) { # too short? - P = c() ; switchoffLD = 1 - if (blocknumber == 1) filequality$filetooshort = TRUE - } else { - if (nrow(P$data) < ((sf * ws * 2) + 1)) { - P = c() ; switchoffLD = 1 - if (blocknumber == 1) filequality$filetooshort = TRUE - } - } - } else { + if (length(P) == 0) { # If data reading is not successful then try following steps to retrieve issue # I am not sure if this is still relevant after all the improvements to GGIRread # but leaving this in just in case it is still needed @@ -229,125 +173,134 @@ g.readaccfile = function(filename, blocksize, blocknumber, filequality, # read the entire block: P = apply_readAxivity(bstart = startpage, bend = endpage) if (length(P) > 1) { # data reading succesful - if (length(P$data) == 0) { # if this still does not work then - P = c() ; switchoffLD = 1 - if (blocknumber == 1) filequality$filetooshort = TRUE - } else { - if (nrow(P$data) < ((sf * ws * 2) + 1)) { - P = c() ; switchoffLD = 1 - if (blocknumber == 1) filequality$filetooshort = TRUE - } else { - filequality$NFilePagesSkipped = NFilePagesSkipped # store number of pages jumped - } - } + filequality$NFilePagesSkipped = NFilePagesSkipped # store number of pages jumped + # Add replications of Ptest to the beginning of P to achieve same data # length as under normal conditions P$data = rbind(do.call("rbind", replicate(NFilePagesSkipped, PtestStartPage$data, simplify = FALSE)), P$data) - } else { # Data reading still not succesful, so classify file as corrupt - P = c() - if (blocknumber == 1) filequality$filecorrupt = TRUE } - } else { - P = c() - if (blocknumber == 1) filequality$filecorrupt = TRUE } } + if ("temp" %in% colnames(P$data)) { + colnames(P$data)[colnames(P$data) == "temp"] = "temperature" + } } else if (mon == MONITOR$AXIVITY && dformat == FORMAT$CSV) { - freadheader = FALSE - headerlength = 0 - startpage = (headerlength + (blocksize * 300 * (blocknumber - 1))) - deltapage = (blocksize*300) - UPI = updatepageindexing(startpage = startpage, deltapage = deltapage, - blocknumber = blocknumber, - PreviousEndPage = PreviousEndPage, mon = mon, dformat = dformat) - startpage = UPI$startpage; endpage = UPI$endpage try(expr = { - P = as.data.frame( - data.table::fread(filename, nrows = deltapage, - skip = startpage, - dec = decn, showProgress = FALSE, header = freadheader), - stringsAsFactors = TRUE) + rawData = data.table::fread(filename, nrows = blocksize, + skip = startpage, + dec = decn, showProgress = FALSE, header = FALSE, + data.table=FALSE, stringsAsFactors=FALSE) }, silent = TRUE) - if (length(P) > 1) { - if (nrow(P) < ((sf * ws * 2) + 1) & blocknumber == 1) { - P = c() ; switchoffLD = 1 #added 30-6-2012 - filequality$filetooshort = TRUE + if (length(rawData) > 0) { + if (nrow(rawData) < blocksize) { + isLastBlock = TRUE } - if (nrow(P) < (deltapage)) { #last block - switchoffLD = 1 + + rawTime = rawData[,1] + + if(class(rawTime)[1] == "character") { + # If timestamps in the csv file are formatted (Y-M-D h:m:s.f), but some of the dates are formatted poorly, + # data.table::fread() might have trouble parsing them as timestamps, and will instead return them as strings. + # as.POSIXct() is less brittle and might still be able to parse such timestamps, but it will be a very slow process. + # For instance, one omGUI export contained timestamps like "2023-11-11 15:22:60.000" (which should normally be "2023-11-11 15:23:00.000"). + # data.table::fread() couldn't parse this but as.POSIXct() could, albeit very slowly. + + rawTime = as.POSIXct(rawTime, tz=configtz, origin = "1970-01-01") + + # if as.POSIXct() also failed to parse this data as timestamps, there's nothing else we can do. + if(class(rawTime)[1] != "POSIXct") { + stop(paste0("Corrupt timestamp data in ", filename), call. = FALSE) + } else { + warning(paste0("Corrupt timestamp data in ", filename, + ". This will greatly slow down processing. To avoid this, use the original .cwa file, ", + "or export your data with Unix timestamps instead."), call. = FALSE) + } + } else { + # If timestamps in the csv file are formatted (Y-M-D h:m:s.f), data.table::fread assumes them to be + # in the timezone of the current device (i.e. as if configtz == ""). + # If this is not the case, we need to force the correct timezone (force, as opposed to converting to that timezone). + if (!is.numeric(rawTime) && configtz != "") { + rawTime = lubridate::force_tz(rawTime, configtz) + } + + # A similar thing needs to be done if the csv file contains Unix timestamps as well. + # OmGui converts device timestamps to Unix timestamps as if the timestamps were originally in UTC. + # So we need to convert the Unix timestamp into a hh:mm:ss format, then force that timestamp + # from UTC into configtz timzone. + if (is.numeric(rawTime)) { + rawTime = as.POSIXct(rawTime, tz="UTC", origin = "1970-01-01") + rawTime = lubridate::force_tz(rawTime, configtz) + } } + + rawTime = as.numeric(rawTime) + # resample the acceleration data, because AX3 data is stored at irregular time points - rawTime = vector(mode = "numeric", nrow(P)) - if (length(params_general[["desiredtz"]]) == 0 & blocknumber == 1) { - cat("Forgot to specify argument desiredtz? Now Europe/London assumed") - params_general[["desiredtz"]] = "Europe/London" - } - rawTime = as.numeric(as.POSIXlt(P[,1],tz = params_general[["desiredtz"]])) - rawAccel = as.matrix(P[,2:4]) + rawAccel = as.matrix(rawData[,2:4]) step = 1/sf - start = rawTime[1] - end = rawTime[length(rawTime)] - timeRes = seq(start, end, step) - nr = length(timeRes) - 1 - timeRes = as.vector(timeRes[1:nr]) - accelRes = matrix(0,nrow = nr, ncol = 3, dimnames = list(NULL, c("x", "y", "z"))) - # at the moment the function is designed for reading the r3 acceleration channels only, + timeRes = seq(rawTime[1], rawTime[length(rawTime)], step) + timeRes = timeRes[1 : (length(timeRes) - 1)] + + # at the moment the function is designed for reading the 3 acceleration channels only, # because that is the situation of the use-case we had. - rawLast = nrow(rawAccel) - accelRes = GGIRread::resample(rawAccel, rawTime, timeRes, rawLast, params_rawdata[["interpolationType"]]) # this is now the resampled acceleration data - P = cbind(timeRes,accelRes) - } else { - P = c() + accelRes = GGIRread::resample(rawAccel, rawTime, timeRes, nrow(rawAccel), params_rawdata[["interpolationType"]]) # this is now the resampled acceleration data + P$data = data.frame(timeRes, accelRes) + colnames(P$data) = c("time", "x", "y", "z") } } else if (mon == MONITOR$MOVISENS && dformat == FORMAT$BIN) { - startpage = blocksize * (blocknumber - 1) + 1 - deltapage = blocksize - UPI = updatepageindexing(startpage = startpage, deltapage = deltapage, - blocknumber = blocknumber, PreviousEndPage = PreviousEndPage, mon = mon, dformat = dformat) - startpage = UPI$startpage; endpage = UPI$endpage file_length = unisensR::getUnisensSignalSampleCount(dirname(filename), "acc.bin") if (endpage > file_length) { endpage = file_length - switchoffLD = 1 - } - P = unisensR::readUnisensSignalEntry(dirname(filename), "acc.bin", - startIndex = startpage, - endIndex = endpage) - P = as.matrix(P) - if (nrow(P) < ((sf * ws * 2) + 1) & blocknumber == 1) { - P = c() - switchoffLD = 1 - filequality$filetooshort = TRUE + isLastBlock = TRUE } + P$data = unisensR::readUnisensSignalEntry(dirname(filename), "acc.bin", + startIndex = startpage, + endIndex = endpage) + if (length(P$data) > 0) { + if (ncol(P$data) < 3) { + P$data = c() + } else { + colnames(P$data) = c("x", "y", "z") + # there may or may not be a temp.bin file containing temperature + try(expr = {P$data$temperature = g.readtemp_movisens(filename, + from = startpage, to = endpage, + acc_sf = sf, acc_length = nrow(P$data), + interpolationType = params_rawdata[["interpolationType"]]) + }, silent = TRUE) + } + } } else if (mon == MONITOR$ACTIGRAPH && dformat == FORMAT$GT3X) { - startpage = blocksize * (blocknumber - 1) + 1 - deltapage = blocksize - UPI = updatepageindexing(startpage = startpage, deltapage = deltapage, - blocknumber = blocknumber, PreviousEndPage = PreviousEndPage, mon = mon, dformat = dformat) - startpage = UPI$startpage; endpage = UPI$endpage - P = try(expr = {as.data.frame(read.gt3x::read.gt3x(path = filename, batch_begin = startpage, - batch_end = endpage,asDataFrame = TRUE))}, silent = TRUE) - if (length(P) == 0 | inherits(P, "try-error") == TRUE) { # too short or not data at all - P = c() ; switchoffLD = 1 - if (blocknumber == 1) filequality$filetooshort = TRUE - if (blocknumber == 1) filequality$filecorrupt = TRUE - } else { - if (nrow(P) < ((sf * ws * 2) + 1)) { - P = c() ; switchoffLD = 1 - if (blocknumber == 1) filequality$filetooshort = TRUE - } # If data passes these checks then it is usefull + P$data = try(expr = {read.gt3x::read.gt3x(path = filename, batch_begin = startpage, + batch_end = endpage, asDataFrame = TRUE)}, silent = TRUE) + if (length(P$data) == 0 || inherits(P$data, "try-error") == TRUE) { # too short or no data at all + P$data = c() + } else { # If data passes these checks then it is usefull + colnames(P$data)[colnames(P$data) == "X"] = "x" + colnames(P$data)[colnames(P$data) == "Y"] = "y" + colnames(P$data)[colnames(P$data) == "Z"] = "z" + + # read.gt3x::read.gt3x returns timestamps as POSIXct with GMT timezone, but they are actally in local time of the device. + # Don't just convert timezones, instead force the correct local timezone of the device (configtz) + # while keeping the same hh:mm:ss time. + P$data$time = lubridate::force_tz(P$data$time, configtz) } } else if (mon == MONITOR$AD_HOC && dformat == FORMAT$AD_HOC_CSV) { # user-specified csv format - startpage = (1 + (blocksize * 300 * (blocknumber - 1))) - deltapage = (blocksize*300) - UPI = updatepageindexing(startpage = startpage,deltapage = deltapage, - blocknumber = blocknumber,PreviousEndPage = PreviousEndPage, - mon = mon, dformat = dformat) - startpage = UPI$startpage; endpage = UPI$endpage + # skip 1 more row only if rmc.firstrow.acc points at a row containing column names. + # This is only relevant for the first chunk of data. + if (blocknumber == 1) { + testheader = data.table::fread(filename, nrows = 2, skip = params_rawdata[["rmc.firstrow.acc"]]-1, + dec = decn, showProgress = FALSE, + header = TRUE, data.table=FALSE, stringsAsFactors=FALSE) + if (suppressWarnings(is.na(as.numeric(colnames(testheader)[1])))) { # first value is *not* a number, so file starts with a header + startpage = startpage + 1 + endpage = endpage + 1 + } + } + try(expr = {P = read.myacc.csv(rmc.file = filename, - rmc.nrow = deltapage, rmc.skip = startpage, + rmc.nrow = blocksize, rmc.skip = startpage, rmc.dec = params_rawdata[["rmc.dec"]], rmc.firstrow.acc = params_rawdata[["rmc.firstrow.acc"]], rmc.firstrow.header = params_rawdata[["rmc.firstrow.header"]], @@ -377,23 +330,63 @@ g.readaccfile = function(filename, blocksize, blocknumber, filequality, interpolationType = params_rawdata[["interpolationType"]], PreviousLastValue = PreviousLastValue, PreviousLastTime = PreviousLastTime, - epochsize = params_general[["windowsizes"]][1:2], - desiredtz = params_general[["desiredtz"]], - configtz = params_general[["configtz"]]) + desiredtz = desiredtz, + configtz = configtz, + header = header) }, silent = TRUE) if (length(sf) == 0) sf = params_rawdata[["rmc.sf"]] - if (length(P) == 4) { # added PreviousLastValue and PreviousLastTime as output of read.myacc.csv - # P = as.matrix(P) # turned off 21-5-2019 - if (nrow(P$data) < ((sf * ws * 2) + 1) & blocknumber == 1) { - P = c() ; switchoffLD = 1 #added 30-6-2012 - filequality$filetooshort = TRUE - } - } else { + } + + # if first block isn't read then the file is probably corrupt + if (length(P$data) <= 1 || nrow(P$data) == 0) { + P = c() + isLastBlock = TRUE + if (blocknumber == 1) { + warning('\nFile empty, possibly corrupt.\n') + filequality$filetooshort = TRUE + filequality$filecorrupt = TRUE + } + } else if (nrow(P$data) < (sf * ws * 2 + 1)) { + # a shorter chunk of data than expected was read + isLastBlock = TRUE + + if (blocknumber == 1) { + # not enough data for analysis P = c() + filequality$filetooshort = TRUE } } + + # remove any columns we don't need/expect + P$data = P$data[,which(colnames(P$data) %in% c("x", "y", "z", "time", "light", "temperature", "wear"))] + + # every column except for time and wear should be numeric. + # If it isn't, some non-numeric input was present, and we don't know how to deal with it. + for (col in c("x", "y", "z", "light", "temperature")) { + if ((col %in% colnames(P$data)) && !is.numeric(P$data[, col])) { + stop(paste0("Corrupt file. ", col, " column contains non-numeric data.")) + } + } + + # the wear column should be logical, but we will coerse it to numeric right away, + # so that later it could be combined with the other numeric columns into a numeric matrix + + if ("wear" %in% colnames(P$data)) { + if (!is.logical(P$data$wear)) { + stop("Corrupt file. The wear column should contail TRUE/FALSE values.") + } + P$data$wear = as.numeric(P$data$wear) + } + + # the time column at this point will be either Unix timestamps or POSIXct objects. + # If POSIXct, we'll convert them to Unix timestamps, so that later this column + # could be combined with the other numeric columns into a numeric matrix + if (("time" %in% colnames(P$data)) && !is.numeric(P$data$time)) { + P$data$time = as.numeric(P$data$time) + } + invisible(list(P = P, filequality = filequality, - switchoffLD = switchoffLD, + isLastBlock = isLastBlock, endpage = endpage, startpage = startpage)) } diff --git a/R/g.readtemp_movisens.R b/R/g.readtemp_movisens.R index 0f04d3fa9..8dbb3165f 100644 --- a/R/g.readtemp_movisens.R +++ b/R/g.readtemp_movisens.R @@ -1,25 +1,33 @@ -g.readtemp_movisens = function(datafile, desiredtz = "", from = c(), to = c(), interpolationType=1) { - temperature = unisensR::readUnisensSignalEntry(dirname(datafile), "temp.bin") - temperature = as.data.frame(temperature) - origin = unisensR::readUnisensStartTime(dirname(datafile)) - temperature$timestamp = seq(origin, origin + nrow(temperature) - 1, by = 1) - rawTime = vector(mode = "numeric", nrow(temperature)) - rawTime = as.numeric(as.POSIXlt(temperature$timestamp,tz = desiredtz)) - rawTemp = as.matrix(temperature[,-c(which(colnames(temperature) == "timestamp"))]) - acc_length = unisensR::getUnisensSignalSampleCount(dirname(datafile), "acc.bin") - step = (nrow(temperature) - 1) / acc_length #ratio of temp sf to acc sf in movisens data - start = rawTime[1] - end = rawTime[length(rawTime)] - timeRes = seq(start, end, step) - nr = length(timeRes) - 1 - timeRes = as.vector(timeRes[1:nr]) - tempRes = matrix(0,nrow = nr, ncol = ncol(rawTemp), dimnames = list(NULL,colnames(rawTemp))) - rawLast = nrow(rawTemp) - tempRes = GGIRread::resample(rawTemp, rawTime, timeRes, rawLast, type=interpolationType) # this is now the resampled temp data - if(length(from) > 0 & length(to) > 0) { - temperature = tempRes[from:to] - } else { - temperature = tempRes +g.readtemp_movisens = function(datafile, from = c(), to = c(), acc_sf, acc_length, interpolationType=1) { + # Acceleration data and temperature were sampled at different rates, + # so we need to resample temperature to get the same sampling rate + # as what we have for acceleration. + + temp_sf = 1 # temperature is most likely sampled at 1Hz, but we'll double-check later + + temp_from = ceiling(from / acc_sf * temp_sf) + temp_to = ceiling(to / acc_sf * temp_sf) + + temperature = unisensR::readUnisensSignalEntry(dirname(datafile), "temp.bin", + startIndex = temp_from, endIndex = temp_to) + new_temp_sf = attr(temperature, "sampleRate") + + # We had guessed that temperature was sampled at 1Hz. Let's check, and if we were wrong, then re-do. + if (temp_sf != new_temp_sf) { + temp_from = ceiling(from / acc_sf * new_temp_sf) + temp_to = ceiling(to / acc_sf * new_temp_sf) + + temperature = unisensR::readUnisensSignalEntry(dirname(datafile), "temp.bin", + startIndex = temp_from, endIndex = temp_to) } + + temperature = temperature$temp + + # we don't care about the exact timestamp values because we'll throw the timestamps away anyway. + rawTime = seq_len(length(temperature)) + timeRes = seq(from = 1, to = rawTime[length(rawTime)], length.out = acc_length) + + temperature = GGIRread::resample(as.matrix(temperature), rawTime, timeRes, length(temperature), type=interpolationType) + invisible(temperature) } diff --git a/R/get_nw_clip_block_params.R b/R/get_nw_clip_block_params.R index e088d96bf..802df6eb7 100644 --- a/R/get_nw_clip_block_params.R +++ b/R/get_nw_clip_block_params.R @@ -1,4 +1,4 @@ -get_nw_clip_block_params = function(chunksize, dynrange, monc, rmc.noise=c(), sf, dformat, +get_nw_clip_block_params = function(chunksize, dynrange, monc, dformat, deviceSerialNumber = "", rmc.noise=c(), sf, rmc.dynamic_range) { blocksize = round(14512 * (sf/50) * chunksize) if (monc == MONITOR$GENEA) blocksize = round(21467 * (sf/80) * chunksize) @@ -17,41 +17,35 @@ get_nw_clip_block_params = function(chunksize, dynrange, monc, rmc.noise=c(), sf if (monc == MONITOR$MOVISENS) blocksize = sf * 60 * 1440 if (monc == MONITOR$VERISENSE && dformat == FORMAT$CSV) blocksize = round(blocksize) - # Clipping threshold: estimate number of data points of clipping based on raw data at about 87 Hz + + if (monc == MONITOR$ACTIGRAPH) { + # If Actigraph then try to specify dynamic range based on Actigraph model + if (length(grep(pattern = "CLE", x = deviceSerialNumber)) == 1) { + dynrange = 6 + } else if (length(grep(pattern = "MOS", x = deviceSerialNumber)) == 1) { + dynrange = 8 + } else if (length(grep(pattern = "NEO", x = deviceSerialNumber)) == 1) { + dynrange = 6 + } + } + + # Clipping threshold if (length(dynrange) > 0) { clipthres = dynrange - 0.5 } else { - if (monc == MONITOR$GENEA) { - clipthres = 5.5 - } else if (monc == MONITOR$GENEACTIV) { - clipthres = 7.5 - } else if (monc == MONITOR$ACTIGRAPH) { - clipthres = 7.5 # hard coded assumption that dynamic range is 8g - } else if (monc == MONITOR$AXIVITY) { - clipthres = 7.5 # hard coded assumption that dynamic range is 8g - } else if (monc == MONITOR$MOVISENS) { + clipthres = 7.5 # hard-coded assumption that dynamic range is 8g + #if (monc == MONITOR$GENEA) clipthres = 5.5 + if (monc == MONITOR$MOVISENS) { clipthres = 15.5 # hard coded assumption that dynamic range is 16g - } else if (monc == MONITOR$VERISENSE) { - clipthres = 7.5 } else if (monc == MONITOR$AD_HOC) { clipthres = rmc.dynamic_range } } # Nonwear threshold: non-wear criteria are monitor-specific racriter = 0.15 # very likely irrelevant parameters, but leave in for consistency - if (monc == MONITOR$GENEA) { - sdcriter = 0.003 - racriter = 0.05 - } else if (monc == MONITOR$GENEACTIV) { - sdcriter = 0.013 - } else if (monc == MONITOR$ACTIGRAPH) { - sdcriter = 0.013 - } else if (monc == MONITOR$AXIVITY) { - sdcriter = 0.013 - } else if (monc == MONITOR$MOVISENS) { - sdcriter = 0.013 - } else if (monc == MONITOR$VERISENSE) { - sdcriter = 0.013 + sdcriter = 0.013 + #if (monc == MONITOR$GENEA) { sdcriter = 0.003; racriter = 0.05 } + if (monc == MONITOR$VERISENSE) { racriter = 0.20 } else if (monc == MONITOR$AD_HOC) { if (length(rmc.noise) == 0) { diff --git a/R/get_starttime_weekday_meantemp_truncdata.R b/R/get_starttime_weekday_meantemp_truncdata.R deleted file mode 100644 index 604c90cbb..000000000 --- a/R/get_starttime_weekday_meantemp_truncdata.R +++ /dev/null @@ -1,154 +0,0 @@ -get_starttime_weekday_meantemp_truncdata = function(temp.available, monc, dformat, data, - P, header, desiredtz, sf, i, datafile, - ws2, starttime, wday, wdayname, configtz = NULL) { - #ensures that first window starts at logical timepoint relative to its size - # (15,30,45 or 60 minutes of each hour) - start_meas = ws2/60 - if (temp.available == TRUE) { - use.temp = TRUE - } else { - use.temp = FALSE - } - meantemp = c() - if (monc == MONITOR$GENEACTIV || (monc == MONITOR$AXIVITY && dformat == FORMAT$CWA) || monc == MONITOR$MOVISENS || (monc == MONITOR$AD_HOC && use.temp == TRUE)) { - if (monc == MONITOR$GENEACTIV) { - if ("temperature" %in% colnames(data)) { - tempcolumn = which(colnames(data) == "temperature") #GGIRread - } else { - tempcolumn = 7 - } - } - if (monc == MONITOR$AXIVITY || monc == MONITOR$AD_HOC) tempcolumn = 5 - if (monc == MONITOR$MOVISENS) tempcolumn = 4 - meantemp = mean(as.numeric(data[, tempcolumn]), na.rm = TRUE) - if (is.na(meantemp) == T) { #mean(as.numeric(data[1:10,7])) - warning("temperature is NA", call. = FALSE) - meantemp = 0 - use.temp = FALSE - } else if (mean(as.numeric(data[1:10,tempcolumn])) > 50) { - warning("temperature value is unreaslistically high (> 50 Celcius)", call. = FALSE) - meantemp = 0 - use.temp = FALSE - } - } - # extraction and modification of starting point of measurement - if (i == 1) { #only do this for first block of data - starttime = g.getstarttime( - datafile = datafile, - P = P, - header = header, - mon = monc, - dformat = dformat, - desiredtz = desiredtz, - configtz = configtz - ) - if (exists("P")) rm(P); gc() - #================================================== - #inspection timezone - timezone = attr(unclass(as.POSIXlt(starttime[1])), which = "tzone") - starttimebefore = as.POSIXlt(starttime) - # assuming that timestamps is good, but that timezone might be lost in conversion from string to POSIXct - if (dformat == FORMAT$BIN) { #not sure whether this is required for csv-format (2) - if (length(which(timezone == "GMT")) > 0) { - if (length(desiredtz) == 0) { - warning("desiredtz not specified, local timezoneused as default", call. = FALSE) - desiredtz = "" - } - starttime = as.POSIXlt(starttime[1], tz = desiredtz) - } - } - #================================================ - #assess weekday - wday = unclass(as.POSIXlt(starttime[1]))$wday #day of the week 0-6 and 0 is Sunday - wday = wday + 1 - weekdays = c("Sunday","Monday","Tuesday","Wednesday","Thursday","Friday","Saturday") - wdayname = weekdays[wday] - #====================================================== - #assess how much data to delete till next 15 minute period - tmp = unlist(strsplit(format(starttime)," ")) - if (length(tmp) > 1) { - starttime2 = as.numeric(unlist(strsplit(tmp[2],":"))) - } else if (length(tmp) == 1 && !grepl("T", tmp)) { - starttime2 = as.numeric(unlist(strsplit(tmp[2],":"))) - } else { - # first get char to POSIX - tmp = iso8601chartime2POSIX(starttime, tz = desiredtz) - tmp = unlist(strsplit(format(tmp)," ")) # to keep it consistent with what we had - starttime2 = as.numeric(unlist(strsplit(format(tmp[2]),":"))) - } - if (length(which(is.na(starttime2) == TRUE)) > 0 | - length(starttime2) == 0) { - starttime2 = c(0,0,0) - } - start_hr = as.numeric(starttime2[1]) - start_min = as.numeric(starttime2[2]) - start_sec = as.numeric(starttime2[3]) - - secshift = 60 - start_sec #shift in seconds needed - if (secshift != 60) { - start_min = start_min + 1 #shift in minutes needed (+1 one to account for seconds comp) - } - if (secshift == 60) secshift = 0 # if starttime is 00:00 then we do not want to remove data - minshift = start_meas - (((start_min/start_meas) - floor(start_min/start_meas)) * start_meas) - if (minshift == start_meas) { - minshift = 0; - - } - sampleshift = ((minshift)*60*sf) + (secshift*sf) #derive sample shift - if (floor(sampleshift) > 1) { - data = data[-c(1:floor(sampleshift)),] #delete data accordingly - } - newmin = start_min + minshift #recalculate first timestamp - newsec = 0 - remem2add24 = FALSE - if (newmin >= 60) { - newmin = newmin - 60 - start_hr = start_hr + 1 - if (start_hr == 24) { #if measurement is started in 15 minutes before midnight - start_hr = 0 - remem2add24 = TRUE #remember to add 24 hours because this is now the wrong day - } - } - starttime3 = paste0(tmp[1], " ", start_hr, ":", newmin, ":", newsec) - #create timestamp from string (now desiredtz is added) - if (length(desiredtz) == 0) { - warning("desiredtz not specified, local timezone used as default", call. = FALSE) - desiredtz = "" - } - starttime_a = as.POSIXct(starttime3, format = "%d/%m/%Y %H:%M:%S", tz = desiredtz) - starttime_b = as.POSIXct(starttime3, format = "%d-%m-%Y %H:%M:%S", tz = desiredtz) - starttime_c = as.POSIXct(starttime3, format = "%Y/%m/%d %H:%M:%S", tz = desiredtz) - starttime_d = as.POSIXct(starttime3, format = "%Y-%m-%d %H:%M:%S", tz = desiredtz) - if (is.na(starttime_a) == FALSE) { - starttime = starttime_a - } else { - if (is.na(starttime_b) == FALSE) { - starttime = starttime_b - } else { - if (is.na(starttime_c) == FALSE) { - starttime = starttime_c - } else { - if (is.na(starttime_d) == FALSE) { - starttime = starttime_d - } else { - warning("date not recognized") - } - } - } - } - if (remem2add24 == TRUE) { - starttime = as.POSIXlt(as.numeric(starttime) + (24 * 3600), origin = "1970-01-01") - } - } - invisible( - list( - starttime = starttime, - meantemp = meantemp, - use.temp = use.temp, - wday = wday, - wdayname = wdayname, - desiredtz = desiredtz, - data = data - ) - ) -} diff --git a/R/get_starttime_weekday_truncdata.R b/R/get_starttime_weekday_truncdata.R new file mode 100644 index 000000000..8f6c3c248 --- /dev/null +++ b/R/get_starttime_weekday_truncdata.R @@ -0,0 +1,54 @@ +get_starttime_weekday_truncdata = function(monc, dformat, data, + header, desiredtz, sf, datafile, + ws2, configtz = NULL) { + # ensures that first window starts at logical timepoint relative to its size + # (15,30,45 or 60 minutes of each hour) + start_meas = ws2/60 + + starttime = g.getstarttime( # this will return a POSIXlt object + datafile = datafile, + data = data, + mon = monc, + dformat = dformat, + desiredtz = desiredtz, + configtz = configtz + ) + + # assess weekday + wday = starttime$wday #day of the week 0-6 and 0 is Sunday + wday = wday + 1 + weekdays = c("Sunday","Monday","Tuesday","Wednesday","Thursday","Friday","Saturday") + wdayname = weekdays[wday] + + # assess how much data to delete till next 15 minute period + secshift = 60 - starttime$sec # shift in seconds needed + if (secshift == 60) { + secshift = 0 # if starttime is 00:00 then we do not want to remove data + } else { + starttime$min = starttime$min + 1 # shift in minutes needed (+1 one to account for seconds comp) + } + + minshift = start_meas - (starttime$min %% start_meas) + if (minshift == start_meas) { + minshift = 0 + } + + sampleshift = (minshift * 60 * sf) + (secshift * sf) # derive sample shift + sampleshift = floor(sampleshift) + if (sampleshift > 1) { + data = data[-c(1 : sampleshift),] # delete data accordingly + } + + # recalculate the timestamp + starttime$min = starttime$min + minshift # if the result is >= 60, hours (and possibly date) will adjust automatically + starttime$sec = 0 + + invisible( + list( + starttime = starttime, + wday = wday, + wdayname = wdayname, + data = data + ) + ) +} diff --git a/R/ismovisens.R b/R/ismovisens.R index 20af3126f..b5c65c27b 100644 --- a/R/ismovisens.R +++ b/R/ismovisens.R @@ -1,15 +1,15 @@ ismovisens = function(data){ # --------------------------- - # Is data a directory? Then use the first file to test if recorded with movisens - directory = dir(data, recursive = T, full.names = T)[1] - isdir = !is.na(directory) - if (isdir == TRUE) data = directory + # Is data a directory? Then use the first file in this directory + # to test if recorded with movisens + first_file = dir(data, recursive = T, full.names = T)[1] + isdir = !is.na(first_file) + if (isdir) data = first_file # --------------------------- - # remove problematic dot at the beginning - data = gsub("^./", "", data) - # Now, I focus on the participant directory - p_dir = strsplit(data, ".", fixed = T) - data_tmp = strsplit(p_dir[[1]], "/") + # At this point, data is guaranteed to be a path to a file. + # Now, focus on the directory containing this file. + data_tmp = strsplit(data, "/") + data_ln = length(data_tmp[[1]]) - 1 data = paste(data_tmp[[1]][1:data_ln], collapse = "/") # --------------------------- diff --git a/R/read.myacc.csv.R b/R/read.myacc.csv.R index acca94ab8..6feb08c4c 100644 --- a/R/read.myacc.csv.R +++ b/R/read.myacc.csv.R @@ -22,145 +22,149 @@ read.myacc.csv = function(rmc.file=c(), rmc.nrow=Inf, rmc.skip=c(), rmc.dec=".", interpolationType = 1, PreviousLastValue = c(0, 0, 1), PreviousLastTime = NULL, - epochsize = NULL, desiredtz = NULL, - configtz = NULL) { + configtz = NULL, + header = NULL) { - if (!is.null(rmc.desiredtz) | !is.null(rmc.configtz)) { + if (length(rmc.col.time) > 0 && !(rmc.unit.time %in% c("POSIX", "character", "UNIXsec", "ActivPAL"))) { + stop("\nUnrecognized rmc.col.time value. The only accepted values are \"POSIX\", \"character\", \"UNIXsec\", and \"ActivPAL\".", call. = FALSE) + } + + if (!is.null(rmc.desiredtz) || !is.null(rmc.configtz)) { generalWarning = paste0("Argument rmc.desiredtz and rmc.configtz are scheduled to be deprecated", " and will be replaced by the existing arguments desiredtz and configtz, respectively.") - showGeneralWarning = TRUE + # Check if both types of tz are provided: if (!is.null(desiredtz) && desiredtz != "" && !is.null(rmc.desiredtz)) { if (rmc.desiredtz != desiredtz) { # if different --> error (don't know which one to use) - showGeneralWarning = FALSE stop(paste0("\n", generalWarning, "Please, specify only desiredtz and set ", "rmc.desiredtz to NULL to ensure it is no longer used.")) } } - if (!is.null(configtz) & !is.null(rmc.configtz)) { # then both provided + if (!is.null(configtz) && !is.null(rmc.configtz)) { # then both provided if (rmc.configtz != configtz) { # if different --> error (don't know which one to use) - showGeneralWarning = FALSE stop(paste0("\n", generalWarning, "Please, specify only configtz and set ", "rmc.configtz to NULL to ensure it is no longer used.")) } } - if (showGeneralWarning == TRUE) { - warning(paste0("\n", generalWarning)) - } + warning(paste0("\n", generalWarning)) # Until deprecation still allow rmc. to be used, - # so us it to overwrite normal tz in this function: + # so use it to overwrite normal tz in this function: if (is.null(desiredtz)) desiredtz = rmc.desiredtz if (desiredtz == "" && !is.null(rmc.desiredtz)) desiredtz = rmc.desiredtz if (is.null(configtz)) configtz = rmc.configtz } # check if none of desiredtz and rmc.desiredtz are provided - if (is.null(desiredtz) & is.null(rmc.desiredtz)) { + if (is.null(desiredtz) && is.null(rmc.desiredtz)) { stop(paste0("Timezone not specified, please provide at least desiredtz", " and consider specifying configtz.")) } - - # bitrate should be or header item name as character, or the actual numeric bit rate - # unit.temp can take C(elsius), F(ahrenheit), and K(elvin) and converts it into Celsius - # Note all argument names start with rmc (read myacc csv) to avoid name clashes when passed on throughout GGIR - if (length(rmc.firstrow.header) == 0) { # no header block - if (rmc.firstrow.acc == 2) { - freadheader = TRUE + if (is.null(rmc.firstrow.acc) || rmc.firstrow.acc < 1) { + rmc.firstrow.acc = 1 + } + skip = rmc.firstrow.acc - 1 + if (!is.null(rmc.skip) && length(rmc.skip) > 0) { + skip = skip + rmc.skip + } + + # only extract the header if it hasn't been extracted for this file before + if (is.null(header)) { + # bitrate should be either header item name as character, or the actual numeric bit rate. + # unit.temp can take C(elsius), F(ahrenheit), and K(elvin) and converts it into Celsius + # Note all argument names start with rmc (read myacc csv) to avoid name clashes when passed on throughout GGIR + if (length(rmc.firstrow.header) == 0) { # no header block + sf = rmc.sf + header = "no header" } else { - freadheader = FALSE - } - skip = rmc.firstrow.acc - sf = rmc.sf - header = "no header" - } else { - # extract header information: - if (length(rmc.header.length) == 0) { - rmc.header.length = rmc.firstrow.acc - 1 - } - - options(warn = -1) # fread complains about quote in first row for some file types - header_tmp = as.data.frame(data.table::fread(file = rmc.file, - nrows = rmc.header.length, - skip = rmc.firstrow.header - 1, - dec = rmc.dec, showProgress = FALSE, header = FALSE, - stringsAsFactors = TRUE, - blank.lines.skip = TRUE)) - options(warn = 0) - if (length(rmc.header.structure) != 0) { # header is stored in 1 column, with strings that need to be split - if (length(header_tmp) == 1) { # one header item - header_tmp = as.matrix(unlist(strsplit(as.character(header_tmp[,1]), rmc.header.structure))) - } else { # multiple header items - mysplit = function(x){ - tmp = strsplit(as.character(x), rmc.header.structure) - tmp = unlist(tmp) - return(tmp) + # extract header information: + if (length(rmc.header.length) == 0) { + rmc.header.length = rmc.firstrow.acc - 1 + } + + options(warn = -1) # fread complains about quote in first row for some file types + header_tmp = data.table::fread(file = rmc.file, + nrows = rmc.header.length, + skip = rmc.firstrow.header - 1, + dec = rmc.dec, showProgress = FALSE, header = FALSE, + blank.lines.skip = TRUE, + data.table=FALSE, stringsAsFactors=FALSE) + options(warn = 0) + if (length(rmc.header.structure) != 0) { # header is stored in 1 column, with strings that need to be split + if (length(header_tmp) == 1) { # one header item + header_tmp = as.matrix(unlist(strsplit(as.character(header_tmp[,1]), rmc.header.structure))) + } else { # multiple header items + mysplit = function(x){ + tmp = strsplit(as.character(x), rmc.header.structure) + tmp = unlist(tmp) + return(tmp) + } + header_tmp0 = header_tmp + header_tmp = unlist(lapply(header_tmp[,1], FUN = mysplit)) + if (length(header_tmp) > 2) { + header_tmp = data.frame(matrix(unlist(header_tmp), nrow = nrow(header_tmp0), byrow = T), stringsAsFactors = FALSE) + colnames(header_tmp) = NULL + } else { + header_tmp = data.frame(matrix(unlist(header_tmp), nrow = 1, byrow = T), stringsAsFactors = FALSE) + colnames(header_tmp) = NULL + } + } + if (ncol(header_tmp) == 1) header_tmp = t(header_tmp) + header_tmp2 = as.data.frame(as.character(unlist(header_tmp[,2])), stringsAsFactors = FALSE) + row.names(header_tmp2) = header_tmp[,1] + colnames(header_tmp2) = NULL + header = header_tmp2 + } else { # column 1 is header name, column 2 is header value + colnames(header_tmp) = NULL + validrows = which(is.na(header_tmp[,1]) == FALSE & header_tmp[,1] != "") + header_tmp = header_tmp[validrows,1:2] + header_tmp2 = as.data.frame(header_tmp[,2], stringsAsFactors = FALSE) + row.names(header_tmp2) = header_tmp[,1] + colnames(header_tmp2) = NULL + header = header_tmp2 + } + # assess whether accelerometer data conversion is needed + if (length(rmc.bitrate) > 0 && length(rmc.dynamic_range) > 0 && rmc.unit.acc == "bit") { + if (is.character(rmc.bitrate[1]) == TRUE) { # extract bitrate if it is in the header + rmc.bitrate = as.numeric(header[which(row.names(header) == rmc.bitrate[1]),1]) } - header_tmp0 = header_tmp - header_tmp = unlist(lapply(header_tmp[,1], FUN = mysplit)) - if (length(header_tmp) > 2) { - header_tmp = data.frame(matrix(unlist(header_tmp), nrow = nrow(header_tmp0), byrow = T), stringsAsFactors = TRUE) - colnames(header_tmp) = NULL - } else { - header_tmp = data.frame(matrix(unlist(header_tmp), nrow = 1, byrow = T), stringsAsFactors = TRUE) - colnames(header_tmp) = NULL + if (is.character(rmc.dynamic_range[1]) == TRUE) { # extract dynamic range if it is in the header + rmc.dynamic_range = as.numeric(header[which(row.names(header) == rmc.dynamic_range[1]),1]) + } + } + # extract sample frequency: + sf = as.numeric(header[which(row.names(header) == rmc.headername.sf[1]),1]) + + if (is.na(sf)) { # sf not retrieved from header + # first see if maybe sf *is* in the header, just not under the rmc.headername.sf name + sf = as.numeric(header[which(row.names(header) == "sample_rate"),1]) + # if sf isn't in the header under the default name either, then use the default value + if (is.na(sf)) { + sf = rmc.sf + header = rbind(header, sf) # add it also to the header + row.names(header)[nrow(header)] = "sample_rate" } } - if (ncol(header_tmp) == 1) header_tmp = t(header_tmp) - header_tmp2 = as.data.frame(as.character(unlist(header_tmp[,2])), stringsAsFactors = FALSE) - row.names(header_tmp2) = header_tmp[,1] - colnames(header_tmp2) = NULL - header = header_tmp2 - } else { # column 1 is header name, column 2 is header value - colnames(header_tmp) = NULL - validrows = which(is.na(header_tmp[,1]) == FALSE & header_tmp[,1] != "") - header_tmp = header_tmp[validrows,1:2] - header_tmp2 = as.data.frame(header_tmp[,2], stringsAsFactors = FALSE) - row.names(header_tmp2) = header_tmp[,1] - colnames(header_tmp2) = NULL - header = header_tmp2 - } - skip = rmc.firstrow.acc - 1 - freadheader = TRUE - # assess whether accelerometer data conversion is needed - if (length(rmc.bitrate) > 0 & length(rmc.dynamic_range) > 0 & rmc.unit.acc == "bit") { - if (is.character(rmc.bitrate[1]) == TRUE) { # extract bitrate if it is in the header - rmc.bitrate = as.numeric(header[which(row.names(header) == rmc.bitrate[1]),1]) + + # standardise key header names to ease use elsewhere in GGIR: + if (length(rmc.headername.sf) > 0) { + row.names(header)[which(row.names(header) == rmc.headername.sf[1])] = "sample_rate" + } + if (length(rmc.headername.sn) > 0) { + row.names(header)[which(row.names(header) == rmc.headername.sn[1])] = "device_serial_number" + } + if (length(rmc.headername.recordingid) > 0) { + row.names(header)[which(row.names(header) == rmc.headername.recordingid[1])] = "recordingID" } - if (is.character(rmc.dynamic_range[1]) == TRUE) { # extract dynamic range if it is in the header - rmc.dynamic_range = as.numeric(header[which(row.names(header) == rmc.dynamic_range[1]),1]) - } - } - # extract sample frequency: - sf = as.numeric(header[which(row.names(header) == rmc.headername.sf[1]),1]) - sn = as.numeric(header[which(row.names(header) == rmc.headername.sn[1]),1]) - ID = as.numeric(header[which(row.names(header) == rmc.headername.recordingid[1]),1]) - - # standardise key header names to ease use elsewhere in GGIR: - if (length(rmc.headername.sf) > 0) { - row.names(header)[which(row.names(header) == rmc.headername.sf[1])] = "sample_rate" - } - if (length(rmc.headername.sn) > 0) { - row.names(header)[which(row.names(header) == rmc.headername.sn[1])] = "device_serial_number" - } - if (length(rmc.headername.recordingid) > 0) { - row.names(header)[which(row.names(header) == rmc.headername.recordingid[1])] = "recordingID" - } - if (length(sf) == 0) { - sf = rmc.sf # if sf not retrieved from header than use default - header = rbind(header,1) # add it also to the header - row.names(header)[nrow(header)] = "sample_rate" } } - if (length(rmc.skip) > 0) { - skip = skip + rmc.skip - } # read data from file - P = as.data.frame(data.table::fread(rmc.file,nrows = rmc.nrow, skip = skip, - dec = rmc.dec, showProgress = FALSE, header = freadheader), - stringsAsFactors = TRUE) + P = data.table::fread(rmc.file, nrows = rmc.nrow, skip = skip, + dec = rmc.dec, showProgress = FALSE, header = "auto", + data.table=FALSE, stringsAsFactors=FALSE) + if (length(configtz) == 0) { configtz = desiredtz } @@ -170,28 +174,28 @@ read.myacc.csv = function(rmc.file=c(), rmc.nrow=Inf, rmc.skip=c(), rmc.dec=".", # select relevant columns, add standard column names P = P[,c(rmc.col.time, rmc.col.acc, rmc.col.temp)] if (length(rmc.col.time) > 0 & length(rmc.col.temp) > 0) { - colnames(P) = c("timestamp","accx","accy","accz","temperature") + colnames(P) = c("time","x","y","z","temperature") } else if (length(rmc.col.time) > 0 & length(rmc.col.temp) == 0) { - colnames(P) = c("timestamp","accx","accy","accz") + colnames(P) = c("time","x","y","z") } else if (length(rmc.col.time) == 0 & length(rmc.col.temp) > 0) { - colnames(P) = c("accx","accy","accz","temperature") + colnames(P) = c("x","y","z","temperature") } else if (length(rmc.col.time) == 0 & length(rmc.col.temp) == 0) { - colnames(P) = c("accx","accy","accz") + colnames(P) = c("x","y","z") } # acceleration and temperature as numeric - P$accx = as.numeric(P$accx) - P$accy = as.numeric(P$accy) - P$accz = as.numeric(P$accz) + P$x = as.numeric(P$x) + P$y = as.numeric(P$y) + P$z = as.numeric(P$z) if (length(rmc.col.temp) > 0) P$temperature = as.numeric(P$temperature) # Convert timestamps if (length(rmc.col.time) > 0) { if (rmc.unit.time == "POSIX") { - P$timestamp = as.POSIXct(format(P$timestamp), origin = rmc.origin, tz = configtz, format = rmc.format.time) + P$time = as.POSIXct(format(P$time), origin = rmc.origin, tz = configtz, format = rmc.format.time) checkdec = function(x) { # function to check whether timestamp has decimal places return(length(unlist(strsplit(as.character(x), "[.]|[,]"))) == 1) } - first_chunk_time = P$timestamp[1:pmin(nrow(P), 1000)] + first_chunk_time = P$time[1:pmin(nrow(P), 1000)] checkMissingDecPlaces = unlist(lapply(first_chunk_time, FUN = checkdec)) if (all(checkMissingDecPlaces) & !is.null(rmc.sf) & @@ -204,23 +208,23 @@ read.myacc.csv = function(rmc.file=c(), rmc.nrow=Inf, rmc.skip=c(), rmc.dec=".", # rmc.sf = 10 # P = data.frame(timestamps = c(rep(ttt - 1, 3), rep(ttt, 10), rep(ttt + 1, 9), rep(ttt + 2, 10), rep(ttt + 3, 4))) #------ - trans = unique(c(1, which(diff(P$timestamp) > 0), nrow(P))) + trans = unique(c(1, which(diff(P$time) > 0), nrow(P))) sf_tmp = diff(trans) timeIncrement = seq(0, 1 - (1/rmc.sf), by = 1/rmc.sf) # expected time increment per second # All seconds with exactly the sample frequency trans_1 = trans[which(sf_tmp == rmc.sf)] indices_1 = sort(unlist(lapply(trans_1, FUN = function(x){x + (1:rmc.sf)}))) - P$timestamp[indices_1] = P$timestamp[indices_1] + rep(timeIncrement, length(trans_1)) + P$time[indices_1] = P$time[indices_1] + rep(timeIncrement, length(trans_1)) # First second if (sf_tmp[1] != rmc.sf) { indices_2 = 1:trans[2] - P$timestamp[indices_2] = P$timestamp[indices_2] + seq(1 - (trans[2]/rmc.sf), 1 - (1/rmc.sf), by = 1/rmc.sf) + P$time[indices_2] = P$time[indices_2] + seq(1 - (trans[2]/rmc.sf), 1 - (1/rmc.sf), by = 1/rmc.sf) } # Last second if (sf_tmp[length(sf_tmp)] != rmc.sf) { indices_3 = (trans[length(trans)-1] + 1):trans[length(trans)] - P$timestamp[indices_3] = P$timestamp[indices_3] + timeIncrement[1:length(indices_3)] + P$time[indices_3] = P$time[indices_3] + timeIncrement[1:length(indices_3)] } # All seconds with other sample frequency and not the first or last second # This code now assumes that most samples in the second were sampled @@ -242,54 +246,53 @@ read.myacc.csv = function(rmc.file=c(), rmc.nrow=Inf, rmc.skip=c(), rmc.dec=".", } else if (length(timeIncrement) < sf2) { timeIncrement2 = c(timeIncrement, rep(timeIncrement[rmc.sf], sf2 - rmc.sf)) } - P$timestamp[indices_4] = P$timestamp[indices_4] + rep(timeIncrement2, length(trans_4)) + P$time[indices_4] = P$time[indices_4] + rep(timeIncrement2, length(trans_4)) } } } } } else if (rmc.unit.time == "character") { - P$timestamp = as.POSIXct(P$timestamp,format = rmc.format.time, tz = configtz) - } else if (rmc.unit.time == "UNIXsec") { - P$timestamp = as.POSIXct(P$timestamp, origin = rmc.origin, tz = configtz) + P$time = as.POSIXct(P$time, format = rmc.format.time, origin = rmc.origin, tz = configtz) + } else if (rmc.unit.time == "UNIXsec" && rmc.origin != "1970-01-01") { + P$time = as.POSIXct(P$time, origin = rmc.origin, tz = desiredtz) } else if (rmc.unit.time == "ActivPAL") { # origin should be specified as: "1899-12-30" rmc.origin = "1899-12-30" - datecode = round(P$timestamp) * 3600*24 - tmp2 = P$timestamp - round(P$timestamp) + datecode = round(P$time) * 3600*24 + tmp2 = P$time - round(P$time) timecode = ((tmp2 * 10^10) * 8.64) / 1000000 numerictime = datecode + timecode - P$timestamp = as.POSIXct(numerictime, origin = rmc.origin, tz = configtz) + P$time = as.POSIXct(numerictime, origin = rmc.origin, tz = desiredtz) } - if (length(which(is.na(P$timestamp) == FALSE)) == 0) { + if (length(which(is.na(P$time) == FALSE)) == 0) { stop("\nExtraction of timestamps unsuccesful, check timestamp format arguments") } - } - if (configtz != desiredtz) { - P$timestamp = as.POSIXct(as.numeric(P$timestamp), - tz = desiredtz, origin = "1970-01-01") + if(!is.numeric(P$time)) { # we'll return Unix timestamps + P$time = as.numeric(P$time) + } } # If acceleration is stored in mg units then convert to gravitational units if (rmc.unit.acc == "mg") { - P$accx = P$accx * 1000 - P$accy = P$accy * 1000 - P$accz = P$accz * 1000 + P$x = P$x * 1000 + P$y = P$y * 1000 + P$z = P$z * 1000 } if (rmc.scalefactor.acc != 1) { - P$accx = P$accx * rmc.scalefactor.acc - P$accy = P$accy * rmc.scalefactor.acc - P$accz = P$accz * rmc.scalefactor.acc + P$x = P$x * rmc.scalefactor.acc + P$y = P$y * rmc.scalefactor.acc + P$z = P$z * rmc.scalefactor.acc } # If acceleration is stored in bit values then convert to gravitational unit if (length(rmc.bitrate) > 0 & length(rmc.dynamic_range) > 0 & rmc.unit.acc == "bit") { if (rmc.unsignedbit == TRUE) { - P$accx = ((P$accx / (2^rmc.bitrate)) - 0.5) * 2 * rmc.dynamic_range - P$accy = ((P$accy / (2^rmc.bitrate)) - 0.5) * 2 * rmc.dynamic_range - P$accz = ((P$accz / (2^rmc.bitrate)) - 0.5) * 2 * rmc.dynamic_range + P$x = ((P$x / (2^rmc.bitrate)) - 0.5) * 2 * rmc.dynamic_range + P$y = ((P$y / (2^rmc.bitrate)) - 0.5) * 2 * rmc.dynamic_range + P$z = ((P$z / (2^rmc.bitrate)) - 0.5) * 2 * rmc.dynamic_range } else if (rmc.unsignedbit == FALSE) { # signed bit - P$accx = (P$accx / ((2^rmc.bitrate)/2)) * rmc.dynamic_range - P$accy = (P$accy / ((2^rmc.bitrate)/2)) * rmc.dynamic_range - P$accz = (P$accz / ((2^rmc.bitrate)/2)) * rmc.dynamic_range + P$x = (P$x / ((2^rmc.bitrate)/2)) * rmc.dynamic_range + P$y = (P$y / ((2^rmc.bitrate)/2)) * rmc.dynamic_range + P$z = (P$z / ((2^rmc.bitrate)/2)) * rmc.dynamic_range } } # Convert temperature units @@ -304,38 +307,26 @@ read.myacc.csv = function(rmc.file=c(), rmc.nrow=Inf, rmc.skip=c(), rmc.dec=".", # check for jumps in time and impute if (rmc.check4timegaps == TRUE) { if (length(sf) == 0) { # estimate sample frequency if not given in header - deltatime = abs(diff(as.numeric(P$timestamp))) + deltatime = abs(diff(as.numeric(P$time))) gapsi = which(deltatime > 0.25) - sf = (P$timestamp[gapsi[1]] - P$timestamp[1]) / (gapsi[1] - 1) + sf = (P$time[gapsi[1]] - P$time[1]) / (gapsi[1] - 1) } - P = g.imputeTimegaps(P, xyzCol = c("accx", "accy", "accz"), timeCol = "timestamp", sf = sf, k = 0.25, + P = g.imputeTimegaps(P, sf = sf, k = 0.25, PreviousLastValue = PreviousLastValue, PreviousLastTime = PreviousLastTime, epochsize = NULL) P = P$x - PreviousLastValue = as.numeric(P[nrow(P), c("accx", "accy", "accz")]) - PreviousLastTime = as.POSIXct(P[nrow(P), "timestamp"]) + PreviousLastValue = P[nrow(P), c("x", "y", "z")] + PreviousLastTime = as.POSIXct(P[nrow(P), "time"], origin = "1970-01-01") } - if (rmc.doresample == TRUE) { #resample - rawTime = vector(mode = "numeric", nrow(P)) - rawTime = as.numeric(as.POSIXct(P$timestamp,tz = configtz)) - rawAccel = as.matrix(P[,-c(which(colnames(P) == "timestamp"))]) - step = 1/sf - start = rawTime[1] - end = rawTime[length(rawTime)] - timeRes = seq(start, end, step) - nr = length(timeRes) - 1 - timeRes = as.vector(timeRes[1:nr]) - accelRes = matrix(0,nrow = nr, ncol = ncol(rawAccel), dimnames = list(NULL,colnames(rawAccel))) - rawLast = nrow(rawAccel) - accelRes = GGIRread::resample(rawAccel, rawTime, timeRes, rawLast, interpolationType) # this is now the resampled acceleration data - colnamesP = colnames(P) - timeRes = as.POSIXct(timeRes, origin = rmc.origin, tz = configtz) - P = as.data.frame(accelRes, stringsAsFactors = TRUE) - P$timestamp = timeRes - P = P[,c(ncol(P),1:(ncol(P) - 1))] + if (rmc.doresample == TRUE && ("time" %in% colnames(P))) { # resample + rawTime = P$time + rawAccel = as.matrix(P[,-c(which(colnames(P) == "time"))]) + timeRes = seq(from = rawTime[1], to = rawTime[length(rawTime)], by = 1/sf) + accelRes = GGIRread::resample(rawAccel, rawTime, timeRes, nrow(rawAccel), interpolationType) # this is now the resampled acceleration data + colnamesP = colnames(P)[-which(colnames(P) == "time")] + P = as.data.frame(accelRes, stringsAsFactors = FALSE) colnames(P) = colnamesP - P$timestamp = as.POSIXct(as.numeric(P$timestamp), - tz = desiredtz, origin = "1970-01-01") + P$time = timeRes } return(list(data = P, header = header, PreviousLastValue = PreviousLastValue, diff --git a/inst/testfiles/ax3_testfile_unix_timestamps.csv b/inst/testfiles/ax3_testfile_unix_timestamps.csv new file mode 100644 index 000000000..ec79b182b --- /dev/null +++ b/inst/testfiles/ax3_testfile_unix_timestamps.csv @@ -0,0 +1,210 @@ +1551178506.0000,0.328125,0.984375,0.203125 +1551178506.0100,0.828125,-0.359375,-0.375000 +1551178506.0200,0.875000,-0.390625,-0.390625 +1551178506.0300,0.890625,-0.390625,-0.390625 +1551178506.0400,0.890625,-0.375000,-0.390625 +1551178506.0500,0.875000,-0.359375,-0.375000 +1551178506.0600,0.875000,-0.359375,-0.390625 +1551178506.0700,0.875000,-0.359375,-0.406250 +1551178506.0800,0.875000,-0.359375,-0.406250 +1551178506.0900,0.890625,-0.343750,-0.406250 +1551178506.1000,0.890625,-0.343750,-0.421875 +1551178506.1100,0.890625,-0.359375,-0.421875 +1551178506.1200,0.890625,-0.359375,-0.421875 +1551178506.1300,0.906250,-0.359375,-0.437500 +1551178506.1400,0.906250,-0.359375,-0.421875 +1551178506.1500,0.906250,-0.343750,-0.421875 +1551178506.1600,0.890625,-0.343750,-0.421875 +1551178506.1700,0.890625,-0.343750,-0.421875 +1551178506.1800,0.875000,-0.343750,-0.421875 +1551178506.1900,0.875000,-0.359375,-0.421875 +1551178506.2000,0.875000,-0.359375,-0.421875 +1551178506.2100,0.875000,-0.359375,-0.406250 +1551178506.2300,0.875000,-0.375000,-0.406250 +1551178506.2400,0.890625,-0.390625,-0.421875 +1551178506.2500,0.906250,-0.375000,-0.421875 +1551178506.2600,0.906250,-0.343750,-0.421875 +1551178506.2700,0.890625,-0.328125,-0.437500 +1551178506.2800,0.875000,-0.312500,-0.437500 +1551178506.2900,0.890625,-0.328125,-0.421875 +1551178506.3000,0.890625,-0.312500,-0.421875 +1551178506.3100,0.906250,-0.328125,-0.421875 +1551178506.3200,0.906250,-0.328125,-0.421875 +1551178506.3300,0.906250,-0.343750,-0.421875 +1551178506.3400,0.906250,-0.343750,-0.421875 +1551178506.3500,0.921875,-0.343750,-0.437500 +1551178506.3600,0.906250,-0.343750,-0.437500 +1551178506.3700,0.906250,-0.343750,-0.437500 +1551178506.3800,0.906250,-0.359375,-0.421875 +1551178506.3900,0.906250,-0.359375,-0.421875 +1551178506.4000,0.890625,-0.343750,-0.406250 +1551178506.4100,0.890625,-0.343750,-0.406250 +1551178506.4200,0.890625,-0.328125,-0.406250 +1551178506.4300,0.890625,-0.328125,-0.421875 +1551178506.4400,0.906250,-0.328125,-0.421875 +1551178506.4500,0.906250,-0.312500,-0.437500 +1551178506.4600,0.906250,-0.328125,-0.437500 +1551178506.4700,0.906250,-0.312500,-0.437500 +1551178506.4800,0.906250,-0.328125,-0.437500 +1551178506.4900,0.906250,-0.328125,-0.437500 +1551178506.5000,0.906250,-0.328125,-0.421875 +1551178506.5100,0.890625,-0.343750,-0.406250 +1551178506.5200,0.890625,-0.343750,-0.421875 +1551178506.5300,0.890625,-0.343750,-0.406250 +1551178506.5400,0.890625,-0.328125,-0.406250 +1551178506.5500,0.890625,-0.328125,-0.406250 +1551178506.5600,0.875000,-0.343750,-0.390625 +1551178506.5700,0.890625,-0.343750,-0.406250 +1551178506.5800,0.890625,-0.343750,-0.421875 +1551178506.5900,0.906250,-0.359375,-0.437500 +1551178506.6000,0.906250,-0.359375,-0.437500 +1551178506.6100,0.906250,-0.343750,-0.421875 +1551178506.6200,0.906250,-0.328125,-0.421875 +1551178506.6300,0.906250,-0.328125,-0.421875 +1551178506.6400,0.906250,-0.328125,-0.421875 +1551178506.6500,0.890625,-0.312500,-0.421875 +1551178506.6600,0.890625,-0.328125,-0.437500 +1551178506.6700,0.890625,-0.312500,-0.437500 +1551178506.6800,0.875000,-0.296875,-0.437500 +1551178506.6900,0.890625,-0.312500,-0.453125 +1551178506.7000,0.875000,-0.312500,-0.468750 +1551178506.7100,0.875000,-0.328125,-0.500000 +1551178506.7200,0.859375,-0.343750,-0.500000 +1551178506.7300,0.843750,-0.343750,-0.500000 +1551178506.7400,0.843750,-0.359375,-0.500000 +1551178506.7500,0.859375,-0.359375,-0.515625 +1551178506.7600,0.843750,-0.343750,-0.515625 +1551178506.7700,0.843750,-0.359375,-0.515625 +1551178506.7800,0.875000,-0.359375,-0.531250 +1551178506.7900,0.890625,-0.375000,-0.500000 +1551178506.8000,0.843750,-0.406250,-0.500000 +1551178506.8100,0.843750,-0.359375,-0.453125 +1551178506.8200,0.765625,-0.359375,-0.453125 +1551178506.8300,0.796875,-0.390625,-0.500000 +1551178506.8400,0.781250,-0.421875,-0.484375 +1551178506.8500,0.750000,-0.421875,-0.500000 +1551178506.8600,0.734375,-0.406250,-0.500000 +1551178506.8700,0.703125,-0.375000,-0.531250 +1551178506.8800,0.734375,-0.375000,-0.546875 +1551178506.8900,0.812500,-0.359375,-0.546875 +1551178506.9000,0.843750,-0.343750,-0.546875 +1551178506.9100,0.812500,-0.343750,-0.531250 +1551178506.9200,0.750000,-0.359375,-0.500000 +1551178506.9300,0.687500,-0.359375,-0.500000 +1551178506.9400,0.671875,-0.343750,-0.500000 +1551178506.9500,0.718750,-0.312500,-0.515625 +1551178506.9600,0.750000,-0.343750,-0.515625 +1551178506.9700,0.781250,-0.359375,-0.546875 +1551178506.9800,0.875000,-0.328125,-0.578125 +1551178506.9900,0.921875,-0.343750,-0.609375 +1551178507.0000,0.906250,-0.359375,-0.656250 +1551178507.0100,0.859375,-0.375000,-0.671875 +1551178507.0200,0.812500,-0.375000,-0.703125 +1551178507.0300,0.781250,-0.359375,-0.734375 +1551178507.0400,0.781250,-0.343750,-0.734375 +1551178507.0500,0.796875,-0.328125,-0.718750 +1551178507.0600,0.812500,-0.328125,-0.703125 +1551178507.0700,0.812500,-0.328125,-0.687500 +1551178507.0800,0.796875,-0.343750,-0.671875 +1551178507.0900,0.781250,-0.343750,-0.656250 +1551178507.1000,0.781250,-0.328125,-0.640625 +1551178507.1100,0.781250,-0.328125,-0.625000 +1551178507.1200,0.765625,-0.312500,-0.625000 +1551178507.1300,0.765625,-0.328125,-0.625000 +1551178507.1400,0.781250,-0.343750,-0.640625 +1551178507.1500,0.796875,-0.359375,-0.656250 +1551178507.1600,0.843750,-0.359375,-0.625000 +1551178507.1700,0.843750,-0.375000,-0.609375 +1551178507.1800,0.828125,-0.359375,-0.609375 +1551178507.1900,0.796875,-0.328125,-0.593750 +1551178507.2000,0.765625,-0.296875,-0.578125 +1551178507.2101,0.734375,-0.296875,-0.578125 +1551178507.2202,0.765625,-0.312500,-0.593750 +1551178507.2303,0.750000,-0.312500,-0.562500 +1551178507.2403,0.703125,-0.312500,-0.593750 +1551178507.2504,0.718750,-0.343750,-0.640625 +1551178507.2605,0.734375,-0.375000,-0.671875 +1551178507.2706,0.750000,-0.375000,-0.671875 +1551178507.2807,0.750000,-0.359375,-0.687500 +1551178507.2908,0.734375,-0.343750,-0.687500 +1551178507.3008,0.734375,-0.343750,-0.687500 +1551178507.3109,0.750000,-0.328125,-0.703125 +1551178507.3210,0.765625,-0.296875,-0.703125 +1551178507.3311,0.781250,-0.296875,-0.703125 +1551178507.3412,0.765625,-0.312500,-0.687500 +1551178507.3513,0.750000,-0.328125,-0.687500 +1551178507.3613,0.765625,-0.328125,-0.687500 +1551178507.3714,0.765625,-0.328125,-0.687500 +1551178507.3815,0.750000,-0.312500,-0.703125 +1551178507.3916,0.750000,-0.312500,-0.703125 +1551178507.4017,0.750000,-0.312500,-0.718750 +1551178507.4118,0.750000,-0.328125,-0.718750 +1551178507.4218,0.765625,-0.343750,-0.734375 +1551178507.4319,0.765625,-0.343750,-0.734375 +1551178507.4420,0.765625,-0.343750,-0.718750 +1551178507.4521,0.750000,-0.328125,-0.703125 +1551178507.4622,0.750000,-0.328125,-0.703125 +1551178507.4722,0.734375,-0.312500,-0.687500 +1551178507.4823,0.734375,-0.312500,-0.687500 +1551178507.4924,0.718750,-0.312500,-0.671875 +1551178507.5025,0.734375,-0.312500,-0.671875 +1551178507.5126,0.734375,-0.312500,-0.671875 +1551178507.5227,0.734375,-0.328125,-0.671875 +1551178507.5328,0.734375,-0.343750,-0.671875 +1551178507.5428,0.750000,-0.328125,-0.656250 +1551178507.5529,0.750000,-0.312500,-0.671875 +1551178507.5630,0.750000,-0.296875,-0.671875 +1551178507.5731,0.750000,-0.296875,-0.687500 +1551178507.5832,0.765625,-0.312500,-0.703125 +1551178507.5933,0.765625,-0.328125,-0.703125 +1551178507.6033,0.765625,-0.343750,-0.703125 +1551178507.6134,0.781250,-0.343750,-0.687500 +1551178507.6235,0.781250,-0.328125,-0.687500 +1551178507.6336,0.765625,-0.312500,-0.671875 +1551178507.6437,0.765625,-0.312500,-0.656250 +1551178507.6537,0.765625,-0.312500,-0.656250 +1551178507.6638,0.765625,-0.312500,-0.656250 +1551178507.6739,0.765625,-0.328125,-0.656250 +1551178507.6840,0.765625,-0.328125,-0.656250 +1551178507.6941,0.765625,-0.343750,-0.671875 +1551178507.7042,0.781250,-0.359375,-0.671875 +1551178507.7143,0.796875,-0.328125,-0.687500 +1551178507.7243,0.781250,-0.312500,-0.687500 +1551178507.7344,0.765625,-0.328125,-0.703125 +1551178507.7445,0.765625,-0.375000,-0.703125 +1551178507.7546,0.781250,-0.453125,-0.703125 +1551178507.7647,0.796875,-0.484375,-0.687500 +1551178507.7747,0.812500,-0.484375,-0.656250 +1551178507.7848,0.828125,-0.437500,-0.625000 +1551178507.7949,0.828125,-0.390625,-0.609375 +1551178507.8050,0.828125,-0.375000,-0.593750 +1551178507.8151,0.828125,-0.390625,-0.578125 +1551178507.8252,0.828125,-0.437500,-0.546875 +1551178507.8353,0.828125,-0.468750,-0.531250 +1551178507.8453,0.828125,-0.437500,-0.546875 +1551178507.8554,0.843750,-0.406250,-0.546875 +1551178507.8655,0.843750,-0.359375,-0.562500 +1551178507.8756,0.828125,-0.359375,-0.562500 +1551178507.8857,0.812500,-0.343750,-0.562500 +1551178507.8958,0.828125,-0.359375,-0.562500 +1551178507.9058,0.828125,-0.359375,-0.546875 +1551178507.9159,0.843750,-0.343750,-0.546875 +1551178507.9260,0.859375,-0.312500,-0.562500 +1551178507.9361,0.843750,-0.312500,-0.593750 +1551178507.9462,0.859375,-0.343750,-0.609375 +1551178507.9563,0.859375,-0.390625,-0.593750 +1551178507.9663,0.875000,-0.390625,-0.562500 +1551178507.9764,0.890625,-0.343750,-0.531250 +1551178507.9865,0.890625,-0.296875,-0.515625 +1551178507.9966,0.906250,-0.296875,-0.531250 +1551178508.0067,0.890625,-0.343750,-0.546875 +1551178508.0168,0.875000,-0.390625,-0.578125 +1551178508.0268,0.843750,-0.390625,-0.593750 +1551178508.0369,0.828125,-0.343750,-0.578125 +1551178508.0470,0.828125,-0.296875,-0.546875 +1551178508.0571,0.859375,-0.281250,-0.531250 +1551178508.0672,0.875000,-0.296875,-0.515625 +1551178508.0773,0.890625,-0.343750,-0.515625 +1551178508.0873,0.890625,-0.375000,-0.500000 +1551178508.0974,0.875000,-0.406250,-0.484375 +1551178508.1075,0.812500,-0.390625,-0.453125 \ No newline at end of file diff --git a/inst/testfiles/ax3test.wav b/inst/testfiles/ax3test.wav deleted file mode 100755 index 7a2e85448..000000000 Binary files a/inst/testfiles/ax3test.wav and /dev/null differ diff --git a/inst/testfiles/ax6_testfile_formatted_timestamps.csv b/inst/testfiles/ax6_testfile_formatted_timestamps.csv new file mode 100755 index 000000000..7acda6793 --- /dev/null +++ b/inst/testfiles/ax6_testfile_formatted_timestamps.csv @@ -0,0 +1,210 @@ +2019-12-23 21:04:06.690,0.007324,0.071289,0.008789,0.274658,-0.503540,15.769958 +2019-12-23 21:04:06.699,0.001953,0.066406,0.007813,0.282288,-0.480652,15.792846 +2019-12-23 21:04:06.710,0.007813,0.068359,0.001953,0.274658,-0.503540,15.769958 +2019-12-23 21:04:06.720,0.007813,0.078125,0.005859,0.305176,-0.488281,15.754699 +2019-12-23 21:04:06.730,0.004883,0.073242,0.009766,0.282288,-0.511169,15.762328 +2019-12-23 21:04:06.739,0.000977,0.068848,0.008301,0.274658,-0.549316,15.762328 +2019-12-23 21:04:06.750,0.001953,0.069336,0.003418,0.312805,-0.495911,15.769958 +2019-12-23 21:04:06.760,0.000977,0.062012,0.004883,0.320435,-0.473022,15.731811 +2019-12-23 21:04:06.770,-0.007324,0.062500,0.008301,0.289917,-0.518799,15.762328 +2019-12-23 21:04:06.779,-0.017090,0.061035,0.009277,0.328064,-0.503540,15.792846 +2019-12-23 21:04:06.789,-0.007324,0.061035,0.009277,0.312805,-0.526428,15.777587 +2019-12-23 21:04:06.800,0.005859,0.064941,0.009766,0.259399,-0.511169,15.739440 +2019-12-23 21:04:06.810,0.010742,0.067383,0.008301,0.289917,-0.518799,15.708922 +2019-12-23 21:04:06.820,0.002930,0.071777,0.007813,0.320435,-0.503540,15.762328 +2019-12-23 21:04:06.829,-0.003418,0.075195,0.006348,0.267029,-0.511169,15.800475 +2019-12-23 21:04:06.840,-0.003418,0.077148,0.004883,0.312805,-0.503540,15.853881 +2019-12-23 21:04:06.850,-0.000488,0.076172,0.004883,0.328064,-0.488281,15.830993 +2019-12-23 21:04:06.860,0.000000,0.068359,0.006348,0.350952,-0.495911,15.792846 +2019-12-23 21:04:06.869,0.000977,0.072754,0.011719,0.320435,-0.549316,15.739440 +2019-12-23 21:04:06.880,0.000000,0.074219,0.007324,0.297546,-0.556946,15.762328 +2019-12-23 21:04:06.890,-0.000488,0.076660,0.002441,0.282288,-0.511169,15.777587 +2019-12-23 21:04:06.900,-0.002441,0.078613,0.008789,0.320435,-0.480652,15.769958 +2019-12-23 21:04:06.909,-0.009766,0.078125,0.011719,0.289917,-0.511169,15.762328 +2019-12-23 21:04:06.920,0.000488,0.074707,0.004883,0.305176,-0.480652,15.823363 +2019-12-23 21:04:06.930,-0.005859,0.067383,0.006348,0.312805,-0.488281,15.724181 +2019-12-23 21:04:06.940,-0.011719,0.063477,0.006836,0.305176,-0.518799,15.739440 +2019-12-23 21:04:06.949,-0.004395,0.066406,0.004395,0.282288,-0.511169,15.777587 +2019-12-23 21:04:06.960,-0.002930,0.070313,0.004883,0.289917,-0.511169,15.777587 +2019-12-23 21:04:06.970,0.008301,0.064453,0.009766,0.320435,-0.511169,15.716552 +2019-12-23 21:04:06.980,0.008301,0.065430,0.008789,0.305176,-0.518799,15.777587 +2019-12-23 21:04:06.990,0.008789,0.070313,0.006348,0.274658,-0.495911,15.815734 +2019-12-23 21:04:07.000,0.004395,0.069336,0.003906,0.259399,-0.473022,15.777587 +2019-12-23 21:04:07.010,-0.000977,0.067383,0.007813,0.274658,-0.518799,15.754699 +2019-12-23 21:04:07.020,0.003906,0.071289,0.010742,0.289917,-0.541687,15.769958 +2019-12-23 21:04:07.030,0.000977,0.070313,0.004395,0.312805,-0.511169,15.762328 +2019-12-23 21:04:07.039,-0.004883,0.066406,0.001465,0.305176,-0.518799,15.762328 +2019-12-23 21:04:07.050,-0.006348,0.071777,0.008789,0.328064,-0.518799,15.747069 +2019-12-23 21:04:07.060,-0.002930,0.072754,0.008789,0.328064,-0.503540,15.754699 +2019-12-23 21:04:07.070,-0.004395,0.073730,0.006348,0.320435,-0.518799,15.747069 +2019-12-23 21:04:07.079,-0.006836,0.072266,0.003418,0.282288,-0.488281,15.724181 +2019-12-23 21:04:07.090,-0.000977,0.070313,0.008301,0.267029,-0.503540,15.769958 +2019-12-23 21:04:07.100,0.002441,0.075195,0.006348,0.305176,-0.526428,15.785216 +2019-12-23 21:04:07.110,-0.001953,0.069824,0.006348,0.320435,-0.473022,15.800475 +2019-12-23 21:04:07.120,-0.005371,0.069336,0.009277,0.305176,-0.503540,15.777587 +2019-12-23 21:04:07.131,-0.006348,0.068848,0.005371,0.320435,-0.518799,15.731811 +2019-12-23 21:04:07.141,-0.008301,0.062012,0.007813,0.259399,-0.511169,15.747069 +2019-12-23 21:04:07.151,-0.000977,0.063965,0.003906,0.305176,-0.518799,15.785216 +2019-12-23 21:04:07.161,0.004395,0.078125,0.003418,0.289917,-0.503540,15.785216 +2019-12-23 21:04:07.171,0.007813,0.081055,0.012207,0.297546,-0.511169,15.762328 +2019-12-23 21:04:07.182,0.005371,0.074219,0.007813,0.312805,-0.495911,15.792846 +2019-12-23 21:04:07.192,-0.013184,0.075195,0.005859,0.366211,-0.495911,15.792846 +2019-12-23 21:04:07.202,-0.020020,0.078125,0.005371,0.335693,-0.495911,15.777587 +2019-12-23 21:04:07.213,-0.027832,0.075195,0.004395,0.282288,-0.503540,15.777587 +2019-12-23 21:04:07.223,-0.027832,0.071777,0.009277,0.305176,-0.511169,15.769958 +2019-12-23 21:04:07.233,-0.023926,0.070313,0.005371,0.297546,-0.534058,15.769958 +2019-12-23 21:04:07.243,-0.020020,0.064941,0.002930,0.274658,-0.503540,15.808105 +2019-12-23 21:04:07.254,-0.020508,0.068359,0.000977,0.305176,-0.518799,15.785216 +2019-12-23 21:04:07.264,-0.012207,0.073242,0.002441,0.305176,-0.495911,15.762328 +2019-12-23 21:04:07.274,-0.003906,0.078613,0.007324,0.328064,-0.518799,15.754699 +2019-12-23 21:04:07.284,0.007324,0.075684,0.007324,0.305176,-0.534058,15.769958 +2019-12-23 21:04:07.295,0.012695,0.079102,0.008301,0.282288,-0.526428,15.777587 +2019-12-23 21:04:07.305,0.019043,0.082520,0.012207,0.320435,-0.495911,15.762328 +2019-12-23 21:04:07.315,0.020508,0.076660,0.008789,0.274658,-0.526428,15.747069 +2019-12-23 21:04:07.325,0.014648,0.071289,0.009277,0.343323,-0.488281,15.708922 +2019-12-23 21:04:07.335,0.007813,0.068848,0.008301,0.312805,-0.480652,15.792846 +2019-12-23 21:04:07.346,-0.005371,0.064941,0.005859,0.289917,-0.511169,15.754699 +2019-12-23 21:04:07.356,-0.001465,0.063965,0.002930,0.305176,-0.541687,15.754699 +2019-12-23 21:04:07.366,0.005859,0.062500,0.008789,0.305176,-0.511169,15.754699 +2019-12-23 21:04:07.376,-0.004883,0.060547,0.002930,0.282288,-0.488281,15.762328 +2019-12-23 21:04:07.387,-0.011230,0.063477,0.000488,0.305176,-0.518799,15.754699 +2019-12-23 21:04:07.397,-0.001953,0.061523,0.004395,0.297546,-0.534058,15.792846 +2019-12-23 21:04:07.407,-0.004883,0.064941,0.003906,0.297546,-0.503540,15.792846 +2019-12-23 21:04:07.418,0.002930,0.064941,0.002930,0.282288,-0.488281,15.762328 +2019-12-23 21:04:07.428,0.006348,0.068848,0.011719,0.274658,-0.503540,15.716552 +2019-12-23 21:04:07.438,0.003906,0.071289,0.011719,0.320435,-0.495911,15.739440 +2019-12-23 21:04:07.448,0.004883,0.074219,0.013184,0.305176,-0.495911,15.785216 +2019-12-23 21:04:07.459,0.009766,0.076172,0.011230,0.305176,-0.495911,15.747069 +2019-12-23 21:04:07.469,0.005371,0.076660,0.005859,0.312805,-0.511169,15.754699 +2019-12-23 21:04:07.479,-0.004395,0.068359,0.004395,0.312805,-0.503540,15.777587 +2019-12-23 21:04:07.489,-0.004883,0.071289,0.008789,0.343323,-0.488281,15.830993 +2019-12-23 21:04:07.500,-0.015137,0.073730,0.008301,0.335693,-0.503540,15.792846 +2019-12-23 21:04:07.510,-0.015137,0.072754,0.010742,0.320435,-0.480652,15.754699 +2019-12-23 21:04:07.520,-0.014648,0.076660,0.011230,0.320435,-0.503540,15.769958 +2019-12-23 21:04:07.529,-0.010742,0.074219,0.002441,0.312805,-0.518799,15.785216 +2019-12-23 21:04:07.539,0.000000,0.077637,0.007324,0.328064,-0.518799,15.754699 +2019-12-23 21:04:07.550,-0.002930,0.076172,0.008789,0.305176,-0.511169,15.777587 +2019-12-23 21:04:07.559,-0.006836,0.071777,0.007813,0.305176,-0.495911,15.800475 +2019-12-23 21:04:07.569,-0.004883,0.073730,0.006348,0.297546,-0.526428,15.762328 +2019-12-23 21:04:07.579,0.000488,0.073730,0.007813,0.297546,-0.526428,15.762328 +2019-12-23 21:04:07.590,-0.002441,0.071289,0.012207,0.267029,-0.503540,15.792846 +2019-12-23 21:04:07.599,-0.003418,0.066895,0.012695,0.289917,-0.518799,15.769958 +2019-12-23 21:04:07.610,-0.005859,0.061035,0.012207,0.274658,-0.480652,15.800475 +2019-12-23 21:04:07.619,-0.012695,0.063477,0.000977,0.274658,-0.488281,15.815734 +2019-12-23 21:04:07.630,-0.009277,0.059570,0.006836,0.297546,-0.511169,15.785216 +2019-12-23 21:04:07.640,-0.004395,0.065918,0.011230,0.328064,-0.526428,15.792846 +2019-12-23 21:04:07.650,-0.005859,0.066406,0.009766,0.289917,-0.549316,15.800475 +2019-12-23 21:04:07.659,-0.004395,0.070801,0.006348,0.282288,-0.526428,15.792846 +2019-12-23 21:04:07.670,-0.002441,0.078125,0.009766,0.289917,-0.518799,15.815734 +2019-12-23 21:04:07.680,0.007324,0.080078,0.011230,0.320435,-0.511169,15.739440 +2019-12-23 21:04:07.690,-0.001953,0.078613,0.011230,0.305176,-0.503540,15.747069 +2019-12-23 21:04:07.699,-0.017578,0.074707,0.007324,0.282288,-0.511169,15.747069 +2019-12-23 21:04:07.710,-0.015137,0.075195,0.003418,0.274658,-0.503540,15.777587 +2019-12-23 21:04:07.720,0.001953,0.074219,0.007813,0.305176,-0.518799,15.792846 +2019-12-23 21:04:07.730,0.005859,0.076172,0.012695,0.350952,-0.503540,15.785216 +2019-12-23 21:04:07.739,0.002930,0.074707,0.010254,0.289917,-0.495911,15.769958 +2019-12-23 21:04:07.750,-0.007813,0.077637,0.011719,0.282288,-0.488281,15.777587 +2019-12-23 21:04:07.760,-0.016113,0.074219,0.004395,0.328064,-0.473022,15.823363 +2019-12-23 21:04:07.770,-0.005859,0.071777,0.002930,0.335693,-0.495911,15.800475 +2019-12-23 21:04:07.779,0.001465,0.071777,0.009766,0.335693,-0.541687,15.754699 +2019-12-23 21:04:07.789,0.006348,0.070801,0.009766,0.267029,-0.518799,15.731811 +2019-12-23 21:04:07.800,-0.001953,0.069336,0.002930,0.297546,-0.526428,15.769958 +2019-12-23 21:04:07.810,-0.001953,0.064453,0.007813,0.251770,-0.526428,15.800475 +2019-12-23 21:04:07.820,-0.000488,0.062988,0.005371,0.274658,-0.503540,15.731811 +2019-12-23 21:04:07.829,-0.003906,0.058105,0.007813,0.297546,-0.511169,15.754699 +2019-12-23 21:04:07.840,0.001953,0.054688,0.005371,0.274658,-0.526428,15.762328 +2019-12-23 21:04:07.850,0.003418,0.054688,0.007813,0.335693,-0.473022,15.716552 +2019-12-23 21:04:07.860,0.008301,0.055664,0.005371,0.350952,-0.511169,15.754699 +2019-12-23 21:04:07.869,0.004395,0.061035,0.007324,0.297546,-0.556946,15.762328 +2019-12-23 21:04:07.880,-0.001953,0.070801,0.001953,0.282288,-0.518799,15.777587 +2019-12-23 21:04:07.890,-0.005859,0.076660,0.005859,0.305176,-0.503540,15.769958 +2019-12-23 21:04:07.900,-0.004883,0.078125,0.006836,0.297546,-0.495911,15.769958 +2019-12-23 21:04:07.910,-0.008301,0.074707,0.005859,0.320435,-0.503540,15.800475 +2019-12-23 21:04:07.920,-0.019043,0.072754,0.005371,0.312805,-0.511169,15.823363 +2019-12-23 21:04:07.930,-0.024414,0.072266,0.014648,0.328064,-0.549316,15.724181 +2019-12-23 21:04:07.941,-0.025391,0.077637,0.012695,0.282288,-0.549316,15.800475 +2019-12-23 21:04:07.951,-0.023926,0.071777,0.007324,0.328064,-0.488281,15.800475 +2019-12-23 21:04:07.961,-0.022949,0.074219,0.005859,0.320435,-0.518799,15.769958 +2019-12-23 21:04:07.971,-0.006348,0.084961,0.008301,0.312805,-0.526428,15.731811 +2019-12-23 21:04:07.982,-0.004395,0.077637,0.007813,0.267029,-0.495911,15.785216 +2019-12-23 21:04:07.992,-0.014648,0.069824,0.011230,0.282288,-0.518799,15.830993 +2019-12-23 21:04:08.002,-0.009766,0.072754,0.007813,0.289917,-0.518799,15.792846 +2019-12-23 21:04:08.012,0.005371,0.075684,0.006348,0.305176,-0.511169,15.754699 +2019-12-23 21:04:08.022,0.004395,0.072754,0.005371,0.343323,-0.526428,15.800475 +2019-12-23 21:04:08.033,0.003906,0.069824,0.004883,0.289917,-0.518799,15.815734 +2019-12-23 21:04:08.043,0.002930,0.066895,0.007324,0.289917,-0.503540,15.747069 +2019-12-23 21:04:08.053,-0.009277,0.063477,0.009277,0.320435,-0.503540,15.769958 +2019-12-23 21:04:08.064,-0.006348,0.064453,0.009277,0.328064,-0.511169,15.823363 +2019-12-23 21:04:08.074,-0.008301,0.064941,0.008301,0.320435,-0.534058,15.800475 +2019-12-23 21:04:08.084,-0.011719,0.063965,0.006348,0.297546,-0.534058,15.762328 +2019-12-23 21:04:08.094,-0.011719,0.062988,0.004883,0.274658,-0.511169,15.769958 +2019-12-23 21:04:08.104,-0.000488,0.071777,0.000488,0.328064,-0.511169,15.762328 +2019-12-23 21:04:08.115,0.004883,0.075684,0.000488,0.328064,-0.518799,15.769958 +2019-12-23 21:04:08.125,0.012207,0.081055,0.007813,0.312805,-0.534058,15.777587 +2019-12-23 21:04:08.135,0.011230,0.077148,0.008789,0.328064,-0.526428,15.800475 +2019-12-23 21:04:08.145,0.008301,0.071777,0.010254,0.305176,-0.518799,15.777587 +2019-12-23 21:04:08.156,0.015625,0.075684,0.007324,0.297546,-0.534058,15.800475 +2019-12-23 21:04:08.166,0.013672,0.073242,0.005371,0.297546,-0.541687,15.800475 +2019-12-23 21:04:08.176,0.004395,0.066406,0.007324,0.282288,-0.511169,15.754699 +2019-12-23 21:04:08.187,0.004395,0.067383,0.004395,0.289917,-0.534058,15.739440 +2019-12-23 21:04:08.197,0.005371,0.068359,0.006348,0.328064,-0.526428,15.724181 +2019-12-23 21:04:08.207,-0.000488,0.066406,0.008789,0.320435,-0.526428,15.777587 +2019-12-23 21:04:08.217,-0.003906,0.065430,0.005859,0.274658,-0.518799,15.769958 +2019-12-23 21:04:08.227,-0.009766,0.062012,0.006836,0.289917,-0.495911,15.792846 +2019-12-23 21:04:08.238,0.006348,0.067383,0.002930,0.297546,-0.511169,15.769958 +2019-12-23 21:04:08.248,0.006836,0.071289,0.002930,0.320435,-0.534058,15.769958 +2019-12-23 21:04:08.258,0.003418,0.066895,0.008301,0.282288,-0.534058,15.747069 +2019-12-23 21:04:08.269,-0.000488,0.069336,0.008301,0.312805,-0.541687,15.747069 +2019-12-23 21:04:08.279,-0.005859,0.074219,0.007813,0.282288,-0.518799,15.808105 +2019-12-23 21:04:08.289,-0.008789,0.076660,0.010742,0.305176,-0.488281,15.838622 +2019-12-23 21:04:08.299,-0.017578,0.075684,0.003418,0.343323,-0.526428,15.808105 +2019-12-23 21:04:08.309,-0.018555,0.072754,0.007813,0.328064,-0.534058,15.769958 +2019-12-23 21:04:08.319,-0.016113,0.070801,0.008301,0.289917,-0.511169,15.769958 +2019-12-23 21:04:08.329,-0.007324,0.076660,0.005859,0.305176,-0.518799,15.754699 +2019-12-23 21:04:08.340,0.001465,0.077637,0.010254,0.282288,-0.549316,15.769958 +2019-12-23 21:04:08.350,0.006348,0.080566,0.008301,0.312805,-0.518799,15.716552 +2019-12-23 21:04:08.359,0.001465,0.075195,0.007324,0.305176,-0.518799,15.823363 +2019-12-23 21:04:08.369,-0.013672,0.067871,0.006836,0.305176,-0.488281,15.830993 +2019-12-23 21:04:08.380,-0.031738,0.062500,0.007813,0.289917,-0.488281,15.769958 +2019-12-23 21:04:08.389,-0.036621,0.065430,0.005371,0.289917,-0.495911,15.769958 +2019-12-23 21:04:08.399,-0.036621,0.059082,0.002930,0.282288,-0.518799,15.769958 +2019-12-23 21:04:08.409,-0.040527,0.053711,0.005371,0.335693,-0.534058,15.823363 +2019-12-23 21:04:08.420,-0.027344,0.061523,0.003418,0.312805,-0.556946,15.800475 +2019-12-23 21:04:08.430,-0.014160,0.059570,0.008301,0.282288,-0.503540,15.800475 +2019-12-23 21:04:08.440,-0.015137,0.066895,0.011719,0.320435,-0.534058,15.808105 +2019-12-23 21:04:08.449,-0.024414,0.065430,0.004395,0.282288,-0.534058,15.777587 +2019-12-23 21:04:08.460,-0.027344,0.062500,0.002441,0.289917,-0.511169,15.747069 +2019-12-23 21:04:08.470,-0.025391,0.069824,0.009766,0.305176,-0.511169,15.747069 +2019-12-23 21:04:08.479,-0.017090,0.068848,0.007813,0.305176,-0.511169,15.762328 +2019-12-23 21:04:08.489,-0.012207,0.073242,0.001953,0.297546,-0.526428,15.838622 +2019-12-23 21:04:08.500,0.006836,0.075684,0.002441,0.328064,-0.564575,15.800475 +2019-12-23 21:04:08.510,0.019043,0.077148,0.005371,0.305176,-0.526428,15.777587 +2019-12-23 21:04:08.520,0.023438,0.077148,0.006348,0.312805,-0.549316,15.792846 +2019-12-23 21:04:08.529,0.024902,0.071289,0.008301,0.328064,-0.526428,15.800475 +2019-12-23 21:04:08.539,0.028320,0.071289,0.008301,0.328064,-0.503540,15.762328 +2019-12-23 21:04:08.550,0.019043,0.068359,0.010254,0.328064,-0.518799,15.747069 +2019-12-23 21:04:08.559,0.009766,0.060059,0.010254,0.282288,-0.518799,15.777587 +2019-12-23 21:04:08.569,0.000488,0.055664,0.010254,0.274658,-0.526428,15.754699 +2019-12-23 21:04:08.579,-0.006348,0.052734,0.004883,0.282288,-0.534058,15.754699 +2019-12-23 21:04:08.590,-0.006348,0.060547,0.004883,0.305176,-0.526428,15.800475 +2019-12-23 21:04:08.600,0.002930,0.066406,0.001953,0.320435,-0.511169,15.846251 +2019-12-23 21:04:08.610,0.002930,0.068848,0.008301,0.335693,-0.526428,15.754699 +2019-12-23 21:04:08.619,-0.003418,0.075684,0.010254,0.335693,-0.503540,15.777587 +2019-12-23 21:04:08.630,-0.006348,0.079590,0.005371,0.305176,-0.534058,15.785216 +2019-12-23 21:04:08.640,-0.005371,0.083008,0.006348,0.320435,-0.526428,15.769958 +2019-12-23 21:04:08.649,-0.006348,0.084473,0.010254,0.305176,-0.495911,15.792846 +2019-12-23 21:04:08.659,-0.003906,0.087891,0.011230,0.259399,-0.495911,15.739440 +2019-12-23 21:04:08.670,-0.000488,0.084473,0.008789,0.320435,-0.495911,15.731811 +2019-12-23 21:04:08.680,0.003418,0.079590,0.005371,0.289917,-0.511169,15.769958 +2019-12-23 21:04:08.690,0.003418,0.079102,0.014160,0.305176,-0.518799,15.792846 +2019-12-23 21:04:08.699,-0.005371,0.069336,0.001465,0.320435,-0.526428,15.785216 +2019-12-23 21:04:08.710,-0.005371,0.066895,0.004395,0.305176,-0.518799,15.769958 +2019-12-23 21:04:08.720,-0.001465,0.066406,0.006348,0.335693,-0.549316,15.785216 +2019-12-23 21:04:08.729,-0.006348,0.064941,0.006348,0.328064,-0.541687,15.800475 +2019-12-23 21:04:08.739,-0.001953,0.061035,0.009766,0.305176,-0.503540,15.777587 +2019-12-23 21:04:08.750,0.005859,0.064941,0.000000,0.312805,-0.518799,15.769958 +2019-12-23 21:04:08.760,0.011230,0.062988,0.002930,0.282288,-0.534058,15.785216 +2019-12-23 21:04:08.770,0.000488,0.069336,0.008301,0.312805,-0.534058,15.785216 +2019-12-23 21:04:08.779,-0.004883,0.068848,0.005859,0.297546,-0.526428,15.777587 +2019-12-23 21:04:08.789,-0.013184,0.071777,0.007324,0.305176,-0.511169,15.815734 +2019-12-23 21:04:08.800,-0.011719,0.076172,0.007324,0.320435,-0.511169,15.792846 \ No newline at end of file diff --git a/man/g.applymetrics.Rd b/man/g.applymetrics.Rd index 0bd4a58cf..e98be9889 100644 --- a/man/g.applymetrics.Rd +++ b/man/g.applymetrics.Rd @@ -64,6 +64,7 @@ Gy = runif(n=10000,min=1,max=3) Gz = runif(n=10000,min=0,max=2) data = cbind(Gx, Gy, Gz) + colnames(data) = c("x", "y", "z") metrics2do = data.frame(do.bfen=TRUE,do.enmo=TRUE,do.lfenmo=FALSE, do.en=FALSE,do.hfen=FALSE,do.hfenplus=FALSE,do.mad=FALSE,do.anglex=FALSE, do.angley=FALSE,do.anglez=FALSE,do.roll_med_acc_x=FALSE, diff --git a/man/g.dotorcomma.Rd b/man/g.dotorcomma.Rd index d9da0de49..e3f1b8812 100644 --- a/man/g.dotorcomma.Rd +++ b/man/g.dotorcomma.Rd @@ -9,7 +9,7 @@ should be interpretted } \usage{ - g.dotorcomma(inputfile,dformat,mon, desiredtz = "", loadGENEActiv = "GGIRread", ...) + g.dotorcomma(inputfile, dformat, mon, ...) } \arguments{ \item{inputfile}{ @@ -23,12 +23,6 @@ Monitor code (accelorometer brand): 0=undefined, 1=GENEA, 2=GENEActiv, 3=Actigraph, 4=Axivity, 5=Movisense, 6=Verisense } - \item{desiredtz}{ - Desired timezone, see documentation \link{g.getmeta} - } - \item{loadGENEActiv}{ - See \link{GGIR} - } \item{...}{ Any input arguments needed for function \link{read.myacc.csv} if you are working with a non-standard csv formatted files. diff --git a/man/g.downsample.Rd b/man/g.downsample.Rd deleted file mode 100644 index 4427d3c6a..000000000 --- a/man/g.downsample.Rd +++ /dev/null @@ -1,41 +0,0 @@ -\name{g.downsample} -\alias{g.downsample} -\title{ - Downsample a vector of numeric values at three time resolutions -} -\description{ - Downsamples a vector of numeric values at three time resolutions: - 1 seconds, ws3 seconds, and ws2 second. Function is not intended - for direct interaction by package end user - -} -\usage{ - g.downsample(sig,fs,ws3,ws2) -} -\arguments{ - \item{sig}{ - Vector of numeric values - } - \item{fs}{ - Sample frequency - } - \item{ws3}{ - ws3 epoch size, e.g. 5 seconds - } - \item{ws2}{ - ws2 epoch size, e.g. 90 seconds - } -} -\value{ -List with three object: var1, var2, and var3 -corresponding to downsample time series at -1 seconds, ws2 seconds, and ws3 seconds resoluton, respectively -} -\examples{ - sig = runif(n=10000,min=1,max=10) - downsampled_sig = g.downsample(sig,fs=20,ws3=5,ws2=15) -} -\keyword{internal} -\author{ -Vincent T van Hees -} diff --git a/man/g.getstarttime.Rd b/man/g.getstarttime.Rd index 32bdfddc3..12ace85bc 100644 --- a/man/g.getstarttime.Rd +++ b/man/g.getstarttime.Rd @@ -9,18 +9,15 @@ for direct use by package user } \usage{ - g.getstarttime(datafile, P, header, mon, dformat, desiredtz, + g.getstarttime(datafile, data, mon, dformat, desiredtz, configtz = NULL) } \arguments{ \item{datafile}{ Full path to data file } - \item{P}{ - Object extracted with \link{g.readaccfile} - } - \item{header}{ - File header extracted with \link{g.inspectfile} + \item{data}{ + Data part of \link{g.readaccfile} output } \item{mon}{ Same as in \link{g.dotorcomma} diff --git a/man/g.imputeTimegaps.Rd b/man/g.imputeTimegaps.Rd index 2f42c4b6c..742b2e7ac 100644 --- a/man/g.imputeTimegaps.Rd +++ b/man/g.imputeTimegaps.Rd @@ -8,7 +8,7 @@ imputes time gaps by the last recorded value per axis normalised to 1 _g_ } \usage{ - g.imputeTimegaps(x, xyzCol, timeCol = c(), sf, k = 0.25, impute = TRUE, + g.imputeTimegaps(x, sf, k = 0.25, impute = TRUE, PreviousLastValue = c(0,0,1), PreviousLastTime = NULL, epochsize = NULL) } @@ -17,12 +17,6 @@ Data.frame with raw accelerometer data, and a timestamp column with millisecond resolution. } - \item{xyzCol}{ - Columnnames or numbers for the x, y and z column - } - \item{timeCol}{ - Column name or number for the timestamp column - } \item{sf}{ Sample frequency in Hertz } diff --git a/man/g.readaccfile.Rd b/man/g.readaccfile.Rd index 8f622863a..56519dce9 100644 --- a/man/g.readaccfile.Rd +++ b/man/g.readaccfile.Rd @@ -12,7 +12,7 @@ g.readaccfile(filename, blocksize, blocknumber, filequality, ws, PreviousEndPage = 1, inspectfileobject = c(), PreviousLastValue = c(0,0,1), PreviousLastTime = NULL, - params_rawdata = c(), params_general = c(), ...) + params_rawdata = c(), params_general = c(), header = NULL, ...) } \arguments{ \item{filename}{ @@ -22,7 +22,7 @@ Size of blocks (in file pages) to be read } \item{blocknumber}{ - Block number relative to start of file + Block number relative to start of file, starting with 1. } \item{filequality}{ Single row dataframe with columns: filetooshort, filecorrupt, @@ -50,6 +50,10 @@ \item{params_general}{ See \link{g.part1} } + \item{header}{ + Header information that was extracted the previous time this file was read, + to be re-used instead of being extracted again. + } \item{...}{ Any input arguments needed for function \link{read.myacc.csv} if you are working with a non-standard csv formatted files. Furter, @@ -62,8 +66,7 @@ \item \code{P} Block object extracted from file with format specific to accelerometer brand \item \code{filequality} Same as in function arguments - \item \code{switchoffLD} Boolean to indicate whether it is worth - continueing to read the next block of data or not + \item \code{isLastBlock} Boolean indicating whether this was the last block to read \item \code{endpage} Page number on which blocked ends, this will be used as input for argument PreviousEndPage when reading the next block. } diff --git a/man/g.readtemp_movisens.Rd b/man/g.readtemp_movisens.Rd index a283f985d..27f0bec76 100644 --- a/man/g.readtemp_movisens.Rd +++ b/man/g.readtemp_movisens.Rd @@ -8,8 +8,8 @@ Reads the temperature from movisens files. it to the matrix where accelerations are stored } \usage{ -g.readtemp_movisens(datafile, desiredtz = "", from = c(), to = c(), - interpolationType=1) +g.readtemp_movisens(datafile, from = c(), to = c(), acc_sf, acc_length, + interpolationType=1) } \arguments{ \item{datafile}{ @@ -17,9 +17,6 @@ g.readtemp_movisens(datafile, desiredtz = "", from = c(), to = c(), movisens store a set of bin file in one folder per recording. GGIR will read the pertinent bin file to access to the temperature data. } - \item{desiredtz}{ - See \link{g.getmeta} - } \item{from}{ Origin point to derive the temperature from movisens files (automatically calculated by GGIR) @@ -28,6 +25,12 @@ g.readtemp_movisens(datafile, desiredtz = "", from = c(), to = c(), End point to derive the temperature from movisens files (automatically calculated by GGIR) } + \item{acc_sf}{ + Sample frequency of acceleration data + } + \item{acc_length}{ + number of acceleration data samples + } \item{interpolationType}{ Integer to indicate type of interpolation to be used when resampling time series (mainly relevant for Axivity sensors), 1=linear, 2=nearest neighbour. } @@ -38,6 +41,6 @@ g.readtemp_movisens(datafile, desiredtz = "", from = c(), to = c(), \keyword{internal} \examples{ \dontrun{ - P = g.readtemp_movisens(datafile, desiredtz = "", from = c(), to = c()) + P = g.readtemp_movisens(datafile, from = c(), to = c(), acc_sf = 64, acc_length = 3000) } } diff --git a/man/get_nw_clip_block_params.Rd b/man/get_nw_clip_block_params.Rd index 52308fd30..bf8046ea6 100644 --- a/man/get_nw_clip_block_params.Rd +++ b/man/get_nw_clip_block_params.Rd @@ -9,8 +9,8 @@ Not designed for direct use by user. } \usage{ - get_nw_clip_block_params(chunksize, dynrange, monc, rmc.noise=c(), - sf, dformat, rmc.dynamic_range) + get_nw_clip_block_params(chunksize, dynrange, monc, dformat, deviceSerialNumber, + rmc.noise=c(), sf, rmc.dynamic_range) } \arguments{ \item{chunksize}{ @@ -22,6 +22,12 @@ \item{monc}{ See \link{g.inspectfile} } + \item{dformat}{ + See \link{g.dotorcomma} + } + \item{deviceSerialNumber}{ + As produced by \link{g.extractheadervars} + } \item{rmc.noise}{ Noise level of acceleration signal in _g_-units, used when working ad-hoc .csv data formats using \link{read.myacc.csv}. The \link{read.myacc.csv} does not take rmc.noise as argument, @@ -31,9 +37,6 @@ \item{sf}{ Numeric, sample frequency in Hertz } - \item{dformat}{ - See \link{g.dotorcomma} - } \item{rmc.dynamic_range}{ Optional, please see \link{read.myacc.csv} } diff --git a/man/get_starttime_weekday_meantemp_truncdata.Rd b/man/get_starttime_weekday_meantemp_truncdata.Rd deleted file mode 100644 index a7e6a0838..000000000 --- a/man/get_starttime_weekday_meantemp_truncdata.Rd +++ /dev/null @@ -1,73 +0,0 @@ -\name{get_starttime_weekday_meantemp_truncdata } -\alias{get_starttime_weekday_meantemp_truncdata} -\title{ - Get starttime (adjusted), weekday, mean temp, and adjust data accordingly. -} -\description{ - Function not intended for direct use by user. - Used inside \link{g.getmeta} as an intermediate step between - loading the raw data and calibrating it. This step includes extracting - the starttime and adjusting it to nearest integer number of long epoch window - lengths in an hour, truncating the data accordingly, extracting the - corresponding weekday and mean temperature (if temperature is available). -} -\usage{ - get_starttime_weekday_meantemp_truncdata(temp.available, monc, - dformat, data, P, header, desiredtz, sf, i, - datafile, ws2, starttime, wday, wdayname, configtz = NULL) -} -\arguments{ - \item{temp.available}{ - Boolean whether temperate is available. - } - \item{monc}{ - See \link{g.inspectfile} - } - \item{dformat}{ - See \link{g.dotorcomma} - } - \item{data}{ - Data part of \link{g.readaccfile} output - } - \item{P}{ - data loaded from accelerometer file with \link{g.readaccfile} - } - \item{header}{ - Header part of \link{g.readaccfile} output - } - \item{desiredtz}{ - See \link{g.getmeta} - } - \item{sf}{ - Numeric, sample frequency in Hertz - } - \item{i}{ - Integer index of passed on from \link{g.getmeta} - to indicate what data block is being read. - } - \item{datafile}{ - See \link{g.getmeta} - } - \item{ws2}{ - Long epoch length - } - \item{starttime}{ - Once calculate it is remembered and fed into this function again, - such that it does not have to be recalulated. - } - \item{wday}{ - Once calculate it is remembered and fed into this function again, - such that it does not have to be recalulated. - } - \item{wdayname}{ - Once calculate it is remembered and fed into this function again, - such that it does not have to be recalulated. - } - \item{configtz}{ - See \link{g.getmeta} - } -} -\keyword{internal} -\author{ -Vincent T van Hees -} diff --git a/man/get_starttime_weekday_truncdata.Rd b/man/get_starttime_weekday_truncdata.Rd new file mode 100644 index 000000000..127356eaf --- /dev/null +++ b/man/get_starttime_weekday_truncdata.Rd @@ -0,0 +1,51 @@ +\name{get_starttime_weekday_truncdata } +\alias{get_starttime_weekday_truncdata} +\title{ + Get starttime (adjusted), weekday, and adjust data accordingly. +} +\description{ + Function not intended for direct use by user. + Used inside \link{g.getmeta} as an intermediate step between + loading the raw data and calibrating it. This step includes extracting + the starttime and adjusting it to nearest integer number of long epoch window + lengths in an hour, truncating the data accordingly, and extracting the + corresponding weekday. +} +\usage{ + get_starttime_weekday_truncdata(monc, + dformat, data, header, desiredtz, sf, + datafile, ws2, configtz = NULL) +} +\arguments{ + \item{monc}{ + See \link{g.inspectfile} + } + \item{dformat}{ + See \link{g.dotorcomma} + } + \item{data}{ + Data part of \link{g.readaccfile} output + } + \item{header}{ + Header part of \link{g.readaccfile} output + } + \item{desiredtz}{ + See \link{g.getmeta} + } + \item{sf}{ + Numeric, sample frequency in Hertz + } + \item{datafile}{ + See \link{g.getmeta} + } + \item{ws2}{ + Long epoch length + } + \item{configtz}{ + See \link{g.getmeta} + } +} +\keyword{internal} +\author{ +Vincent T van Hees +} diff --git a/man/read.myacc.csv.Rd b/man/read.myacc.csv.Rd index 73e6320a6..75c355689 100644 --- a/man/read.myacc.csv.Rd +++ b/man/read.myacc.csv.Rd @@ -35,9 +35,9 @@ interpolationType=1, PreviousLastValue = c(0, 0, 1), PreviousLastTime = NULL, - epochsize = NULL, desiredtz = NULL, - configtz = NULL) + configtz = NULL, + header = NULL) } \arguments{ \item{rmc.file}{ @@ -163,9 +163,6 @@ \item{PreviousLastTime}{ Automatically identified last timestamp in previous chunk of data read. } - \item{epochsize}{ - Numeric vector of length two, with short and long epoch sizes. Only used by GGIR internally. - } \item{desiredtz}{ Timezone in which device was worn. } @@ -173,6 +170,10 @@ Timezone in which device was configured. If equal to desiredtz you can leave this in its default value. } + \item{header}{ + Header information that was extracted the previous time this file was read, + to be re-used instead of being extracted again. + } } \details{ To use this function in the context of GGIR use all arguments from this function, diff --git a/tests/testthat/test_detect_nonwear.R b/tests/testthat/test_detect_nonwear.R index f6dd93c4f..6f3c7d67b 100644 --- a/tests/testthat/test_detect_nonwear.R +++ b/tests/testthat/test_detect_nonwear.R @@ -11,6 +11,7 @@ test_that("detects non wear time", { create_test_acc_csv(Nmin = Ndays * 1440, sf = sf) data = as.matrix(read.csv("123A_testaccfile.csv", skip = 10)) + colnames(data) = c("x", "y", "z") # 2013 algorithm ------ # clipthres to 1.7 to test the clipping detection @@ -40,3 +41,49 @@ test_that("detects non wear time", { # remove generated file ------ if (file.exists("123A_testaccfile.csv")) file.remove("123A_testaccfile.csv") }) + +test_that("the wear column is processed correctly", { + skip_on_cran() + + N = 3000 + sf = 10 + + op = options(digits.secs = 4) + on.exit(options(op)) + + timeseq = ((0:(N - 1))/sf) + + set.seed(100) + accx = rnorm(N) + set.seed(200) + accy = rnorm(N) + set.seed(300) + accz = rnorm(N) + wear = c(rep(TRUE,N/3),rep(FALSE,N/4),rep(TRUE,N/3),rep(FALSE,N/12)) + + # create a csv file with a wear channel + S = data.frame(x = accx, y = accy, z = accz, time = timeseq, wear = wear, stringsAsFactors = TRUE) + testfile = "testcsv.csv" + write.csv(S, file = testfile, row.names = FALSE) + on.exit({if (file.exists(testfile)) file.remove(testfile)}, add = TRUE) + + + I_obj = g.inspectfile(testfile, + rmc.dec = ".", + rmc.firstrow.acc = 1, rmc.firstrow.header = c(), rmc.unit.time = "UNIXsec", + rmc.col.acc = 1:3, rmc.col.temp = c(), rmc.col.time = 4, + rmc.sf = sf, desiredtz = "", + rmc.col.wear = 5) + + filequality = data.frame(filetooshort = FALSE, filecorrupt = FALSE, + filedoesnotholdday = FALSE, NFilePagesSkipped = 0) + + M_obj = g.getmeta(datafile = testfile, desiredtz = "", windowsizes = c(1,60,120), + inspectfileobject = I_obj, + rmc.dec = ".", + rmc.firstrow.acc = 1, rmc.firstrow.header = c(), rmc.unit.time = "UNIXsec", + rmc.col.acc = 1:3, rmc.col.temp = c(), rmc.col.time = 4, + rmc.sf = sf, desiredtz = "", + rmc.col.wear = 5) + expect_equal(M_obj$metalong$nonwearscore, c(3,0,3,3,3)) +}) diff --git a/tests/testthat/test_greadaccfile.R b/tests/testthat/test_greadaccfile.R index ffd457488..cfc752cca 100644 --- a/tests/testthat/test_greadaccfile.R +++ b/tests/testthat/test_greadaccfile.R @@ -1,24 +1,47 @@ library(GGIR) context("g.readaccfile") -test_that("g.readaccfile and g.inspectfile can read gt3x, cwa, and actigraph csv files correctly", { +test_that("g.readaccfile and g.inspectfile can read movisens, gt3x, cwa, Axivity csv, and actigraph csv files correctly", { skip_on_cran() - cat("\nActigraph .csv") - for (filename in c("ActiGraph13.csv", "ActiGraph61.csv")) { - csvfile = system.file(paste0("testfiles/", filename), package = "GGIR")[1] + desiredtz = "Pacific/Auckland" + configtz = "Europe/Berlin" + params = extract_params(input = list(frequency_tol = 0.1, interpolationType = 1, + desiredtz = desiredtz, configtz = configtz)) + params_rawdata = params$params_rawdata + params_general = params$params_general - Icsv = g.inspectfile(csvfile, desiredtz = "Europe/London") - expect_equal(Icsv$monc, MONITOR$ACTIGRAPH) - expect_equal(Icsv$dformc, FORMAT$CSV) - } + filequality = list(filecorrupt = FALSE, filetooshort = FALSE) + dayborder = 0 + + # For Axivity csv files, we'll be able to read files with both unix and formatted (Y-M-D h:m:s) timestamps + Ax3CsvFile = system.file("testfiles/ax3_testfile_unix_timestamps.csv", package = "GGIR")[1] + Ax6CsvFile = system.file("testfiles/ax6_testfile_formatted_timestamps.csv", package = "GGIR")[1] + + Ax3CwaFile = system.file("testfiles/ax3_testfile.cwa", package = "GGIRread")[1] + Ax6CwaFile = system.file("testfiles/ax6_testfile.cwa", package = "GGIRread")[1] - cwafile = system.file("testfiles/ax3_testfile.cwa", package = "GGIRread")[1] GAfile = system.file("testfiles/GENEActiv_testfile.bin", package = "GGIRread")[1] gt3xfile = system.file("testfiles/actigraph_testfile.gt3x", package = "GGIR")[1] - desiredtz = "Europe/London" - params_rawdata = list(frequency_tol = 0.1, interpolationType = 1, - desiredtz = "Europe/London", configtz = "Europe/London") + cat("\nActigraph .csv") + + create_test_acc_csv() + filename = "123A_testaccfile.csv" + on.exit({if (file.exists(filename)) file.remove(filename)}, add = TRUE) + + Icsv = g.inspectfile(filename, desiredtz = desiredtz) + expect_equal(Icsv$monc, MONITOR$ACTIGRAPH) + expect_equal(Icsv$dformc, FORMAT$CSV) + + csv_read = g.readaccfile(filename, blocksize = 10, blocknumber = 1, filequality = filequality, + dayborder = dayborder, ws = 3, + PreviousEndPage = 1, inspectfileobject = Icsv, + params_rawdata = params_rawdata, params_general = params_general) + expect_equal(nrow(csv_read$P$data), 3000) + expect_false(csv_read$filequality$filecorrupt) + expect_false(csv_read$filequality$filetooshort) + expect_equal(sum(csv_read$P$data), 3151.11, tolerance = .01, scale = 1) + cat("\nActigraph .gt3x") # actigraph .gt3x Igt3x = g.inspectfile(gt3xfile, desiredtz = desiredtz) @@ -27,23 +50,84 @@ test_that("g.readaccfile and g.inspectfile can read gt3x, cwa, and actigraph csv expect_equal(Igt3x$sf, 30) EHV = g.extractheadervars(Igt3x) expect_equal(EHV$deviceSerialNumber, "MOS2E39180594_firmware_1.9.2") + + gt3x_read = g.readaccfile(gt3xfile, blocksize = 3000, blocknumber = 1, filequality = filequality, + dayborder = dayborder, ws = 3, + PreviousEndPage = 1, inspectfileobject = Igt3x, + params_rawdata = params_rawdata, params_general = params_general) + expect_equal(nrow(gt3x_read$P$data), 17640) + expect_false(gt3x_read$filequality$filecorrupt) + expect_false(gt3x_read$filequality$filetooshort) + expect_equal(sum(gt3x_read$P$data[c("x","y","z")]), 2732.35, tolerance = .01, scale = 1) + Mgt3x = g.getmeta(datafile = gt3xfile, desiredtz = desiredtz, windowsize = c(1,300,300), inspectfileobject = Igt3x) expect_true(Mgt3x$filetooshort) expect_false(Mgt3x$filecorrupt) + cat("\nAxivity .cwa") - # axivity .cwa - Icwa = g.inspectfile(cwafile, desiredtz = desiredtz, params_rawdata = params_rawdata) + + Icwa = g.inspectfile(Ax3CwaFile, params_rawdata = params_rawdata, params_general = params_general) expect_equal(Icwa$monc, MONITOR$AXIVITY) expect_equal(Icwa$dformc, FORMAT$CWA) expect_equal(Icwa$sf, 100) EHV = g.extractheadervars(Icwa) expect_equal(EHV$deviceSerialNumber,"39434") - Mcwa = g.getmeta(cwafile, desiredtz = desiredtz, windowsize = c(1,300,300), + + cwa_read = g.readaccfile(Ax3CwaFile, blocksize = 10, blocknumber = 1, filequality = filequality, + dayborder = dayborder, ws = 2, + PreviousEndPage = 1, inspectfileobject = Icwa, + params_rawdata = params_rawdata, params_general = params_general) + expect_equal(cwa_read$P$header$blocks, 145) + expect_equal(sum(cwa_read$P$data[c("x","y","z")]), 280.53, tolerance = .01, scale = 1) + + Mcwa = g.getmeta(Ax3CwaFile, desiredtz = desiredtz, windowsize = c(1,300,300), inspectfileobject = Icwa) expect_true(Mcwa$filetooshort) expect_false(Mcwa$filecorrupt) - + + ax3_start_timestamp = cwa_read$P$data$time[1] + + Icwa = g.inspectfile(Ax6CwaFile, params_rawdata = params_rawdata, params_general = params_general) + cwa_read = g.readaccfile(Ax6CwaFile, blocksize = 10, blocknumber = 1, filequality = filequality, + dayborder = dayborder, ws = 2, + PreviousEndPage = 1, inspectfileobject = Icwa, + params_rawdata = params_rawdata, params_general = params_general) + ax6_start_timestamp = cwa_read$P$data$time[1] + + cat("\nAxivity .csv") + + for (csvData in list(list(Ax3CsvFile, 200, -11.80, ax3_start_timestamp), + list(Ax6CsvFile, 200, 14.84, ax6_start_timestamp))) { + IAxivityCsv = g.inspectfile(csvData[[1]], params_rawdata = params_rawdata, params_general = params_general) + expect_equal(IAxivityCsv$monc, MONITOR$AXIVITY) + expect_equal(IAxivityCsv$dformc, FORMAT$CSV) + + csv_read = g.readaccfile(csvData[[1]], blocksize = 1, blocknumber = 1, filequality = filequality, + dayborder = dayborder, ws = 1, + PreviousEndPage = 1, inspectfileobject = IAxivityCsv, + params_rawdata = params_rawdata, params_general = params_general) + + # For both ax3 and ax6 files, we expect 4 columns: timestamp and XYZ. + # All gyro data in ax6 files gets ignored. + expect_equal(ncol(csv_read$P$data), 4) + + expect_equal(nrow(csv_read$P$data), csvData[[2]]) + expect_false(csv_read$filequality$filecorrupt) + expect_false(csv_read$filequality$filetooshort) + expect_equal(sum(csv_read$P$data[c("x","y","z")]), csvData[[3]], tolerance = .01, scale = 1) + + # check that the timestamps for the Axivity csv look the same as they did for + # the original cwa version of the same file (this verifies that timestamp conversion + # worked the same for both formats) + expect_equal(csv_read$P$data$time[1], csvData[[4]], tolerance = .01, scale = 1) + + MAxCsv = g.getmeta(datafile = Ax3CsvFile, desiredtz = desiredtz, windowsize = c(1,300,300), + inspectfileobject = IAxivityCsv) + expect_true(MAxCsv$filetooshort) + expect_false(MAxCsv$filecorrupt) + } + cat("\nGENEActiv .bin") # GENEActiv .bin IGA = g.inspectfile(GAfile, desiredtz = desiredtz) @@ -53,38 +137,167 @@ test_that("g.readaccfile and g.inspectfile can read gt3x, cwa, and actigraph csv EHV = g.extractheadervars(IGA) expect_equal(EHV$deviceSerialNumber,"012967") + + GA_num_blocks = 2 + GA_read = g.readaccfile(GAfile, blocksize = GA_num_blocks, blocknumber = 1, filequality = filequality, + dayborder = dayborder, ws = 3, + desiredtz = desiredtz, PreviousEndPage = 1, inspectfileobject = IGA) + + # As of R 4.0, an extra header row is extracted, which affects the positioning of the values. + # expect_equal(as.numeric(as.character(wav_read$P$header$hvalues[7])),17) + expect_equal(round(sum(GA_read$P$data[, 2:4]), digits = 2), -271.97) + expect_equal(GA_read$endpage, GA_num_blocks) + + # print(GA_read$P$header) + # expect_equal(as.character(unlist(GA_read$P$header[3, 1])), "216 Hours") + MGA = g.getmeta(GAfile, desiredtz = desiredtz, windowsize = c(1,300,300), verbose = FALSE, inspectfileobject = IGA) expect_true(MGA$filetooshort) + + cat("\n Movisens") + + output_dir = "output_unisensExample" + on.exit({if (file.exists(output_dir)) unlink(output_dir, recursive = TRUE)}, add = TRUE) + if (file.exists(output_dir)) unlink(output_dir, recursive = TRUE) + + zip_file = "0.3.4.zip" + on.exit({if (file.exists(zip_file)) unlink(zip_file)}, add = TRUE) + if (!file.exists(zip_file)) { + # link to a tagged release of Unisens/unisensR github repo + movisens_url = "https://github.com/Unisens/unisensR/archive/refs/tags/0.3.4.zip" + download.file(url = movisens_url, destfile = zip_file, quiet = TRUE) + } + movisens_dir = "unisensR-0.3.4" + on.exit({if (file.exists(movisens_dir)) unlink(movisens_dir, recursive = TRUE)}, add = TRUE) + if (file.exists(movisens_dir)) { + unlink(movisens_dir, recursive = TRUE) + } + unzip(zipfile = zip_file, exdir = ".") + movisensFile = file.path(getwd(), "unisensR-0.3.4/tests/unisensExample/acc.bin") + + Mcsv = g.inspectfile(movisensFile, desiredtz = desiredtz) + expect_equal(Mcsv$monc, MONITOR$MOVISENS) + expect_equal(Mcsv$dformc, FORMAT$BIN) + expect_equal(Mcsv$sf, 64) + + movisens_blocksize = 3000 + movisens_read = g.readaccfile(movisensFile, blocksize = movisens_blocksize, blocknumber = 1, filequality = filequality, + dayborder = dayborder, ws = 3, + PreviousEndPage = 1, inspectfileobject = Mcsv, + params_rawdata = params_rawdata, params_general = params_general) + expect_equal(nrow(movisens_read$P$data), movisens_blocksize) + expect_false(movisens_read$filequality$filecorrupt) + expect_false(movisens_read$filequality$filetooshort) + expect_equal(sum(movisens_read$P$data[c("x","y","z")]), 4383.67, tolerance = .01, scale = 1) + expect_equal(movisens_read$endpage, movisens_blocksize) + + # read the next block (set PreviousEndPage to movisens_read$endpage) + movisens_read2 = g.readaccfile(movisensFile, blocksize = movisens_blocksize, blocknumber = 2, filequality = filequality, + dayborder = dayborder, ws = 3, + PreviousEndPage = movisens_read$endpage, inspectfileobject = Mcsv, + params_rawdata = params_rawdata, params_general = params_general) + expect_equal(nrow(movisens_read2$P$data), movisens_blocksize) + expect_equal(movisens_read2$endpage, movisens_blocksize * 2) + # if the 1st sample of 2nd block is identical to the last sample of the 1st block, + # this means that we calculated the startpage of the 2nd block incorrectly. + expect_false(any(movisens_read2$P$data[1,] == movisens_read$P$data[nrow(movisens_read$P$data),])) + + # ad-hoc csv file + + # Create test file: No header, with temperature, with time. + # The first row of the file will contain column names. + # Have this file contain 6000 samples. We'll read it in 2 sets of 3000 lines, + # and if there are any lines unnecessarily skipped, then the second attempt to + # read a block of 3000 will return fewer than 3000 lines. + N = 6000 + sf = 30 + x = Sys.time()+((0:(N-1))/sf) + timestamps = as.POSIXlt(x, origin="1970-1-1", tz = configtz) + mydata = data.frame(Xcol = rnorm(N), timecol = timestamps, Ycol = rnorm(N), Zcol = rnorm(N), + tempcol = rnorm(N) + 20) + testfile = "testcsv1.csv" + on.exit({if (file.exists(testfile)) file.remove(testfile)}, add = TRUE) + + write.csv(mydata, file = testfile, row.names = FALSE) + + AHcsv = g.inspectfile(testfile, + rmc.dec=".", rmc.sf=30, rmc.unit.time="POSIX", + rmc.firstrow.acc = 1, rmc.firstrow.header=c(), + rmc.col.acc = c(1,3,4), rmc.col.temp = 5, rmc.col.time=2, + rmc.unit.acc = "g", rmc.unit.temp = "C", rmc.origin = "1970-01-01") + + expect_equal(AHcsv$monc, MONITOR$AD_HOC) + expect_equal(AHcsv$dformc, FORMAT$AD_HOC_CSV) + expect_equal(AHcsv$sf, 30) + + # Read the file starting with row 1 (rmc.firstrow.acc = 1); this row contains column names. + # Verify that full 3000 rows are still read. + csv_read = g.readaccfile(testfile, blocksize = 10, blocknumber = 1, filequality = filequality, # blocksize is # of pages of 300 samples + dayborder = dayborder, ws = 3, desiredtz = desiredtz, configtz = configtz, + PreviousEndPage = c(), inspectfileobject = AHcsv, + rmc.dec=".", rmc.sf=30, rmc.unit.time="POSIX", + rmc.firstrow.acc = 1, rmc.firstrow.header=c(), + rmc.col.acc = c(1,3,4), rmc.col.temp = 5, rmc.col.time=2, + rmc.unit.acc = "g", rmc.unit.temp = "C", rmc.origin = "1970-01-01") + + expect_equal(nrow(csv_read$P$data), 3000) + expect_false(csv_read$filequality$filecorrupt) + expect_false(csv_read$filequality$filetooshort) + expect_equal(sum(csv_read$P$data[c("x","y","z")]), -38.90, tolerance = .01, scale = 1) + # endpage should be (# of rows read + 1) because for the next block we'll need to skip + # not just the rows read but also the row containing column names. + expect_equal(csv_read$endpage, 3001) + + # since the 1st row of the file contains column names, pointing rmc.firstrow.acc 2 + # should lead to the same eaxt 3000 lines being read (the lines after the column names). + csv_read2 = g.readaccfile(testfile, blocksize = 10, blocknumber = 1, filequality = filequality, + dayborder = dayborder, ws = 3, desiredtz = desiredtz, configtz = configtz, + PreviousEndPage = c(), inspectfileobject = AHcsv, + rmc.dec=".", rmc.sf=30, rmc.unit.time="POSIX", + rmc.firstrow.acc = 2, rmc.firstrow.header=c(), + rmc.col.acc = c(1,3,4), rmc.col.temp = 5, rmc.col.time=2, + rmc.unit.acc = "g", rmc.unit.temp = "C", rmc.origin = "1970-01-01") + + expect_equal(nrow(csv_read2$P$data), 3000) + expect_false(csv_read2$filequality$filecorrupt) + expect_false(csv_read2$filequality$filetooshort) + expect_equal(sum(csv_read2$P$data[c("x","y","z")]), sum(csv_read$P$data[c("x","y","z")]), tolerance = .01, scale = 1) + expect_equal(csv_read2$endpage, 3000) + + # reading the next 3000 lines should also give the same result for rmc.firstrow.acc == 1 or 2. + csv_read3 = g.readaccfile(testfile, blocksize = 10, blocknumber = 2, filequality = filequality, + dayborder = dayborder, ws = 3, desiredtz = desiredtz, configtz = configtz, + PreviousEndPage = csv_read$endpage, inspectfileobject = AHcsv, + rmc.dec=".", rmc.sf=30, rmc.unit.time="POSIX", + rmc.firstrow.acc = 1, rmc.firstrow.header=c(), + rmc.col.acc = c(1,3,4), rmc.col.temp = 5, rmc.col.time=2, + rmc.unit.acc = "g", rmc.unit.temp = "C", rmc.origin = "1970-01-01") + + expect_equal(nrow(csv_read3$P$data), 3000) + + csv_read4 = g.readaccfile(testfile, blocksize = 10, blocknumber = 2, filequality = filequality, + dayborder = dayborder, ws = 3, desiredtz = desiredtz, configtz = configtz, + PreviousEndPage = csv_read2$endpage, inspectfileobject = AHcsv, + rmc.dec=".", rmc.sf=30, rmc.unit.time="POSIX", + rmc.firstrow.acc = 2, rmc.firstrow.header=c(), + rmc.col.acc = c(1,3,4), rmc.col.temp = 5, rmc.col.time=2, + rmc.unit.acc = "g", rmc.unit.temp = "C", rmc.origin = "1970-01-01") + + expect_equal(nrow(csv_read4$P$data), 3000) + expect_equal(sum(csv_read3$P$data[c("x","y","z")]), sum(csv_read4$P$data[c("x","y","z")]), tolerance = .01, scale = 1) + # test decimal separator recognition extraction - decn = g.dotorcomma(cwafile,dformat = FORMAT$CWA, mon = MONITOR$AXIVITY, desiredtz = desiredtz) + decn = g.dotorcomma(Ax3CwaFile,dformat = FORMAT$CWA, mon = MONITOR$AXIVITY, desiredtz = desiredtz) expect_equal(decn,".") decn = g.dotorcomma(GAfile, dformat = FORMAT$BIN, mon = MONITOR$GENEACTIV, desiredtz = desiredtz) expect_equal(decn,".") - filequality = list(filecorrupt = FALSE, filetooshort = FALSE) - dayborder = 0 - - cwa_read = g.readaccfile(cwafile, blocksize = 10, blocknumber = 1, filequality = filequality, - dayborder,ws = 3, desiredtz = desiredtz, - PreviousEndPage = 1, inspectfileobject = Icwa, - params_rawdata = params_rawdata) - GA_read = g.readaccfile(GAfile, blocksize = 2, blocknumber = 1, filequality = filequality, - dayborder = dayborder, ws = 3, - desiredtz = desiredtz, PreviousEndPage = 1, inspectfileobject = IGA) - expect_equal(cwa_read$P$header$blocks, 145) - expect_equal(round(cwa_read$P$data[200, 6], digits = 4), 0) - - # As of R 4.0, an extra header row is extracted, which affects the positioning of the values. - # expect_equal(as.numeric(as.character(wav_read$P$header$hvalues[7])),17) - expect_equal(round(sum(GA_read$P$data.out[, 2:4]), digits = 2), -467.59) - # print(GA_read$P$header) - # expect_equal(as.character(unlist(GA_read$P$header[3, 1])), "216 Hours") - + #also test one small other function: datadir = system.file("testfiles", package = "GGIR")[1] fnames = datadir2fnames(datadir = datadir, filelist = FALSE) - expect_equal(length(fnames$fnames), 6) - expect_equal(length(fnames$fnamesfull), 6) + expect_equal(length(fnames$fnames), 7) + expect_equal(length(fnames$fnamesfull), 7) }) diff --git a/tests/testthat/test_imputeTimegaps.R b/tests/testthat/test_imputeTimegaps.R index 86b44f269..3813b4d8a 100644 --- a/tests/testthat/test_imputeTimegaps.R +++ b/tests/testthat/test_imputeTimegaps.R @@ -4,37 +4,46 @@ test_that("timegaps are correctly imputed", { N = 10000 sf = 20 x = data.frame(time = as.POSIXct(x = (1:N)/sf, tz = "", origin = "1970/1/1"), - X = 1:N, Y = 1:N, Z = 1:N) - xyzCol = c("X", "Y", "Z") + x = 1:N, y = 1:N, z = 1:N) + xyzCol = c("x", "y", "z") + x_without_time = data.frame(x = 1:N, y = 1:N, z = 1:N) + xyzCol = c("x", "y", "z") - x_without_time = data.frame(X = 1:N, Y = 1:N, Z = 1:N) - xyzCol = c("X", "Y", "Z") - + # Duration of each consecutive gap is equal to the distance netween + # the sample right before and the sample right after the zeros that got removed to create this gap. + # So the duration of each gap is equal to the number of zeros + 1. + ngaps = 4 zeros = c(5:200, 6000:6500, 7000:7500, 8000:9500) + gaps_duration = length(zeros) + ngaps + gaps_duration = gaps_duration/sf/60 # TEST THAT SAME FILE WITH DIFFERENT FORMATS IS IMPUTED IN THE SAME WAY ---- # Format 1: with timestamp & with timegaps (no zeroes, incomplete dataset) x1 = x[-zeros,] - x1_imputed = g.imputeTimegaps(x1, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1)) + x1_imputed = g.imputeTimegaps(x1, sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1)) x1_imputed_QClog = x1_imputed$QClog; x1_imputed = x1_imputed$x - x1_removed = g.imputeTimegaps(x1, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1)) + x1_removed = g.imputeTimegaps(x1, sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1)) x1_removed_QClog = x1_removed$QClog; x1_removed = x1_removed$x - + + # make sure the timestamps got imputed correctly + # (here we are checking that the last imputed timestamp is correct relative to the first timestamp after the gap) + expect_equal(as.numeric(x1_imputed$time[201] - x1_imputed$time[200]), 1/sf, tolerance = .01, scale = 1) + # Format 2: with timestamp & with zeros (complete dataset) x2 = x x2[zeros, xyzCol] = 0 - x2_imputed = g.imputeTimegaps(x2, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1)) + x2_imputed = g.imputeTimegaps(x2, sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1)) x2_imputed_QClog = x2_imputed$QClog; x2_imputed = x2_imputed$x - x2_removed = g.imputeTimegaps(x2, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1)) + x2_removed = g.imputeTimegaps(x2, sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1)) x2_removed_QClog = x2_removed$QClog; x2_removed = x2_removed$x # Format 3: without timestamp & with zeros (complete dataset) x3 = x_without_time x3[zeros, xyzCol] = 0 - x3_imputed = g.imputeTimegaps(x3, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1)) + x3_imputed = g.imputeTimegaps(x3, sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1)) x3_imputed_QClog = x3_imputed$QClog; x3_imputed = x3_imputed$x - x3_removed = g.imputeTimegaps(x3, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1)) + x3_removed = g.imputeTimegaps(x3, sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1)) x3_removed_QClog = x3_removed$QClog; x3_removed = x3_removed$x # tests number of rows @@ -58,19 +67,19 @@ test_that("timegaps are correctly imputed", { expect_equal(x2_imputed_QClog$timegaps_n, 4) expect_equal(x3_imputed_QClog$timegaps_n, 4) - expect_equal(x1_imputed_QClog$timegaps_min, length(zeros)/sf/60) - expect_equal(x2_imputed_QClog$timegaps_min, length(zeros)/sf/60) - expect_equal(x3_imputed_QClog$timegaps_min, length(zeros)/sf/60) + expect_equal(x1_imputed_QClog$timegaps_min, gaps_duration) + expect_equal(x2_imputed_QClog$timegaps_min, gaps_duration) + expect_equal(x3_imputed_QClog$timegaps_min, gaps_duration) # TEST IMPUTATION WHEN FIRST ROW IS NOT CONSECUTIVE TO PREVIOUS CHUNK ---- # Format 4: with timestamp & with timegaps (no zeroes, incomplete dataset) x4 = x[-zeros,] PreviousLastTime = x[1,"time"] - 30 # dummy gap of 30 seconds between chunks suppressWarnings({ # warning arising from made up PreviousLastTime - x4_imputed = g.imputeTimegaps(x4, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = TRUE, + x4_imputed = g.imputeTimegaps(x4, sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1), PreviousLastTime = PreviousLastTime) x4_imputed_QClog = x4_imputed$QClog; x4_imputed = x4_imputed$x - x4_removed = g.imputeTimegaps(x4, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = FALSE, + x4_removed = g.imputeTimegaps(x4, sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1), PreviousLastTime = PreviousLastTime) x4_removed_QClog = x4_removed$QClog; x4_removed = x4_removed$x }) @@ -83,9 +92,9 @@ test_that("timegaps are correctly imputed", { # Format 5: with timestamp & with zeros (complete dataset) x5 = x x5[zeros, xyzCol] = 0 - x5_imputed = g.imputeTimegaps(x5, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1)) + x5_imputed = g.imputeTimegaps(x5, sf = sf, k = 2/sf, impute = TRUE, PreviousLastValue = c(0,0,1)) x5_imputed_QClog = x5_imputed$QClog; x5_imputed = x5_imputed$x - x5_removed = g.imputeTimegaps(x5, xyzCol, timeCol = "time", sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1)) + x5_removed = g.imputeTimegaps(x5, sf = sf, k = 2/sf, impute = FALSE, PreviousLastValue = c(0,0,1)) x5_removed_QClog = x5_removed$QClog; x5_removed = x5_removed$x expect_equal(nrow(x5_imputed), N) diff --git a/tests/testthat/test_read.myacc.csv.R b/tests/testthat/test_read.myacc.csv.R index de7fde65d..567d6501a 100644 --- a/tests/testthat/test_read.myacc.csv.R +++ b/tests/testthat/test_read.myacc.csv.R @@ -115,8 +115,8 @@ test_that("read.myacc.csv can handle files without header, no decimal places in # Evaluate with decimal places in seconds expect_equal(nrow(D1$data), 20) expect_equal(ncol(D1$data), 5) - expect_equal(strftime(D1$data$time[1:5], format = '%Y-%m-%d %H:%M:%OS2', - tz = "Europe/London"), + expect_equal(strftime(as.POSIXct(D1$data$time[1:5], tz = "Europe/London", origin = "1970-01-01"), + format = '%Y-%m-%d %H:%M:%OS2', tz = "Europe/London"), c("2022-11-02 13:01:16.00", "2022-11-02 13:01:16.03", "2022-11-02 13:01:16.06", @@ -138,8 +138,8 @@ test_that("read.myacc.csv can handle files without header, no decimal places in rmc.headername.recordingid = "ID") expect_equal(nrow(D2$data), 20) expect_equal(ncol(D2$data), 5) - expect_equal(strftime(D2$data$time[1:5], format = '%Y-%m-%d %H:%M:%OS2', - tz = "Europe/London"), + expect_equal(strftime(as.POSIXct(D2$data$time[1:5], tz = "Europe/London", origin = "1970-01-01"), + format = '%Y-%m-%d %H:%M:%OS2', tz = "Europe/London"), c("2022-11-02 18:01:16.50", "2022-11-02 18:01:16.53", "2022-11-02 18:01:16.56", "2022-11-02 18:01:16.59", "2022-11-02 18:01:16.63")) @@ -397,7 +397,7 @@ test_that("read.myacc.csv can handle gaps in time and irregular sample rate", { rmc.headername.sn = "serial_number", rmc.headername.recordingid = "ID", rmc.check4timegaps = TRUE, rmc.doresample = TRUE) - expect_that(nrow(D2$data), equals(169)) # because data expands with 5 seconds that are now imputed + expect_that(nrow(D2$data), equals(170)) # because data expands with 5 seconds that are now imputed expect_that(ncol(D2$data), equals(5))