Opened 8 years ago

Closed 3 years ago

#3837 closed defect (fixed)

AREA_OR_POINT=Point is ignored in georeferencing

Reported by: rhijmans Owned by: warmerdam
Priority: high Milestone:
Component: GDAL_Raster Version: unspecified
Severity: major Keywords: AREA_OR_POINT gtiff
Cc:

Description

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 Changed 8 years ago by warmerdam

Status: newassigned

To my immense consternation section 2.5.2.2 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.

comment:2 Changed 8 years ago by warmerdam

Keywords: gtiff added

comment:3 Changed 8 years ago by warmerdam

Milestone: 1.6.4
Resolution: fixed
Status: assignedclosed

Fix addressed in:

http://trac.osgeo.org/gdal/wiki/rfc33_gtiff_pixelispoint

Applied in trunk (r21158), 1.7 (r21159), 1.6 (r21160) and 1.6-esri (r21161).

comment:4 Changed 5 years ago by Even Rouault

Milestone: 1.6.4

Milestone 1.6.4 deleted

comment:5 Changed 3 years ago by feurbano

Resolution: fixed
Status: closedreopened

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)

[...] Metadata:

AREA_OR_POINT=Point TIFFTAG_SOFTWARE=MODIS Reprojection Tool v4.1 March 2009

Image Structure Metadata:

INTERLEAVE=BAND

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.

comment:6 Changed 3 years ago by Even Rouault

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 in reply to:  6 Changed 3 years ago by feurbano

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 Changed 3 years ago by Even Rouault

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 Changed 3 years ago by feurbano

Resolution: fixed
Status: reopenedclosed

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).

Note: See TracTickets for help on using tickets.