У меня есть большие ежедневные наборы данных в растровом формате. Я хочу рассчитать общее количество пикселей на основе значений в разных многоугольниках в одном шейп-файле. Шейп-файл представляет собой 90-метровую модель STRM DEM, классифицированную по 24 высотным зонам. Эти 24 зоны высот представляют 24 полигона в одном шейп-файле. Я хочу проверить, сколько пикселей занимает каждый многоугольник.
Прежде всего, мне нужно проверить общее количество пикселей на основе следующих значений (200210240250) в каждом многоугольнике и, наконец, сохранить его в CSV.
Я уже разработал код: но столкнулся с проблемой в конце.
library(sp)
library(rgdal)
library(raster)
mod = raster("MOYDGL06_Maximum_Snow_Extent_2004097.tif")
shp= readOGR("Gilgit_DEM_24.shp")
mod_ext = extract(mod,shp,df=T,na.rm=T)
mod_mask = mask(mod,shp)
plot(r2,axes = TRUE,ext = extent(shp))
r3_200 = rasterToPoints(mod_mask,function(x){ x ==200 },spatial = TRUE)
r3_210 = rasterToPoints(mod_mask,function(x){ x ==210 },spatial = TRUE)
r3_240 = rasterToPoints(mod_mask,function(x){ x ==240 },spatial = TRUE)
r3_250 = rasterToPoints(mod_mask,function(x){ x ==250 },spatial = TRUE)
r3_200_1 = raster::intersect(shp,r3_200)
write.csv(r3_200_1,file = 'r2_extract_gilgit.csv')
Изображение и код R доступны в по этой ссылке