14

A reverse clip saves only the part of your spatial object that is outside the bounds of another object, as opposed to a regular clip which saves the parts that are inside the other object.

Performing reverse clip in ArcMap? shows how to do it in ArcMap.

How do I do this in R?

Reproducible example (on Linux machines):

system("wget 'https://github.com/Robinlovelace/Creating-maps-in-R/archive/master.zip' -P /tmp/")
unzip("/tmp/master.zip", exdir = "/tmp/master")
uk <- readOGR("/tmp/master/Creating-maps-in-R-master/data/", "ukbord")
lnd <- readOGR("/tmp/master/Creating-maps-in-R-master/data/", "LondonBoroughs")
plot(uk)
plot(lnd, add = T, col = "black")

What I want here to do is to save all of the UK except for London. Visually, I want the black shape in the resulting image to be a hole.

enter image description here

4

Answer for Simple Features:

sf package draws on Geometry Engine Open Source, and as such can access the list of commands such as st_within etc.

One such command, st_difference, will do the job:

require(sf)

# make a square simple feature
s <- rbind(c(1,1), c(1,5), c(5,5), c(5,1), c(1,1))
s.sf <-st_sfc(st_polygon(list(s)))
s.pol = st_sf(ID = "sq", s.sf)

# make a smaller triangle simple feature
s2 <- rbind(c(2,2), c(3,4), c(4,2), c(2,2))
s2.sf <-st_sfc(st_polygon(list(s2)))
s2.pol = st_sf(ID = "tr", s2.sf)

# find the 'difference', i.e. reverse of st_intersection
t <- st_difference(s.pol,s2.pol)

plot(t)

# have a look at the new geometry, a well known text format with exterior followed by hole
st_geometry(t)[[1]]
POLYGON((1 1, 1 5, 5 5, 5 1, 1 1), (2 2, 4 2, 3 4, 2 2))

also see towards the bottom of this article

can also be done by coercing Sp to sf with st_as_sf. Heed the warnings as attributes can be tricky to manage!

12

Seems to be a simple application of gDifference from the rgeos package:

> require(rgeos)
> ukhole = gDifference(uk, lnd)
Warning message:
In RGEOSBinTopoFunc(spgeom1, spgeom2, byid, id, "rgeos_difference") :
  spgeom1 and spgeom2 have different proj4 strings
> plot(ukhole)

The projection warning is because the LondonBoroughs shapefile doesn't have a .prj file.

Just to make sure that's a hole and not an outline or another solid polygon:

> gArea(lnd) + gArea(ukhole) - gArea(uk)
[1] 0
5

Bit late to the party, but there is a simple way to do this with mask using the 'inverse' argument;

ukhole <- mask(uk, lnd, inverse = TRUE)
  • From the raster package. And from sf any ideas? – RobinLovelace Aug 3 '17 at 10:39

Your Answer

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

Not the answer you're looking for? Browse other questions tagged or ask your own question.