Reversed raster after readGDAL()
On Wed, 10 Feb 2010, Agustin Lobo wrote:
I think I tripped on this same stone, but cannot find anything in my records, perhaps it's on the lost computer.
Sorry, can't find such a posting - maybe my search keys were wrong - you could try the nabble archive.
The file test2.tif was made by plugin Interpolation in QGIS out of a vector layer of points with nitrate concentration values, just a simple inverse distance method. I think that the author f that pluginis Marco Hugentobler (Marco, sorry if it's not you) to whom I'm forwarding. Relevant info from the metadata in QGIS for this raster: Layer Spatial Reference System: +proj=utm +zone=31 +ellps=intl +units=m +no_defs Origin 424389,4.68549e+06 Pixel Size: 279.887,-279.887
No, it is 345.0078, 194.0195.
NW corner in QGIS (interactive):424388.820,4685490.931 SW corner in QGIS (interactive):424389.122,4635951.740
SW in R is 424389, 4635822
or SW cell centre 424561.5, 4635919.0 which is the SW corner plus 0.5* the
resolution.
Try:
GDALinfo(system.file("pictures/cea.tif", package = "rgdal")[1])
system(paste("gdalinfo", system.file("pictures/cea.tif", package =
"rgdal")[1]))
to see the oddness of gdalinfo on your raster. The SW (LL) corner is north
of the NW (UL) one. The provided raster is OK, so I guess one might think
that yours isn't?
Image with grid: http://sites.google.com/site/eospansite/dummy/test2gridqgis.jpg The weird thing is that QGIS uses gdal to read the geotif files.
How is the plugin writing the GTiff? Roger
Agus 2010/2/10 Roger Bivand <Roger.Bivand at nhh.no>
On Wed, 10 Feb 2010, Agustin Lobo wrote:
Hi! After
r3 <- readGDAL("test2.tif")
test2.tif has GDAL driver GTiff and has 256 rows and 256 columns
image(r3,col=rainbow(128))
I get the image N-S reversed. I've put an screenshot and the actual geotif image here: http://sites.google.com/site/eospansite/dummy/test2_screenshot2.jpeg http://sites.google.com/site/eospansite/dummy/test2.tif
Collected from: http://sites.google.com/site/eospansite/dummy Yes, it is flipped on the Y axis. What software made it? Could you please show the images with axes, as I think that the GTiff has a wrong sign on its y-step? What are the coordinates of the NW and SW corners of the image in QGIS? GDAL expects images to be read by row N to S, so a y-step with the wrong sign (+ instead of -) flips the image northward. Roger PS: gdalinfo externally has lower left north of upper left:
system("gdalinfo test2.tif")
Driver: GTiff/GeoTIFF Files: test2.tif Size is 256, 256 Coordinate System is: PROJCS["ED50 / UTM zone 31N", ? ?GEOGCS["ED50", ? ? ? ?DATUM["European_Datum_1950", ? ? ? ? ? ?SPHEROID["International 1924",6378388,297.0000000000014, ? ? ? ? ? ? ? ?AUTHORITY["EPSG","7022"]], ? ? ? ? ? ?AUTHORITY["EPSG","6230"]], ? ? ? ?PRIMEM["Greenwich",0], ? ? ? ?UNIT["degree",0.0174532925199433], ? ? ? ?AUTHORITY["EPSG","4230"]], ? ?PROJECTION["Transverse_Mercator"], ? ?PARAMETER["latitude_of_origin",0], ? ?PARAMETER["central_meridian",3], ? ?PARAMETER["scale_factor",0.9996], ? ?PARAMETER["false_easting",500000], ? ?PARAMETER["false_northing",0], ? ?UNIT["metre",1, ? ? ? ?AUTHORITY["EPSG","9001"]], ? ?AUTHORITY["EPSG","23031"]] Origin = (424389.000000000000000,4635822.000000000000000) Pixel Size = (345.007812500000000,194.019531250000000) Metadata: ?AREA_OR_POINT=Area Image Structure Metadata: ?INTERLEAVE=BAND Corner Coordinates: Upper Left ?( ?424389.000, 4635822.000) ( ?2d 5'20.10"E, 41d52'11.88"N) Lower Left ?( ?424389.000, 4685491.000) ( ?2d 4'56.97"E, 42d19'2.06"N) Upper Right ( ?512711.000, 4635822.000) ( ?3d 9'11.42"E, 41d52'24.52"N) Lower Right ( ?512711.000, 4685491.000) ( ?3d 9'15.31"E, 42d19'14.90"N) Center ? ? ?( ?468550.000, 4660656.500) ( ?2d37'10.89"E, 42d 5'47.83"N) Band 1 Block=256x4 Type=Float64, ColorInterp=Gray and:
GDALinfo("test2.tif")
rows ? ? ? ?256 columns ? ? 256 bands ? ? ? 1 origin.x ? ? ? ?424389 origin.y ? ? ? ?4636016 res.x ? ? ? 345.0078 res.y ? ? ? 194.0195 oblique.x ? 0 oblique.y ? 0 driver ? ? ?GTiff projection ?+proj=utm +zone=31 +ellps=intl +units=m +no_defs file ? ? ? ?test2.tif apparent band summary: ? GDType ? ? ? ?Bmin ? ? ? Bmax 1 Float64 -4294967295 4294967295 Metadata: AREA_OR_POINT=Area
Agus ? ? ? ?[[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
-- 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
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