20

I am working with a shapefile that has many polygons.

How do I add one more field named "area_sqkm" and calculate area for each polygon in the shapefile?

PolyGeo
  • 65,136
  • 29
  • 109
  • 338
HLe
  • 321
  • 1
  • 2
  • 4
  • 1
    What libraries do you have installed and what have you tried? – GISHuman Jun 29 '16 at 15:57
  • 5
    rgeos::gArea(x,byid=TRUE) – mdsumner Jun 29 '16 at 17:21
  • 1
    Welcome to gis.stackexchange! Please note that a good question on this site is expected to show some degree of research on your part, i.e. what you have tried and - if applicable - code so far. For more info, you can check our [faq]. – underdark Jun 29 '16 at 19:53

1 Answers1

30

You can do

library(raster)
x <- shapefile('file.shp')
crs(x)
x$area_sqkm <- area(x) / 1000000

Or with terra

library(terra)
x <- vect('file.shp')
x$area_sqkm <- expanse(x) / 1000000

That works for both an angular CRS (longitude/latitude) and for planar CRSs. It is generally best to use a longitude/latitude CRS to compute area a planar CRS can distort size with a lot (How much? It depends on the CRS and on the location and extent of your spatial data).

Robert Hijmans
  • 10,683
  • 25
  • 35
  • 1
    is raster::area() with unprojected data more accurate than rgeos::gArea() with projected data? – Richard DiSalvo Apr 30 '19 at 22:45
  • 2
    It can be. Because projections distort. If the inter-node distance of the polygon is high, projection can distort the area further (see the geosphere vignette). – Robert Hijmans May 01 '19 at 00:14
  • 2
    For other readers: This returns a numeric vector equal to the number of polygons (it looked to me as if it might be just for one polygon at a time) – geneorama Aug 27 '21 at 21:23