Hoping to get some help - likley more user error than a problem...
I am working with 7 large (5GB) raster stacks (1984-2021) that cover my area of interest. I am trying to crop/mask a long list of polygons to the raster that they fall within. The desired output would be a spatraster stack of just the polygon extent (or rectangle including the polygon) that is more manageble for analysis.
library(terra)
##simplified example
r_1 <- rast("file_name_raster_west.tif")
r_2 <-rast("file_name_raster_east.tif")
vect <- vect("file_name_polygon.tif")
#logic
crop(r_1, vect) # returns spatraster - vect is in west
crop(r_2, vect) # returns error - not in extent
#
clip_to_rast <- function (vect, r1, r2){
test <-crop(r1, aoi_v) #tried to use terra::extract to
test2 <-crop(r2, aoi_v)
if (nrow(test)>0)
return(test)
else (nrow(test2)>0) #tried = F,
return(test2)
ifelse("no match")
}
output <- clip_to_rast(vect, r_1, r_2)
##Error: [crop] extents do not overlap #dont understand why it doesnt return the first match in this case.
Will need to work over 100 polygons and check 7 rasters.
I have tried to (1) merge the rasters; but it crashes R (2) moasic; same issue (3) if/else statement; unsuccessful - I am not sure how to cycle across the 7 rasters and intigrate crop into a boolean; probaly simple I just cant get it to work... (4) try-catch (and a few other SO with DE-9IM); but I dont want to know if it throws an error; I want it to just move on to the next raster the excute the operation. (5) I tried to use terra::extract above (so that I could use !is.na = F); but it retured TRUE for each extracted pixel Would need to find a way to crop based on that result.
I cant imagine this is a unique operation; so throwing it out to hopfully learn a flow/function that helps. Thank you for any advise!
Update: Tried 2 other approaches
##Try Catch##
raster_files <- list(r_1, r_2)
clipped_polygons <- list()
# Loop over each raster file
for (i in seq_along(raster_files)) {
# Attempt to clip the polygon to the raster
tryCatch({
clipped_polygon <- crop(i, polygon, mask = T )
clipped_polygons[i] <- clipped_polygon
}, error = function(e) {
print("not in this raster")
})
}
###this returns 2x["not in this raster"] and a empty list (though I know it falls within r_1)
##if else approach
raster_files <- list(r_1, r_2)
clipped_polygons <- list()
for (t in raster_files){
clipped_polygon <- crop(t, polygon, mask = T )
if(nrow(clipped_polygon)>0){
clipped_polygons_3[i] <- clipped_polygon
} else {
print("nope")
}
}
## this works ONLY is r_1 is first in the list...otherwise, I get the Error: [crop] extents do not overlap and a empty list (it does not move to the second list)
Any thoughts?