R/tangler.R

In tangles: Anonymisation of Spatial Point Patterns and Grids

## Function for tangling original coordinates using an existing transformation plan
##
## Applies an existing transformation sequence (from a tanglerInfo object) to a new dataset.
## Useful for tangling related data to match previously Anonymised datasets.
##
## Inputs:
## - data: Either a 2-column matrix, data.frame, sf POINT object, or terra raster
## - tanglerInfo: Detangler object previously generated by tangles()
## - raster_object: Logical; if TRUE, the data is treated as a raster object
## - stub: Character; optional identifier for file naming
## - saveTangles: Logical; if TRUE, output is saved as RDS (and GeoTIFF for raster)
## - exportShapefile: Logical; if TRUE, exports tangler output as shapefile (no CRS)
## - path: Directory for saving output

tangler <- function(
    data = NULL,
    tanglerInfo = NULL,
    raster_object = FALSE,
    stub = NULL,
    saveTangles = FALSE,
    exportShapefile = FALSE,
    path = NULL
) {
  # If the user wants to save outputs but didn't supply a path,
  # default into tempdir()
  if ((saveTangles || exportShapefile) && is.null(path)) {
    path <- tempdir()}
  
  ## Check for raster-unsafe transformations
  if (raster_object && any(!(tanglerInfo$unpicker$degree %in% c(0, 90, 180, 270, NA)))) {
    stop("Rotation angles in detangler object are not suitable for raster alignment (must be multiples of 90 degrees).")
  }
  
  ## Extract coordinates
  if (raster_object) {
    tempD <- data.frame(cellNos = seq_len(terra::ncell(data)))
    vals <- as.data.frame(terra::values(data))
    tempD <- cbind(tempD, vals)
    cellNos <- tempD$cellNos
    gXY <- data.frame(terra::xyFromCell(data, cellNos))
    xyData <- as.matrix(gXY)
  } else if (inherits(data, "sf")) {
    xyData <- sf::st_coordinates(data)
  } else {
    xyData <- data
  }
  
  ###### Internal step functions (reapply transformations)
  
  leap_Xr <- function(xyData, r.num) {
    xyData[, 1] <- xyData[, 1] + r.num
    xyData
  }
  
  leap_Yr <- function(xyData, r.num) {
    xyData[, 2] <- xyData[, 2] + r.num
    xyData
  }
  
  rotate_XYr <- function(xyData, deg, origin.point) {
    x <- t(xyData[, 1])
    y <- t(xyData[, 2])
    v <- rbind(x, y)
    
    center <- v
    center[1, ] <- origin.point[1]
    center[2, ] <- origin.point[2]
    
    theta <- (deg * pi) / 180
    R <- matrix(c(cos(theta), -sin(theta), sin(theta), cos(theta)), nrow = 2)
    
    s <- v - center
    so <- R %*% s
    vo <- so + center
    
    xyData <- cbind(vo[1, ], vo[2, ])
    xyData
  }
  
  ###### Apply transformations in forward order
  
  for (i in seq_len(nrow(tanglerInfo$unpicker))) {
    seq.step <- tanglerInfo$unpicker$step[i]
    
    if (seq.step == 1) {
      xyData <- leap_Xr(xyData, tanglerInfo$unpicker$leap_dist[i])
    } else if (seq.step == 2) {
      xyData <- leap_Yr(xyData, tanglerInfo$unpicker$leap_dist[i])
    } else if (seq.step == 3) {
      xyData <- rotate_XYr(
        xyData,
        deg = tanglerInfo$unpicker$degree[i],
        origin.point = c(tanglerInfo$unpicker$origin_X[i],
                         tanglerInfo$unpicker$origin_Y[i])
      )
    }
  }
  
  xyData <- as.data.frame(xyData)
  names(xyData) <- c("X", "Y")
  hash.out <- tanglerInfo$hash
  
  ## Reconstruct raster if needed
  if (raster_object) {
    tDat <- cbind(xyData, tempD)
    
    if (ncol(tDat) > 4) {
      value_cols <- names(tDat)[-(1:3)]
      r_list <- lapply(value_cols, function(col) {
        xyz <- tDat[, c("X", "Y", col)]
        names(xyz) <- c("x", "y", "value")
        terra::rast(x = xyz, type = "xyz")
      })
      rasterOuts <- do.call(c, r_list)
      names(rasterOuts) <- value_cols
    } else {
      xyz <- tDat[, c(1, 2, 4)]
      names(xyz) <- c("x", "y", "value")
      rasterOuts <- terra::rast(x = xyz, type = "xyz")
    }
    
    if (saveTangles) {
      saveRDS(rasterOuts, file = file.path(path, paste0("tangledXY_raster_", hash.out, ".rds")))
      for (i in seq_len(terra::nlyr(rasterOuts))) {
        rz <- rasterOuts[[i]]
        out.name <- file.path(path, paste0("tangledXY_raster_", names(rasterOuts[[i]]), "_", hash.out, ".tif"))
        terra::writeRaster(rz, filename = out.name, overwrite = TRUE)
      }
    }
    
    return(list(rasterOuts))
    
  } else {
    ## Save output as RDS
    if (saveTangles) {
      rds_out <- file.path(path, paste0("tanglerXY_", stub, "_", hash.out, ".rds"))
      saveRDS(xyData, file = rds_out)
    }
    
    ## Optionally export to shapefile
    if (exportShapefile) {
      pts_sf <- sf::st_as_sf(xyData, coords = c("X", "Y"))
      sf::st_crs(pts_sf) <- NA
      shp_out <- file.path(path, paste0("tanglerXY_", stub, "_", hash.out, ".shp"))
      suppressWarnings(suppressMessages(
        sf::st_write(pts_sf, shp_out, delete_layer = TRUE, quiet = TRUE)
      ))
    }
    
    return(xyData)
  }
}