If we have a raster, say of integer elevation data for a country, and one polygon shapefile, say of 300 river basins in that country, with a unique name for each, how would we most easily get output like this for them all?
basinID, gridcellelev a, 320 a, 321 a, 320 b, 17 b, 18 b, 19
The most burdensome way seems to be filtering/converting the single shapefile into 300 shapefiles, clipping the raster 300 times into 300 uniqueID rasters, reading them back in, generating individual tables for each basin, then combining them all together.
On the other hand the ideal way seems to be skipping the file generation, not saving the xy data, and creating the same table from just using one raster and one shapefile - by iteratively selecting the cells within a basin, stamping them with the basinID, creating a table, losing the coordinates, and keep repeating and appending that table until the 300th basin.
I'm not looking for any statistics, just the raw data listing of the grid cell elevations that would have been part of some standard clip method. I believe the attribute table that comes with the clipped raster from ArcMap is the counts/frequency of the cells. That output table format works for me too.
I don't know to minimally reproduce a raster and a polygon shapefile, so I'd just be grateful for any tips/libraries/functions/examples. If starting with R, here's a starting point:
library(tidyverse) library(raster) library(rgdal) library(sf) elev_raster <- raster("spain_elev_meters.tif") #integer raster basins <- readOGR("spainbasins.shp", "spainbasins") %>% st_as_sf() #unique basin ID column: `basinID`
I'd prefer to do everything in R, but happy to try any batch methods in ArcGIS as well (I have 10.6 with Spatial Analyst, but not Pro).