library(sp)
library(dismo)
# spatial data (political districts of Austria)
gadm <- getData('GADM', country = "AT", level = 2)
# view
plot(gadm)
# some addresses
pts <- geocode(c("Aldrans, Grubenweg", "Wien, Stephansdom", "Salzburg, Mozartplatz"))
# make pts spatial
coords <- SpatialPoints(pts[, c("longitude", "latitude")])
spdf_pts <- SpatialPointsDataFrame(coords, pts)
# assign CRS/projection (which is WGS 1984)
crs <- CRS(" +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0")
proj4string(spdf_pts) <- crs
# check data
str(spdf_pts)
str(gadm)
# plot pts on top
plot(spdf_pts, cex = 2, col = 2, add = T)
# do an intersection (points in polygon)
# yielding the polygon's attribute data
over(spdf_pts, gadm)
20 Jun 2013
Spatial Overlays with R - Retrieving Polygon Attributes for a Set of Points
A short tutorial for spatial overlays using R-GIS..
Subscribe to:
Post Comments
(
Atom
)
FYI I need to change the code to
ReplyDeletecoords <- SpatialPoints(pts[, c("lon", "lat")])
to make it work. Very cool though!