Using variable window ("kernel") size for raster processing in QGIS?
Input square window size raster:
Script:
import arcpy, os, traceback, sys, numpy
from arcpy import env
env.overwriteOutput = True
dem=r'D:\Scratch\demG'
wsize=r'D:\Scratch\wsize'
d=arcpy.Describe(dem)
try:
def showPyMessage():
arcpy.AddMessage(str(time.ctime()) + " - " + message)
demArray = arcpy.RasterToNumPyArray(dem,"","","",-9999)
sizeArray=arcpy.RasterToNumPyArray(wsize,"","","",-9999)
blankArray=arcpy.RasterToNumPyArray(dem,"","","",-9999)
d=arcpy.Describe(wsize)
origin=d.extent.lowerLeft
cSize=arcpy.Raster(wsize).meanCellHeight
nRows,nCols=demArray.shape
cellsTotal=nCols*nRows
# main loop
arcpy.SetProgressor("step", "", 0, nRows)
for nRow in range (nRows):
for nCol in range (nCols):
S=sizeArray[nRow,nCol]
shift=S/2
top=max(0,nRow-shift)
bottom=min(nRows-1,nRow+shift)
rows=demArray[top:bottom+1]
left=max(0,nCol-shift)
right=min(nCols-1,nCol+shift)
kernel=rows[:,range(left,right+1)]
nR,nC=kernel.shape
blankArray[nRow,nCol]=float(numpy.sum(kernel))/nR/nC
arcpy.SetProgressorPosition()
myRaster = arcpy.NumPyArrayToRaster(blankArray,origin,cSize,cSize)
oneGrid=Con(myRaster<>-9999,myRaster)
oneGrid.save(r'D:\Aerials\nztm\toriver')
del demArray,blankArray,sizeArray
except:
message = "\n*** PYTHON ERRORS *** "; showPyMessage()
message = "Python Traceback Info: " + traceback.format_tb(sys.exc_info()[2])[0]; showPyMessage()
message = "Python Error Info: " + str(sys.exc_type)+ ": " + str(sys.exc_value) + "\n"; showPyMessage()
Hillshade before:
and after