14

後で空間分析に使用するために、ボロノイポリゴンをマップと組み合わせたいと思います。結合してシェープファイル/空間ポリゴンとして保存したいポイントとシェープファイルの数があります。ボロノイポリゴンを取得するには、このトピックの関数を使用します。

私のコードは次のとおりです。

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))

私の地図はここから入手できます:https ://docs.google.com/file/d/0B-ZJyVlQBsqlSURiN284dF9YNUk/edit

library(rgdal)
voronoipolygons <- function(x) {
  require(deldir)
  if (.hasSlot(x, 'coords')) {
    crds <- x@coords  
  } 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'))))
}

そして、voronoipolygonsを取得するための私のコード:

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

結果: ここに画像の説明を入力してください

このボロノイポリゴンを新しい空間ポリゴンとしてマップ内に配置したいと思います。

私は答えに感謝します!

編集:

明確にするために、私はこのようなことを達成したいと思います(この線はボロノイポリゴンから作成する必要があります):

ここに画像の説明を入力してください

4

1 に答える 1

12

少し変更された関数は、追加の空間ポリゴン引数を取り、そのボックスに拡張します:

voronoipolygons <- function(x,poly) {
  require(deldir)
  if (.hasSlot(x, 'coords')) {
    crds <- x@coords  
  } 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)

}

次に、次のようにします。

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

ggは SpatialPolygons オブジェクトであり、proj4 文字列の不一致に関する警告が表示される場合があることに注意してください。proj4 文字列をいずれかのオブジェクトに割り当てる必要がある場合があります。

于 2012-08-28T13:02:49.143 に答える