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 acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Browse other questions tagged or ask your own question.