三、Arcpy基础【ArcGIS Python系列】

时间：2023年8月24日 阅读时间：15~20分钟 在阅读本篇和后续篇章时，推荐要有一点点Python基础。 代码练习notebook为4.2.3-arcpy基础(代码练习).ipynb

1.地理数据处理基础知识

ArcGIS中的地理处理允许您执行空间分析和建模以及自动执行GIS任务。典型的地理处理工具获取输入数据（要素类、栅格或表），执行地理处理任务，然后生成输出数据作为结果。ArcGIS包含数百种地理处理工具。地理处理工具的示例包括用于创建缓冲区、用于向表添加字段以及用于对地址表进行地理编码的工具。

地理处理通过创建组合了一系列工具的序列来支持工作流的自动化。一个工具的输出实际上成为下一个工具的输入。通过使用模型（model builder）和脚本，可以在ArcGIS中地理处理工具的自动化工作流。

2.从导入ArcPy开始

ArcPy 包含许多模块、类和函数，这使得可以在 Python 脚本中使用 ArcGIS Pro 中的所有地理处理工具。

import arcpy

执行上述语句后，就可以运行随 ArcGIS Pro 安装的工具箱中的所有地理处理工具。包括用于处理数据的模块 (arcpy.da)、地图脚本模块 (arcpy.mp)、用于图像分析和解释的模块 (arcpy.ia) 以及用于地图代数和栅格分析的模块（arcpy.sa）。导入 ArcPy 后，您就可以开始使用其模块、函数和类。

在脚本中导入 ArcPy 不仅会导入 ArcPy 的功能，还会执行两项重要检查：ArcPy 的可用性和许可证的可用性。如果输出RuntimeError: NotInitialized错误消息，请参照4.2.1-arcpy介绍和安装.md安装ArcGIS Pro。

3.设置工作区

首先得理解Python中绝对路径和相对路径的区别，简单提示一下：

• 绝对路径（Absolute Path）是从文件系统的根目录开始的完整路径。它包含了从根目录到目标文件或目录的所有目录层级。在不同的操作系统中，根目录的表示方式可能不同。例如，在Windows系统中，绝对路径可以以盘符（如C:\）开始，而在Linux或Mac系统中，绝对路径以斜杠（/）开始。在代码中如果是反斜杠 ""，应该改为 “/”(正斜杠)或''\'(两个反斜杠)。或者写成r"C:/data1"。
• 相对路径（Relative Path）是相对于当前工作目录的路径。当前工作目录是指运行Python程序时所在的目录。相对路径指定了从当前工作目录到目标文件或目录的路径。相对路径可以是简单的文件名或目录名，也可以是包含目录层级关系的路径。

ArcPy中的工作空间指定的就是工作目录，对应的可以使用相对路径引用。独立的 Python 脚本默认情况有一个当前工作目录，默认情况下该目录是脚本的位置。当设置 arcpy.env.workspace 时，ArcGIS Pro 将会在该路径下查找和操作数据。

tip: 您可以使用 os.getcwd() 获取当前工作目录，并且可以使用 os.chdir("/path") 更改当前工作目录。这样我们就能够在工作目录中使用相对路径指定路径了，保证了代码的可移植性。

arcpy.env.workspace 本质是一个Python类

注意理解env 是一个Python类（class），workspace 是该类的一个属性（property）。arcpy.env.workspace对应arcpy.<class>.<property>，所以arcpy.<classname>.<property> = <value>就是工作空间的属性值。

例如，你有一个名为 "C:\Data" 的文件夹，其中包含了你要使用的地理数据，你可以通过以下方式将它设置为工作空间：

import arcpy

arcpy.env.workspace = r"C:\Data"

# 创建地理数据库
arcpy.CreateFileGDB_management(arcpy.env.workspace, "myGDB.gdb") # 在工作空间下创建名为myGDB.gdb的地理数据库

在这个例子中，arcpy.env.workspace 被设置为 "C:\Data"，这意味着在执行地理处理脚本时，ArcGIS Pro 将会在该文件夹下查找和操作数据。

使用 arcpy.env.workspace 的好处是，它可以确保地理处理脚本在不同的环境中都能正常工作，无论是在 Windows 还是其他操作系统上。它提供了一种统一的方式来设置工作空间，使得脚本可以在不同的计算机上或不同的工作目录中运行，而不需要手动更改路径。（此方法和python的相对路径的作用相同）例如你可以这样指定工作空间：

import os

# 在整个脚本前指定一次绝对路径
data_dir = r'C:\Users\<用户名>\Documents\Python_\Github\arcgis-notebooks-tutorial\hurricane_analysis\data'
# 以后路径都是用相对路径 利用os.path.join处理路径能避免许多问题
hurricanes_raw_dir = os.path.join(data_dir,'hurricanes_raw')

# 利用mkdir创建检查和创建目录 
if not os.path.exists(hurricanes_raw_dir):
    os.mkdir(os.path.join(data_dir,'hurricanes_raw'))

总而言之，arcpy.env.workspace 是一个用于设置地理处理脚本工作空间的变量，它确保脚本能够在不同的环境中正确访问和操作数据。

4.使用地理处理工具

ArcPy 使您可以访问 ArcGIS Pro 中的所有地理处理工具。打开软件能看到有很多地理处理工具：

image-20230824120611163

（1）调用工具的方法两种方法

• 作为python函数调用：arcpy.<toolname_toolboxalias>(<parameters>)

例如调用裁剪工具：

import arcpy
arcpy.env.workspace = "C:/Data"
arcpy.Clip_analysis("streams.shp", "study.shp", "result.shp")

• 还可以通过使用与工具箱别名匹配的模块来使用工具。首先将工具箱作为模块调用，然后将工具作为该模块的函数调用，最后是工具的参数。 arcpy.<toolboxalias>.<toolname>(<parameters>)

import arcpy
arcpy.env.workspace = "C:/Data"
arcpy.analysis.Clip("streams.shp", "study.shp", "result.shp")

两种方法都是正确的。

小tips：

1. Python区分大小写，Clip不等于clip
2. 在代码行中空格对执行没有影响，但是对可读性有影响👀。
3. 函数和变量之间不要有空格，() 正确， <toolname> (<parameters>)不正确。
4. arcpy .analysis. Clip()不产生错误但是不是正确格式。

（2）示例：使用缓冲区buffer

而具体函数的使用可以参照帮助文档。我们从缓冲区buffer的帮助文档工具中举例说明：

缓冲工具图示

1. 程序中通过搜索找到buffer工具，可以看到通过此地图处理工具的可视化操作的参数：带星号的是必填此参数，分别是输入要素、输出要素类和距离。

image-20230809141458109

1. 点击右上角的问号❔，进入帮助页面，找到Python代码:

image-20230809122137740

1. 可以看到缓冲工具有四个参数，"{}"中的参数是可选的，可以不填或者用None和""表示, 其余的是必选参数。缓冲工具的语法是:

arcpy.analysis.Buffer(in_features, out_feature_class, buffer_distance_or_field, {line_side}, {line_end_type}, {dissolve_option}, {dissolve_field}, {method})

Python函数中的参数和软件中的参数基本是一一对应的，前三个参数为必选参数。我们简单浏览整个表格后然后一一说明：

1）必选参数

2）可选参数：

看起来很复杂，如果你了解ArcGIS Pro中使用缓冲区的方法，其实只需要将相应参数填入函数内就可以了，例如：

import arcpy, os

# 用os.getcwd()设置文件夹
data1_dir = os.path.join(os.getcwd(), "resource/data1")

arcpy.env.workspace = data1_dir

# 将shp导入数据库
arcpy.CopyFeatures_management("streets.shp", os.path.join("demo.gdb", "streets"))

#  可以更改工作空间
arcpy.env.workspace = os.path.join(data1_dir, "demo.gdb") # 改成你的data1文件夹位置

# 缓冲区分析
arcpy.analysis.Buffer("streets", "streets_Buffered", "20 Meters", "FULL", "ROUND", "LIST", "LABEL_CLAS")

我们将输出的文件streets_Buffered拖入地图中：

image-20230812071525712

以上也可以改写成以下形式方便阅读：

arcpy.analysis.Buffer(in_features="streets",
                      out_feature_class= "streets_Buffered_1", 
                      buffer_distance_or_field="20 Meters", 
                      line_side="FULL", 
                      line_end_type="ROUND", 
                      dissolve_option="LIST", 
                      dissolve_field="LABEL_CLAS")

你也可以单独定义变量，方便代码复用和制作脚本：

in_features="streets"
out_feature_class="streets_Buffered_2"
buffer_distance_or_field="20 Meters"
line_side="FULL"
line_end_type="ROUND"
dissolve_option="LIST"
dissolve_field="LABEL_CLAS"

arcpy.analysis.Buffer(in_features,out_feature_class, buffer_distance_or_field, line_side, line_end_type, dissolve_option, dissolve_field)

（3）调用工具的小技巧

在软件中打开某个地理处理工具的界面，你可以点击右上角的问号进入帮助页面，在填好参数之后就可以将此操作复制为python命令，甚至你可以在运行完地理处理工具确保没有错误之后再复制为python命令。如图所示：

你也可以打开历史记录，复制你运行过的命令：

5.空间参考

在 ArcGIS 中，为什么统一地图标准，通常使用两种坐标系：

1. 地理坐标系（Geographic Coordinate System）：地理坐标系使用经度（Longitude）和纬度（Latitude）来表示地球上的位置。经度表示东西方向上的位置，纬度表示南北方向上的位置。常用的地理坐标系包括经度-纬度坐标系，如 WGS84（世界大地测量系统1984）， GCJ02（火星坐标系），BD09（百度坐标系）等。前者是目前GPS使用的坐标系，后两者是国内使用常使用的坐标系，被加密，WGS84转后者可以使用百度或高德提供的地图转换服务，反过来转为WGS84需要用单独的方法。此处有吐槽。
2. 投影坐标系（Projected Coordinate System）：投影坐标系是将地理坐标系映射到平面上的二维坐标系。它使用笛卡尔坐标系，其中位置由 X 和 Y 坐标表示。投影坐标系通常用于地图制作和空间分析。常见的投影坐标系包括通用横轴墨卡托投影（Universal Transverse Mercator，UTM）。

（1）理解空间参考类

我们通过空间参考类（SpatialReference）来指定和引用空间参考。一般在创建空白要素类的时候以及投影转换的时候使用。

此类具有多个属性，包括坐标系参数。但是，若要使用这些属性，必须实例化类 （instantiated），需要为此类创建一个对象。类就像一个此对象的蓝图，你可以通过实例化类在此蓝图的基础上创建一个对象。例如我们查看shpfile系列文件中.prj文件来看看其空间参考属性：

import arcpy, os
prjfile = os.path.join(data1_dir, "streets.prj")
spatialref = arcpy.SpatialReference(prjfile)
spatialref

spatialref对象

上图可以spatialref的全部属性。我们也可以单独获取其属性，比如

spatialref.name # 获取空间参考文件的名称
# >>> 'NAD_1983_HARN_Adj_MN_Clay_Feet'

以下演示定义一个地理坐标系空间参考对象：

sr1 = arcpy.SpatialReference("GCS_WGS_1984")
sr2 = arcpy.SpatialReference(4326) # ！！！数字类型不是字符串

# 判断两个参考系是否相等
sr1 == sr2 
>>> True # 证明相等

可以同时对空间参考对象定义地理坐标系和投影坐标系。

（2）投影的概念

投影是一种将地球表面上的三维地理坐标（经度、纬度和高程）映射到二维平面上的方法。由于地球是一个三维椭球体，将其映射到平面上会引入形状、距离和方向的变形。因此，选择适当的投影方法对于特定的地理区域和任务非常重要。ArcGIS 提供了各种投影方法，包括等距圆柱投影、等距圆锥投影、等面积投影、等角投影等。每种投影方法都有其特定的优势和适用范围。ArcGIS中使用投影和投影栅格工具进行投影变换，对应的Arcpy方法是arcpy.management.Project和arcpy.management.ProjectRaster，如果还未定义投影需要用定义投影工具。

img

（3）哪些情况需要采用投影坐标系

选择投影坐标系解决地球表面的曲面到平面的映射问题，避免地球曲面产生的误差。以下情况需要使用投影坐标系：

1. 地图制作：当需要制作地图时，通常需要将地球表面的曲面映射到平面上。由于地球是一个三维椭球体，直接在平面上表示地球上的地理坐标会引入形状、距离和方向的变形。通过采用适当的投影坐标系，可以将地理坐标转换为平面坐标，以在地图上准确地表示地理特征、距离和方向。
2. 空间分析：在进行空间分析时，需要进行地理数据的测量、叠加和分析。在地理坐标系下，直接进行距离、面积和方向的计算可能不准确，因为地球的曲面会引入误差。通过将数据转换到适当的投影坐标系，可以进行准确的空间分析，确保测量和计算的精度。
3. 数据叠加：当需要将来自不同数据源的地理数据进行叠加时，这些数据可能使用不同的地理坐标系。为了进行准确的叠加，您需要将数据转换到相同的投影坐标系，以确保它们在平面上的位置和几何关系正确匹配。
4. 可视化和展示：在将地理数据可视化和展示时，使用投影坐标系可以确保地图的形状和比例符合实际。通过选择适当的投影坐标系，可以在地图上准确地显示地理特征和空间分布，使观众能够更好地理解和解读地理信息。

（4）WKID代号的查询

WKID 空间参考代号查询方式：

通过查询空间参考代号，可以找到对应的空间参考文件，从而确定空间参考。空间参考代号可以通过以下两个文件查询：

1. geographic_coordinate_systems.pdf
2. projected_coordinate_systems.pdf

比如我们要查询WGS1984的WKID和标准命名，打开文件1搜索：

image-20230812074808315

也可以在arcgis pro中点开任何一个图层或者使用打开投影工具，进行可视化操作：

image-20230824132330857

image-20230824132354731

（5）投影相关函数的使用

简单来说，对于矢量数据采用投影arcpy.management.Project，对于栅格数据采用投影栅格arcpy.management.ProjectRaster，如果没有数据空间参考采用定义投影arcpy.management.DefineProjection。他们都可以传入空间参考类的实例化对象作为参数传入，拿定义投影举例：

import arcpy
infc = r"C:\data\demo.shp"
sr = arcpy.SpatialReference(4326) # 建议输入WKID
arcpy.DefineProjection_management(infc, sr)

6.示例

我们看个简单的Python例子，如果要将文件夹中所有的shapefile复制到地理数据库，我们应该怎么办：

# 1导入包
import arcpy, os

# 2定义相关参数
gdb = "demo.gdb"
workspace = r"Z:\Sync\Urban-Spatial-Data-Analysis-For-Beginners\4-空间数据分析\4.2-arcpy\resource\data1"
arcpy.env.workspace = workspace

# 3创建地理数据库
arcpy.CreateFileGDB_management(workspace, gdb) 

# 列出文件夹中的要素类 shp文件
fc_list = arcpy.ListFeatureClasses()

for fc in fc_list:
    fc_desc = arcpy.Describe(fc)
    if fc_desc.shapeType == "Polyline":
        newfc = os.path.join(gdb, "Polyline",fc_desc.basename)
        
prj = arcpy.SpatialReference(os.path.join(workspace, "streets.prj")) # 读取shp文件的投影信息
arcpy.CreateFeatureDataset_management(gdb, "Polyline", prj) # 在数据库中创建名叫Polyline的空白要素类 指定其空间参考为 "streets.prj"的空间参考

arcpy.CopyFeatures_management(fc, newfc)

你可以将以下代码按照我的分块形式复制到jupyter notebook中执行，以便更清楚的了解每行代码发生了什么：

在第3步代码运行之后，你会发现data1文件夹下多了一个空的gdb数据库：

第4步我们想把data1文件夹里所有(其实只有一个)多段线的要素导入到此数据库，首先列出当前工作空间的要素类：

fc_list = arcpy.ListFeatureClasses()
fc_list
"""
>>> (">>>"我用来表示notebook单元格的输出，下文不再另做说明)
>>> ['parcels.shp', 'streets.shp']
"""

然后把这两个要素类通过数据arcpy.Describe返回的对象中的数据类型shapeType进行判定，如果是多段线则构建输出文件名。

for fc in fc_list:
    fc_desc = arcpy.Describe(fc)
    if fc_desc.shapeType == "Polyline":
        newfc = os.path.join(gdb, "Polyline",fc_desc.basename)
newfc
"""
>>> 'demo.gdb\\Polyline\\streets'
"""

我们复制此要素到数据库：

prj = arcpy.SpatialReference(os.path.join(workspace, "streets.prj")) # 读取shp文件的投影信息
arcpy.CreateFeatureDataset_management(gdb, "Polyline", prj) # 在数据库中创建名叫Polyline的空白要素类

arcpy.CopyFeatures_management(fc, newfc)

创建成功后此要素会被自动加载到arcgis的地图中：

streets文件创建成功

以上示例通过遍历列表的方式可以将重复性的导入工作自动化，省时省力~