Opened 13 years ago
Closed 8 years ago
#3837 closed defect (fixed)
AREA_OR_POINT=Point is ignored in georeferencing
|Reported by:||rhijmans||Owned by:||warmerdam|
GDAL docs ( http://www.gdal.org/gdal_datamodel.html ) say:
The following metadata items have well defined semantics in the default domain: AREA_OR_POINT: May be either "Area" (the default) or "Point". Indicates whether a pixel value should be assumed to represent a sampling over the region of the pixel or a point sample at the center of the pixel. This is not intended to influence interpretation of georeferencing which remains area oriented.
AREA_OR_POINT=Point appears in geotiff files produced by the HEG-Modis Reprojection Tool. Presumably these files use the flag "PixelIsPoint" mentioned in the geoTiff docs ( http://www.remotesensing.org/geotiff/spec/geotiff2.5.html ).
These geotiff docs say that "PixelIsArea" and "PixelIsPoint" _do_ have implications for georeferencing (that's the point).
It appears that gdalinfo indeed ignores this flag. And it seems to me that this is wrong. Or does gdal do adjustments "under the hood"?
In contrast, it seems that ArcMap 9.3 does not ignore it (although it seems to make a mistake in how it corrects for it)
Change History (9)
comment:1 by , 13 years ago
|Status:||new → assigned|
comment:2 by , 13 years ago
comment:3 by , 13 years ago
|Status:||assigned → closed|
Fix addressed in:
Applied in trunk (r21158), 1.7 (r21159), 1.6 (r21160) and 1.6-esri (r21161).
comment:4 by , 10 years ago
Milestone 1.6.4 deleted
comment:5 by , 8 years ago
|Status:||closed → reopened|
I generated a time series of MODIS NDVI using HEG-Modis Reprojection Tool. The output is a (set of) geotiff with this metainformation (from gdalinfo)
AREA_OR_POINT=Point TIFFTAG_SOFTWARE=MODIS Reprojection Tool v4.1 March 2009
Image Structure Metadata:
Corner Coordinates: Upper Left ( -10.0000000, 65.0000000) ( 10d 0' 0.00"W, 65d 0' 0.00"N) [...]
This is formally correct, but when I visualize the image in QGIS (2.6) and when I import it into Postgis (2.1.3) the image is shifted half-pixel, which means that "AREA_OR_POINT=Point" is ignored. Instead, Arcgis (10.something) correctly visualize it. maybe I am missing something, but the problem seems to be the same reported in this old ticket so I re-open it.
follow-up: 7 comment:6 by , 8 years ago
The behaviour of georeferencing can be changed by setting the GTIFF_POINT_GEO_IGNORE configuration option to NO to disable the AREA_OR_POINT=Point correction. Attaching or providing a link to the above mentionned image would be helpful too. An upper let corner with exact coordinates is suspicious. It should be shifted by a half pixel to take account the GDAL convention, so I suspect that the file itself might be incorrect.
comment:7 by , 8 years ago
The entire image is quite big (it covers all europe, around 200GB)): https://www.dropbox.com/s/zw8ybiz2f8djq91/MOD13Q1_2014_12_03.250m_16_days_NDVI.zip?dl=0
here I have cut a small piece (few KB): https://www.dropbox.com/s/fk2k624qfwwjolw/MOD13Q1_2014_12_03.250m_16_days_NDVI_cut.zip?dl=0
I am probably misunderstanding the behavior of the default value of GTIFF_POINT_GEO_IGNORE: once "AREA_OR_POINT=Point", I expect GDAL to consider the UL coordinates as the centroid of the corner pixel and thus shift the image, but this doesn't happen.
Or at least this is what I understand from http://trac.osgeo.org/gdal/wiki/rfc33_gtiff_pixelispoint ("Interpretation of the raster space from the GeoTIFF tie points will be offset by half a pixel in the PixelIsPoint case").
comment:8 by , 8 years ago
The raw registration of the file as indicated by listgeo is
Geotiff_Information: Version: 1 Key_Revision: 1.0 Tagged_Information: ModelTiepointTag (2,3): 0 0 0 -9.99895835 64.99895835 0 ModelPixelScaleTag (1,3): 0.0020833 0.0020833 0 End_Of_Tags. Keyed_Information: GTModelTypeGeoKey (Short,1): ModelTypeGeographic GTRasterTypeGeoKey (Short,1): RasterPixelIsPoint GeographicTypeGeoKey (Short,1): GCS_WGS_84 GeogCitationGeoKey (Ascii,7): "WGS 84" GeogAngularUnitsGeoKey (Short,1): Angular_Degree GeogSemiMajorAxisGeoKey (Double,1): 6378137 GeogInvFlatteningGeoKey (Double,1): 298.257223563 End_Of_Keys. End_Of_Geotiff.
So in pixel-is-point model, the center of the top-left pixel is -9.99895835 64.99895835. But in GDAL data model, geotransform is alwas reported with the pixel is area convention, so in order for the center of the GDAL pixel to match -9.99895835 64.99895835, GDAL advertizes a top-left cordinate of the top-left pixel with a half-pixel shift towards the north-west, hence -10.0000000, 65.0000000. So I'd say GDAL behaviour is correct.
comment:9 by , 8 years ago
|Status:||reopened → closed|
Thanks a lot for the quick answer, and sorry for reporting this "false" alarm. Now I have a picture of the whole thing and I confirm that GDAL behaves correctly (while Arcgis makes an error in the visualization).
Just for documentation, here below all the information that I gathered:
In the log-file of the HEG-Modis Reprojection Tool, the coordinates are (-10,65), but: "The upper left (UL) corner specified for output in GeoTIFF refers to the center of the corner pixel" while "Any output corner coordinates specified by the GUI, in the status box, or by the command line, in standard output or to the log file, represent the outer extent of the pixel"
The geotiff is correctly registered with reference to the center of the corner pixel and tagged as "GTRasterTypeGeoKey (Short,1): RasterPixelIsPoint"
The standard output of GDALINFO reports that "AREA_OR_POINT=Point" and shows the coordinates related to the outer extent of the file (with the correction of half pixel).
To my immense consternation section 126.96.36.199 of the geotiff specification does make it clear that the georeferencing semantics are altered by half a pixel in the pixel-as-area case refuting 15 years of stubborn behavior on my part.
I shall attempt to prepare a GDAL RFC to address how this transition will be handled.
Details to follow.