[R-sig-Geo] raster: croping and aggregating

[Barry Rowlingson]

On Fri, Jul 15, 2011 at 11:10 AM, Agustin Lobo <alobolistas at gmail.com> wrote:
> I have 2 raster objects
>
>> show(CODEQP)
> class       : RasterLayer
> dimensions  : 10, 10, 100  (nrow, ncol, ncell)
> resolution  : 100, 100  (x, y)
> extent      : 438048.4, 439048.4, 4633415, 4634415  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=utm +zone=31 +ellps=intl +units=m +no_defs
> values      : in memory
> min value   : 1
> max value   : 1
>
>> dem = raster("mde10mny.tif")
>> show(dem)
> class       : RasterLayer
> dimensions  : 2333, 2788, 6504404  (nrow, ncol, ncell)
> resolution  : 10, 10  (x, y)
> extent      : 435926.9, 463806.9, 4613063, 4636393  (xmin, xmax, ymin, ymax)
> coord. ref. : NA
> values      : /media/Iomega_HDD/DIPU2011/RVISITANTS/mde10mny.tif
> min         : ?
> max         : ?
>
> and want to make a thrid object "dem100" perfectly aligned to CODEQP.
>
> I've tried with crop() but, as noted in the help page, the input
> extent is actually modified to get aligned to the input
> raster pixels and I actually would like the opposite:
>
>> dem100 = crop(dem,extent(CODEQP))
>> show(dem100)
> class       : RasterLayer
> dimensions  : 100, 100, 10000  (nrow, ncol, ncell)
> resolution  : 10, 10  (x, y)
> extent      : 438046.9, 439046.9, 4633413, 4634413  (xmin, xmax, ymin, ymax)
> coord. ref. : NA
> values      : in memory
> min value   : 549.4215
> max value   : 627.9377
>
>> show(CODEQP)
> class       : RasterLayer
> dimensions  : 10, 10, 100  (nrow, ncol, ncell)
> resolution  : 100, 100  (x, y)
> extent      : 438048.4, 439048.4, 4633415, 4634415  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=utm +zone=31 +ellps=intl +units=m +no_defs
> values      : in memory
> min value   : 1
> max value   : 1
>
> Note they are not identical.
>
> Is there any way (besides using grass from within R) to get dem100 to
> have the same geometry (extent and pixel size) than CODEQP?

 I think that doing a crop() can only return cells that are directly
taken from the input raster. Your CODEQP raster doesn't align to the
dem raster, so what you want is not a crop but will involve some kind
of interpolation, as your target extent cells overlap more than one
input raster cells.

 Possibly extract() will do it?

Barry