(三十五)arcpy开发&计算polyline的起点和结束点
今天我们来学习一下,使用arcpy来寻找polyline要素类的起始端点和结束端点。首先,程序开始会对要素类polyline进行遍历,然后利用半正矢公式对两个端点进行计算。这里主要利用到半正矢公式,以及在arcpy中创建要素类字段,并添加相应的值。我们来看一下具体的实现代码。
import arcpy
import sys
import math
def haversine(point1, point2):
lat1, long1 = point1
lat2, long2 = point2
diffLat = math.radians(lat2 - lat1)
diffLong = math.radians(long2 - long1)
radius = 6378137
a = (math.sin(diffLat / 2) * math.sin(diffLat / 2) +
math.cos(math.radians(lat1)) * math.cos(math.radians(lat2)) *
math.sin(diffLong / 2) * math.sin(diffLong / 2))
c = 2 * math.atan2(math.sqrt(a), math.sqrt(1 - a))
d = c * radius
return d
def find_distant_pts(allpts, beginpts, endpts, second_check=False):
beg = beginpts[0]
end = endpts[0]
biggest_distance = 0
if len(allpts) > 2:
for pt in allpts:
if allpts.count(pt) == 1:
for pt2 in allpts:
if allpts.count(pt2) == 1 or second_check is True:
test_distance = haversine(pt, pt2)
if test_distance > biggest_distance:
biggest_distance = test_distance
if pt in beginpts:
beg = pt
end = pt2
else:
beg = pt2
end = pt
return beg, end
def calc_endpoint_coord(in_polyline):
beginlat = "BegLat"
beginlong = "BegLong"
endlat = "EndLat"
endlong = "EndLong"
field_list = [beginlat, beginlong, endlat, endlong]
existing = [field.name for field in arcpy.ListFields(in_polyline)]
field_error = False
for f in field_list:
if f not in existing:
arcpy.AddField_management(in_polyline, f, "DOUBLE")
else:
arcpy.AddError("ERROR - field named {} already exists".format(f))
field_error = True
if field_error:
sys.exit()
with arcpy.da.UpdateCursor(in_polyline, ['SHAPE@', beginlong, beginlat, endlong, endlat]) as cur:
for row in cur:
beginpts = [(part.getObject(0).Y, part.getObject(0).X) for part in row[0]]
endpts = [(part.getObject(part.count - 1).Y, part.getObject(part.count - 1).X) for part in row[0]]
allpts = beginpts + endpts
beg, end = find_distant_pts(allpts, beginpts, endpts)
if beg == end:
beg, end = find_distant_pts(allpts, beginpts, endpts, True)
row[1] = beg[1]
row[2] = beg[0]
row[3] = end[1]
row[4] = end[0]
cur.updateRow(row)
in_polyline=r"C:\Users\qin\Desktop\zxy\shp_lb.shp"
calc_endpoint_coord(in_polyline)
在代码的中,注意使用for循环的便捷形式。
beginpts = [(part.getObject(0).Y, part.getObject(0).X) for part in row[0]]
endpts = [(part.getObject(part.count - 1).Y, part.getObject(part.count - 1).X) for part in row[0]]
最后,我们来看一下具体的实现结果。
更多内容,请微信扫二维码关注公众号,或者加入arcpy开发qq学习群:487352121
转载自:https://blog.csdn.net/u010608964/article/details/89412976