CODEWITHSUNDEEP
rprojectionr-raster
Matifou
Matifou

Reputation: 8880

r raster to polygon shows a regular grid even if cells had different area?

I am converting a raster to polygons. In my initial raster, cells have different area, but when I convert to polygon using rasterToPolygons and plot it, I see an equi-distant grid, even though initial cells had different area?

Is it an issue with my R code, or an issue with my bad understanding of projections? How could I do to represent the polygons with area proportional to the cells area?

library(raster)
r <- raster(nrow=18, ncol=36)
r$value <- 1:ncell(r)
r$area <- as.data.frame(area(r))$layer

head(as.data.frame(r))

ras_to_pol <- rasterToPolygons(r)
spplot(ras_to_pol, "area", main="Area of initial cells converted to polygons")

enter image description here

Thanks!

Upvotes: 1

Views: 506

Answers (1)

lbusett
lbusett

Reputation: 5932

In your code, you are not specifying the resolution of the cells nor the extent of the raster or the projection. Therefore raster defaults to:

extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0

and since you are using nrow=18, ncol=36 you will get cells of 10 degrees both on latitude and longitude, and therefore "square" cells when you plot them in "cartesian coordinates" on a plane.

However, in geographic coordinates the area of the surface of the "spheroid" corresponding to a 10x10 degrees cell actually changes with latitude (See here https://gis.stackexchange.com/a/29743 and here https://badc.nerc.ac.uk/help/coordinates/cell-surf-area.html for an explanation), as shown by the "colours" on your plot.

If you want to have polygons with variable "area" you have to reproject to an equal-area projection. There are many: here I'm using a Conical Equal Area as an example: 

ras_to_pol_new  <- spTransform(ras_to_pol, CRS("+init=epsg:3410"))
spplot(ras_to_pol_new, "area", main="Area of initial cells converted to polygons")

enter image description here

You can see how the area of the polygons decreases with latitude. (For a very good primer on projections, you can see for example: https://kartoweb.itc.nl/geometrics/Map%20projections/mappro.html)

Upvotes: 1

Related Questions