Histograma zonal para cada región (fila) de un polígono (tabla) con python script

Question

Tengo un DEM.tif (la simbología utiliza un intervalo igual - 100m) y un archivo de polígonos (regions.shp) con muchas regiones (los nombres están definidos en la columna 'reg_name'). Me gustaría utilizar el Histograma Zonal para producir tablas (y posteriormente archivos xlxs) con dos columnas: la primera debería ser la clase (200-300,300-400 etc.) y la segunda la frecuencia. Tengo DEM.tif y regions.shp en mi documento de mapa y si uso el código en la ventana de python para una región (sin usar un bucle for) funciona. Pero cuando quiero ejecutar el código para generar tablas para cada región no funciona. El problema principal es probablemente porque después de la línea extracted_raster=arcpy.sa.ExtractByMask(in_raster, lyr) el 'extracted_raster' no aparece en el índice o en la vista de datos. Así, arcpy.ApplySymbologyFromLayer_management("extracted_raster",in_raster) no se puede aplicar. Me he dado cuenta de que si ejecuto el código sólo hasta arcpy.ApplySymbologyFromLayer_management("extracted_raster",in_raster) el 'extracted_raster' se muestra en mi TOC sólo después de la última iteración (la última región de la lista). Todo 'extracted_raster' Los rastreos antes de la última región no se muestran... ¿Cómo puedo hacer que aparezcan? Estoy utilizando arcpy.RefreshActiveView() pero no muestra la trama... Y sé que el Histograma Zonal produce una tabla con columnas para todas las regiones a la vez, pero me gustaría obtener tablas con los mínimos y máximos reales en una tabla para cada región (sin clases con frecuencia cero antes y después). Aquí está mi script:

#Set MapDocument to "CURRENT"...
mxd=arcpy.mapping.MapDocument("CURRENT")

#Defining parameters...
in_raster='DEM.tif' #Raster to be extracted and used for histogram data production...
zone_layer='regions.shp' #This is a layer with regions...
zone_field='reg_name' #The name of the field with names of the regions...

#Code for listing names of regions...
def getValueList (inputTable, field):
    values = set()
    rows = arcpy.SearchCursor(inputTable)
    for row in rows:
            values.add(row.getValue(field))
    return sorted(values)
arcpy.AddMessage("Koda getValuesList je definirana!")

#Listing names of regions...
list=getValueList(zone_layer, zone_field)
lyr=arcpy.mapping.Layer(zone_layer)

where_clause=(str(zone_field))+"="

#Loop for
for zone in list:
    lyr.definitionQuery=where_clause+"'"+str(zone)+"'"
    arcpy.RefreshActiveView()
    extracted_raster=arcpy.sa.ExtractByMask(in_raster, lyr)
    arcpy.RefreshActiveView()
    arcpy.ApplySymbologyFromLayer_management("extracted_raster",in_raster)
    output_table="histogram_data"+str(zone)+".dbf"
    arcpy.sa.ZonalHistogram(lyr,zone_field,"extracted_raster",output_table)
print "Done!

"

Answer 1

Un par de pensamientos... Creo que deberías seguir el ejemplo aquí y guarda los datos extraídos si quieres que se conserven. No creo que se guarde de otra manera, por lo que se sobrescribe cuando se itera a través del bucle for. Podrías utilizar el mismo método que tienes para nombrar la tabla añadiendo tu zona al nombre de la ruta del raster. A continuación de la ayuda de Arc.

outExtractByMask = ExtractByMask("elevation", "mask.shp")
outExtractByMask.save("C:/sapyexamples/output/maskextract")

Puede utilizar MakeRasterLayer_management a partir de los rásteres extraídos, lo que creará capas temporales en su documento cartográfico. Entonces podrá aplicar su simbología.

Answer 2

Aquí está mi código que funciona. El resultado no se muestra en mi vista TOC o Datos (pensé que esto es necesario para la herramienta ApplySimbology..), pero funciona. Hice un script que puede ser utilizado como una herramienta.

import arcpy
from string import maketrans
arcpy.CheckOutExtension("Spatial")

in_raster=arcpy.GetParameterAsText(0)
zone_layer=arcpy.GetParameterAsText(1)
zone_field=arcpy.GetParameterAsText(2)
files_folder=arcpy.GetParameterAsText(3)
excel_files_folder=arcpy.GetParameterAsText(4)
tmpraster_folder=arcpy.GetParameterAsText(5)
tmpshp_buffer=arcpy.GetParameterAsText(6)

def getValueList (inputTable, field):
    values = set()
    rows = arcpy.SearchCursor(inputTable)
    for row in rows:
            values.add(row.getValue(field))
    return sorted(values)

list=getValueList(zone_layer, zone_field)
lyr=arcpy.mapping.Layer(zone_layer)

where_clause=(str(zone_field))+"="

inchar = " "
outchar = "_"
change = maketrans(inchar, outchar)

x=1
for zone in list:

    name_raster=str(x)+str(zone)[:7]
    strzone=str(zone)

    zone_name_corr=strzone.translate(change)
    lyr.definitionQuery=where_clause+"'"+str(zone)+"'"

    name_buffer=tmpshp_buffer+"\\"+name_raster+"_b50.shp"
    arcpy.Buffer_analysis(lyr, name_buffer, 50)
    extracted_raster=arcpy.sa.ExtractByMask(in_raster, name_buffer)
    extracted_raster.save(tmpraster_folder+"\\"+name_raster)
    raster_lyr=str(zone)+"ly"
    arcpy.MakeRasterLayer_management(tmpraster_folder+"\\"+name_raster, raster_lyr)
    arcpy.RefreshActiveView()
    arcpy.ApplySymbologyFromLayer_management(raster_lyr,in_raster)

    output_table=files_folder+"\\histogram_data"+zone_name_corr+".dbf"

    arcpy.sa.ZonalHistogram(lyr,zone_field,raster_lyr,output_table)

    x=x+1