CODEWITHSUNDEEP
rgisintersectionr-sf
M.Teich
M.Teich

Reputation: 607

Executing sf operations on grouped objects

I am trying to find intersecting sfc_LINESTRINGs by group, where both intersecting sfc_LINESTRING objects are grouped. This question specifically relates to st_intersects, but I think the workflow will be the same for other operations within sf. Similar questions have been answered here and here, but in both cases only one of the intersecting objects was grouped and in any case I haven't been smart enough to make it work for my data set.

library(tidyverse)
library(sf)

line1 <- data.frame(lon = c(1, 1.5, 2, 2.5, 2.0, 2.5),
                lat = c(1, 1.5, 2, 2.5, 2.1, 2.6),
                grp = c("A", "A", "B", "B", "C", "C")
                ) %>% 
  st_as_sf(coords = c("lon","lat"), dim = "XY") %>% 
  group_by(grp) %>% 
  summarise(do_union = FALSE) %>% 
  st_cast("LINESTRING") %>% 
  mutate(grp = c("A","A","B"))

line2 <- data.frame(lon = c(1.2, 1.3, 2.2, 2.3),
                lat = c(1, 1.5, 2, 2.5),
                grp = c("A", "A", "B", "B")
                ) %>% 
  st_as_sf(coords = c("lon","lat"), dim = "XY") %>% 
  group_by(grp) %>% 
  summarise(do_union = FALSE) %>% 
  st_cast("LINESTRING")

plot(line1, lwd = 2, key.pos = 4, reset = FALSE)
plot(line2, lwd = 2, lty = 2, add = TRUE)

enter image description here

I would now like to find the intersecting sfc_LINESTRINGs by the relevant group. Unfortunately, the simple solution I had hoped for using group_by throws an error:

int <- line1 %>% 
  group_by(grp) %>% 
  filter(lengths(st_intersects(line1, line2)) >0)

Using my actual data set, I can get it to run, but the group_by command seems to be ignored, so that I get all intersecting objects, not just those which match by group. To be specific, in this case the result should only select one sfc_LINESTRING for each group, i.e. matched continuous and dashed green line and matched continuous and dashed red line. The dashed red line should NOT create an intersection with the green line. Any suggestions will be much appreciated.

Upvotes: 1

Views: 242

Answers (1)

lovalery
lovalery

Reputation: 4652

Not completely sure if I understood your request correctly, but I am giving it a try! So, please find below one possible solution using the st_join() function of the sf library, and then some dplyr functions to get the final result.

Reprex

  • Code
library(dplyr)
library(sf)

results <- line1 %>% 
  st_join(., line2, join = st_intersects) %>% 
  filter(., grp.x == grp.y) %>% 
  rename(grp = grp.x) %>% 
  select(-grp.y)
  • Output
results
#> Simple feature collection with 2 features and 1 field
#> Geometry type: LINESTRING
#> Dimension:     XY
#> Bounding box:  xmin: 1 ymin: 1 xmax: 2.5 ymax: 2.6
#> CRS:           NA
#> # A tibble: 2 x 2
#>   grp           geometry
#>   <chr>     <LINESTRING>
#> 1 A       (1 1, 1.5 1.5)
#> 2 B     (2 2.1, 2.5 2.6)
  • Visualization
plot(results)

PS: In case you want to get the two intersecting lines in each group within the same sf object, you can extend the code as follows:

  • Complementary code
results2 <- results %>% 
  st_union(.,line2, by_feature = TRUE) %>% 
  filter(., grp == grp.1) %>% 
  select(-grp.1)
#> Warning: attribute variables are assumed to be spatially constant throughout all
#> geometries
  • Output 2
results2
#> Simple feature collection with 2 features and 1 field
#> Geometry type: MULTILINESTRING
#> Dimension:     XY
#> Bounding box:  xmin: 1 ymin: 1 xmax: 2.5 ymax: 2.6
#> CRS:           NA
#> # A tibble: 2 x 2
#>   grp                                                                   geometry
#>   <chr>                                                        <MULTILINESTRING>
#> 1 A     ((1 1, 1.25 1.25), (1.25 1.25, 1.5 1.5), (1.2 1, 1.25 1.25), (1.25 1.25~
#> 2 B     ((2 2.1, 2.275 2.375), (2.275 2.375, 2.5 2.6), (2.2 2, 2.275 2.375), (2~
  • Visualization
plot(results2)

Created on 2022-01-04 by the reprex package (v2.0.1)

Upvotes: 1

Related Questions