Split A polygon into multiple where two or more polylines intersect each other
Hi All,
I have thousands of polygons, each polygon contains polylines inside.
as per requirement from customer, need make separate polygon where two polylines intersects each other, when you extend polylines programmatically using .net .
please suggest me how to complete this task.
Accepted Solutions
Thank you Xander,
i have checked the script and executed on original data, result output is very good.
Thanks once again for your swift response 
Hi Xander!
Thanks for your query and below are my response for your query .
Question 1: how are you going to define which of the lines to extend first?
Response: Yes it is difficult to say which one has to extend first,
- so first both polylines has to extend upto intersection of polygon
- once polygons split into multiple user has to check visually which one has to ignore
- what i have shown earlier that is final output
- please look into below image
Question 2: what the data represents ?
actually the data belongs to agriculture (polygons) and fences (polylines (property boundaries))
Question 3: could do this type of things with python too and that could result in a lot less code!
yes, i agree, but my all earlier tools for this project has developed on.net,
if python works on the same feature class instead of creating another output
,it's good and requesting you to provide solution and sample code for python too.
Another thing to keep in mind is how far you will extend the line:
Should each line divide the polygon only in two or may there be more resulting parts?
As a side note; for these type of scenarios I would always write to a new featureclass. Is there a special reason why you want (or need) to write to the same featureclass?
OK, so let's visually show what I did (with python). Below the test data I used for writing the code. There are two lines per polygon (can be more or less) and some are partly outside the polygons, and some are completely inside the polygons:
To match the lines with the polygons, and cut off the parts outside the polygons, I did an intersect (manually) to create this:
Not only do I now have the parts of the lines inside the polygons, I also (more importantly) have the attributes of the polygons for each line (an identifier that I can use to match) the lines to the polygons, without the need to spatially select the lines for each polygon.
The results is this:
Each polygon is divided by the lines. In an additional step the multi part polygons are converted to single parts. All attributes are preserved as well, but I didn't spend any time in reconstructing any relationshipclasses in the code.
Below the code use:
import arcpy
def main():
import os
arcpy.env.overwriteOutput = True
fc_pol = r'C:\GeoNet\SplitPolygons\gdb\data.gdb\Polygons'
fc_lines = r'C:\GeoNet\SplitPolygons\gdb\data.gdb\Polylines2'
fc_out = r'C:\GeoNet\SplitPolygons\gdb\data.gdb\Polygons_out3'
fld_id = 'Name'
# distance to extent each lines in each direction
extend_dist = 50000
# create ouput polygon fc
print "create ouput polygon fc"
sr = arcpy.Describe(fc_pol).spatialReference
ws, fc_name = os.path.split(fc_out)
arcpy.CreateFeatureclass_management(ws, fc_name, "POLYGON", fc_pol, None, None, sr)
# create list of lines per id
print "create list of lines per id"
dct_lines = {}
flds = ('SHAPE@', fld_id)
with arcpy.da.SearchCursor(fc_lines, flds) as curs:
for row in curs:
line = row[0]
name = row[1]
if name in dct_lines:
lines = dct_lines[name]
lines.append(line)
dct_lines[name] = lines
else:
dct_lines[name] = [line]
# create list of fields
print "create list of fields"
flds = getValidFields(fc_pol)
index_id = flds.index(fld_id)
# start insert cursor
print "start insert cursor"
with arcpy.da.InsertCursor(fc_out, flds) as curs_out:
with arcpy.da.SearchCursor(fc_pol, flds) as curs:
for row in curs:
polygon = row[0]
name = row[index_id]
print "Processing polygon", name
if name in dct_lines:
# polygon has lines
lst_lines = dct_lines[name]
print " - Lines found:", len(lst_lines)
lst_polygons = [polygon]
for polyline in lst_lines:
polyline2 = ExtendLine(polyline, extend_dist)
lst_polygons = getSplitPolygons(lst_polygons, polyline2)
lst_polygons = createSingleParts(lst_polygons)
print " - Polygons created:", len(lst_polygons)
for polygon in lst_polygons:
lst_row = list(row)
lst_row[0] = polygon
row_out = tuple(lst_row)
curs_out.insertRow(row_out)
else:
# polygon has no lines, write the original polygon
curs_out.insertRow(row)
def getSplitPolygons(lst_polygons, polyline):
lst_new = []
lst_i = range(len(lst_polygons))
lst_i.reverse()
for i in lst_i:
polygon = lst_polygons.pop(i)
try:
lst_pols = polygon.cut(polyline)
except:
lst_pols = [polygon]
lst_new.extend(lst_pols)
return lst_new
def createSingleParts(lst_polygons):
lst_new = []
for polygon in lst_polygons:
if polygon.partCount > 1:
for i in range(polygon.partCount):
part = polygon.getPart(i)
polygon_part = arcpy.Polygon(arcpy.Array(part), polygon.spatialReference)
lst_new.append(polygon_part)
else:
lst_new.append(polygon)
return lst_new
def getValidFields(fc):
lst = []
for fld in arcpy.ListFields(fc):
if not fld.type in ['OID', 'Geometry']:
if not fld.name in ['SHAPE_Length', 'SHAPE_Area']:
lst.append(fld.name)
lst.insert(0, 'SHAPE@')
return lst
def ExtendLine(polyline, extend_dist):
sr = polyline.spatialReference
start_pnt = polyline.firstPoint
end_pnt = polyline.lastPoint
pntg1 = arcpy.PointGeometry(start_pnt, sr)
pntg2 = arcpy.PointGeometry(end_pnt, sr)
angle = getAngle(pntg1, pntg2, sr)
angle180 = angle + 180
pntg1a = pntg1.pointFromAngleAndDistance(angle180, extend_dist, 'PLANAR')
pntg2a = pntg2.pointFromAngleAndDistance(angle, extend_dist, 'PLANAR')
return arcpy.Polyline(arcpy.Array([pntg1a.firstPoint, pntg2a.firstPoint]), sr)
def getAngle(pntg1, pntg2, sr):
return pntg1.angleAndDistanceTo(pntg2, method='PLANAR')[0]
if __name__ == '__main__':
main()Some explanation:
- line 18, the output polygon featureclass is created using the schema of the input polygon featureclass (same fields)
- lines 22 - 33, create a dictionary were each key is the id (name) of the polygon, and the value is a list of lines for that polygon name
- For each name of the polygon (line 46) on line 49 it is checked if the polygon has any lines (using the dictionary)
- For each polyline (line 54), the line is extended (line 55), the divided polygons are returned (line 56) and converted to singe parts (line 57)
- Each created polygon (lines 60 - 64) are written to the output featureclass with the attributes of the original polygon
- starting on line 70 the supporting functions are defined.
Kind regards, Xander
Thank you Xander,
i have checked the script and executed on original data, result output is very good.
Thanks once again for your swift response
