一. 使用 Python(ArcPy) 绘制shp1.什么是ArcPyArcPy 是一个安装 ArcGIS 会附带的站点包，通过 Python 实现。简言之，能通过 Python 直接调用 arcpy 执行地理数据分析、数据转换、数据管理和地图自动化以及多种客制化的需求。 ArcGIS 帮助链接：https://desktop./zh-cn/arcmap/10.3/analyze/arcpy/what-is-arcpy-.htm 需要注意的是 ArcPy 是通过 Python2 实现的，不兼容 Python3！ 需要注意的是 ArcPy 是通过 Python2 实现的，不兼容 Python3！ 需要注意的是 ArcPy 是通过 Python2 实现的，不兼容 Python3！ 除非使用 ArcGIS pro，要花钱哦！ 2.如何构造shp(面)主要涉及类Array 类 ArcPy 中自带的数组对象，可包含点和数组，用于构造几何对象。 官方文档：https://desktop./zh-cn/arcmap/10.3/analyze/arcpy-classes/array.htm Polygon 类 可以使用 Polygon 类从头开始创建几何对象。接收 Array 类中数据然后使用 Polygon 类完成要素类的创建。 官方文档：https://desktop./zh-cn/arcmap/10.3/analyze/arcpy-classes/polygon.htm 构造面点连成线，线组成了一个面。只需在 arcpy.Polygon() 类中传入一个点集就能形成一个面，点集（数组）通过 arcpy.Array() 构造传入（最后一个点不必与首个点相同，就算没有也会自动闭合）。 [[20.0, 20.0],[30.0, 20.0],[30.0, 10.0],[20.0, 10.0]] 多个面就是多个点集。两个点集就是两个面、两个图形。 [[[5.0, 3.0],[3.0, 3.0], [3.0, 5.0],[5.0, 3.0]], [[7.0, 5.0],[5.0, 5.0], [5.0, 7.0],[7.0, 5.0]] ] 那么如果是环岛、孔洞等图形呢？ 不过是多个点集，当某个点集在另一个点集范围内时，arcpy.Polygon() 类可以自动把重叠部分完全扣除掉，实现环岛、孔洞等。所以使用该方法绘制图形是非常方便的。 [[[40.0, 40.0], [50.0, 40.0], [50.0, 30.0], [40.0, 30.0]], [[45.0, 35.0], [48.0, 35.0], [48.0, 36.0], [45.0, 36.0]] ] 以下是 使用 arcpy.Array() 和 arcpy.Polygon() 实现简单的点集转 shp 代码示例： # -*- coding:utf-8 -*-# --------------------------------------------------------------------------- # Author: Lcc hygnic # Created on: 2020/10/16 15:55 # Reference: """ Description:Python2.7 Usage: """ # --------------------------------------------------------------------------- import arcpy import os def draw_poly(coord_list, sr, y, x): """ 创建多边形 coord_list(List)：多个点组成的坐标 sr: 投影系 y(Int): y坐标列 x(Int): x坐标列 """ parts = arcpy.Array() yuans = arcpy.Array() yuan = arcpy.Array() for part in coord_list: for pnt in part: if pnt: yuan.add(arcpy.Point(pnt[y], pnt[x])) else: # null point - we are at the start of a new ring yuans.add(yuan) yuan.removeAll() # we have our last ring, add it yuans.add(yuan) yuan.removeAll() # if we only have one ring: remove nesting if len(yuans) == 1: yuans = yuans.getObject(0) parts.add(yuans) yuans.removeAll() # 只有一个，单个图形 if len(parts) == 1: parts = parts.getObject(0) return arcpy.Polygon(parts, sr) if __name__ == '__main__': # 点集，将之前的示例的三个图形放到一起 poly = [ [[20.0, 20.0], [30.0, 20.0], [30.0, 10.0], [20.0, 10.0]], [[5.0, 3.0], [3.0, 3.0], [3.0, 5.0], [5.0, 3.0]], [[7.0, 5.0], [5.0, 5.0], [5.0, 7.0], [7.0, 5.0]], [[40.0, 40.0], [50.0, 40.0], [50.0, 30.0], [40.0, 30.0]], [[45.0, 35.0], [48.0, 35.0], [48.0, 36.0], [45.0, 36.0]] ] arcpy.env.overwriteOutput = True output_folder = os.getcwd() # 创建空白 shp name = "NewShp" # ▶注释1◀ blank_shp = arcpy.CreateFeatureclass_management( output_folder, name, "Polygon", spatial_reference=None) # 创建、写入面 Rows = arcpy.da.InsertCursor(blank_shp, "SHAPE@") f = poly p = draw_poly(f, sr=None, y=0, x=1) Rows.insertRow([p]) write_sth = "Output path: " + os.path.join(output_folder, name) print write_sth del Rows 运行该程序会创建一个名为 NewShp.shp 的文件。 ▶注释1◀： arcpy.CreateFeatureclass_management() 创建一个空要素类，该方法可以传入参考系参数 spatial_reference。在这里没有设置参考系，如果之后需要的话，需要使用定义投影工具定义投影系或者提前就设置好投影系。 二. 土地报备坐标txt文件（坐标交换数据）数据主要依据标准《勘测定界界址点坐标交换格式》📚我们拿到的 TXT 文件中的数据基本按照以下标准： 附录一：《勘测定界界址点坐标交换格式》 ...坐标交换格式具有两种格式，分别如下： 文本格式 [属性描述] 格式版本号= 数据产生单位= 数据产生日期= 坐标系= 几度分带= 投影类型= 计量单位= 带号= 精度= 转换参数=X平移,Y平移,Z平移,X旋转,Y旋转,Z旋转,尺度参数 [地块坐标] 界址点数,地块面积,地块编号,地块名称,记录图形属性(点、线、面),图幅号,地块用途,地类编码,@ {点号,地块圈号,X坐标,Y坐标 ... ... 点号,地块圈号,X坐标,Y坐标} 界址点数,地块面积,地块编号,地块名称,记录图形属性(点、线、面),图幅号,地块用途,地类编码,@ {点号,地块圈号,X坐标,Y坐标 ... ... 点号,地块圈号,X坐标,Y坐标} 注意： 所有的逗号分隔符都必须是英文输入法状态下的逗号；地块圈号不能小于零；数据产生日期的格式为：2000-12-12；坐标系为54北京坐标系或80国家大地坐标系；投影类型为高斯克吕格或等角多圆锥；几度分带为3或6；带号、精度、转换参数、界址点数、地块面积、地块圈号,X坐标,Y坐标必须为数字型；且不能用该（9999,000,000）方式表示；地块编号、地块名称、记录图形属性(点、线、面)、图幅号、地块用途、地类编码、点号的每项里不能含有“,” 、“@”符号。 点号中不能有字母，如J01等，建议为1。（预审报盘坐标格式要求） ... 例子：（该例子中数据经过处理，如果使用该数据生成shp,结果可能不符合预期） [属性描述] 格式版本号= 数据产生单位=国土资源部 数据产生日期=2020-9-10 坐标系=2000国家大地坐标系 几度分带=6 投影类型=高斯克吕格 计量单位=米 带号=35 精度=0.001 转换参数=0,0,0,0,0,0,0 [地块坐标] 8137,559.4952,1,火星乡地球村改造项目,面,,,,@ J1,1,3433398.731421,35572460.011417 J2,1,3433413.686436,35572361.307118 J3,1,3433422.659303,35572214.746216 J4,1,3433425.650143,35572095.104671 J5,1,3433431.631890,35571897.696121 J6,1,3433446.586915,35571804.9739 J7,1,3433437.613643,35571673.368316 J8,1,3433377.792817,35571493.906440 J1,1,3433398.731421,35572460.011417 18400,998.7312,1,火星乡地球村改造项目,面,,,,@ J1,1,3437490.7824,35577811.1468 J2,1,3437488.2225,35577838.46 J3,1,3437488.2225,35577864.9191 J4,1,3437503.5857,35577884.5502 J5,1,3437514.6817,35577899.9133 J6,1,3437514.6818,35577928.0804 J7,1,3437506.9998,35577940.0295 J8,1,3437506.5118,35577944.1285 J9,1,3437503.9028,35577960.8496 J1,1,3437490.7824,35577811.1468 数据结构说明： 三.转换过程中的坑那么把土地报备坐标 txt 文件转换成 shp，主要做的就是 把TXT文件中的点正确的分成一个个点集（数组），然后传入 arcpy.Polygon() 类就行了。 然而在实际操作中，由于各种原因，有很多坑，故可能出现下面这些情况： 这种情况出现的原因通常只有一个，就是没有正确的把一个个点集分开，导致本应该闭合的一个图形没有结束，然后最后一个点又连接到了下一个图形，出现这种线“乱飞乱连”的错误情况。 而为什么没有正确把点集分开，是因为对数据的理解不全面（或者说数据有几种样式，txt 文本自己本人从相关官网下载，不存在其他人私自修改的可能），上文所说的标准是我在网上找到的，感觉说的不甚清楚且不全。 实际在实践中，出现了以下几种样式（情况）： 情况一无地块描述，如 18400,998.7312,1,火星乡地球村改造项目,面,,,,@ 。 这种用 @ 符号结尾的段落，是每一个单独地块都有的描述符，本可以作为地块之间区分的标志物。然而实际上部分 txt 文件中直接就没有地块描述符 @ ，就如下图两个单独地块之间没有用地块描述符 @ 分开。 情况二序号混乱。 下图的这种情况可以理解为多部件（一个地块实际是由两个单独地块组成的），就像这样 如果第二个部件在第一个部件范围内，就会变成孔洞、环岛，就像这样 从 txt 文本可以看到部件从1变成2，序号是连续的，从J1019 到 J1020（中间的 J1 是部件1的起始点，可有可无，在这里不影响）。
然后奇怪的事情来了，或者说，“另一种情况吧”。 下图也是多部件，部件从1变成2，但是在这个 txt 文本中， 序号没有紧接上面从 J6077 到 J6078，而是变成了 J1，然后 J2，重新开始了。
情况三更多的奇怪情况。 从上面来看，由于序号规则不详、不完善或者说数据错误，非常容易混淆多部件和地块之间的编号。 更糟糕的是，我还遇到过以下几种情况，正如下面三个不同 txt 文本的截图所示: 从左到右，第一个是直接开头序号都不为1了，直接从 J9677 开始； 第二个是部件开头部不为1，为0； 第三个是部件计数跳跃，直接从1变成174，完全不知道为什么会这样！
四.解决思路 🤔
那么如何处理这些情况呢？咋一看上去非常复杂。 实际上非常简单，那就是不去纠结部件和地块序号，直接全部当成不同的地块处理（相当于有多部件的地块被打散了，这完全不影响出图，反正在一张shp文件上）。 到了这里，我已经逐渐理解一切，关键在于序号列和部件列，一旦序号列从1重新开始或者部件发生改变，数据就进行分割。 分好的数据（一个地块）是一个列表，然后一个 txt 文件中的所有数据组成一个大列表进行处理。 五.实现代码
完整实现代码（脚本形式）: ✌️好耶 points_genarator 方法将输入的 txt 文件正确分割成多个点集（数组）；然后使用 draw_poly 方法将该点集绘制成图形。 ▶注释1◀： 由于我们不再详细区分多部件等，所以不再需要描述符 @，况且有的 txt 文件直接都没有该描述符。所以直接将存在 @ 的行删除。 ▶注释2◀： 这个循环中的代码看上去比较复杂，其作用是判断每一行的部件或者序号的顺序是否发生改变，从来将不同点集直接分开。 ▶注释3◀： 文本中的第三行和第四行对应的不是 x y 坐标，而是 y x 坐标。所以 y=3, x=2。 ▶注释4◀： 输入 txt 文件所在地址和输出的 shp 文件的地址，然后运行即可。
正确分割数据得到的正确 shp 文件与错误 shp 文件的对比：
结束语
参考文章 关于构造面：https://gis./questions/14952/arcgis-10-0-arcpy-how-to-create-a-polygon-geometry-from-inner-and-outer-ring-po 使用版本 ArcGIS10.3 Python2.7 源代码及其PDF文档下载地址： 链接：https://pan.baidu.com/s/1ihBKcufLuRkAmXlJFFsDDA 提取码：3t9e 如果下载地址失效的话，可以回复 d1 获取最新下载地址。 很难的文章系列：
... 不难的有趣文章系列： ... 更多文章可以使用搜索哦 |
|
|
