¿Obtener la extensión de cada polígono en shapefile usando ArcPy?

18

En ArcGIS 10 y Python quiero obtener la información de extensión (xmax, ymax, xmin, ymin) de cada uno de los polígonos en un shapefile.

Puedo obtener la extensión de todo el shapefile usando

file=r"D:\SCRATCH\ARCGIS0k_trc_tiles_TVM.shp"
desc=arcpy.Describe(file)
print desc.extent.Xmax
  

394551.52085039532

Pero no puedo encontrar la manera de obtener la misma información para cada fila en el conjunto de datos.

rows = arcpy.SearchCursor("100k_trc_tiles_TVM")
for row in rows:
 print row

imprime las 31 filas en el conjunto de datos pero

for row in rows:
 desc=arcpy.Describe(row)
 print desc.extent.Xmax

da un error.

  

Error de tiempo de ejecución: Objeto: Describe entrada   el valor no es un tipo válido

Estaba pensando en agregar los valores de extensión a la tabla usando "calcular geometría" pero esto solo da el centroide. Entonces supongo que podemos usar algo como row.GetValue ("xmax").

Dicho esto, sé que podemos crear X / Y, max / min utilizando la función de enlace pero sería mejor si pudiéramos evitar tener que agregar campos, especialmente si ArcPy puede obtener estos valores.

Básicamente, necesito que las extensiones se incorporen a la herramienta de recorte para recortar 30 áreas de datos (según las hojas de mapa de 1: 100,000) para el geoprocesamiento, ya que la herramienta Dividir falla debido al gran tamaño del conjunto de datos (ver ¿Por qué Intersect da ERROR 999999: Error al ejecutar la función Topología no válida [demasiados puntos finales de líneas)? ). Quiero automatizar esto ya que se repite en una serie de conjuntos de datos.

=== script de trabajo ===

# Emulates Arc Info SPLIT tool by using Clip but
# Requires a FC from which each row is used as the input clip feature.
# Each row must be rectangular.
# Used on 12GB FGDB with 100 million records.


#Licence: Creative Commons
#Created by: George Corea; [email protected], [email protected]
import arcpy, string

#inFrame=arcpy.GetParameterAsText(0) # Input dataframe FC
#inFile=arcpy.GetParameterAsText(1) # Input FC for splitting
#outDir=arcpy.GetParameterAsText(2) # Output FGDB

inFrame=r"D:\SCRATCH\ARCGIS0k_trc_tiles_TVM.shp"
inFile=r"c:\junk6\data_Merge.gdb\FullRez_m2b"
outDir=r"D:\SCRATCH\Projects6\datasplit\test_slaasp.gdb"
#NameField="Name_1"

#arcpy.env.workspace = r"C:/Workspace"
arcpy.env.overwriteOutput = True

rows = arcpy.SearchCursor(inFrame)
shapeName = arcpy.Describe(inFrame).shapeFieldName
for row in rows:
    feat = row.getValue(shapeName)
    Name = row.Name_1
    print "Executing clip on: "+str(Name)
    extent = feat.extent
    #print extent.XMin,extent.YMin,extent.XMax,extent.YMax
# Create an in_memory polygon
    XMAX = extent.XMax
    XMIN = extent.XMin
    YMAX = extent.YMax
    YMIN = extent.YMin
    pnt1 = arcpy.Point(XMIN, YMIN)
    pnt2 = arcpy.Point(XMIN, YMAX)
    pnt3 = arcpy.Point(XMAX, YMAX)
    pnt4 = arcpy.Point(XMAX, YMIN)
    array = arcpy.Array()
    array.add(pnt1)
    array.add(pnt2)
    array.add(pnt3)
    array.add(pnt4)
    array.add(pnt1)
    polygon = arcpy.Polygon(array)
    ShapeFile = outDir+"\temp_poly"
    arcpy.CopyFeatures_management(polygon, ShapeFile)

    #print Name
### Set local variables
    in_features = inFile
    clip_features = ShapeFile
    out_feature_class = outDir+"\"+Name
    xy_tolerance = "0.22"

    # Execute Clip

    try:
        arcpy.Clip_analysis(in_features, clip_features, out_feature_class, xy_tolerance)
        print "Completed: "+str(Name)
    except:
        error = arcpy.GetMessages()
        print "Failed on: "+str(Name)+" due to "+str(error)
    
pregunta GeorgeC 04.11.2011 - 07:11

6 respuestas

25

Obtenga el objeto de forma en su cursor y acceda a su propiedad de extensión. Consulte la Ayuda de ArcGIS Trabajo con geometría en Python :

shapeName = arcpy.Describe(inFeatures).shapeFieldName
for row in rows:
    feat = row.getValue(shapeName)
    extent = feat.extent
    print extent.XMin,extent.YMin,extent.XMax,extent.YMax
    
respondido por el Luke 04.11.2011 - 20:35
7

El conjunto de herramientas Contenedor de límites hace exactamente lo que desea. Si solo desea fragmentos de código, examine las funciones dentro de los scripts, uno trata explícitamente la extensión.

EDIT

Debo agregar que la secuencia de comandos agregará valores a los campos Izquierda, Derecha, Superior e Inferior en el archivo de salida creado que se pueden usar para el procesamiento posterior

    
respondido por el user681 04.11.2011 - 10:16
2

Acabo de probar Geometría de delimitación mínima (sobre) (en Gestión de datos) en ArcGIS 10 y parece que hace exactamente lo mismo, para todos los campos.

    
respondido por el Sid 22.10.2012 - 12:55
1

Como se describe en ¿Extracción de coordenadas de vértices de polígono en ArcMap?  puede obtener los vértices de un polígono y luego agregar las coordenadas x e y de cada vértice como campos en la tabla de atributos. Esto tiene la limitación de no adjuntar las coordenadas max / min directamente a cada polígono, pero esto se puede lograr de varias maneras.

El método con el que estoy más familiarizado es leer los campos x e y en las listas de python utilizando el módulo pyshp , que luego puede ordenarse para encontrar los valores máximo y mínimo para cada polígono. Pyshp se puede usar para abrir una clase de escritor para agregar nuevos campos a los polígonos originales y escribir estos valores máximos y mínimos en el polígono correcto.

Creo que esto se puede hacer usando el arcpy, pero tuve muchos problemas al escribir en shapefiles en 9.3 usando el geoprocesador, así que prefiero el método pyshp, pero no estoy seguro si el módulo arcpy ha resuelto estos problemas.

    
respondido por el sgrieve 18.09.2012 - 09:30
0

¿Intentó usar mayúsculas en la "M" en "XMax"? Creo que se supone que es:

print desc.extent.XMax

en lugar de

print desc.extent.Xmax

de acuerdo con la documentación . Por supuesto, eso me hace preguntarme cómo funcionó tu primer fragmento de código. De cualquier manera, ¡pruébalo!

    
respondido por el dmahr 04.11.2011 - 14:10
0

Otra forma sería hacer un SearchCursor () en el shapefile, luego puedes usar row.shape.extent:

rows = arcpy.SearchCursor(shapefileName)

for row in rows:
   extent = row.shape.extent
   ...
   ...
    
respondido por el Ruth 20.01.2012 - 15:45

Lea otras preguntas en las etiquetas