I have problems while using the ZonalHistogram function in Python. It seems to give strange results, while, with the same code run in arcgis' python window everything works fine. I am trying to study the snow cover per elevation ranges, thus including in the function a MODIS image (showing snow cover) as the input raster, and a reclissified DEM as the input value raster. I tried to modify the reclassification many times, but still not solving the problem... Moreover, inverting the input raster with the input value raster (meaning, inverting the function parameters), everything works fine, but it's a mess since I'm working on thousands of images and having the MODIS not always the same number of columns (while being the DEM always the same). Hope someone could help me. Here is my script, so that you can have a look at it in order to better understand my request:
import arcpy, os, sys, shutil
newpath= r'G:\\Chile_2012\\MODIS10A2.V005\\provvisoria\\Zonal_Histogram\\temporary'
if not os.path.exists (newpath):
os.makedirs (newpath)
arcpy.env.workspace= 'G:\\Chile_2012\\MODIS10A2.V005\\Cobertura_nival\\Cuencas_Copia\\campione'
arcpy.env.overwriteOutput= True
from arcpy.sa import *
listRaster= arcpy.ListRasters ()
elev_ranges= 'G:\\Chile_2012\\MODIS10A2.V005\\DEM\\elev_ranges'
out_dir= 'G:\\Chile_2012\\MODIS10A2.V005\\provvisoria\\Zonal_Histogram\\'
try:
for raster in listRaster:
arcpy.CheckOutExtension ('Spatial')
ZonalHistogram (raster, 'Value', elev_ranges, newpath+ '\\'+ raster)
arcpy.ClearEnvironment ('workspace')
arcpy.env.workspace= newpath
listTable = arcpy.ListTables ()
print listTable
arcpy.Merge_management (listTable, out_dir+ 'Rio_Aconcagua_Merge.dbf')
shutil.rmtree (newpath)
except:
print arcpy.GetMessages ()
I would add that my idea would be to work with layers (by means of functions that could transform rasters in layers) since the table returned in arcgis has a "LABEL" field (coming from the classification in symbology) which is the thing I would like to include in the final table and which is not included while working out of arcgis (in pywin), in order to know on what I am working on.