Skip to content

Instantly share code, notes, and snippets.

@stefanocudini
Created March 20, 2013 01:45
Show Gist options
  • Star 11 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save stefanocudini/5201689 to your computer and use it in GitHub Desktop.
Save stefanocudini/5201689 to your computer and use it in GitHub Desktop.
read pixel value from GeoTiff raster file by lat lon
#!/usr/bin/env python
#http://geoinformaticstutorial.blogspot.it/2012/09/reading-raster-data-with-python-and-gdal.html
#http://www.gis.usu.edu/~chrisg/python/2009/lectures/ospy_slides4.pdf
from osgeo import gdal,ogr
from osgeo.gdalconst import *
import struct
import sys
lon = 12.502742
lat = 42.243713
lat = float(sys.argv[2])
lon = float(sys.argv[3])
def pt2fmt(pt):
fmttypes = {
GDT_Byte: 'B',
GDT_Int16: 'h',
GDT_UInt16: 'H',
GDT_Int32: 'i',
GDT_UInt32: 'I',
GDT_Float32: 'f',
GDT_Float64: 'f'
}
return fmttypes.get(pt, 'x')
ds = gdal.Open(sys.argv[1], GA_ReadOnly)
if ds is None:
print 'Failed open file'
sys.exit(1)
transf = ds.GetGeoTransform()
cols = ds.RasterXSize
rows = ds.RasterYSize
bands = ds.RasterCount #1
band = ds.GetRasterBand(1)
bandtype = gdal.GetDataTypeName(band.DataType) #Int16
driver = ds.GetDriver().LongName #'GeoTIFF'
success, transfInv = gdal.InvGeoTransform(transf)
if not success:
print "Failed InvGeoTransform()"
sys.exit(1)
px, py = gdal.ApplyGeoTransform(transfInv, lon, lat)
structval = band.ReadRaster(int(px), int(py), 1,1, buf_type = band.DataType )
fmt = pt2fmt(band.DataType)
intval = struct.unpack(fmt , structval)
print round(intval[0],2) #intval is a tuple, length=1 as we only asked for 1 pixel value
@PARPedraza
Copy link

My friend, I have a question. I try your code. I'd read a raster, but I have an error.
What lat lon need your code?

raster='C:...\IMG_DATA\T13QFD_20171020T172341_B01.jp2' # input raster
lon = -103.538
lat = 21.144
success, transfInv = gdal.InvGeoTransform(transf)
ValueError: too many values to unpack

@logi-26
Copy link

logi-26 commented Jul 31, 2018

I was going to use this script, but I discovered that GDAL has a built in function for this very purpose (maybe it did not when this script was developed). It is called gdallocationinfo: https://www.gdal.org/gdallocationinfo.html
To read a pixel's value using lat/lon co-ordinates, just specify the '-wgs84' flag.

@weagle08
Copy link

needed this in nodejs, and this solution got me there. thanks!

@stefanocudini
Copy link
Author

@weagle08 In nodejs I read a tiff using a spawn of command gdallocationinfo

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment