Hello everyone, I've been exploring ways to "clip a grid/raster with a polygon/shape file" in R but what I have seen so far seems more suited for relatively fine grids and large polygons (e.g. many grids fall within each polygon). I'm confronted with a relatively coarse data grid for the US (0.3 deg ~ 32 km) and would like to clip these grids using county boundaries and calculate mean values for each polygon/county. Basically, many counties will fall "in between" several grids and so, ideally, I need to calculate a weighted mean for each polygon by weighting each grid value by the fraction of the polygon it is covering. I found a way to do this with the "extract" function in the raster package using the "weight" option (see example below and ?extract), but it takes about 60 minutes for one grid/raster layer for all US counties (~ 1 second/county/layer). My problem is that I need to do this for about 30 years * 365 days * 8 daily snapshots = 87,600 layers / variable. Enough to keep me 8+ years waiting for each variable, and I have about 25 variables to process... Do you know a (hopefully way) more efficient way to do this? Thank you in advance for any guidance! Ariel PS: I'm running on a 2.66 GHz core 2 duo iMac / 8 GB RAM / Mac OS X 10.7 (Lion) / R 2.12.2 FILES: ------- shape file: http://terpconnect.umd.edu/~ortiz/R/US48_rr-fixed.grb (~500KB) GRIB grid file: http://terpconnect.umd.edu/~ortiz/R/US48_rr-fixed.grb (~12MB) CODE: ------- #used libraries library(rgdal) library(raster) #setting up data grid <- readGDAL("US48_rr-fixed.grb") proj <- CRS(proj4string(fixed)) countyshape0 <- readShapePoly("co99_d00.shp", proj4string = CRS("+proj=longlat")) countyshape <- spTransform(countyshape0, proj) removestate <- countyshape$STATE %in% c("02","15","72") # leave out Alaska, Hawai and Puerto Rico countyshape48 <- countyshape[!removestate,] rasteredgrid <-raster(grid, 3) # chooses one layer #visual glimpse: image(rasteredgrid) plot(countyshape48, add=TRUE) #cookie cutting/raster clipping timer <- proc.time() # start the clock test <- extract(flxraster, countyshape48, weights=TRUE, fun=mean) proc.time() - timer # stop the clock ----- Ariel Ortiz-Bobea PhD Candidate in Agricultural & Resource Economics University of Maryland - College Park -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/How-to-efficiently-clip-large-grids-raster-with-polygons-tp6804923p6804923.html Sent from the R-sig-geo mailing list archive at Nabble.com.
How to efficiently clip large grids/raster with polygons?
4 messages · Robert J. Hijmans, Ariel Ortiz-Bobea
An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20110917/845bb65b/attachment.pl>
Thanks Robert, It works great. It took the same time to clip the RasterStack than the single RasterLayer. I read your "Introduction to the ?raster? package (version 1.9-11)" and I found very useful too to understand the different objects. Thanks for developing this resource! I wonder about a way to combine several variables (23) and time periods (~90,000) to efficiently do the polygon clipping on a large RasterBrick/Stacks. I feel that I need to choose between: 1- performing the clipping on 23 RasterStacks (one for each variable, e.g. temperature) in which each band corresponds to a time period (90,000 of them) 2- performing the clipping on 90,000 RasterStacks (one for each time period, e.g. 1980/1/1 09:00) in which each band corresponds to a variable (23 of them). ... and option 1 would be expectedly faster. But, is there a way to combine both options to do the clipping "ALL at once" of all the time periods and variables? I read about the "setZ" function for time-tagging layers and the ability to perform functions for each raster cell over time... but I'm not sure I can make use of it in this context (stacking 2 time-period layers doubles my band number... basically, I feel it does not lead to an additional dimension in the structure of the object). Any thoughts would be greatly appreciated! Regards, Ariel ----- Ariel Ortiz-Bobea PhD Candidate in Agricultural & Resource Economics University of Maryland - College Park -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/How-to-efficiently-clip-large-grids-raster-with-polygons-tp6804923p6806329.html Sent from the R-sig-geo mailing list archive at Nabble.com.
An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20110918/8a0d90e5/attachment.pl>