Quiero una imagen de varianza local con un 3x3 de una imagen raster geoespacial usando python. Mi enfoque hasta ahora era leer en la banda raster como una matriz, a continuación, utilizando la notación matricial para ejecutar una ventana móvil y escribir la matriz en una nueva imagen raster. Este enfoque funcionó bien para un filtro de paso alto como se describe en este tutorial: http://www.gis.usu.edu/~chrisg/python/2009/lectures/ospy_slides6.pdf
Luego intenté calcular la varianza con varias aproximaciones, la última usando numpy (como np), pero sólo obtengo una imagen gris con el mismo valor en todas partes. Estoy abierto a cualquier tipo de solución. Si al final me da la varianza local media, sería aún mejor Gracias, Mace
rows = srcDS.RasterYSize
#read in as array
data = srcBand.ReadAsArray(0,0, cols, rows).astype(np.int)
#calculate the variance for a 3x3 window
outVariance = np.zeros((rows, cols), np.float)
outVariance[1:rows-1,1:cols-1] = np.var([(data[0:rows-2,0:cols-2]),
(data[0:rows-2,1:cols-1]),
(data[0:rows-2,2:cols] ),
(data[1:rows-1,0:cols-2]),
(data[1:rows-1,1:cols-1]),
(data[1:rows-1,2:cols] ),
(data[2:rows,0:cols-2] ),
(data[2:rows,1:cols-1] ),
(data[2:rows,2:cols] )])
#output
outDS = driver.Create(outFN, cols, rows, 1, GDT_Float32)
outDS.SetGeoTransform(srcDS.GetGeoTransform())
outDS.SetProjection(srcDS.GetProjection())
outBand = outDS.GetRasterBand(1)
outBand.WriteArray(outVariance,0,0)
...
