# Clipping lines from intersecting polygons with sf in R

I am attempting to clip sf lines where they overlap with polygons in R - the opposite of the clip in the image below:

I've created a reprex below, where the blue lines outside of the circle should be what is left after the clip.

``````library(sf)
#> Linking to GEOS 3.10.2, GDAL 3.4.2, PROJ 8.2.1; sf_use_s2() is TRUE
library(tmap)

begin.coords = data.frame(lon = -78,
lat = 43.45)
end.coords = data.frame(lon = -78,
lat = 43.55)

l_sf2 <- vector("list", nrow(begin.coords))
for (i in seq_along(l_sf2)){
l_sf2[[i]] <- st_linestring(as.matrix(rbind(begin.coords[i, ], end.coords[i,])))
}

line_df = st_sfc(l_sf2) %>%
st_as_sf() %>%
st_set_crs(.,4326) %>%
st_transform(.,crs=32618)

poly_df = data.frame(lon = -78,
lat = 43.5) %>%
st_as_sf(coords=c('lon','lat')) %>%
st_set_crs(.,4326) %>%
st_transform(.,crs=32618) %>%
st_buffer(.,dist=4000)

tm_shape(line_df)+
tm_lines(col='blue')+
tm_shape(poly_df)+
tm_polygons(col='red',alpha=.5)
``````

Created on 2023-07-12 with reprex v2.0.2

I understand this process is likely not an easy one, given that line geometry is defined by starting and ending points, thus, such a clip would have to create separate lines with new starting and ending points.

EDIT:

After incorporating Spacedman's answer, it seems that the real issue occurs when there are multiple linestrings that change direction. In the reprex below, it appears as though `st_difference()` looks to clip lines only when they are overlapping locations where all of the polygons overlap:

``````library(sf)
library(tmap)
library(dplyr)

df1 = data.frame(lon=c(-78,-77.95,-77.9,-77.85,-77.8),
lat=c(43.45,43.5,43.55,43.5,43.45),
id = c(1:5)) %>%
st_as_sf(coords=c('lon','lat'),crs=4326) %>%
st_transform(crs=32618)

df2 = data.frame(lon=c(-78,-77.95,-77.9,-77.85,-77.8),
lat=c(43.5,43.45,43.5,43.45,43.5),
id=c(1:5)) %>%
st_as_sf(coords=c('lon','lat'),crs=4326) %>%
st_transform(crs=32618)

df1 %>%
summarise(do_union = F) %>%
st_cast("LINESTRING") -> sf_line

df2 %>%
summarise(do_union = F) %>%
st_cast("LINESTRING") -> sf_line2

sf_bothlines = rbind(sf_line,sf_line2)

poly_df = data.frame(lon = c(-77.9,-77.95),
lat = c(43.5,43.5)) %>%
st_as_sf(coords=c('lon','lat')) %>%
st_set_crs(.,4326) %>%
st_transform(.,crs=32618) %>%
st_buffer(.,dist=4000)

tm_shape(sf_bothlines)+
tm_lines(col='blue')+
tm_shape(poly_df)+
tm_polygons(col='red',alpha=.5)
``````

``````
clipped=st_difference(sf_bothlines,poly_df)

tm_shape(clipped)+
tm_lines(col='blue')+
tm_shape(poly_df)+
tm_polygons(col='red',alpha=.5)
``````

Created on 2023-07-12 with reprex v2.0.2

## 1 Answer

Use `st_difference`:

``````> dd = st_difference(line_df, poly_df)
> plot(dd)
``````

Note that for this to work as a "clip" to the polygon you need the polygon as a single feature, otherwise the lines will get clipped to each feature in turn, possibly resulting in many lines inside feature A because they were the bits resulting from the difference with feature B. To fix this, unite all the polygons into a single feature using `st_union`.

So given your object with two features:

``````> poly_df
Simple feature collection with 2 features and 0 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 257503.6 ymin: 4816427 xmax: 269546.1 ymax: 4824570
Projected CRS: WGS 84 / UTM zone 18N
geometry
1 POLYGON ((269546.1 4820427,...
2 POLYGON ((265503.6 4820570,...
``````

The union creates a single combined feature:

``````> st_union(poly_df)
Geometry set for 1 feature
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 257503.6 ymin: 4816427 xmax: 269546.1 ymax: 4824570
Projected CRS: WGS 84 / UTM zone 18N
POLYGON ((269546.1 4820427, 269540.6 4820218, 2...
``````

Then the difference of the lines with the union is:

``````> plot(st_difference(sf_bothlines, st_union(poly_df)))
> plot(poly_df, add=TRUE, col="#80808080")
``````

• Weird. I swear this wasn't working when I tried earlier. ¯_(ツ)_/¯ Commented Jul 12, 2023 at 19:56
• I do get some weird behavior when doing this for curved lines though - 3,130 lines becomes over 466,000 lines. Is there a reason for this or a way around it? Commented Jul 12, 2023 at 20:04
• It could happen, depends how often your lines cut the polygon, they could could get chopped up into hundreds of pieces. See what happens to a small sample of those 3130. Commented Jul 12, 2023 at 21:02
• It's behaving very strangely. It's still creating large numbers - 31 lines is turning into 4,619, but they are still overlapping and aren't being clipped. I'll try to create another reprex. Commented Jul 13, 2023 at 0:13
• I made a major edit after figuring out what the problem is. Since this answered the original question, I will leave it as accepted until another answer that addresses the edit is posted. Commented Jul 13, 2023 at 0:48