Форумы GIS-Lab.info

Код: Выделить всё

import pyproj
from rdp import rdp
from openpyxl import load_workbook

xlssourse = 'coord.xlsx'
wbsource = load_workbook(xlssourse)
sheetsource = wbsource['Лист1']

PointList = []
geodesic = pyproj.Geod(ellps='WGS84')

row = 2
fileId = sheetsource.cell(column=1, row=row).value

while fileId is not None:
    long = sheetsource.cell(column=1, row=row).value
    lat = sheetsource.cell(column=2, row=row).value
    if (type(long) == float) and (type(lat) == float):
        if [long, lat] not in PointList:
            PointList.append([long, lat])
    row += 1
    fileId = sheetsource.cell(column=1, row=row).value

PointListRdp = rdp(PointList, epsilon=0.0000005)

row = 1
for element in PointListRdp:
    pntOneX = element[0]
    pntOneY = element[1]
    if ('pntTwoX' in locals()) and ('pntTwoY' in locals()):
        cur_azimuth, back_azimuth, distance = geodesic.inv(pntOneX, pntOneY, pntTwoX, pntTwoY)
        if 'new_azimuth' in locals():
            dif_azimuth = new_azimuth - cur_azimuth
            if (abs(dif_azimuth) > 120) and (abs(dif_azimuth) < 240):
                print(pntTwoX, pntTwoY)
        new_azimuth = cur_azimuth
    pntTwoX = pntOneX
    pntTwoY = pntOneY
    row += 1