I have a legacy program that attempts to read USGS 1/3 arcsec data. When I compare results to GDAL, I find that the lat/lon to pixel conversion is off by 1 line. The reason (I think) is because the legacy code calculates the latitude of line 0 as:

lat0 = yllcorner + (nrows - 1) * cellsize

where yllcorner, nrows and cellsize come from the .HDR file.

I see in the GDAL code in frmts/raw/ehdrdataset.cpp that there is no "subtract 1" term:

 dfULYMap = dfYLLCorner + nRows * dfYDim;

I'm assuming that GDAL is right. Can anyone explain why so that a noob can understand?

  • 1
    It depends on how you look at it: the GDAL code gives you the Y-value for the top edge of the top row (line 0), the legacy code gives you the value for the bottom edge of line 0 (which also happens to be the top edge of line 1). – Hermann Aug 24 '13 at 7:51
  • Hermann, that makes total sense and solves my problem! Thanks. – Reilly Aug 24 '13 at 20:49

There is actually an alternative explanation which has to do with the concepts that raster pixels could either represent points or areas. So it again depends on how you look at it. Below the example of a 3 x 3 pixel raster.

| * | * | * |
| * | * | * |
| * | * | * |

Pixel is Point: you have to go n-1 rows from the lower left point (asterisk) to the top left one

Pixel is Area: you have to go n rows from the lower left corner of the lower left area (box) to the upper left corner of the upper left area.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.