Search code examples
pythonmatplotlibscipygdalgeotiff

calculate latitude and longitude from a geotiff image

How can I calculate mean latitude and mean longitude of each block from a geotiff image separating it into regular blocks (say, 50 by 50 pixels).

The input data is just the imagery,for example as downloaded from: http://eoimages.gsfc.nasa.gov/images/imagerecords/57000/57752/land_shallow_topo_2048.tif

This can be opened into python using gdal as folows:

import gdal
geotiff = gdal.Open ('land_shallow_topo_2048.tif')
colum_numbers,row_numbers,band_numbers=geotiff.RasterXSize,
                                       geotiff.RasterYSize,geotiff.RasterCount
print (colum_numbers,row_numbers,band_numbers)
2048 1024 3
Solution

  • Take a look at this question: Obtain Latitude and Longitude from a GeoTIFF File

    Since your image is already in lat-lon coordinates, the following values should be in lat-lon:

    gt = geotiff.GetGeoTransform()
minx = gt[0]
miny = gt[3] + width*gt[4] + height*gt[5] 
maxx = gt[0] + width*gt[1] + height*gt[2]
maxy = gt[3]

    Now the size of a pixel is:

    latPxSz = (maxy - miny) / row_numbers
lonPxSz = (maxx - minx) / column_numbers

    The center of a 50x50 box that is i box rows and j box columns away from the corner is:

    boxCenterLat = (i + 0.5) * 50 * latPxSz + miny
boxCenterLon = (j + 0.5) * 50 * lonPxSz + minx

    If you are using some other coordinate system, you can achieve a similar result by doing the transformation from the other, similar, question.