# Optimising Min/Max temporal raster search in python/gdal

I need to search across ~100 1200x1200 raster tiles and produce a min and a max raster. Along the way i also want to exclude values that are over 1 or less than 0. What would be the quickest method to do this in python with gdal? I was previously loading all the rasters into memory and then looping over the pixels, but this seems the slowest possible method.

Thanks!

• try to use the matrix (extract it from the raster) and look for optimization methods to find the max and min value of big matrix in python. – Gago-Silva Aug 13 '13 at 11:42

## 2 Answers

You could use numpy max/min methods and a masked array.

For example:

``````from osgeo import gdal
import numpy
import glob

#Loop through and open all rasters and stack them into a 3d array
rasterpaths=glob.glob(r'*.tif')
for ras in rasterpaths:
ds=gdal.Open(ras)
dat=ds.GetRasterBand(1).ReadAsArray()
dat=dat[None,:,:]               #turn 2d array into 3d

try:stk=numpy.vstack((stk,dat)) #Do we already have a 3d stack?
except NameError:stk=dat        #Nope, this is the first time through the loop

#Create masked array where valid values are >=0 and <=1
stk=numpy.ma.masked_outside(stk,0,1)

#Get max and min
max=stk.max(axis=0)
min=stk.min(axis=0)

#Write out to new rasters with gdal if you like...
``````

For anyone in the future: I solved this by essentially only loading two bands into memory (i and i+1) and then comparing just these two. If i+1 met was greater/smaller it was accepted as the new maximum/minimum.

The theory is in very similar to https://stackoverflow.com/questions/11964450/python-order-a-list-of-numbers-without-built-in-sort-min-max-function