Hi everyone,
I've some set of points which I'm joining to make a line. Refer the Picture below.
The problem I'm facing is there are no points(elevation points) on the polygon (River Bank). Is there any way of extending the line to the polygon?
At present m using Extend line to layer, but it taking too much of time and inconsistent in the lines. Refer the Picture below.
Accepted Solutions
I just downloaded the sample data, and ran my code from the Draw a Line to the nearest Polygon thread, and it ran without any errors and did exactly what you wanted. I did have to increase the dist value since you are using State Plane, and the gaps between the lines and river banks are greater than 100 feet.
Here is the code that worked for me, using your 2 shape files (load the shape files into ArcMap first and run from interactive Python window):
>>> import math
>>> dist = 500
>>>
>>> fc_line = "Pointoline"
>>> fc_poly = "RB_23"
>>> with arcpy.da.SearchCursor(fc_poly, "SHAPE@") as scur:
... for poly, in scur:
... boundary = poly.boundary()
... arcpy.SelectLayerByLocation_management(fc_line, "WITHIN", poly)
... with arcpy.da.UpdateCursor(fc_line, "SHAPE@") as ucur:
... for line, in ucur:
... arr, = line.getPart()
... SR = line.spatialReference
...
... p1, p2 = arr[:2]
... angle = math.atan2(p2.Y - p1.Y, p2.X - p1.X)
... p = arcpy.Point(p1.X - dist * math.cos(angle),
... p1.Y - dist * math.sin(angle))
... arr.insert(0, p)
...
... pn1, pn = arr[-2:]
... angle = math.atan2(pn.Y - pn1.Y, pn.X - pn1.X)
... p = arcpy.Point(pn.X + dist * math.cos(angle),
... pn.Y + dist * math.sin(angle))
... arr.append(p)
...
... line = arcpy.Polyline(arr, SR)
... line = line.cut(boundary)[1]
... ucur.updateRow([line])
...
>>>
I just downloaded the sample data, and ran my code from the Draw a Line to the nearest Polygon thread, and it ran without any errors and did exactly what you wanted. I did have to increase the dist value since you are using State Plane, and the gaps between the lines and river banks are greater than 100 feet.
Here is the code that worked for me, using your 2 shape files (load the shape files into ArcMap first and run from interactive Python window):
>>> import math
>>> dist = 500
>>>
>>> fc_line = "Pointoline"
>>> fc_poly = "RB_23"
>>> with arcpy.da.SearchCursor(fc_poly, "SHAPE@") as scur:
... for poly, in scur:
... boundary = poly.boundary()
... arcpy.SelectLayerByLocation_management(fc_line, "WITHIN", poly)
... with arcpy.da.UpdateCursor(fc_line, "SHAPE@") as ucur:
... for line, in ucur:
... arr, = line.getPart()
... SR = line.spatialReference
...
... p1, p2 = arr[:2]
... angle = math.atan2(p2.Y - p1.Y, p2.X - p1.X)
... p = arcpy.Point(p1.X - dist * math.cos(angle),
... p1.Y - dist * math.sin(angle))
... arr.insert(0, p)
...
... pn1, pn = arr[-2:]
... angle = math.atan2(pn.Y - pn1.Y, pn.X - pn1.X)
... p = arcpy.Point(pn.X + dist * math.cos(angle),
... pn.Y + dist * math.sin(angle))
... arr.append(p)
...
... line = arcpy.Polyline(arr, SR)
... line = line.cut(boundary)[1]
... ucur.updateRow([line])
...
>>>
The option provided by Joshua Bixby is a very good one. It takes the direction of the first and last part of the line and extends it to the polygon boundary. As an alternative you can simply snap the first and last point to polygon boundary (see code below).
#-------------------------------------------------------------------------------
# Name: snapline2polygon.py
# Purpose:
#
# Author: xbakker
#
# Created: 17/04/2017
#-------------------------------------------------------------------------------
def main():
print "Start process..."
import arcpy
fc_line = r'C:\GeoNet\SnapLine2Polygon\Pointoline.shp'
fc_pol = r'C:\GeoNet\SnapLine2Polygon\RB_23.shp'
fc_out = r'C:\GeoNet\SnapLine2Polygon\snapped_lines.shp'
# fech polygon and get perimeter
print "Get polygon and boundary"
polygon = arcpy.da.SearchCursor(fc_pol, ('SHAPE@')).next()[0]
boundary = polygon.boundary()
print "start loop..."
lst_feats = []
cnt = 0
with arcpy.da.SearchCursor(fc_line, ('SHAPE@')) as curs:
for row in curs:
cnt += 1
if cnt % 10 == 0:
print " - processing line", cnt
polyline = row[0]
pnt1 = polyline.firstPoint
pnt2 = polyline.lastPoint
pntg1_snap = boundary.snapToLine(pnt1)
pntg2_snap = boundary.snapToLine(pnt2)
polyline_out = ExtendPolyline(polyline, pntg1_snap.firstPoint, pntg2_snap.firstPoint)
lst_feats.append(polyline_out)
# store results
print "store results..."
arcpy.CopyFeatures_management(lst_feats, fc_out)
def ExtendPolyline(polyline, pnt1, pnt2):
# will snap first and last line to nearest part of given polyline
# does not account for small segments!
sr = polyline.spatialReference
lst_pnts = []
for part in polyline:
for pnt in part:
lst_pnts.append(pnt)
lst_pnts.insert(0, pnt1)
lst_pnts.append(pnt2)
return arcpy.Polyline(arcpy.Array(lst_pnts), sr)
if __name__ == '__main__':
main()
This will produce errors when you have small segments and the first of last point is closer to the wrong side of the polygon boundary:
To avoid this you would have to edit the line data and create the correct lines to be extended or use the logic suggested by Joshua.
