The context: I have a border and routes followed by vehicles. I want to
- decompose the routes, cutting them everytime they cross a border, and
- for each part, know what are the distances performed.
Simple reproducible example, part 1: input data: as input, one square as polygon for the border, three routes as linestrings:
library(sf) library(lwgeom) # Defining the "border", a polygon four_corners <- matrix(c(0, 0, 1, 0, 1, 1, 0, 1, 0, 0), ncol=2, byrow=TRUE) pts <- list(four_corners) my_polygon <- st_polygon(pts) plot(my_polygon) # Defining the first "route", a linestring three_points <- matrix(c(-1, 0, 0.5, 0.5, 2, 0), ncol=2, byrow=TRUE) my_linestring1 <- st_linestring(three_points) print(my_linestring1) plot(my_linestring1, add=TRUE) # Defining the second "route", a linestring two_points <- matrix(c(0.5, -0.5, 0.5, 0.25), ncol=2, byrow=TRUE) my_linestring2 <- st_linestring(two_points) print(my_linestring2) plot(my_linestring2, add=TRUE) # Defining the third "route", a linestring four_points <- matrix(c(0.25, 0.75, 0.25, 0.6, 0.75, 0.6, 0.75, 0.75), ncol=2, byrow=TRUE) my_linestring3 <- st_linestring(four_points) print(my_linestring3) plot(my_linestring3, add=TRUE) # Group the three routes in a spatial dataframe with ID g <- st_sfc(my_linestring1, my_linestring2, my_linestring3) df_routes <- data.frame(id_column = 1:3) st_geometry(df_routes) <- g
df_routes looks like this:
id_column geometry 1 1 LINESTRING (-1 0, 0.5 0.5, ... 2 2 LINESTRING (0.5 -0.5, 0.5 0... 3 3 LINESTRING (0.25 0.75, 0.25...
With one internal route, one route crossing the border once, and one route crossing the border twice:
Simple reproducible example, part 2: the solution
# Split the route in different subroutes cut by the border splitted_routes <- st_split(df_routes, my_polygon) # Decompose the GEOMETRYCOLLECTION of LINESTRINGs in a list of LINESTRINGs splitted_routes_as_single_linestring <- st_collection_extract(splitted_routes, "LINESTRING") # Get the length of the LINESTRING list_of_length = st_length(splitted_routes_as_single_linestring) # Remove the geometry, in order to get a dataframe as output splitted_routes_as_single_linestring$geometry <- NULL # Just merge the length in a column distance in the original dataframe ID_and_length <- cbind(splitted_routes_as_single_linestring, distance = list_of_length)
This approach works on this simple example and generates the distances of each part on each side of the border:
id_column distance 1 1 1.054093 1.1 1 1.054093 1.2 1 1.054093 2 2 0.500000 2.1 2 0.250000 3 3 0.800000
In particular, the 3rd route is internal and the output is simply one distance corresponding to the complete distance of the original route.
Now, when using this approach on a real life example, some internal routes are decomposed in small pieces. They should be decomposed as one, but somehow, my approach is not working. When I am extracting LINESTRINGS out of the GEOMETRYCOLLECTIONS coming from st_split, I got several objects, which was not my goal/expectation.