5

SpatialPolygonsDataFrameオブジェクトとしてreadOGR（）を使用してRに読み込んだ2つのシェープファイルがあります。どちらも、内部の境界が異なるニュージーランドの地図です。1つには、領土の権限の境界を表す約70のポリゴンがあります。もう1つは、面積単位を表す約1900です。

私の目的は、より大きなプロジェクトの厄介な基本部分ですが、これらのマップを使用して、面積単位を検索し、それが主に属する地域の権限を返すことができる参照テーブルを作成することです。over（）を使用して、どのポリゴンを見つけることができます。重複しますが、多くの場合、エリアユニットは、少なくとも一部は複数の領土当局内にあるようです。個々のケースを見ると、通常、エリアユニットの90％以上が単一の領土当局にあることがわかります。

over（）が行うことを実行するが、すべての重複するポリゴンだけでなく、いくつかの重複するポリゴンのどれがそれぞれの場合に最も重複しているのかを識別できる、すぐに使える手段はありますか？

4

2 に答える 2

7

これが@Silverfishの答えを利用して仕事をしたコードです

library(sp)
library(rgeos)
library(rgdal)

###
# Read in Area Unit (AU) boundaries
au <- readOGR("C:/Users/Peter Ellis/Documents/NZ", layer="AU12")

# Read in Territorial Authority (TA) boundaries
ta <- readOGR("C:/Users/Peter Ellis/Documents/NZ", layer="TA12")

###
# First cut - works ok when only one TA per area unit
x1 <- over(au, ta)
au_to_ta <- data.frame(au@data, TAid = x1)

###
# Second cut - find those with multiple intersections
# and replace TAid with that with the greatest area.

x2 <- over(au, ta, returnList=TRUE)

# This next loop takes around 10 minutes to run:
for (i in 1:nrow(au_to_ta)){
    tmp <- length(x2[[i]])
    if (tmp>1){
        areas <- numeric(tmp)
        for (j in 1:tmp){
            areas[j] <- gArea(gIntersection(au[i,], ta[x2[[i]][j],]))
            }
#       Next line adds some tiny random jittering because
#       there is a case (Kawerau) which is an exact tie
#       in intersection area with two TAs - Rotorua and Whakatane

        areas <- areas * rnorm(tmp,1,0.0001)

        au_to_ta[i, "TAid"] <- x2[[i]][which(areas==max(areas))]
    }

}


# Add names of TAs
au_to_ta$TA <- ta@data[au_to_ta$TAid, "NAME"]

####
# Draw map to check came out ok
png("check NZ maps for TAs.png", 900, 600, res=80)
par(mfrow=c(1,2), fg="grey")
plot(ta, col=ta@data$NAME)

title(main="Original TA boundaries")
par(fg=NA)
plot(au, col=au_to_ta$TAid)
title(main="TA boundaries from aggregated\nArea Unit boundaries")
dev.off()

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

于 2013-01-09T05:58:13.573 に答える
5

私はあなたが望むものが地理情報システムSEでカバーされていると思います：

https://gis.stackexchange.com/questions/40517/using-r-to-calculate-the-area-of-multiple-polygons-on-a-map-that-intersect-with?rq=1

特に、テリトリーポリゴンがT1、T2、T3などで、分類しようとしているポリゴンがAである場合は、おそらくAとT1、次にAとT2、次にAとT3などで使用gAreaし、どちらを選択しますか。gIntersection最大面積があります。rgeos（パッケージが必要になります。）

于 2013-01-08T04:32:23.550 に答える