Is there an way to perform the equivalent of the identity tool from ArcGIS where all boundaries of spatial data x are retained and areas intersecting with spatial data y are updated with attributes and border from y.
Answering in case someone else running into this issue. One of my coworkers helped me get to the solution of binding the results of the intersection of the two layers, with the difference of the focal layer and intersecting layer.
arc.ident <- function(layer_a, layer_b){
int_a_b <- st_intersection(layer_a, layer_b)
rest_of_a <- st_difference(layer_a, st_union(layer_b))
output <- bind_rows(int_a_b, rest_of_a)
return(st_as_sf(output))
}
or as a tidyverse pipe
arc.ident.output <- st_intersection(layer_a, layer_b) %>%
bind_rows(st_difference(layer_a, st_union(layer_b)))