[R-sig-Geo] Creating polygons out of an irregular data set to represent specific classes
Sébastien Durand
v8extra at gmail.com
Fri Dec 3 06:03:29 CET 2010
Ok here is my latest development on the topic! But I still have a question! At the end I still get segmented polygons!
library(akima)
library(maptools)
library(raster)
data(akima)
# define de dimension of grid
dimgrid=50
# Interpolate to regular grid
ak.li=interp(akima$x, akima$y, akima$z, xo=seq(min(akima$x), max(akima$x), length = dimgrid),yo=seq(min(akima$y), max(akima$y), length = dimgrid),linear = TRUE, extrap=FALSE, duplicate = "mean")
# Show interpolation
image(ak.li)
points(akima)
with(akima, text(x, y, formatC(z,dig=2), adj = -0.1))
# Create a raster with the data
r <- raster(nrow=length(ak.li$y), ncol=length(ak.li$x), crs=NA)
xmin(r) <- min(ak.li$x, na.rm=TRUE)
xmax(r) <- max(ak.li$x, na.rm=TRUE)
ymin(r) <- min(ak.li$y, na.rm=TRUE)
ymax(r) <- max(ak.li$y, na.rm=TRUE)
tmp=ak.li$z
tmp=apply(tmp, 1, rev)
r[] <- tmp
# Show the raster
plot(r)
# Convert the data found between ]20 and 25[ into polygons
pol <- rasterToPolygons(r, fun=function(x){x>20 & x<25})
# Show the polygons
plot(pol, add=T, col='red')
# Attempt to merge similar adjacent into one polygon
union = unionSpatialPolygons(pol, ID=rep(1,
times=length(pol at polygons)))
# Show
plot(union, add=TRUE, col="red")
All polygons were not merged into a single one and I don't know how to fix that!
Any help?
More information about the R-sig-Geo
mailing list