Skip to content

Instantly share code, notes, and snippets.

@robinkraft
Created February 10, 2012 23:48
Show Gist options
  • Save robinkraft/1794188 to your computer and use it in GitHub Desktop.
Save robinkraft/1794188 to your computer and use it in GitHub Desktop.
Load hdf file with GDAL and Python, get NDVI
from osgeo import gdal
import os
layer_dict = {"reli":"reliability", "qual":"Quality", "ndvi":"NDVI", "evi":"EVI"}
def print_data(layer, data):
print data
print "data:", layer
print "type:", data.dtype
print "mean:", data.mean()
print "max :", data.max()
print "min :", data.min()
return
def load_hdf(fname, layer="ndvi"):
print "\nLoading %s" % os.path.split(fname)[1]
# convenience lookup so we can use short names
layer_key = layer_dict.get(layer)
hdf = gdal.Open(fname)
sdsdict = hdf.GetMetadata('SUBDATASETS')
sdslist =[sdsdict[k] for k in sdsdict.keys() if '_NAME' in k]
sds = []
for n in sdslist:
sds.append(gdal.Open(n))
# make sure the layer we want is in the file
if layer_key:
full_layer = [i for i in sdslist if layer_key in i] # returns empty if not found
idx = sdslist.index(full_layer[0])
data = sds[idx].ReadAsArray()
print_data(layer, data)
return data
else:
print "Layer %s not found" % layer
@tkorting
Copy link

tkorting commented Nov 1, 2017

Thanks a lot for this source code.

@iqbalhabibie
Copy link

what type of MODIS do this source code use?

because i have error in line 7
Syntax Error: Missing parentheses in call to 'print': E:\load hdf file\load_hdf.py, line 7, pos 14
print data

@kannes
Copy link

kannes commented Sep 12, 2019

@iqbalhabibie That's a Python 2 vs Python 3 issue. The code above is in Python 2.

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