I have a SpatialLinesDataFrame representing streets. I need to find the intersection of 2 streets.
I'll copy a toy example from How to get the intersection point of two vector?
library(rgdal)
library(rgeos)
a = c(1,5,2,6,3,6,3,5,7)
b = c(5,3,5,7,2,6,9,3,6)
SL1 <- SpatialLines(list(Lines(Line(cbind(seq_along(a),a)), "A")))
SL2 <- SpatialLines(list(Lines(Line(cbind(seq_along(b),b)), "B")))
# Build my own SpatialLines with both objects, couldn't find a way
# to build it from SL1 and SL2
SL <- SpatialLines(list(Lines(Line(cbind(seq_along(a),a)), "A"), Lines(Line(cbind(seq_along(b),b)), "B")))
SL.df <- SpatialLinesDataFrame(SL, data=data.frame(OBJECTID=paste(1:2), NAME=c("st 1", "av B"), row.names=row.names(SL)))
From the answer to the question I copied, I know I can get the intersection as
# Find intersections
coords <- coordinates(gIntersection(SL1, SL2))
But suppose I have read SL.df with readOGR and don't have access to the original components. How should I extract individual streets?
After quite some wandering, I've come up with
st.1 <- SpatialLines(SL.df@lines[which(row.names(SL.df) == row.names(SL.df)[which(SL.df@data$NAME == 'st 1')])])
av.B <- SpatialLines(SL.df@lines[which(row.names(SL.df) == row.names(SL.df)[which(SL.df@data$NAME == 'av B')])])
coords <- coordinates(gIntersection(st.1, av.B))
It works, but the extraction of the individual streets doesn't seem very elegant. Is there a better way to do it?
Thanks!
A couple of options:
coordinates(gIntersection(subset(SL.df, NAME=="st 1"),
subset(SL.df, NAME=="av B")))
# x y
# 1 1.666667 3.666667
# 1 2.400000 3.800000
# 1 4.500000 4.500000
# 1 6.000000 6.000000
# 1 7.750000 4.500000
coordinates(gIntersection(SL.df[SL.df$NAME=="st 1",],
SL.df[SL.df$NAME=="av B",]))
# x y
# 1 1.666667 3.666667
# 1 2.400000 3.800000
# 1 4.500000 4.500000
# 1 6.000000 6.000000
# 1 7.750000 4.500000