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]]
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
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.biv...@nhh.no
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
