# -*- coding: utf-8 -*-
"""
ArcGIS自相交修复工具
功能:
1. 查找自相交点
2. 创建0.01矩形面擦除原始数据,解决自相交问题
"""
import arcpy
import os
import uuid
arcpy.env.overwriteOutput = True
def find_self_intersections(input_fc):
"""
使用检查几何工具查找自相交问题
"""
scratch_gdb = arcpy.env.scratchGDB
check_table = os.path.join(scratch_gdb, "check_geom_" + uuid.uuid4().hex[:8])
# 检查几何
arcpy.management.CheckGeometry(
input_fc,
check_table,
"SELF_INTERSECTION"
)
# 统计自相交数量
count = 0
self_intersect_oids = []
with arcpy.da.SearchCursor(check_table, ["FEATURE_OID", "PROBLEM"]) as cursor:
for oid, problem in cursor:
if problem == "Self intersection":
count += 1
self_intersect_oids.append(oid)
# 清理临时表
arcpy.management.Delete(check_table)
return count, self_intersect_oids
def create_small_rectangle(center_x, center_y, size, sr):
"""
在指定位置创建一个小的矩形面
size: 矩形的边长(单位与坐标系一致)
"""
half_size = size / 2
# 创建矩形四个角点
points = [
arcpy.Point(center_x - half_size, center_y - half_size),
arcpy.Point(center_x + half_size, center_y - half_size),
arcpy.Point(center_x
+ half_size, center_y + half_size),
arcpy.Point(center_x - half_size, center_y + half_size),
arcpy.Point(center_x - half_size, center_y - half_size) # 闭合
]
array = arcpy.Array(points)
polygon = arcpy.Polygon(array, sr)
return polygon
def fix_self_intersection_with_erase(input_fc, output_fc, rect_size=0.01):
"""
通过创建小矩形面擦除来修复自相交
参数:
input_fc: 输入面要素类
output_fc: 输出要素类
rect_size: 矩形边长(默认0.01)
"""
desc = arcpy.Describe(input_fc)
sr = desc.spatialReference
scratch_gdb = arcpy.env.scratchGDB
# 第一步:检查自相交
print("正在检查自相交问题...")
# 使用更直接的方法:检查每个要素的自相交
# 对于面要素,自相交通常发生在环相交的地方
# 创建临时要素类用于存储小矩形
rect_fc = os.path.join(scratch_gdb, "temp_rect_" + uuid.uuid4().hex[:8])
arcpy.management.CreateFeatureclass(
scratch_gdb,
os.path.basename(rect_fc),
"POLYGON",
spatial_reference=sr
)
# 查找所有自相交点
self_intersect_points = []
with arcpy.da.SearchCursor(input_fc, ["OID@", "SHAPE@"]) as cursor:
for oid, shape in cursor:
if shape is None:
continue
# 获取面的所有环
for part in shape:
points = [
p for p in part if p is not None]
# 检查线段相交
for i in range(len(points) - 1):
for j in range(i + 2, len(points) - 1):
# 避免相邻线段(它们共享一个端点)
if j == i + 1:
continue
p1, p2 = points[i], points[i + 1]
p3, p4 = points[j], points[j + 1]
# 计算两条线段的交点
intersect_point = line_intersection(p1, p2, p3, p4)
if intersect_point:
self_intersect_points.append((intersect_point.X, intersect_point.Y, oid))
print(f"发现 {len(self_intersect_points)} 个自相交点")
# 创建小矩形面用于擦除
with arcpy.da.InsertCursor(rect_fc, ["SHAPE@"]) as icur:
for x, y, oid in self_intersect_points:
rect = create_small_rectangle(x, y, rect_size, sr)
icur.insertRow([rect])
print(f"创建了 {len(self_intersect_points)} 个擦除矩形")
# 使用擦除工具
if len(self_intersect_points) > 0:
print("正在执行擦除操作...")
arcpy.analysis.Erase(
input_fc,
rect_fc,
output_fc
)
else:
# 没有自相交,直接复制
print("没有发现自相交问题,直接复制...")
arcpy.management.CopyFeatures(input_fc, output_fc)
# 清理临时数据
try:
arcpy.management.Delete(rect_fc)
except:
pass
print(f"\n处理完成!")
print(f"原始数据: {input_fc}")
print(f"输出数据: {output_fc}")
print(f"修复了 {len(self_intersect_points)} 个自相交点")
def line_intersection(p1, p2, p3, p4):
"""
计算两条线段的交点
返回交点坐标,如果不相交则返回None
"""
x1, y1 = p1.X, p1.Y
x2, y2 = p2.X, p2.Y
x3, y3 = p3.X, p3.Y
x4, y4 = p4.X, p4.Y
# 计算方向向量
dx1 = x2 - x1
dy1 = y2 - y1
dx2 = x4 - x3
dy2 = y4 - y3
# 计算叉积
cross = dx1 * dy2 - dy1 * dx2
if abs(cross) 1e-10:
# 平行或共线
return None
# 计算参数
t = ((x3 - x1) * dy2 - (y3 - y1) * dx2) / cross
u = ((x3 - x1) * dy1 - (y3 - y1) * dx1) / cross
# 检查交点是否在两条线段内
if 0 <= t <= 1 and 0 <= u <= 1:
# 计算交点坐标
intersect_x = x1 + t * dx1
intersect_y = y1 + t * dy1
return arcpy.Point(intersect_x, intersect_y)
return None
def main():
input_fc = r"C:\Users\yl\Documents\ArcGIS\Projects\MyProject17\MyProject17.gdb\m"
output_fc = r"C:\Users\yl\Documents\ArcGIS\Projects\MyProject17\MyProject17.gdb\密码"
rect_size = 0.01
fix_self_intersection_with_erase(input_fc, output_fc, rect_size)
if __name__ == "__main__":
main()
