How to extract coordinates values from a shapefile?
On Wednesday 09 June 2010, Hadley Wickham wrote:
On Wed, Jun 9, 2010 at 3:24 PM, Edzer Pebesma <edzer.pebesma at uni-muenster.de> wrote:
The example provided by Matt assumes that each polygon consists of a
single ring, and doesn't have islands, lakes etc. The function below
pasts all coordinates to a single 2-column matrix. For clarity's sake, I
avoided nested sapply's.
library(maptools)
xx <- readShapePoly(system.file("shapes/sids.shp",
package="maptools")[1], IDvar="FIPSNO", proj4string=CRS("+proj=longlat
+ellps=clrk66")) allcoordinates = function(x) {
? ?ret = NULL
? ?polys = x at polygons
? ?for(i in 1:length(polys)) {
? ? ? ?pp = polys[[i]]@Polygons
? ? ? ?for (j in 1:length(pp))
? ? ? ? ? ?ret = rbind(ret, coordinates(pp[[j]]))
? ?}
? ?ret
}
q = allcoordinates(xx)
# check:
plot(xx)
lines(q, col='blue')
And that's basically what fortify does, except it covers a few more cases. Hadley
Not to be nit-picky, but wouldn't it be safer to pre-allocated "ret" instead of dynamically appending? This would allow the suggested function to scale to very large geometries. Cheers, Dylan
Dylan Beaudette Soil Resource Laboratory http://casoilresource.lawr.ucdavis.edu/ University of California at Davis 530.754.7341