2012-02-22 17 views
7

Tengo varios conjuntos de puntos (para diferentes años ~ 20)¿Cómo crear polígonos Thiessen a partir de puntos utilizando paquetes R?

Quiero generar polígonos de Thiessen para cada conjunto de puntos usando paquetes espaciales r.

Sé que esto puede hacerse utilizando SIG, pero como yo quiero un proceso por lotes algo en I habría

útil.

+2

'install.packages ("SOS"); biblioteca ("sos"); findFn ("thiessen") ' –

+0

Estoy usando ArcGIS para lo mismo ahora ... – user2760

Respuesta

16

No nos ha dado acceso a sus datos, pero aquí hay un ejemplo de los puntos que representan las ciudades del mundo, utilizando un enfoque descrito por Carson Farmer en his blog. De esperar que va a comenzar ...

# Carson's Voronoi polygons function 
voronoipolygons <- function(x) { 
    require(deldir) 
    require(sp) 
    if (.hasSlot(x, 'coords')) { 
    crds <- [email protected] 
    } else crds <- x 
    z <- deldir(crds[,1], crds[,2]) 
    w <- tile.list(z) 
    polys <- vector(mode='list', length=length(w)) 
    for (i in seq(along=polys)) { 
    pcrds <- cbind(w[[i]]$x, w[[i]]$y) 
    pcrds <- rbind(pcrds, pcrds[1,]) 
    polys[[i]] <- Polygons(list(Polygon(pcrds)), ID=as.character(i)) 
    } 
    SP <- SpatialPolygons(polys) 
    voronoi <- SpatialPolygonsDataFrame(SP, data=data.frame(x=crds[,1], 
    y=crds[,2], row.names=sapply(slot(SP, 'polygons'), 
    function(x) slot(x, 'ID')))) 
} 

Ejemplo 1: Entrada es un SpatialPointsDataFrame:

# Read in a point shapefile to be converted to a Voronoi diagram 
library(rgdal) 
dsn <- system.file("vectors", package = "rgdal")[1] 
cities <- readOGR(dsn=dsn, layer="cities") 

v <- voronoipolygons(cities) 

plot(v) 

Voronoi diagram of cities

Ejemplo 2: Entrada es vectores de x, y coordenadas:

dat <- data.frame(x=runif(100), y=runif(100)) 
v2 <- voronoipolygons(dat) 
plot(v2) 

Another voronoi diagram

+0

He modificado la función para que acepte vectores de coordenadas (esperados en orden x, y) así como un' SpatialPointsDataFrame'. – jbaums

+0

Esto funcionó..TY! – user2760

+0

Me alegra oírlo. – jbaums

1

El mismo principio que se muestra por jbaums, pero el código más simple:

library(dismo) 
library(rgdal) 
cities <- shapefile(file.path(system.file("vectors", package = "rgdal")[1], "cities")) 

v <- voronoi(cities) 
plot(v) 
Cuestiones relacionadas