I have two polygon layers, a set of aquifer polygons and a set of district boundaries. I'm trying to clip the aquifer layers to be in only the districts that I want.

From examining the other posts on the forum I saw that st_crop almost achieves what I want except that it only takes a bbox, and I found that st_intersection is again nearly the exact solution except that SF doesn't like one of the two multipolygons because the boundaries between the polygons touch so I get the Topology Exception. Does anyone know how to handle this? This seems like a pretty straightforward thing, but despite my efforts, I'm unable to overcome it.

I also thought to loop through each element of the aquifer polygons, but I wasn't sure if this would work without dissolving the adjoining district boundaries first. Either way it feels like a clumsy solution because eventually I'd have to rejoin the finished product anyways to one single SF object.

for (i in length(aquifers)){
  aquiferslist[i] <- st_intersection(aquifers[i], dist.bound)
  • hard to test without your data but perhaps you could try rmapshaper::ms_clip – sebdalgarno Jan 10 at 15:56

Your Answer

By clicking "Post Your Answer", you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.