何度も使用geosphere::areaPolygon()
して成功しましたが、今問題が発生しました。
別のポリゴン sp1 を含むポリゴン sp2 があります。したがって、sp2 は sp1 よりも大きな領域を持つ必要があります。で各面積を計算するとareaPolygon()
、逆の結果が得られます。
areasp1 = 10133977
areasp2 = 9858811
gSymdifference を使用して対称的な異なるポリゴン sp3 を見つけました。
areasp3 = 275165.4
sp1 : 赤の破線
sp2 : 黒の線
sp3 : 緑の点線
私が使用したコード:
library(sp)
library(geosphere)
library(rgeos)
library(raster)
s1 <- data.frame(lon= c(130.4327, 129.85127, 121.00775, 84.9176, 60.40123,
58.97929, 58.55841, 94.95237),lat= c(30.25074, 29.83075, 23.7992, 28.25964,
36.89905, 37.72305, 52.58793, 43.47459))
s2 <- data.frame(lon= c(130.4327, 129.85127, 121.00775, 54.35377, 58.55841,
94.95237),lat= c(30.25074, 29.83075, 23.7992, 31.61798, 52.58793, 43.47459))
sp1 <- SpatialPolygons(list(Polygons(list(Polygon(s1)),1)))
crs(sp1) <- CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0")
sp2 <- SpatialPolygons(list(Polygons(list(Polygon(s2)),1)))
crs(sp2) <-CRS ("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0")
areasp1 <- areaPolygon(sp1)/10^6 # to get area in square km
areasp2 <- areaPolygon(sp2)/10^6
sp3 <- gSymdifference(sp1,sp2,checkvalidity = TRUE)
areasp3 <- areaPolygon(sp3)/10^6
すべてのポリゴンは、 の自己交差またはその他の問題についてテストされているため、ここでgIsValid()
言及されている問題とは関係ありません。
sp1がsp2よりも面積が大きい理由について何か考えはありますか?
注:areasp2
+を合計areasp3
すると、ほぼ に等しくなりareasp1
ます。これが偶然かどうかはわかりません。