GridTopology, Gluing SpatialGridDataFrames together from adjacent Rasters
Chris,
Something along these lines might work:
library(raster)
drg85<- raster("~path/drg85.tif")
drg97<- raster("~path/drg97.tif")
drg85_p<- crop(drg85, extent(xmin, xmax, ymin, ymax))
drg97_p<- crop(drg97, extent(xmin2, xmax2, ymin2, ymax2))
ext <- unionExtent(drg85_p, drg97_p)
r1 <- expand(drg85_p, ext)
x <- crop(r1, drg97_p)
x <- resample(drg97_p, x, method='bilinear')
r <- merge(r1, y, filename='test.tif')
Robert
On Thu, Jul 22, 2010 at 10:03 AM, Chris English <sglish at hotmail.com> wrote:
Hi List,
I am new to R spatial. ?I have group of islands in a lake that I kayak around and observe birds, animals, and vegetation. ?I wanted to create a map of my trips and points of observation.
My study area is essentially at the edges of two 7.5 quads, south edge of top and north edge of bottom. ?My goal has been to join parts of these two quad maps and derive a map showing my study area.
To Start:
drg85<- raster("~path/drg85.tif")
drg87<- raster("~path/drg97.tif")
drg85_p<- crop(drg85, extent(xmin, xmax, ymin, ymax))
if(require(rgdal)){drg85_pw<- writeRaster(drg85_p, filename="drg85_p",
format="GTiff", overwrite=TRUE)}
(do the same crop and write for drg97_p)
the cropping also allowed for plotting in a small RAM laptop
drg85_p<- readGDAL("~path/drg85_p.tif")#returns SpatialGridDataFrame
drg97_p<- readGDAL("~path/drg97_p.tif")#returns SpatialGridDataFrame
a more direct approach
drg85<- raster("~path/drg85.tif")
drg97<- raster("~path/drg97.tif")
drg85_p<- crop(drg85, extent(xmin, xmax, ymin, ymax) , filename =
"drg85_p.tif") # filename can be omitted
drg97_p<- crop(drg97, extent(xmin, xmax, ymin, ymax) , filename = "drg85_p.tif")
drg85_p<- as(drg85, 'SpatialGridDataFrame')
drg97_p<- as(drg97, 'SpatialGridDataFrame')
So I tried to cbind: ?drg8597_ps<- cbind(drg85_p, drg97_p) Error in stop.ifnot.equal(x, grds[[1]]) : topology is not equal
getGridTopology(drg85_p)
? ? ? ? ? ? ? ? ? ? ? ? ? ? x ? ? ? ? ? ? ? ? ? ?y cellcentre.offset 1.011904e+06 ? 4.556105e+05 cellsize ? ? ? ? ? ? ?8.252375e+00 ? 8.252375e+00 cells.dim ? ? ? ? ? ?2.001000e+03 ? 8.960000e+02
getGridTopology(drg97_p)
? ? ? ? ? ? ? ? ? ? ? ? ? ? x ? ? ? ? ? ? ? ? ?y cellcentre.offset 1.011901e+06 ? 4.54005e+05 cellsize ? ? ? ? ? ? ?8.251850e+00 ?8.25185e+00 cells.dim ? ? ? ? ? ?2.013000e+03 ? 2.00000e+02 And indeed they're different, most especially between y cells.dim. Checked to see if the grids held data:
fullgrid(drg85_p)
[1] TRUE
fullgrid(drg97_p)
[1] TRUE So I thought, well most of the area is in drg85_p, why not make the topology of drg97_p the same as drg85_p and tried:
drg97_p<- GridTopology(cellcentre.offset(x=1.011904e+06,y=4.556105e+05),cellsize(x=8.252375e+00, y=8.253275e+00), cells.dim(x=2.001000e+03 y=8.960000e+02))
Error in initialize(value, ...) : ?could not find function "cellcentre.offset" OK- get the form right:
drg97_p<- GridTopology(c(1.011904e+06, 4.556105e+05),c(8.252375e+00, 9.252375e+00),c(2.0010000e+03,8.960000e+02)) fullgrid(drg97_p)
[1] FALSE Oops, the data. ?Well, I've made a grid. Go readGDAL to drg97_pd again. And now I can't figure out how to read out the data and write to my equal topology. Any help would be appreciated. Chris
_________________________________________________________________ The New Busy is not the old busy. Search, chat and e-mail from your inbox. N:WL:en-US:WM_HMP:042010_3 ? ? ? ?[[alternative HTML version deleted]] _______________________________________________ R-sig-Geo mailing list R-sig-Geo at stat.math.ethz.ch https://stat.ethz.ch/mailman/listinfo/r-sig-geo