1

I have a raster with integer values graded from 0 to 32. I also have a shapefile of this area subdivided into three basins A, B and C.

I would like to calculate the area (or the number of pixels, then transform it to area) of each raster class based on the shapefile. For example, for Basin "A", what is the area (or number of pixels) with a value of 0? and with a value of 1? and with value 2......?

I don't have a script for that, but I have here what would be the initial to load the files in R.

library(sf)
library(stars)
library(maptools)

rastermap<-read_stars('C:/0mestrado/SRTM_SAGA/duvida/raster_2005_UTM_cut.tif')
shapefilemap<-readShapePoly('C:/0mestrado/SRTM_SAGA/duvida/01vetor_UTM.shp')

data-> https://drive.google.com/file/d/1ylE49TWp22c71Io09MUL33QXRc5ccNcn/view?usp=sharing

2

I suggest using package terra for this purpose.

library(terra)
# open files in terra native format
rastermap<-rast('C:/0mestrado/SRTM_SAGA/duvida/raster_2005_UTM_cut.tif')
shapefilemap<-vect('C:/0mestrado/SRTM_SAGA/duvida/01vetor_UTM.shp')

# perform extraction
ext<-extract(rastermap, shapefilemap)

Now, you have a dataframe with two columns, the first (ID) is the name of your polygons, the second are the pixel values. You can know the number of pixels of each value by using table.

> table(ext$ID, ext$classification_2011)
   
          3       4       9      12      15      20      24      25      33      39      41
  1  856040  219788    5041   22516 2294996      93    1021     418     148  241474   31824
  2  874728  122898    3053   26962 2637298      86    3865    2320    1698   46133   23336
  3 2192903  422501   58432   75469 5433626   44820    1884    1902    7542  122088   30779
3
  • hello friend, I didn't understand two things about your answer: 1) Does terra calculate the number of pixels inside each polygon? 2) why the term ext$classification_2011? – wesleysc352 19 mins ago Jun 4 at 19:34
  • It extracts one value per pixel within each polygon. In other words, you will have as much values as pixels; therefore, you have the count of pixels. With regard to your second question, I have no idea. I used the raster file you provided and I got that name. I presume it was the original name of the raster, but I do not know.
    – HMSP
    Jun 7 at 7:33
  • Hello I found out, the origin of the term classification_2011, when I apply View(ext), it returns a table with a column with that name. Jun 10 at 15:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.