Вычислить площадь в скрипте Python в ArcMap

14

Я пытаюсь вычислить площадь многоугольника в моем скрипте Python. Я создаю новый многоугольник из слияния двух вместе, и я хотел бы добавить область полученного многоугольника к полю в выходном файле. Многоугольник хранится в обычном шейп-файле и проецируется. Площадь желательно в единицах карты.

Я бы подумал, что это довольно распространенная и простая задача, но, несмотря на большое количество Google, я до сих пор не смог найти работающие решения.

Я планировал использовать arcpy.updateCursorдля вставки значения после его вычисления (на этом этапе в FC есть только одна особенность), поэтому проще всего будет вернуть его в виде переменной. Любое альтернативное решение, которое выполняет ту же задачу (получение значения площади в правильном поле), также будет работать.

Я также попробовал Калькулятор поля от Python. Модифицированный со страниц справки я думал, что следующее будет работать, но пока не повезло.

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

Запуск ArcGIS Basic 10.1 SP1 с Python 2.7 в Windows 7.

Соответствующие части моего текущего кода выглядят так:

#/.../
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 ...
Мартин
источник
Какой у вас тип выхода? Shapefile, файловая база геоданных, что-то еще? Ваш выходной файл проецируется или не проецируется?
blord-castillo
Кроме того, не могли бы вы опубликовать немного больше примера кода, в частности курсор, который вы используете для обновления? Скорее всего, вы можете достичь того, что вы хотите, используя SHAPE@AREAкак часть вашего курсора, чтобы прочитать область; но структура кода зависит от того, находится ли ваша область в тех единицах, которые вы хотите выписать.
blord-castillo

Ответы:

29

Существует три различных способа найти и сохранить область многоугольника в классе пространственных объектов с помощью arcpy: 1) калькулятор поля, 2) «классические» курсоры arcpy и 3) arcpy.daкурсоры. Часть этого заимствована из моего предыдущего ответа об использовании SearchCursor .


1. Полевой калькулятор

  • При использовании полевого калькулятора существует три различных типа выражений, которые используют разные анализаторы выражений. Это указано в третьем параметре инструмента геообработки Calculate Field . При доступе к свойствам объекта Geometry с использованием like in !shape.area!, вы должны использовать анализатор Python 9.3.

  • Выражение, которое вы имели до того, split() команду для результата !SHAPE.AREA!. Это возвращает listобъект Python , который нельзя преобразовать в floatобъект.

  • В своем выражении вы можете указать единицу возвращаемой области, используя @SQUAREKILOMETERS флаг, заменяя SQUAREKILOMETERSединицы на странице справки « Рассчитать поле» .

Вот код Python, который я бы использовал для этого метода:

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 - «Классические» курсоры

  • При использовании классических курсоров (то есть arcpy.UpdateCursor) объект курсора является итеративным, содержащим rowобъекты. Вы должны использовать getValueи setValueметоды , чтобы получить геометрию из ряда (как объекта геометрии и установите значение площади в rowкачестве поплавка.

  • Ваша выходная строка хранится во временном пустом месте, пока вы не вызовете updateRow метод на курсоре. Это сохраняет новые данные в фактический набор данных.

Вот код Python, который я бы использовал для этого метода:

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 - курсоры arcpy.da

  • При использовании новых курсоров в модуле доступа к данным (т. Е. arcpy.da.UpdateCursor) Необходимо передать список имен полей в качестве второго параметра в конструкторе курсора. Это требует дополнительной работы заранее, но получающиеся rowобъекты представляют собой списки Python, что упрощает чтение и запись данных при переборе строк курсора. arcpy.da.UpdateCursorтакже имеет лучшую производительность, чем arcpy.UpdateCursorчастично, потому что пропускает неважные поля, особенно геометрию.

  • При чтении геометрии, вы можете выбрать один из нескольких геометрических лексем, например SHAPE@TRUECENTROID, SHAPE@AREA, или SHAPE@. Использование «более простого» токена значительно повышает производительность по сравнению с тем SHAPE@, который содержит всю информацию о геометрии. Полный список токенов находится на arcpy.da.UpdateCursorстранице справки.

  • Как и прежде, ваша выходная строка хранится во временном временном пространстве, пока вы не вызовете updateRow метод для курсора. Это сохраняет новые данные в фактический набор данных.

Вот код Python, который я бы использовал для этого метода:

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
источник
5
Прекрасный ответ. Просто хотел сказать, что с 10.2 вы бы просто сделали, так row[1] = row[0]как areaатрибута больше нет . Вы также можете использовать курсор в качестве диспетчера контекста в withвыражении и не беспокоиться об удалении чего-либо.
Пол Н