Latitude/Longitude of each pixel using python and gdal

I needed lat/long of each pixel of a GeoTiff file. I searched through internet to find a solution and ,thanks to Stack Overflow, found a piece of code that I modified. As usual, python and gdal were used. Just follow the instructions below.

First get the number of rows and columns of you image by using gdalinfo (although it can be automated as well but lets just go with it)

gdalinfo raster.tif

Now edit the python code below and set rows and columns according to your image (don’t forget to change the image name).

from osgeo import gdal

# Open tif file
ds = gdal.Open('raster.tif')

# GDAL affine transform parameters, According to gdal documentation xoff/yoff are image left corner, a/e are pixel wight/height and b/d is rotation and is zero if image is north up. 
xoff, a, b, yoff, d, e = ds.GetGeoTransform()

def pixel2coord(x, y):
 """Returns global coordinates from pixel x, y coords"""
 xp = a * x + b * y + xoff
 yp = d * x + e * y + yoff
 return(xp, yp)

# get columns and rows of your image from gdalinfo
rows = 36+1
colms = 34+1

if __name__ == "__main__":
 for row in  range(0,rows):
  for col in  range(0,colms): 
   print pixel2coord(col,row)

Finally, execute this python script using following command (provided that you have named your file script.py)

python script.py

Here is the output of my raster.tif

(29.25, 6.0)
(29.5, 6.0)
(29.75, 6.0)
(30.0, 6.0)
(30.25, 6.0)
(30.5, 6.0)
(30.75, 6.0)
(31.0, 6.0)
(31.25, 6.0)
(31.5, 6.0)
.
.
.

In my case, I converted output to csv and opened it in qGIS along with raster.tif (shown below).

for_fb2

and a closer look

for_fb

 

20 thoughts on “Latitude/Longitude of each pixel using python and gdal

    • Modify the following code

      # csv library
      library csv

      # writer for csv
      csvfile = open(‘/home/geoinformatics/myscripts/py/csv/sequence_fmnr.csv’, ‘wb’)
      csv_writer = csv.writer(csvfile, delimiter=’,’,quotechar=’|’, quoting=csv.QUOTE_MINIMAL)

      #now in the final loop use the writer instead of print command
      csv_writer.writerow(str_temp)

      • Thanks for the response. Your code works.
        We get points of leftmost corner. Can i get Lat/Long for centre of pixels(centroid). Do you know any way in python if that can be done?

        Thanks.

  1. Thank you. I got the centre coordinates by modifying the code. The code works awesome on Sentinel 2A pixels.
    But the code doesnt work out on Landsat pixels. In one pixel in Landsat i am getting many number of points through the code.
    Can you advise on this. Thanks

  2. Nice article. Thank You.

    Had a question here :
    How can I get inverse of this.
    def pixel2coord(x, y):
    “””Returns global coordinates from pixel x, y coords”””
    xp = a * x + b * y + xoff
    yp = d * x + e * y + yoff
    return(xp, yp)

    like from global coordinates to pixel ?

  3. Hi, This is really helpful. I’m relatively new to working with rasters in Python. Can you point me to some instruction about how to find rows and column with gdalinfo? Having trouble with this step. Thank you!

    • Thanks.

      first you need to install gdal. Then goto the folder where your raster is (or provide the path in the gdalinfo command). The output will have a section called “Size is” and that is the size of your raster (in terms of pixels) first value is column and the second is rows.

  4. Hi,
    I’m very new to python and am using a jupyter notebook. I’m getting this error:

    NameError Traceback (most recent call last)
    in ()
    2
    3 “””Returns global coordinates from pixel x, y coords”””
    —-> 4 xp = a * x + b * b + xoff
    5 yp = d * e + y * b + yoff
    6 return(xp, yp)

    NameError: name ‘x’ is not defined

Leave a reply to mnahmad Cancel reply