2012-08-28 28 views
12

Daha sonra mekansal analiz için kullanmak üzere Voronoi çokgenlerini harita ile birleştirmek istiyorum. Ben birleştirmek ve daha sonra şekil dosya/uzaysal çokgenler olarak kaydetmek istiyorum nokta ve şekil dosya sayısı var. Voronoi çokgenler almak için this topic işlevini kullanıyorum.Voronoi çokgenlerini ve haritalarını birleştirin

coords<-data.frame(LONG=c(16.9252,16.9363,16.9408,16.8720,16.9167,16.9461,16.9093,16.9457,16.9171,16.8506,16.9471,16.8723,16.9444,16.9212,16.8809,16.9191,16.8968,16.8719,16.9669,16.8845), 
LAT=c(52.4064,52.4266,52.3836,52.3959,52.4496,52.3924,52.4012,52.3924,52.3777,52.4368,52.4574,52.3945,52.4572,52.3962,52.3816,52.3809,52.3956,52.3761,52.4236,52.4539)) 

Haritamın burada mevcuttur: https://docs.google.com/file/d/0B-ZJyVlQBsqlSURiN284dF9YNUk/edit

library(rgdal) 
voronoipolygons <- function(x) { 
    require(deldir) 
    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)) 
    require(sp) 
    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')))) 
} 

Ve kod voronoipolygons almak için:

pzn.coords<-voronoipolygons(coords) 
plot(pznall) 
plot(pzn.coords,add=T) 
points(coords$LONG,coords$LAT) 

sonucu: enter image description here

şöyle

Benim kodudur

Haritanımın içinde yeni bir uzaysal poligon olarak bu voronoi poligonuna sahip olmak istedim.

Anwsers'a minnettar olurum!

Düzenleme:

Açık olmak gerekirse, ben böyle bir şey (bu hatlar voronoi Polygonlardan oluşturulmalıdır) ulaşmak istiyoruz:

enter image description here

+2

Ne istediğinizi net değil. Pzn.coords nesnesinde bir spatialpolygonsdataframe olarak voronoi çokgenleri var. – Spacedman

+0

Sağladığınız bağlantı bir '.rdata' dosyasıdır. 'Save()' özelliğini kullanarak yeni bir 'myvoronietc.rdata' dosyası oluşturabilirsiniz. 'kaydet ('pzn.coords', 'lotsa_other_data', dosya = 'myvoronietc.rdata')'. –

+0

Her iki dosyayı da birleştirmek istiyorum, haritanımdaki çokgenli poligonlara sahip olmak, bu sınırların çokluğunu voronoi çokgenlerine sahip olmamak ama harita sınırlarımda olmak istiyorum. – Maciej

cevap

12

Biraz modifiye işlevi alır Ek bir uzamsal çokgenler argümanı ve bu kutuya genişler:

voronoipolygons <- function(x,poly) { 
    require(deldir) 
    if (.hasSlot(x, 'coords')) { 
    crds <- [email protected] 
    } else crds <- x 
    bb = bbox(poly) 
    rw = as.numeric(t(bb)) 
    z <- deldir(crds[,1], crds[,2],rw=rw) 
    w <- tile.list(z) 
    polys <- vector(mode='list', length=length(w)) 
    require(sp) 
    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')))) 

    return(voronoi) 

} 

Sonra yapın:

pzn.coords<-voronoipolygons(coords,pznall) 
library(rgeos) 
gg = gIntersection(pznall,pzn.coords,byid=TRUE) 
plot(gg) 

Not gg bir SpatialPolygons nesne olduğunu ve uyumsuz proj4 dizeleri hakkında bir uyarı alabilirsiniz. Proj4 dizesini nesnelerin birine veya başkalarına atamanız gerekebilir.

+0

Teşekkür ederim çok çok, istediğim buydu! – Maciej

+0

7 numaralı satırda 'pznall' ("çoktan" olmalıdır) – jbaums

+0

Bende, rw = as.numeric (t (bbox (pznall))) 'rw = as.numeric (t) ile değiştirilmelidir. (bb)) '(bbox (poli) olarak' '' '' '' '' '' '' '' '' '' '' '' '' '' '' '' '' '' '' '' 'ile kaydedilmişti. –