R/mod07_add_hedgerows.R

In ecoservR/ecoserv_tool: Natural Capital Mapping Toolkit (UK)

#######################################
### Worflow functions               ###
### for EcoservR                    ###
### Sandra Angers-Blondin           ###
### 22 April 2021                   ###
#######################################

#' Add Hedgerows
#'
#' This function adds hedgerow data to the basemap. As opposed to the other basemap processing functions that extract information into existing polygons, this function creates new geometries within the basemap. Data input for hedgerows must be linear features.

#' @param mm The mm object loaded in the environment, can be at various stages of updating.
#' @param studyAreaBuffer The buffered study area generated during mod01 or reloaded when resuming a session.
#' @param projectLog The RDS project log file generated by the wizard app and containing all file paths to data inputs and model parameters
#' @return Saves a project_title_MM_07.RDS file to project folder
#' @export

add_hedgerows <- function(mm = parent.frame()$mm,
                          studyAreaBuffer = parent.frame()$studyAreaBuffer,
                          projectLog = parent.frame()$projectLog){

   timeA <- Sys.time() # start time


   ## Extract the file paths and other info from project log ----------------------

   output_temp <- projectLog$output_temp
   title <- projectLog$title
   scratch_path <- file.path(output_temp, "ecoservR_scratch")
   if (!dir.exists(scratch_path)) dir.create(scratch_path)

   # Get path
   hedgepath <- projectLog$df[projectLog$df$dataset == "hedgerows", ][["path"]]  # path to hedge data, if available
   hedgetype <- guessFiletype(hedgepath)

   if (!is.na(hedgepath) & !is.null(hedgepath)){
      message("Preparing to update baseline with hedgerow data...")


      # DATA IMPORT ---------------------------------------------------------------------------------

      # Import hedge data

      hedge <- loadSpatial(hedgepath,
                           filetype = hedgetype)  # using the loading into list function in case there are many layers

      # Tidy up the files so we can combine them together

      hedge <- lapply(hedge, function(x) sf::st_zm(x, drop = TRUE) %>% # remove Z dimension if present
                         checkcrs(., studyAreaBuffer) %>% # check proj
                         sf::st_geometry() %>% # drop attributes
                         sf::st_union() %>%
                         sf::st_as_sf() %>%  #make sure format is sf df
                         sf::st_cast(to="MULTILINESTRING") %>% sf::st_cast(to = "LINESTRING")) #

      hedge <- do.call(rbind, hedge) %>% sf::st_as_sf()  # putting back into one single sf object

      ## make sure the geometry column has the same name in both hedges and mm
      hedge <- ecoservR::rename_geometry(hedge, "geometry")
      mm <- lapply(mm, function(x) ecoservR::rename_geometry(x, "geometry"))


      # DATA PREP -----------------------------------------------------------------------------------

      message("Creating hedgerows from linear features")
      # Keep only the geometry column
      # and buffer to create polygons from linear features
      hedge <- sf::st_geometry(hedge) %>%
         sf::st_buffer(1.25, endCapStyle = "ROUND", joinStyle = "MITRE")
      # rounded ends give best results when unioning, even if add an extra meter to the end of the hedge


      ## Grid the polygonized hedgerows so we have named tiles that match mm

      SAgrid <- ecoservR::grid_study(studyAreaBuffer) # create grid
      SAgrid <- SAgrid[names(SAgrid) %in% names(mm)]
      hedge <- lapply(SAgrid, function(x)  # tile hedges
         sf::st_intersection(hedge, sf::st_geometry(x)) %>%
            ecoservR::checkgeometry(.) %>%
            sf::st_as_sf())

      hedge <- hedge[sapply(hedge, function(x) nrow(x) > 0)] # remove empty tiles

      rm(SAgrid)

      # The union is necessary to remove overlap (takes a while)

      message("Unioning hedgerows to remove overlaps")


      hedge <- lapply(hedge, function(x){

         sf::st_union(x) %>%
            sf::st_as_sf() %>%
            ecoservR::checkgeometry(.) # union and convert to single-part poly
      })



      # Trim hedges -------------------------------------------------------------

      # Trim hedges so that the buffer doesn't encroach on structures (buildings and roads)

      for (i in 1:length(mm)){

         tilename <- names(mm)[[i]]

         if (!tilename %in% names(hedge)){next}


         message("Reshaping hedgerows, tile ", tilename)

         # Get index of manmade features for that tile
         index <- (mm[[i]][["Make"]] == "Manmade" | mm[[i]][["Theme"]] == "Water")

         if (any(index)){  # if we detected some of those features

            # mask them from the hedge tile
            hedge[[tilename]] <- rmapshaper::ms_erase(hedge[[tilename]], mm[[i]][index,]) %>%
               ecoservR::checkgeometry(.)
         }
      }



      # Remove land from mm -----------------------------------------------------

      ## In each tile, identify which mastermap polygons intersect the hedges.
      ## Remove the shape of the hedges from them

      for (i in 1:length(mm)){

         tilename <- names(mm)[[i]]

         if (!tilename %in% names(hedge)){next}

         message("Preparing basemap to receive hedgerows, tile ", tilename)

         # creating indices so we subset out the smaller amount of data possible
         indexmm <- which(lengths(sf::st_intersects(mm[[i]], hedge[[tilename]]))>0)  # to subset mm poly intersecting


         if (length(indexmm) > 0){ # only erase if there are features

            indexh <- which(lengths(sf::st_intersects(hedge[[tilename]], mm[[i]][indexmm,]))>0)

            # create a subset with the revised polygons
            mm_erased <- rmapshaper::ms_erase(mm[[i]][indexmm,], hedge[[tilename]][indexh,]) %>%
               checkgeometry(.)

            # remove the original polys from df
            mm[[i]] <- mm[[i]][-indexmm,]

            # and merge the new set in
            mm[[i]] <- rbind(mm[[i]], mm_erased)
            rm(mm_erased)

         }
      }


      # Burn in hedges ----------------------------------------------------------
      message("Burning hedges into basemap")

      ## Bind hedges into corresponding mm list element

      for (i in 1:length(mm)){

         tilename <- names(mm)[i]

         if (!tilename %in% names(hedge)){next}

         # Make extra sure the geometry cols have the same name (caused problems before)

         if (attr(hedge[[tilename]], "sf_column") != attr(mm[[tilename]], "sf_column")){

            hedge[[tilename]] <- ecoservR::rename_geometry(hedge[[tilename]],
                                                           attr(mm[[tilename]], "sf_column"))
         }


         ## Add columns to hedges to correspond with mm attributes

         mmcols <- setdiff(names(mm[[i]]), names(hedge[[tilename]]))  # mm cols not present in hedge col

         hedge[[tilename]]$rmapshaperid <- NULL # drop rmapshaper column
         hedge[[tilename]][mmcols] <- NA  # ready for binding

         # Populate attributes
         hedge[[tilename]] <- hedge[[tilename]] %>%
            dplyr::mutate(Term = "Hedgerow",
                          Group = "Natural Environment",
                          Make = "Natural",
                          Theme = "Land"
            )





         mm[[i]] <- rbind(mm[[i]], hedge[[tilename]])
      }



      suppressWarnings(rm(indexmm, indexh, i, mmcols, hedge, tilename) )


      # Validate geometries one last time

      mm <- checkgeometry(mm, "POLYGON")

      # SAVE UPDATED MASTER MAP ---------------------------------------------------------------------

      saveRDS(mm, file.path(output_temp, paste0(title, "_MM_07.RDS")))

      # Update the project log with the information that map was updated

      projectLog$last_success <- "MM_07.RDS"

      timeB <- Sys.time() # stop time

      # add performance to log
      projectLog$performance[["add_hedges"]] <- as.numeric(difftime(
         timeB, timeA, units="mins"
      ))

      updateProjectLog(projectLog) # save revised log

      # and delete contents of scratch folder
      cleanUp(scratch_path)


      message(paste0("Finished updating with hedgerow data. Process took ",
                     round(difftime(timeB, timeA, units = "mins"), digits = 1),
                     " minutes."))



   } else {message("No hedgerow data input specified.")} # end of running condition


   # Return mm to environment, whether it has been updated or not.
   return({
      invisible({
         mm <<- mm
         projectLog <<- projectLog
      })
   })

} # end of function