diff --git a/NAMESPACE b/NAMESPACE index 873739a..d0f0416 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -31,6 +31,7 @@ exportMethods("[") exportMethods("[[") exportMethods(as.list) exportMethods(bands) +exportMethods(coordinates) exportMethods(coverages) exportMethods(createPatterns) exportMethods(crop) diff --git a/R/methods.R b/R/methods.R index 8216b4a..53b9e7a 100644 --- a/R/methods.R +++ b/R/methods.R @@ -402,6 +402,17 @@ setMethod("crop", } ) +#' @inheritParams twdtwRaster-class +#' @rdname twdtwRaster-class +#' @param obj object of cless twdtwRaster +#' @export +setMethod("coordinates", + signature = signature("twdtwRaster"), + definition = function(obj, ...){ + coordinates(obj@timeseries[[1]], ...) + } +) + #' @inheritParams twdtwRaster-class #' @rdname twdtwRaster-class #' @export diff --git a/R/twdtwApply.R b/R/twdtwApply.R index ae7769b..9b4a1bf 100644 --- a/R/twdtwApply.R +++ b/R/twdtwApply.R @@ -75,6 +75,8 @@ #' #' @param chunk.size An integer. Set the number of cells for each block, #' see \code{\link[raster]{blockSize}} for details. +#' +#' @param silent if set to TRUE then hide chunk processing message. Default is FALSE. #' #' @references #' [1] Maus V, Camara G, Cartaxo R, Sanchez A, Ramos FM, de Queiroz, GR. @@ -186,7 +188,7 @@ twdtwApply.twdtwTimeSeries = function(x, y, weight.fun, dist.method, step.matrix #' @export setMethod(f = "twdtwApply", "twdtwRaster", def = function(x, y, resample, length, weight.fun, dist.method, step.matrix, n, span, min.length, theta, - breaks=NULL, from=NULL, to=NULL, by=NULL, overlap=0.5, chunk.size=1000, filepath=NULL, ...){ + breaks=NULL, from=NULL, to=NULL, by=NULL, overlap=0.5, chunk.size=1000, filepath=NULL, silent=FALSE, ...){ if(!is(step.matrix, "stepPattern")) stop("step.matrix is not of class stepPattern") if(is.null(weight.fun)) @@ -215,14 +217,12 @@ setMethod(f = "twdtwApply", "twdtwRaster", if(resample) y = resampleTimeSeries(object=y, length=length) twdtwApply.twdtwRaster(x, y, weight.fun, dist.method, step.matrix, n, span, min.length, theta, - breaks, overlap, chunk.size, filepath, ...) + breaks, overlap, chunk.size, filepath, silent, ...) }) twdtwApply.twdtwRaster = function(x, y, weight.fun, dist.method, step.matrix, n, span, min.length, theta, - breaks, overlap, chunk.size, filepath, - mc.preschedule = TRUE, mc.set.seed = TRUE, mc.silent = FALSE, - mc.cores = getOption("mc.cores", 1L), mc.cleanup = TRUE, ...){ + breaks, overlap, chunk.size, filepath, silent, ...){ # Set blocks Multi-thread parameters minblocks = round(nrow(x)*ncol(x) / chunk.size) @@ -235,6 +235,8 @@ twdtwApply.twdtwRaster = function(x, y, weight.fun, dist.method, step.matrix, n, matching_bands = pattern_names[pattern_names %in% raster_bands] if(length(matching_bands)<1) stop(paste0("Attributes (bands) of the raster and patterns do not match")) + if("doy"%in%coverages(x) & !"doy"%in%matching_bands) + matching_bands = c("doy", matching_bands) x = subset(x, layers=matching_bands) raster_bands = coverages(x) @@ -254,42 +256,34 @@ twdtwApply.twdtwRaster = function(x, y, weight.fun, dist.method, step.matrix, n, # Get time line timeline = as.Date(index(x)) - get_aligs = function(x){ - # twdtwApply(x, y=y, weight.fun=weight.fun) - twdtwApply(x, y=y, weight.fun=weight.fun, dist.method=dist.method, step.matrix=step.matrix, n=n, span=span, min.length=min.length, theta=theta, keep=FALSE) - } + # Raster info + proj_str = projection(x) + coord = data.frame(coordinates(x)) + names(coord) = c("longitude","latitude") fun = function(i){ - if(!mc.silent) print(paste0("Procesing chunk ",i,"/",threads[length(threads)])) - - # Get time series from raster - ts_list = lapply(as.list(x), FUN=getValuesBlock, row=blocks$row[i], nrows=blocks$nrows[i]) + if(!silent) print(paste0("Procesing chunk ",i,"/",threads[length(threads)])) - # Create a dummy array - nts = seq(1, nrow(ts_list[[1]])) - m = length(levels) - n = length(breaks)-1 - - # Create zoo time series - ts_zoo = lapply(nts, FUN=.bulidZoo, x=ts_list, timeline=timeline) + # Get twdtwTimeSeries from raster + cells = cellFromRow(x@timeseries$ndvi, blocks$row[i]:blocks$nrows[i]) + ts = getTimeSeries(x, y = coord[cells,], proj4string = proj_str) - # Create twdtwTimeSeries object - ts = try(twdtwTimeSeries(ts_zoo), silent = TRUE) - if(is(ts, "try-error")) - return(lapply(levels, function(l) writeValues(b_files[[l]], matrix(9999, nrow=length(nts), ncol=n), blocks$row[i]))) - # Apply TWDTW analysis - twdtw_results = lapply(as.list(ts), FUN=get_aligs) + twdtw_results = twdtwApply(ts, y=y, weight.fun=weight.fun, dist.method=dist.method, + step.matrix=step.matrix, n=n, span=span, + min.length=min.length, theta=theta, keep=FALSE) # Get best mathces for each point, period, and pattern - A = lapply(twdtw_results, FUN=.lowestDistances, m=m, n=n, levels=levels, breaks=breaks, overlap=overlap, fill=9999) + m = length(levels) + h = length(breaks)-1 + A = lapply(as.list(twdtw_results), FUN=.lowestDistances, m=m, n=h, levels=levels, breaks=breaks, overlap=overlap, fill=9999) # Reshape list to array - A = sapply(A, matrix, nrow=n, ncol=m, simplify = 'array') + A = sapply(A, matrix, nrow=h, ncol=m, simplify = 'array') # Write raster files - lapply(seq_along(levels), function(l) writeValues(b_files[[levels[l]]], matrix(t(A[,l,]),ncol=n), blocks$row[i])) + lapply(seq_along(levels), function(l) writeValues(b_files[[levels[l]]], matrix(t(A[,l,]),ncol=h), blocks$row[i])) } # Apply TWDTW analysis diff --git a/man/twdtwApply.Rd b/man/twdtwApply.Rd index 57641e8..e048c15 100644 --- a/man/twdtwApply.Rd +++ b/man/twdtwApply.Rd @@ -19,7 +19,7 @@ twdtwApply(x, y, resample = TRUE, length = NULL, weight.fun = NULL, \S4method{twdtwApply}{twdtwRaster}(x, y, resample, length, weight.fun, dist.method, step.matrix, n, span, min.length, theta, breaks = NULL, from = NULL, to = NULL, by = NULL, overlap = 0.5, chunk.size = 1000, - filepath = NULL, ...) + filepath = NULL, silent = FALSE, ...) } \arguments{ \item{x}{an object of class twdtw*. This is the target time series. @@ -87,6 +87,8 @@ see \code{\link[raster]{blockSize}} for details.} \item{filepath}{A character. The path to save the raster with results. If not informed the function saves in the current work directory.} + +\item{silent}{if set to TRUE then hide chunk processing message. Default is FALSE.} } \value{ An object of class twdtw*. @@ -137,8 +139,7 @@ time_interval = seq(from=as.Date("2007-09-01"), to=as.Date("2013-09-01"), log_fun = weight.fun=logisticWeight(-0.1,50) r_twdtw = twdtwApply(x=rts, y=patt, weight.fun=log_fun, breaks=time_interval, - filepath="~/test_twdtw", overwrite=TRUE, format="GTiff", mc.cores=3, - chunk.size=1000) + filepath="~/test_twdtw", overwrite=TRUE, format="GTiff") plot(r_twdtw, type="distance") diff --git a/man/twdtwAssess.Rd b/man/twdtwAssess.Rd index 6533a52..7a3c9e7 100644 --- a/man/twdtwAssess.Rd +++ b/man/twdtwAssess.Rd @@ -113,7 +113,7 @@ plot(twdtw_assess, type="accuracy") plot(twdtw_assess, type="area") plot(twdtw_assess, type="map", samples = "all") plot(twdtw_assess, type="map", samples = "incorrect") - +plot(twdtw_assess, type="map", samples = "correct") # Create latex tables twdtwXtable(twdtw_assess, table.type="matrix") diff --git a/man/twdtwRaster-class.Rd b/man/twdtwRaster-class.Rd index ab15896..3169f36 100644 --- a/man/twdtwRaster-class.Rd +++ b/man/twdtwRaster-class.Rd @@ -27,6 +27,7 @@ \alias{[[,twdtwRaster,ANY,ANY-method} \alias{labels,twdtwRaster-method} \alias{crop,twdtwRaster-method} +\alias{coordinates,twdtwRaster-method} \alias{extent,twdtwRaster-method} \alias{show,twdtwRaster-method} \alias{is.twdtwRaster,ANY-method} @@ -74,6 +75,8 @@ \S4method{crop}{twdtwRaster}(x, y, ...) +\S4method{coordinates}{twdtwRaster}(obj, ...) + \S4method{extent}{twdtwRaster}(x, y, ...) \S4method{show}{twdtwRaster}(object) @@ -120,6 +123,8 @@ originally stores in separated files. See details.} \item{object}{an object of class twdtwRaster.} +\item{obj}{object of cless twdtwRaster} + \item{crs}{character or object of class 'CRS'. PROJ.4 description of the coordinate reference system. For other arguments and more details see \code{\link[raster]{projectRaster}}.}