(三十五)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

                                                                     

  • 0
    点赞
  • 5
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 1
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

yGIS

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值