Calcular área dentro de la secuencia de comandos de Python en ArcMap

14

Estoy tratando de calcular el área de un polígono dentro de mi script Python. Creo un nuevo polígono fusionando dos juntos, y me gustaría agregar el área del polígono resultante a un campo en el archivo de salida. El polígono se almacena en un archivo de forma regular y se proyecta. Área preferiblemente en unidades de mapa.

Hubiera pensado que se trataba de una tarea bastante común y simple, pero a pesar de mucho Googleing, hasta ahora no he podido encontrar soluciones que funcionen.

Estaba planeando usar arcpy.updateCursorpara insertar el valor una vez que se calcula (solo hay una característica en el FC en esta etapa), por lo que es más fácil si se puede devolver como una variable. Cualquier solución alternativa que realice la misma tarea (obtener el valor del área en el campo correcto) también funcionará.

También probé la calculadora de campo de Python. Modificado de las páginas de ayuda, pensé que lo siguiente funcionaría, pero hasta ahora no tuve suerte.

arcpy.AddField_management(tempPgs, "Shape_area", 'DOUBLE')
exp = "float(!SHAPE.AREA!.split())"
arcpy.CalculateField_management(tempPgs, "Shape_area", exp)

Ejecutando ArcGIS Basic 10.1 SP1 con Python 2.7 en Windows 7.

Las partes relevantes de mi código actual se ven así:

#/.../
arcpy.Copy_management(inpgs, outpgs)
arcpy.AddField_management(outpgs, 'Shape_area', 'LONG')
fields = AM.FieldLst(outpgs)

#/.../

# Identify and search for shapes smaller than minimum area
where1 = '"' + 'Shape_Area' + '" < ' + str(msz)
polyrows = arcpy.SearchCursor(inpgs, where1)

for prow in polyrows:
    grd1 = prow.GridID   # GridID on the current polygon
    grd2 = nDD.get(grd1) # GridID on the polygon downstream

    # Update features
    if grd2
        geometry1 = prow.Shape
        geometry2 = geometryDictionary[grd2]

        # Update temporary features
        arcpy.Merge_management([geometry1, geometry2], tempMerged)
        arcpy.Dissolve_management(tempMerged, tempPgs)

        fds = AM.FieldLst(tempPgs)

        for field in fields[2:]:
            arcpy.AddField_management(tempPgs, field, 'DOUBLE')

        for fd in fds[2:]:
            arcpy.DeleteField_management(tempPgs, fd)

        exp = "float(!SHAPE.AREA!.split())"
        arcpy.CalculateField_management(tempPgs, "Shape_area", exp)

        # Append them to output FC
        try:
            arcpy.Append_management(tempPgs, outpgs, "TEST")
        except arcgisscripting.ExecuteError:
            arcpy.Append_management(tempPgs, outpgs, "NO_TEST")

    elif ...

    else ...
Martín
fuente
¿Cuál es su tipo de salida? Shapefile, geodatabase de archivos, ¿algo más? ¿Su archivo de salida es proyectado o no proyectado?
blord-castillo
Además, ¿podría publicar un poco más de la muestra de código, en particular el cursor que está utilizando para realizar la actualización? Lo más probable es que pueda lograr lo que desea usando SHAPE@AREAcomo parte de su cursor para leer el área; pero la estructura del código depende de si su área está en las mismas unidades que lo que desea escribir.
blord-castillo

Respuestas:

29

Hay tres formas diferentes de encontrar y almacenar el área de polígonos en una clase de entidad con arcpy: 1) calculadora de campo, 2) cursores de arco "clásicos" y 3) arcpy.dacursores. Algo de esto se tomó prestado de mi respuesta anterior sobre el uso de SearchCursor .


1. Calculadora de campo

  • Cuando se usa la calculadora de campo, hay tres tipos de expresiones diferentes que usan analizadores de expresiones diferentes. Esto se especifica en el tercer parámetro de la herramienta de geoprocesamiento Calcular campo . Al acceder a las propiedades del objeto Geometry usando like in !shape.area!, debe usar el analizador Python 9.3.

  • La expresión que tenía antes realizaba un split()comando sobre el resultado de !SHAPE.AREA!. Esto devuelve un listobjeto Python , que no se puede convertir en un floatobjeto.

  • En su expresión, puede especificar la unidad del área devuelta utilizando la @SQUAREKILOMETERSbandera, reemplazando SQUAREKILOMETERScon las unidades en la página de ayuda Calcular campo .

Aquí está el código de Python que usaría para este método:

tempPgs = "LayerName"
arcpy.AddField_management(tempPgs, "Shape_area", "DOUBLE")
exp = "!SHAPE.AREA@SQUAREKILOMETERS!"
arcpy.CalculateField_management(tempPgs, "Shape_area", exp, "PYTHON_9.3")

2. Arc 10.0 - Cursores "clásicos"

  • Cuando se usan cursores clásicos (es decir arcpy.UpdateCursor), el objeto cursor es un objeto iterativo que contiene rowobjetos. Debe usar los métodos getValuey setValuepara obtener la geometría de la fila (como un objeto de geometría y establecer el valor del área rowcomo flotante).

  • Su fila de salida se almacena en un espacio temporal temporal hasta que llame al updateRowmétodo en el cursor. Esto guarda los nuevos datos en el conjunto de datos real.

Aquí está el código de Python que usaría para este método:

tempPgs = "LayerName"
arcpy.AddField_management(tempPgs, "Shape_area", "DOUBLE")
geometryField = arcpy.Describe(tempPgs).shapeFieldName #Get name of geometry field
cursor = arcpy.UpdateCursor(tempPgs)
for row in cursor:
    AreaValue = row.getValue(geometryField).area #Read area value as double
    row.setValue("Shape_area",AreaValue) #Write area value to field
    cursor.updateRow(row)
del row, cursor #Clean up cursor objects

3. Arc 10.1 - cursores arcpy.da

  • Al usar los nuevos cursores en el módulo de acceso a datos (es decir arcpy.da.UpdateCursor), debe pasar una lista de nombres de campo como el segundo parámetro en el constructor del cursor. Esto requiere un poco más de trabajo por adelantado, pero los rowobjetos resultantes son listas de Python, lo que hace que sea más fácil leer y escribir datos al recorrer las filas del cursor. arcpy.da.UpdateCursortambién tiene un mejor rendimiento que arcpy.UpdateCursor, en parte porque omite campos sin importancia, especialmente la geometría.

  • Al leer la geometría, puede elegir una de una serie de fichas de geometría, por ejemplo SHAPE@TRUECENTROID, SHAPE@AREAo SHAPE@. El uso de un token "más simple" mejora en gran medida el rendimiento en comparación con el SHAPE@que contiene toda la información de geometría. La lista completa de tokens se encuentra en la arcpy.da.UpdateCursorpágina de ayuda.

  • Como antes, su fila de salida se almacena en un espacio temporal temporal hasta que llame al updateRowmétodo en el cursor. Esto guarda los nuevos datos en el conjunto de datos real.

Aquí está el código de Python que usaría para este método:

tempPgs = "LayerName"
arcpy.AddField_management(tempPgs, "Shape_area", "DOUBLE")
CursorFieldNames = ["SHAPE@AREA","Shape_area"]
cursor = arcpy.da.UpdateCursor(tempPgs,CursorFieldNames)
for row in cursor:
    AreaValue = row[0].area #Read area value as double
    row[1] = AreaValue #Write area value to field
    cursor.updateRow(row)
del row, cursor #Clean up cursor objects
dmahr
fuente
55
Maravillosa respuesta. Solo quería decir que a partir de 10.2, simplemente haría row[1] = row[0]ya que ya no hay un areaatributo. También puede usar el cursor como administrador de contexto en una withdeclaración y no tener que preocuparse por eliminar nada.
Paul H