Creating polygon connecting endpoints of multiple lines using ArcPy?

Question

I'm trying to figure out how to create a polygon that connects all the endpoints of a shapefile containing a set of polyilnes with pythonscript in ArcGIS, I'm having trouble doing this as the order of the nodes in the polygon is important. I want to achieve the grey polygon in the picture from the green lines

I want to connect the endpoints of the green lines to create the grey polygon without having to do it manually

• Ian Turton: do your lines have some attribute to give the ordering?
• Matej: first, you need the ordering defined as @iant asked, then you need the rule whether connecting endpoint to the next startpoint or doing it any other way
• Ian Turton: failing that maybe some sort of alpha hull on the end points?
• Amanda: The line does to some degree have attributes to give them order. They have a ID number, but for the example above the rightbranch have ID 1-7, the left 15- 21 and after they're connected the IDs are 22-27
• FelixIP: You can get very close by a) creating TIN, using lines, b) converting TIN to triangles c) selecting triangles that share a boundary with lines. You'll have only 1 polygon to delete at the top

Answer 1

STEPS:

Compute sections centre points: enter image description here

Built their Euclidean minimum spanning tree, dissolve it and compute buffer, distance equal half of shortest section length: enter image description here

Create section end points and compute their chainage (distance along line) on the boundary of the buffer (closed polyline version of buffer): enter image description here

Sort end points in ascending order using chainage field. Points below labelled by their FID:

enter image description here

Create polygon from ordered set of points: enter image description here

Script:

import arcpy, traceback, os, sys,time
from heapq import *
from math import sqrt
import itertools as itt
from collections import defaultdict
try:
 def showPyMessage():
 arcpy.AddMessage(str(time.ctime()) + " - " + message)
 # MST by PRIM's
 def prim( nodes, edges ):
 conn = defaultdict( list )
 for n1,n2,c in edges:
 conn[ n1 ].append( (c, n1, n2) )
 conn[ n2 ].append( (c, n2, n1) )
 mst = []
 used = set( nodes[ 0 ] )
 usable_edges = conn[ nodes[0] ][:]
 heapify( usable_edges )
 while usable_edges:
 cost, n1, n2 = heappop( usable_edges )
 if n2 not in used:
 used.add( n2 )
 mst.append( ( n1, n2, cost ) )
 for e in conn[ n2 ]:
 if e[ 2 ] not in used:
 heappush( usable_edges, e )
 return mst 
 mxd = arcpy.mapping.MapDocument("CURRENT")
 SECTIONS=arcpy.mapping.ListLayers(mxd,"SECTION")[0]
 PGONS=arcpy.mapping.ListLayers(mxd,"RESULT")[0]
 d=arcpy.Describe(SECTIONS)
 SR=d.spatialReference
 cPoints,endPoints,lMin=[],[],1000000
 with arcpy.da.SearchCursor(SECTIONS, "Shape@") as cursor:
 # create centre and end points
 for row in cursor:
 feat=row[0]
 l=feat.length
 lMin=min(lMin,feat.length)
 theP=feat.positionAlongLine (l/2).firstPoint
 cPoints.append(theP)
 theP=feat.firstPoint
 endPoints.append(theP)
 theP=feat.lastPoint
 endPoints.append(theP)
 arcpy.AddMessage('Computing minimum spanning tree')
 m=len(cPoints)
 nodes=[str(i) for i in range(m)]
 p=list(itt.combinations(range(m), 2))
 edges=[]
 for f,t in p:
 p1=cPoints[f]
 p2=cPoints[t]
 dX=p2.X-p1.X;dY=p2.Y-p1.Y
 lenV=sqrt(dX*dX+dY*dY)
 edges.append((str(f),str(t),lenV))
 MST=prim(nodes,edges)
 mLine=[]
 for edge in MST:
 p1=cPoints[int(edge[0])]
 p2=cPoints[int(edge[1])]
 mLine.append([p1,p2])
 pLine=arcpy.Polyline(arcpy.Array(mLine),SR)
 # create buffer and compute chainage
 buf=pLine.buffer(lMin/2)
 outLine=buf.boundary()
 chainage=[]
 for p in endPoints:
 measure=outLine.measureOnLine(p)
 chainage.append([measure,p])
 chainage.sort(key=lambda x: x[0])
 # built polygon
 pGon=arcpy.Array()
 for pair in chainage:
 pGon.add(pair[1])
 pGon=arcpy.Polygon(pGon,SR)
 curT = arcpy.da.InsertCursor(PGONS,"SHAPE@")
 curT.insertRow((pGon,))
 del curT
except:
 message = "\n*** PYTHON ERRORS *** "; showPyMessage()
 message = "Python Traceback Info: " + traceback.format_tb(sys.exc_info()[2])[0]; showPyMessage()
 message = "Python Error Info: " + str(sys.exc_type)+ ": " + str(sys.exc_value) + "\n"; showPyMessage()

I know it is a bicycle, but it’s my own and I like it

Answer 2

I post this solution for QGIS here because it is free software and easy to implement. I considered only the right "branch" of polyline vector layer; as it can be observed at the next image (12 features at attributes table):

enter image description here

The code (algorithm in a one line python list comprehension), for running at the Python Console of QGIS, is:

layer = iface.activeLayer()
features = layer.getFeatures()
features = [feature for feature in features]
n = len(features)
geom = [feature.geometry().asPolyline() for feature in features ]
#multi lines as closed shapes
multi_lines = [[geom[i][0], geom[i][1], geom[i+1][1], geom[i+1][0], geom[i][0]]
 for i in range(n-1)]
#multi polygons
mult_pol = [[] for i in range(n-1)]
for i in range(n-1):
 mult_pol[i].append(multi_lines[i])
#creating a memory layer for multi polygon
crs = layer.crs()
epsg = crs.postgisSrid()
uri = "Polygon?crs=epsg:" + str(epsg) + "&field=id:integer""&index=yes"
mem_layer = QgsVectorLayer(uri,
 "polygon",
 "memory")
QgsMapLayerRegistry.instance().addMapLayer(mem_layer)
mem_layer.startEditing()
#Set features
feature = [QgsFeature() for i in range(n-1)]
for i in range(n-1):
 #set geometry
 feature[i].setGeometry(QgsGeometry.fromPolygon(mult_pol[i]))
 #set attributes values
 feature[i].setAttributes([i])
 mem_layer.addFeature(feature[i], True)
#stop editing and save changes
mem_layer.commitChanges()

After running the code:

enter image description here

it was produced a polygon memory layer (with 11 features at its attributes table). It works nicely.

Answer 3

You could select the endpoints that will participate in a polygon, create a TIN from only those points. Convert the TIN to polygons, dissolve the polygons. The trick to automating this process is deciding which points to contribute to each polygon. If you have lines with valid directions, and those lines all share some common attribute you could write a query to export say the end vertices using line vertices to points, then select by attribute those points that have the common attribute value.
Better would be to extract/select the points, read the x , y values using a cursor, use the x, y values to write a new polygon. I cannot see an attached picture in your post but if point order matters then once you have the x, y values stored in a Python list, sort them. http://resources.arcgis.com/EN/HELP/MAIN/10.1/index.html#//002z0000001v000000

Answer 4

Expanding on @iant comment, the closest geometry to your snapshot is the alpha shape (alpha hull) of the endpoints. Fortunately many well received threads have been already answered on GIS SE. For example:

• Create Complex Polygon from Point Layer using only Boundary Points in ArcGIS

• What are Definition, Algorithms and Practical Solutions for Concave Hull?

To solve your problem, first use Feature To Point to extract the end points. Then use the python tool from this Link to calculate the concave hull.

• PolyGeo: Your first link seems to be broken.