Skip to content

Commit

Permalink
Fix 27 ditch disagg for scaling and resampling. Speed gain
Browse files Browse the repository at this point in the history
  • Loading branch information
mzarowka committed Mar 20, 2024
1 parent 55af691 commit 36e9ecd
Show file tree
Hide file tree
Showing 12 changed files with 523 additions and 72 deletions.
10 changes: 4 additions & 6 deletions DESCRIPTION
Original file line number Diff line number Diff line change
Expand Up @@ -16,19 +16,18 @@ Roxygen: list(markdown = TRUE)
RoxygenNote: 7.3.1
Suggests:
testthat (>= 3.0.0),
arrow,
patchwork,
sf,
knitr,
rmarkdown
rmarkdown,
scales,
prospectr
Config/testthat/edition: 3
URL: https://github.com/mzarowka/HSItools
BugReports: https://github.com/mzarowka/HSItools/issues
Imports:
dplyr,
DT,
ggplot2,
prospectr,
shiny,
shinycssloaders,
shinyFiles,
Expand All @@ -41,11 +40,10 @@ Imports:
rlang,
data.table,
shinyalert,
stringr,
tibble,
readr,
fs,
scales
sf
Depends:
R (>= 4.1.0)
LazyData: true
Expand Down
4 changes: 4 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -28,6 +28,10 @@ export(run_core)
export(spectra_position)
export(spectra_sub)
export(stretch_raster_full)
import(DT)
import(shiny)
import(shinyFiles)
importFrom(ggplot2,theme)
importFrom(rlang,.data)
importFrom(shinyalert,shinyalert)
importFrom(stats,median)
2 changes: 1 addition & 1 deletion R/extractors.R
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ extract_spectral_series <- function(raster, index = NULL, calibration = NULL, ex
} else if (filename == TRUE) {
# Check source
if (terra::sources(raster) == "") {
rlang::warn(message = "In memory object. Using working directory instead.")
rlang::warn(message = "In memory object. Using working directory.")

filename <- paste0(getwd(), "/spectral_profile.csv")

Expand Down
90 changes: 79 additions & 11 deletions R/normalization.R
Original file line number Diff line number Diff line change
Expand Up @@ -179,6 +179,69 @@ raster_crop <- function(raster, type, dir = NULL, roi) {
return(raster)
}

#' #' Create reference SpatRaster
#' #'
#' #' @param raster terra SpatRaster of the captured reference.
#' #' @param roi Region Of Interest: extent to match data raster.
#' #' @param ref_type type of reference, one of "whiteref" or "darkref".
#' #' @param ... additional arguments.
#' #'
#' #' @return a terra SpatRaster of reference matching the data raster extent.
#' #' @export
#' #'
#' #' @description Creating reference SpatRaster covering core extent
#' #' Create one mean reference row SpatRaster by averaging data every column by aggregation
#' #' Create reference SpatRaster matching capture SpatRaster extent by disaggregation.
#' create_reference_raster_disagg <- function(raster, roi, ref_type, ...) {
#' # Store additional parameters
#' params <- rlang::list2(...)
#'
#' # Check if correct class is supplied.
#' if (!inherits(raster, what = "SpatRaster")) {
#' rlang::abort(message = "Supplied data is not a terra SpatRaster.")
#' }
#'
#' # Check number of rois and prepare ids.
#' if (length(roi) == 1) {
#' roi_id <- NULL
#' } else {
#' roi_id <- paste0("ROI_", seq(1:length(roi)))
#' }
#'
#' if (ref_type == "whiteref") {
#' name <- "WHITEREF"
#'
#' } else {
#' name <- "DARKREF"
#' }
#'
#' # Aggregate data into one row SpatRaster, divide by number of rows
#' cli::cli_alert("Aggregate { name }")
#'
#' raster <- terra::aggregate(
#' raster,
#' fact = c(terra::nrow(raster), 1),
#' fun = "mean",
#' overwrite = TRUE,
#' steps = terra::ncell(raster) * terra::nlyr(raster))
#'
#' # Set new extent to match extent of capture SpatRaster
#' terra::ext(raster) <- roi
#'
#' # Disaggregate data over entire extent to match capture SpatRaster extent, multiply by ymax
#' cli::cli_alert("Disaggregate { name }")
#'
#' raster <- terra::disagg(
#' raster,
#' fact = c(terra::ymax(raster), 1),
#' filename = paste0(params$path, "/products/", name, "_", basename(params$path), "_disaggregated.tif"),
#' overwrite = TRUE,
#' steps = terra::ncell(raster) * terra::nlyr(raster))
#'
#' # Return raster
#' return(raster)
#' }

#' Create reference SpatRaster
#'
#' @param raster terra SpatRaster of the captured reference.
Expand Down Expand Up @@ -215,28 +278,33 @@ create_reference_raster <- function(raster, roi, ref_type, ...) {
name <- "DARKREF"
}

# Create 1 row template
template_row <- terra::rast(terra::ext(c(terra::xmin(roi), terra::xmax(roi), 0, 1)), resolution = c(1, 1))

# Create full core template
template_core <- terra::rast(terra::ext(roi), resolution = c(1, 1), nlyrs = terra::nlyr(raster))

# Aggregate data into one row SpatRaster, divide by number of rows
cli::cli_alert("Aggregate { name }")

raster <- terra::aggregate(
raster,
fact = c(terra::nrow(raster), 1),
fun = "mean",
overwrite = TRUE,
steps = terra::ncell(raster) * terra::nlyr(raster))
fun = "mean") |>
terra::resample(
template_row)

# Set new extent to match extent of capture SpatRaster
terra::ext(raster) <- roi

# Disaggregate data over entire extent to match capture SpatRaster extent, multiply by ymax
cli::cli_alert("Disaggregate { name }")
# Scale and resample data over entire extent to match capture SpatRaster extent, multiply by ymax
cli::cli_alert("Scaling and resampling { name }")

raster <- terra::disagg(
raster,
fact = c(terra::ymax(raster), 1),
filename = paste0(params$path, "/products/", name, "_", basename(params$path), "_disaggregated.tif"),
overwrite = TRUE,
steps = terra::ncell(raster) * terra::nlyr(raster))
raster <- raster |>
terra::rescale(fx = 1, fy = terra::nrow(template_core) * 2) |>
terra::resample(
template_core,
filename = paste0(params$path, "/products/", name, "_", basename(params$path), "_disaggregated.tif"))

# Return raster
return(raster)
Expand Down
2 changes: 2 additions & 0 deletions R/utils.R
Original file line number Diff line number Diff line change
Expand Up @@ -155,8 +155,10 @@ pixel_to_distance <- function(core, pixel_ratio, ymax, ymin = 0, sample_start, s
#' @export
#'
#' @examples
#' if (interactive() == TRUE) {
#' a1 <- readRDS(file.path(system.file(package = "HSItools"),"extdata/HSItools_core.rds"))
#' change_output_dir(a1)
#' }
#'
change_output_dir = function(run_core_output){
#currentRoot <- rprojroot::find_root(rprojroot::criteria$is_rstudio_project)
Expand Down
Loading

0 comments on commit 36e9ecd

Please sign in to comment.