R/habitat.analysis.R

In mottensmann/DBChecks: Checking Database

#' Habitat analysis
#'
#' @description Crops vector layers using overlay and recalculates area
#' (if polygon) and length (if line) of objects
#'
#' @details vectorised if length overlay > 1
#'
#' @details
#' 1.) Clip to overlay extent (native:clip)
#' 2.) Recalculate geometry of objects
#' 3.) Output data.frame
#'
#' @param features list of ESRI shapes
#' @param overlay overlay ESRI shape file
#' @inheritParams meancoordinates
#' @import magrittr
#' @export
#'
habitat.analysis <- function(features = NULL,
                             overlay =  NULL,
                             root =  'C:/OSGeo4W64') {

## checks:
## ###########################################################################
if (class(features) != "list") stop("features must be a list of shapes")

## set up RQGIS
## ###########################################################################
if ("RQGIS3" %in% rownames(installed.packages())) {
  library(RQGIS3, quietly = T)
  ## Where to find RQGIS
  if (!is.null(root)) {
    RQGIS3::set_env(root = root)
  } else {
    RQGIS3::set_env()
  }
} else {
  stop("Install RQGIS3 to proceed")
}
## ###########################################################################

## define output file
.RData <- stringr::str_replace(overlay, ".shp", ".RData")

## keep track of temp files to remove them once not needed anymore
unlink_file <- file.path(temp.folder, "unlink.txt")
if (!file.exists(unlink_file)) {
  cat("Create", unlink_file, "\n")
  sink(unlink_file, append = F)
  sink()
}

## Loop through features
## ###########################################################################
results <- lapply(features, function(one_feature) {

  temp.file <- temp.shp()
  ## Clip to overlay extent
  ## ######################
  x <- RQGIS3::run_qgis(
    alg = "native:clip",
    INPUT = one_feature,
    OVERLAY = overlay,
    OUTPUT = temp.file,
    load_output = T,
    show_output_paths = F)
  ## ######################

  ## Recalculate geometry
  ## ####################
  temp.file2 <- temp.shp()
  if (any(names(x) == "Area")) {
    x <- RQGIS3::run_qgis(
      alg = "qgis:fieldcalculator",
      INPUT = temp.file,
      FIELD_NAME = "NEW_AREA",
      FIELD_TYPE = 0,
      FIELD_PRECISION = 3,
      FORMULA = "$area/10000", #ha
      OUTPUT = temp.file2,
      load_output = T,
      show_output_paths = F)
    out <- sum(x[["NEW_AREA"]])
  } else if (any(names(x) %in% c("length", "Length"))) {
    x <- RQGIS3::run_qgis(
      alg = "qgis:fieldcalculator",
      INPUT = temp.file,
      FIELD_NAME = "NEW_LENGTH",
      FIELD_TYPE = 0,
      FIELD_PRECISION = 3,
      FORMULA = "$length", #ha
      OUTPUT = temp.file2,
      load_output = T,
      show_output_paths = F)
    out <- sum(x[["NEW_LENGTH"]])
  } else {
    stop(paste("Unknown geometry in", one_feature))
  }

  ## discard temp files
  ## ##################
  silent <- lapply(c(temp.file, temp.file2), unlink.shp) %>%
    do.call("rbind",.) %>%
    set_rownames(1:nrow(.))

  ## if some can not be removed retry, keep them
  to_unlink <- character()
  to_unlink <- dplyr::filter(silent, success == "No")$file
  #####################

  return(out)
}) %>%
  set_names(names(features))

save(results, file = .RData)

## add to file and try to delete
if (length(to_unlink) > 0) {
  ## write to file
  sink(unlink_file, append = TRUE)
  cat(to_unlink, sep = "\n")
  sink()
  ## read file content
  to_unlink <- readr::read_delim(unlink_file, delim = "\t",
                                 skip_empty_rows = T,
                                 col_names = F,
                                 col_types = "c")[[1]]
  silent <- sapply(to_unlink, unlink)
  to_unlink <- data.frame(file = to_unlink,
                          success = ifelse(silent == 0, "Yes", "No"))
  to_unlink <- dplyr::filter(to_unlink, success == "No")$file
  ## Update
  sink(unlink_file, append = F)
  cat(to_unlink, sep = "\n")
  sink()

}

return(results)
}