[R-sig-Geo] computing square kilometers in shape file polygons

[Dan Putler]

Hi Joe,

I assume you want to append the area to the attribute data frame. If yes, the
function pasted below should help you out, assuming your study area
isn't too large. It is based on re-projecting the data into an area
preserving projection, and then calculating the areas in that
projection. I've mostly worked with it in British Columbia, so use a BC
Albers equal area project (EPSG:3005). The map units of your projection
will likely be in meters, so the calculated areas will be in m^2, which
will be easy to convert.

Dan

The code:

# The function areaPolygons() returns the area of the polygons in a
# SpatialPolygons or SpatialPolygonsDataFrame object. If a proj4string is
# provided then the areas are calculated in the map units of this
projection.
# If it is not provided, then the areas are calculated in the map units
of the
# original SpatialPolygon or SpatialPolygonsDataFrame object. If the
projection
# of the sp object is geographic (i.e., lon/lat) the areas are in
decimal degrees^2.

# Requires the sp and rgdal packages

# Programmer: Dan Putler
# Created: 06Nov09
# Modified: 06Nov09

areaPolygons<- function(spPoly, proj4string = NULL) {
      if(class(spPoly)[[1]] != "SpatialPolygonsDataFrame"&
class(spPoly)[[1]] != "SpatialPolygons") {
          stop("spPoly must be a SpatialPolygonsDataFrame or a
SpatialPolygons object.")
      }
      require(sp)
      require(rgdal)
      if(!is.null(proj4string)) {
          if(class(proj4string)[[1]] != "CRS") {
              stop("The proj4string must be of class CRS")
          }
          spP<- spTransform(spPoly, CRS = proj4string)
      }
      else {
          spP<- spPoly
      }
      areas<- unlist(lapply(spP at polygons, function(x) a<- x at area))
      return(areas)
}


On 09/03/2010 03:30 PM, brwin338 at aol.com wrote:
>  Good Afternoon
>  Does any one know of a package/function that could be used to compute the square kilometers for polygons in a long-lat projected shapefile?
>  I can access the values in the "area" slots in the object resulting from using maptool's  "readShapePoly" command.  I have approximated the km^2 from these values by using one of the polygons's (such as a given state's) known area in km^2 and rescaling the polygon's "area" relative to the state's known area.  I then applied this "scaling factor" to the other polygons in my shapefile object to approximate the other polygons' are in km^2.
>
>  However, if I understand things correctly, my approximate "area to km^2" scaling  will become increasing distorted as I move north or south relative to my original "normalizing" state.  Is there a function that can compensate for changes in latitude when computing the polygon areas associated with long-lat polygons?  I hope this question is somewhat coherent
>  Thanks
>  Joe
>
>
>
>
>  	[[alternative HTML version deleted]]
>
>  _______________________________________________
>  R-sig-Geo mailing list
>  R-sig-Geo at stat.math.ethz.ch
>  https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>