sync stuff done for now

NAMESPACE

@@ -3,38 +3,54 @@
 export(alignBurstSeq)
 export(applyLinCor)
 export(applySync)
+export(applySync_old)
 export(checkInp)
 export(checkInpSync)
 export(checkInpSyncData)
 export(downsampleToaList_random)
 export(downsampleToaList_selective)
 export(fineTuneSyncModel)
 export(getBbox)
+export(getDatSyncCov)
+export(getDatTmbSync)
+export(getDistMat)
 export(getDownsampledToaList)
+export(getEpsLong)
+export(getFixedHydrosVec)
 export(getInp)
 export(getInpSync)
+export(getLinCorr)
 export(getOffsetVals)
 export(getParamsTmbSync)
 export(getRandomTmbSync)
+export(getSmoothOffsets)
+export(getSsVec)
+export(getSyncCheckDat)
 export(getSyncCoverage)
 export(getSyncModel)
 export(getSyncToa)
 export(getToaYaps)
 export(plotBbox)
+export(plotSyncCheck)
+export(plotSyncCov)
+export(plotSyncModel)
 export(plotSyncModelCheck)
 export(plotSyncModelHydros)
 export(plotSyncModelResids)
+export(plotSyncNetwork)
 export(plotYaps)
 export(prepDetections)
+export(prepSyncParams)
 export(runTmb)
 export(runYaps)
 export(simHydros)
 export(simTelemetryTrack)
 export(simToa)
 export(simTrueTrack)
+export(stopSilent)
 export(tempToSs)
 export(testYaps)
-export(yapsifyHydros)
+export(yapsify)
 importFrom(Rcpp,sourceCpp)
 importFrom(data.table,"%between%")
 importFrom(data.table,"%like%")
@@ -50,16 +66,19 @@ importFrom(ggplot2,geom_bar)
 importFrom(ggplot2,geom_boxplot)
 importFrom(ggplot2,geom_histogram)
 importFrom(ggplot2,geom_hline)
+importFrom(ggplot2,geom_label)
 importFrom(ggplot2,geom_line)
 importFrom(ggplot2,geom_point)
 importFrom(ggplot2,geom_ribbon)
+importFrom(ggplot2,geom_segment)
 importFrom(ggplot2,geom_smooth)
 importFrom(ggplot2,geom_tile)
 importFrom(ggplot2,geom_violin)
 importFrom(ggplot2,geom_vline)
 importFrom(ggplot2,ggtitle)
 importFrom(ggplot2,labs)
 importFrom(ggplot2,scale_fill_gradientn)
+importFrom(ggplot2,scale_y_continuous)
 importFrom(ggplot2,scale_y_log10)
 importFrom(ggplot2,theme)
 importFrom(ggplot2,xlab)

NEWS.md

@@ -6,7 +6,7 @@
 * Will break current code 
 * Improve functions to obtain sync models
 * Attempt to improve workflow
-* Remove dependence of deprecated package splusTimeSeries 
+* Remove dependence on deprecated package splusTimeSeries 
 
 
 # yaps v.1.2.5.9000

R/applySync.R → R/applySync_old.R

@@ -7,7 +7,7 @@
 #' @export
 #' @return A `data.table` with the now synchronized time-of-arrivals in column `eposync`. 
 #' @example man/examples/example-yaps_ssu1.R
-applySync <- function(toa, hydros="", sync_model){
+applySync_old <- function(toa, hydros="", sync_model){
 	if(is.matrix(toa)) {type <- "toa_matrix"
 	} else if(data.table::is.data.table(toa)) {type <- "detections_table"}
 

R/getDistMat.R

@@ -0,0 +1,13 @@
+#' Get matrix of hydro-to-hydro distances
+#' @export
+getDistMat <- function(hydros){
+	dist_mat <- matrix(ncol=nrow(hydros), nrow=nrow(hydros))
+	rownames(dist_mat) <- hydros$h_sn
+	colnames(dist_mat) <- hydros$h_sn
+	for(hx in 1:nrow(hydros)){
+		for(hy in 1:nrow(hydros)){
+			dist_mat[hx, hy] <- sqrt((hydros[hx, x] - hydros[hy, x] )^2 + (hydros[hx, y] - hydros[hy, y] )^2 + (hydros[hx, z] - hydros[hy, z] )^2)
+		}
+	}
+	return(dist_mat)
+}

R/getSSVec.R

@@ -0,0 +1,12 @@
+#' Internal function. Extract vector of speed of sounds for each timestamp from supplied data.
+#' Input data must be a data.table() containing columns ts (POSIXct timestamp) and ss (numerical speed of sound)
+#' @inheritParams getInpSync
+#' @noRd
+#' @export
+getSsVec <- function(inp_toa_list, ss_data){
+	roll <- data.table::data.table(ts = as.POSIXct(inp_toa_list$epo_self_vec, origin="1970-01-01", tz="UTC"))
+	data.table::setkey(ss_data, ts)
+	data.table::setkey(roll, ts)
+	ss_data_vec <- ss_data[roll, roll="nearest"]$ss
+	return(ss_data_vec)
+}

R/getSyncCheckDat.R → R/getSyncCheckDat_old.R

@@ -3,7 +3,7 @@
 #' @param extreme_threshold Ignore delta values larger than this threshold. 
 #' @inheritParams getInpSync
 #' @noRd
-getSyncCheckDat <- function(sync_model, extreme_threshold=1000){
+getSyncCheckDat_old <- function(sync_model, extreme_threshold=1000){
 	toa <- sync_model$inp_synced$inp_params$toa
 	toa_sync <- applySync(toa, sync_model=sync_model)
 

R/stopSilent.R

@@ -0,0 +1,8 @@
+#' Helper function to stop script silent
+#' @noRd
+#' @export
+stopSilent <- function() {
+  opt <- options(show.error.messages = FALSE)
+  on.exit(options(opt))
+  stop()
+}

R/syncGetters.R

@@ -1,13 +1,3 @@
-#' Internal function. Extract speed of sounds for each timestamp used in sync-process from supplied data.
-#' @inheritParams getInpSync
-#' @noRd
-getSsDataVec <- function(inp_toa_list, ss_data){
-	roll <- data.table::data.table(ts = as.POSIXct(inp_toa_list$epo_self_vec, origin="1970-01-01", tz="UTC"))
-	data.table::setkey(ss_data, ts)
-	data.table::setkey(roll, ts)
-	ss_data_vec <- ss_data[roll, roll="nearest"]$ss
-	return(ss_data_vec)
-}
 
 # #' Internal function. Get toa for sync from sync_dat
 # #' @inheritParams getInpSync
@@ -40,38 +30,3 @@ getSsVals <- function(inp_toa_list, n_ss_day){
 
 	return(list(n_ss_idx=n_ss_idx, ss_idx=ss_idx, ss_levels=ss_levels))
 }
-
-#' Internal function to get residuals from sync_model in long format
-#' @inheritParams getInpSync
-#' @noRd
-getEpsLong <- function(report, pl, inp_sync){
-
-	if(inp_sync$dat_tmb_sync$ss_data_what == "est"){
-		ss_vec <- pl$SS[inp_sync$dat_tmb_sync$ss_idx]
-	} else {
-		ss_vec <- inp_sync$dat_tmb_sync$ss_data_vec
-	}
-
-
-	eps <- report$eps
-	eps[which(eps==0)] <- NA
-	eps_long <- data.table::data.table(reshape2::melt(eps))
-	if(inp_sync$sync_type == 'top'){
-		colnames(eps_long) <- c('ping', 'hydro_idx', 'E')
-		eps_long[, sync_tag_idx:=rep(inp_sync$dat_tmb_sync$sync_tag_idx_vec, times=ncol(eps))]
-		eps_long[, ss:=rep(ss_vec, times=ncol(eps))]
-	} else {
-		colnames(eps_long) <- c('E')
-		eps_long[ , ping := inp_sync$dat_tmb_sync$toa_delta[,'ping_idx']]
-		eps_long[ , h1_idx := inp_sync$dat_tmb_sync$toa_delta[,'h1']]
-		eps_long[ , h2_idx := inp_sync$dat_tmb_sync$toa_delta[,'h2']]
-		eps_long[ , sync_tag_idx := inp_sync$dat_tmb_sync$toa_delta[,'sync_tag_idx']]
-		eps_long[ , ss := ss_vec[inp_sync$dat_tmb_sync$toa_delta[,'ping_idx']]]
-	}
-	eps_long[, E_m:=E*ss]
-
-	eps_long <- eps_long[!is.na(E)]
-
-	return(eps_long)
-}
-

R/sync_applySync.R

@@ -0,0 +1,86 @@
+#' Apply sync_model to detection data
+#' @export
+applySync <- function(dat_sync, sync_model){
+	dat_temp <- copy(dat_sync)
+
+	lin_corr 	<- sync_model$inp_synced$sync_params$lin_corr
+	gams 		<- sync_model$gams
+	hydros 		<- sync_model$inp_synced$inp_params$hydros
+	epo_range 	<- range(sync_model$inp_synced$inp_params$offset_levels)
+
+	tk <- sync_model$inp_synced$sync_params$time_keeper
+	dat_sync_cov <- getDatSyncCov(sync_model$inp_synced)
+
+
+	dat_synced <- data.table()
+	cat("Applying sync_model to hydro data...\n")
+
+	for(h in 1:nrow(hydros)){
+		gams[[h]]$newdata <- data.table(epofrac=dat_temp[h_sn == hydros[h, h_sn], epofrac])
+	}
+
+	n_cores <- parallel::detectCores()-1
+	cat("...running parallel using ",n_cores," cores\n")
+	cl <- parallel::makeCluster(n_cores)
+	pred_list <- parallel::parLapplyLB(cl, gams,  fun=function(k){
+		mgcv::predict.gam(k, newdata=k$newdata)
+	})
+	parallel::stopCluster(cl)
+
+	for(h in 1:nrow(hydros)){
+		cat(".",h,".")
+		dat_temp_h <- dat_temp[h_sn == hydros[h_idx == h, h_sn]]
+		dat_temp_h[, lin_corr := lin_corr[h, 1] + epofrac * lin_corr[h, 2]]
+		# dat_temp_h[, pred := pred_h <- mgcv::predict.gam(gams[[h]], newdata=dat_temp_h)]
+		dat_temp_h[, pred := pred_list[[h]]]
+		dat_temp_h[, eposync := epofrac - lin_corr - pred]
+
+		dat_temp_h[, lin_corr := NULL]
+		dat_temp_h[, pred := NULL]
+
+		dat_synced <- rbind(dat_synced, dat_temp_h)
+
+		if(h %% 10 == 0){ cat("\n")}
+
+	}
+	cat("\n")
+	setorder(dat_synced, eposync)
+
+	# # # NA'ing data outside sync period
+	dat_synced[!epofrac %between% epo_range, eposync := NA]
+	n_outside_sync_range <- nrow(dat_synced[is.na(eposync)])
+	if(n_outside_sync_range > 0){
+		cat("NOTE: ",n_outside_sync_range," rows are outside the period synced by the model - these are returned as eposync=NA. \n")
+	}
+
+	# # # NA'ing data where time keeper had no data 
+	nas_tk <- dat_sync_cov[h_sn == tk & N == 0]
+	if(nrow(nas_tk) > 0){
+		for(i in 1:nrow(nas_tk)){
+			na_range_i <- sync_model$inp_synced$inp_params$offset_levels[nas_tk[i, offset_idx], ]
+			dat_synced[epofrac %between% na_range_i, eposync := NA]
+		}
+	}
+	n_tk_nas <- nrow(dat_synced[is.na(eposync)]) - n_outside_sync_range
+	if(n_tk_nas > 0){
+		cat("NOTE: ",n_tk_nas," rows are not synced because the time keeper had too few data - these are returned as eposync=NA. This affects all hydros. \n")
+	}
+
+	# # # NAing data where individual hydros had no sync data...
+	nas_non_tk <- dat_sync_cov[h_sn != tk & N == 0]
+	if(nrow(nas_non_tk) > 0){
+		for(i in 1:nrow(nas_non_tk)){
+			na_range_i <- sync_model$inp_synced$inp_params$offset_levels[nas_non_tk[i, offset_idx], ]
+			dat_synced[epofrac %between% na_range_i, eposync := NA]
+		}
+	}
+	n_non_tk_nas <- nrow(dat_synced[is.na(eposync)]) - n_outside_sync_range - n_tk_nas
+	if(n_non_tk_nas > 0){
+		cat("NOTE: ",n_non_tk_nas," rows are not synced because the specific hydro(s) had too few data - these are returned as eposync=NA. \n")
+		cat(paste0("... These hydros are affected: ", paste0(unique(nas_non_tk$h_sn), collapse=","),"\n"))
+	}
+
+	dat_temp <- NULL
+
+	return(dat_synced)
+}