English

How To: Calculate Zonal Statistics for Overlapping Zones

Summary

Instructions provided describe a sample Python script that can be used to calculate zonal statistics for overlapping zones.

When using polygon features as the input zones in the Zonal Statistics as Table tool, the tool rasterizes the vector features. This results in data loss if overlapping polygons are present in the zone feature.

Using Python, the tool's limitations can be overcome by iteratively processing each feature to overcome such a problem.

Procedure

Code:
import arcpy, os, sys, string
from arcpy import env
arcpy.CheckOutExtension("spatial")
def CreateDirectory(DBF_dir):
if not os.path.exists(DBF_dir):
os.mkdir(DBF_dir)
print "created directory {0}".format(DBF_dir)

def ZonalStasAsTable(fc,DBF_dir,raster,zoneField):

for row in arcpy.SearchCursor(fc):
lyr = "Zone_{0}_lyr".format(row.OBJECTID)
tempTable = DBF_dir + os.sep + "zone_{0}.dbf".format(row.OBJECTID)
arcpy.MakeFeatureLayer_management(fc, lyr, "\"OBJECTID\" = {0}".format(row.OBJECTID))
print "Creating layer {0}".format(lyr)
out_layer = DBF_dir + os.sep + lyr + ".lyr"
arcpy.SaveToLayerFile_management(lyr, out_layer, "ABSOLUTE")
print "Saved layer file"
arcpy.gp.ZonalStatisticsAsTable_sa(out_layer, zoneField, raster, tempTable, "DATA", "ALL")
print "Populating zonal stats for {0}".format(lyr)
del row, lyr

def MergeTables(DBF_dir,zstat_table):
arcpy.env.workspace = DBF_dir
tableList = arcpy.ListTables()
arcpy.Merge_management(tableList,zstat_table)
print "Merged tables. Final zonalstat table {0} created. Located at {1}".format(zstat_table,DBF_dir)
del tableList
if __name__ == "__main__":
ws = "C:/TEMP"
DBF_dir = ws + os.sep + "DBFile"
fc = "C:/TEMP/zone/zone_data/zones.gdb/zones"
raster = r"C:/TEMP/zone/raster/zonalstat_DEM.img"
zoneField = "ZoneID"
zstat_table = DBF_dir + os.sep + "Zonalstat.dbf"
CreateDirectory(DBF_dir)
ZonalStasAsTable(fc,DBF_dir,raster,zoneField)
MergeTables(DBF_dir,zstat_table)