前言

最近需要用ArcGIS统计面内点的数量，搜索了一下，官方提供的工具为空间连接，挺好用，但是有一个缺点，一次只能计算一个点要素类，多个点要素类需要重复操作，比较繁琐。搜索一番网上没有较完美的解决方案，刚好最近在研究Arcpy，开始手搓轮子。

实现思路

Arcpy实现判定相交我想到的有两个方法，一个是官方的按位置选择图层，一个是遍历图层用Arcpy的Geometry类中的contains判定相交，此处使用的是前者，Arcpy调用官方工具也十分简单，在具体工具的运行键右边有个向下拉的部分，点击后选择复制成Python代码即可。

踩过的坑

新添加的图层无法使用SelectLayerByLocation

根据官方介绍，需要使用MakeFeatureLayer创建要素图层

MakeFeatureLayer后编辑字段会影响源要素属性

先CopyFeatures_management将要素类复制到数据库，再使用MakeFeatureLayer创建要素图层，操作完成之后不需要重新导出

组内要素的处理问题

.split("\\")[-1]从路径提取中的文件名

要素类行数获取

int(arcpy.management.GetCount("目标要素")[0])

成品代码

注释很详细，看不懂翻翻官方手册就行

# -*- coding: utf-8 -*-
"""
Tool:               统计
Source Name:        <File name>
Version:            Arcgis Pro 2.8
Author:             Chao
Usage:              <Command syntax>
Required Arguments: <point_parameter>
                    <surface_parameter>
Optional Arguments: <point_parameter2>
                    <point_parameter3>
Description:        统计面内多个点要素类的点数量
"""
import arcpy


def ScriptTool(point_parameter, surface_parameter, result_parameter):
    true_result_parameter = result_parameter
    # 根据绝对路径提取文件名
    result_parameter = result_parameter.split("\\")[-1]

    # 复制图层并添加到图层列表，不使用MakeFeatureLayer将无法使用SelectLayerByLocation，不copy直接创建图层要素修改时源要素会变化
    arcpy.CopyFeatures_management(surface_parameter, true_result_parameter)
    arcpy.management.MakeFeatureLayer(true_result_parameter, result_parameter)

    # 获取面文件行数，用于计算进度
    task_total_num = int(arcpy.management.GetCount(surface_parameter)[0])

    # 添加进度条
    arcpy.SetProgressor("step", f"即将开始任务,共计{task_total_num}个面",0, task_total_num , 1)

    # 将点文件文本切割为列表
    point_list = point_parameter.split(";")
    # 获取面要素已有字段
    fieldList = arcpy.ListFields(result_parameter)
    field_name_List = []
    for field in fieldList:
        field_name_List.append(field.name)

    for point_name in point_list:

        # 处理需要添加的字段
        if "\\" in point_name:
            point_name_str = point_name.split("\\")[-1]
        else:
            point_name_str = point_name
        # 判断是否包含该字段，没有则添加
        if point_name_str not in field_name_List:
            # 添加字段
            arcpy.management.AddField(result_parameter, point_name_str, "LONG", None, None, None, point_name_str,
                                      "NULLABLE","NON_REQUIRED",'')
    # 添加游标，用于遍历面要素
    updateCursor = arcpy.UpdateCursor(result_parameter)
    task_num = 1
    for row in updateCursor:
        OBJECTID = row.getValue("OBJECTID")
        # 选择面要素的一个面，通过按位置选择工具选出与其相交的点要素
        arcpy.management.SelectLayerByAttribute(result_parameter, "NEW_SELECTION", f"OBJECTID = {OBJECTID}", None)
        arcpy.management.SelectLayerByLocation(point_parameter, "INTERSECT", result_parameter, None, "NEW_SELECTION",
                                               "NOT_INVERT")
        # 修改进度条显示文字
        arcpy.SetProgressorLabel(f"开始统计面,当前第{task_num}/{task_total_num}个面")
        # 修改进度条进度，向后跳动一个间隔数
        arcpy.SetProgressorPosition()
        for point_name in point_list:
            
            
            num = int(arcpy.management.GetCount(point_name)[0])

            if "\\" in point_name:
                point_name_str = point_name.split("\\")[-1]
            else:
                point_name_str = point_name
            # 为所选字段赋值
            arcpy.management.CalculateField(result_parameter, point_name_str, str(num), "PYTHON3", '', "TEXT",
                                            "NO_ENFORCE_DOMAINS")
        task_num = task_num + 1

    # 清除选择的要素
    arcpy.management.SelectLayerByAttribute(result_parameter, "CLEAR_SELECTION")
    for point_name in point_list:
        arcpy.management.SelectLayerByAttribute(point_name, "CLEAR_SELECTION")

    arcpy.SetProgressorLabel(f"统计完成")
    return


if __name__ == '__main__':
    # 获取脚本输入信息
    # Getpoint_parameterAsText-将输入要素转换为文本，多选要素中间";"相隔
    surface_parameter = arcpy.GetParameterAsText(0)
    point_parameter = arcpy.GetParameterAsText(1)
    result_parameter = arcpy.GetParameterAsText(2)
    ScriptTool(point_parameter, surface_parameter,result_parameter)

脚本参数设置，如图

最后

工具实现的面内多点要素类的统计，统计信息以字段的形式储存在输出要素类中，另外实现了进度条信息的可视化。第二种实现方法在面和点的数量较少的时候运行估计会快些，暂时还没造轮子，留个坑。

补充第二种(contains)方法

脚本属性设置同第一种方法

# -*- coding: utf-8 -*-
"""
Tool:               统计
Source Name:        <File name>
Version:            Arcgis Pro 2.8
Author:             Chao
Usage:              <Command syntax>
Required Arguments: <point_parameter>
                    <surface_parameter>
Optional Arguments: <point_parameter2>
                    <point_parameter3>
Description:        统计面内多个点要素类的点数量
"""
import arcpy


def ScriptTool(point_parameter, surface_parameter, result_parameter):
    true_result_parameter = result_parameter
    # 根据绝对路径提取文件名
    result_parameter = result_parameter.split("\\")[-1]

    # 复制图层并添加到图层列表，不使用MakeFeatureLayer将无法使用SelectLayerByLocation，不copy直接创建图层要素修改时源要素会变化
    arcpy.CopyFeatures_management(surface_parameter, true_result_parameter)
    arcpy.management.MakeFeatureLayer(true_result_parameter, result_parameter)

    # 获取面文件行数，用于计算进度
    task_total_num = int(arcpy.management.GetCount(surface_parameter)[0])

    # 添加进度条
    arcpy.SetProgressor("step", f"即将开始任务,共计{task_total_num}个面",0, task_total_num , 1)

    # 将点文件文本切割为列表
    point_list = point_parameter.split(";")
    # 获取面要素已有字段
    fieldList = arcpy.ListFields(result_parameter)
    field_name_List = []
    for field in fieldList:
        field_name_List.append(field.name)

    for point_name in point_list:

        # 处理需要添加的字段
        if "\\" in point_name:
            point_name_str = point_name.split("\\")[-1]
        else:
            point_name_str = point_name
        # 判断是否包含该字段，没有则添加
        if point_name_str not in field_name_List:
            # 添加字段
            arcpy.management.AddField(result_parameter, point_name_str, "LONG", None, None, None, point_name_str,
                                      "NULLABLE","NON_REQUIRED",'')
    # 添加游标，用于遍历面要素
    updateCursor = arcpy.UpdateCursor(result_parameter)
    task_num = 1
    for row in updateCursor:
        surface_geometry = row.shape
        # 修改进度条显示文字
        arcpy.SetProgressorLabel(f"开始统计面,当前第{task_num}/{task_total_num}个面")
        # 修改进度条进度，向后跳动一个间隔数
        arcpy.SetProgressorPosition()
        for point_name in point_list:
            pointCursor = arcpy.UpdateCursor(point_name)
            num = 0
            for point_row in pointCursor:
                point_geometry =point_row.shape
                ## 判断点是否在面内
                if (surface_geometry.contains(point_geometry)):
                    num = num + 1

            if "\\" in point_name:
                point_name_str = point_name.split("\\")[-1]
            else:
                point_name_str = point_name
            # 为所选字段赋值
            row.setValue(point_name_str, str(num))
            updateCursor.updateRow(row)

            del pointCursor
        task_num = task_num + 1



    arcpy.SetProgressorLabel(f"统计完成")
    return


if __name__ == '__main__':
    # 获取脚本输入信息
    # Getpoint_parameterAsText-将输入要素转换为文本，多选要素中间";"相隔
    surface_parameter = arcpy.GetParameterAsText(0)
    point_parameter = arcpy.GetParameterAsText(1)
    result_parameter = arcpy.GetParameterAsText(2)
    ScriptTool(point_parameter, surface_parameter,result_parameter)
最后修改：2023 年 11 月 26 日
© 允许规范转载
如果觉得我的文章对你有用，请随意赞赏