R/GIS: How to subset a shapefile by a lat-long bounding box?

Question

I want to subset out a shapefile (the .shp and associated files are here) into another one bounded by a set of coordinates, say between longs [80,90] and lats [20,30], and then write this out as another shapefile. If I use the maptools package:

df = readShapeLines("/path/asia_rivers.shp")

and then look at the structure of the file with as.data.frame(df), I can't find any obvious way of subsetting by coordinates. I can use the PBSmapping package to subset:

df = importShapefile("/path/asia_rivers.shp")
df_sub = subset(df, X>=80 & X<=90 & Y >=20 & Y <=30)

but then I can't seem to be able to coerce this into a SpatialLines data frame which can be exported via writeSpatialShape() in maptools. I keep getting this error: Error in PolySet2SpatialLines(df_sub) : unknown coordinate reference system. Surely I am missing something very basic and there should be an easy way of subsetting geo-data by geo-coordinates?

Answer 1

You could try the following:

library(rgeos)
rivers <- readWKT("MULTILINESTRING((15 5, 1 20, 200 25), (-5 -8,-10 -8,-15 -4), (0 10,100 5,20 230))")
bbx <- readWKT("POLYGON((0 40, 20 40, 20 0, 0 0, 0 40))") 

rivers.cut <- gIntersection(rivers, bbx)

plot(rivers, col="grey")
plot(bbx, add=T, lty=2)
plot(rivers.cut, add=T, col="blue")