areaPolygon: Area of a longitude/latitude polygon

Description

Compute the area of a polygon in angular coordinates (longitude/latitude) on an ellipsoid.

Usage

# S4 method for matrix
areaPolygon(x, a=6378137, f=1/298.257223563, ...)
# S4 method for SpatialPolygons
areaPolygon(x, a=6378137, f=1/298.257223563, ...)

Arguments

longitude/latitude of the points forming a polygon; Must be a matrix or data.frame of 2 columns (first one is longitude, second is latitude) or a SpatialPolygons* object

major (equatorial) radius of the ellipsoid

ellipsoid flattening. The default value is for WGS84

...

Additional arguments. None implemented

Value

area in square meters

References

C.F.F. Karney, 2013. Algorithms for geodesics, J. Geodesy 87: 43-55. https://dx.doi.org/10.1007/s00190-012-0578-z. Addenda: http://geographiclib.sf.net/geod-addenda.html. Also see http://geographiclib.sourceforge.net/

Examples

Run this code

# NOT RUN {
p <- rbind(c(-180,-20), c(-140,55), c(10, 0), c(-140,-60), c(-180,-20))
areaPolygon(p)

# Be careful with very large polygons, as they may not be what they seem!
# For example, if you wanted a polygon to compute the area equal to about 1/4 of the ellipsoid
# this won't work:
b <- matrix(c(-180, 0, 90, 90, 0, 0, -180, 0), ncol=2, byrow=TRUE)
areaPolygon(b)
# Becausee the shortest path between (-180,0) and (0,0) is 
# over one of the poles, not along the equator!
# Inserting a point along the equator fixes that
b <- matrix(c(-180, 0, 0, 0, -90,0, -180, 0), ncol=2, byrow=TRUE)
areaPolygon(b)

# }