doubt with readGDAL
On Sun, 4 Oct 2009, "Jos? M. Blanco Moreno" wrote:
Dear R-users community, I have been trying to figure out on my own, but no way... I am trying to read and decimate a geotiff: library(maptools) library(spatial) library(sp) library(rgdal) download.file(url='http://srtm.geog.kcl.ac.uk/portal/srtm41/srtm_data_arcascii/srtm_37_04.zip',destfile='srtm_37_04.zip') shell('unzip -o srtm_37_04.zip') cata <- readGDAL('srtm_37_04.tif',offset=c(2900,0),region.dim=c(2600,5000),output.dim=c(260,500)) cata2 <- cata cata2 at data$band1 <- ifelse(is.na(cata2 at data$band1),-1,cata2 at data$band1) cata2 at data$band1 <- c(50,100,250,500,1000,1500,2000,3000)[cut(cata at data$band1,c(0,50,100,250,500,1000,1500,2000,3000))] image(cata2,col=gray(seq(.9,0.6,length=9))) contour(cata2,add=T,levels=c(0,50,100,250,500,1000,1500,2000,3000),col=gray(.7),drawlabels=F) axis(1) axis(2) But the coordinates don't come out: they are translated in the y-axis *a lot*; bbox(cata) returns: min max x -4.166183e-04 5.000417 y -1.577894e+01 -10.778109 While in the non-decimated version they are: min max x -0.0004166183 4.16625 y 40.4170835754 42.58375
Thanks for a clear report - the y range and the cell sizes were wrong when all three of offset=, region.dim= and output.dim= were used together. With current rgdal 0.6-18, using the Arc ASCII file (not GTiff), I see for the undecimated subscene:
gridparameters(catb)
cellcentre.offset cellsize cells.dim x 4.839851e-08 0.0008333333 5000 y 4.041750e+01 0.0008333333 2600
bbox(catb)
min max x -0.0004166183 4.16625 y 40.4170835754 42.58375 but for the decimated subscene:
gridparameters(cata)
cellcentre.offset cellsize cells.dim x 0.004584215 0.01000167 500 y -15.769325078 0.01923397 260
bbox(cata)
min max x -4.166183e-04 5.000417 y -1.577894e+01 -10.778109 In forthcoming 0.6-19 (committed to the sourceforge rgdal repository) and the undecimated subscene:
gridparameters(catb)
cellcentre.offset cellsize cells.dim x 4.839851e-08 0.0008333333 5000 y 4.041750e+01 0.0008333333 2600
bbox(catb)
min max x -0.0004166183 4.16625 y 40.4170835754 42.58375 with the decimated subscene now corresponding:
gridparameters(cata)
cellcentre.offset cellsize cells.dim x 0.003750048 0.008333333 500 y 40.421250242 0.008333333 260
bbox(cata)
min max x -0.0004166183 4.16625 y 40.4170835754 42.58375 This is with this subscene on this data set - I would be very grateful if anyone could check from a source install after checking out from sourceforge for other, including Southern hemisphere decimated subscenes (there is a sign wrinkle). I'll wait for feedback before submitting to CRAN. Again, thanks for a helpful report! Roger
Could anyone indicate a workaround to preserve the coordinates of the grid? Thank you for any information on this issue (which I suppose is generated by my own incompetence). Jos? M. sessionInfo() R version 2.8.1 (2008-12-22) i386-pc-mingw32 locale: LC_COLLATE=Spanish_Spain.1252;LC_CTYPE=Spanish_Spain.1252;LC_MONETARY=Spanish_Spain.1252;LC_NUMERIC=C;LC_TIME=Spanish_Spain.1252 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] rgdal_0.6-8 spatial_7.2-47 maptools_0.7-25 sp_0.9-37 foreign_0.8-36 loaded via a namespace (and not attached): [1] grid_2.8.1 lattice_0.17-25
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