3

从这个问题中，我获得了一个很棒的功能areaPolygon()，它给了我一个坐标多边形内的区域。但是，当我尝试使用该函数时，计算似乎很奇怪：

我首先创建一组两个点

 require(fields)
 coords <- c(11.3697193956209, 47.233380520521, 11.3723606043791, 
 47.235179479479)
 coords <- matrix(coords, nrow=2, ncol=2, byrow=TRUE)

然后我检查这两者之间的距离：

 rdist.earth(coords,coords,miles=FALSE)[1,2]

获得：0.2827821 公里（这将是矩形的对角线）

我继续创建一个矩形

 polygon <- matrix(coords, nrow=2, ncol=2)
 polygon <- rbind(polygon, polygon)
 polygon[4,2] <- polygon[1,2]
 polygon[4,1] <- polygon[2,1]
 polygon[3,2] <- polygon[2,2]
 polygon[3,1] <- polygon[1,1]
 polygon <- rbind(polygon, polygon[1,])

看看这看起来是否不错：plot(polygon)

第四步：我计算多边形内的面积。

geosphere::areaPolygon(polygon)
[1] 31.99288 #from the help file I know this ought to be square metres.

但是，我会预料200*200=40000 m²到，因为我的斜角的边是 200 x 200 米。这可以通过

rdist.earth(polygon,coords,miles=FALSE)

           [,1]         [,2]
 [1,] 0.0000000 2.827821e-01
 [2,] 0.2827821 9.504539e-05
 [3,] 0.2002671 1.996434e-01
 [4,] 0.1996501 2.002671e-01

所以现在来回答我的问题（最后）我做错了什么？非常感谢您的帮助！

4

1 回答 1

5

你做了一个无效的多边形！如果你用它来绘制它，type="l"你会看到一个蝴蝶结：

> plot(polygon,type="l")

一半的领结将有一个负面积，另一半是正面积，所以你的结果是两半尺寸的差异。它们不会完全一样，因为地球是球形的...

您只需要在 中重新排序您的积分polygon。

于 2012-09-28T11:11:33.663 回答