Ariel,

It appears that you are extracting values layer by layer, in a loop. If you
do that, the weights are re-calculated each iteration. Better to use a
RasteStack or (better yet) a RasterBrick  and extract all values at once.
Also make sure you have the most recent version of the raster package.

If that does not do it, you could disaggregate your raster data (again, in
one step) such that you do not need to use the weights.

Robert


On Sat, Sep 17, 2011 at 6:10 PM, Ariel Ortiz-Bobea <aortizbobea@arec.umd.edu
> wrote:

> 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.
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

	[[alternative HTML version deleted]]

