Clip a raster

On Tue, 3 Feb 2009, Etienne Bellemare Racine wrote:

Is it possible to clip a raster after having it read with readGDAL ?
Instead of using region.dim() and offset(), I would like to clip it "on
the fly" to visualize statistics associated with objects of different
size. The clip is square -specified with upper left and lower right
corners- and fit with the pixels size. If you have an interactive way,
I'm even more interested.

For an interactive method for a SpatialGrid* object, try:

library(sp)
data(meuse.grid)
coordinates(meuse.grid) <- c("x", "y")
gridded(meuse.grid) <- TRUE
fullgrid(meuse.grid) <- TRUE
bbox(meuse.grid)
image(meuse.grid, "dist", breaks=seq(0,1,0.1), col=heat.colors(10))
# interactive choice of rectangle of interest, ll and ur
pts <- locator(2)
pts
points(pts)
grd <- getGridTopology(meuse.grid)
gi <- getGridIndex(coordinates(pts), grd)
cd <- slot(grd, "cells.dim")
ys <- sort(ceiling(gi / cd[1]))
xs <- sort(gi %% cd[1])
meuse_sub <- meuse.grid[ys[1]:ys[2], xs[1]:xs[2]]
bbox(meuse_sub)
image(meuse_sub, "dist", breaks=seq(0,1,0.1), col=grey.colors(10),
 add=TRUE)

What you call clipping is the "[" operator here.

Hope this helps,

Roger

Etienne


       [[alternative HTML version deleted]]

--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no