75、QGIS学习记录

最新推荐文章于 2024-03-21 00:07:09 发布

chaoyuzhang

最新推荐文章于 2024-03-21 00:07:09 发布

阅读量188

点赞数

分类专栏：图形学

本文链接：https://blog.csdn.net/chaoyuzhang/article/details/115485046

版权

图形学专栏收录该内容

3 篇文章 0 订阅

订阅专栏

在console画多边形

wkt2 = 'POLYGON ((11717491.32085317 4192105.64176356, 11717299.80128882 4193264.57267995, 11717299.8012888 4193264.5726799, 11717078.05634601 4192351.42419002, 11715386.86564284 4193357.23269148, 11715382.15507762 4193327.5158613, 11715382.1550776 4193327.5158613, 11715377.1297748 4193340.6659426, 11715390.95950833 4193383.05904149, 11715390.95950837 4193383.05904173, 11715325.6489058 4193440.2176517, 11715339.93261539 4193493.72901487, 11715339.9326154 4193493.7290149, 11715339.93261541 4193493.72901487, 11715575.92693334 4193535.72223333, 11717466.67044189 4193535.72223333, 11717515.67084537 4192091.15999688, 11717491.32085317 4192105.64176356), ())'
temp = QgsVectorLayer('POLYGON?crs=epsg:4326', 'tmpgeo', 'memory')
QgsProject.instance().addMapLayer(temp)
temp.startEditing()
geom = QgsGeometry().fromWkt(wkt2)
feature = QgsFeature()
feature.setGeometry(geom)
temp.dataProvider().addFeatures([feature])
# 提交修改
temp.commitChanges()

剔除尖刺

import math

#清洗成需要的格式
def data_clean(wkt):
    wkt = wkt[10:-3]
    P1, P2 = wkt.split("), (")

    L1 = P1.split(", ")
    L1 = [[float(s) for s in c.split(" ")] for c in L1]

    L2 = P2.split(", ")
    L2 = [[float(s) for s in c.split(" ")] for c in L2]

    return L1, L2

#数据标准化
def CoordinateStandardization(L1, L2):
    x_min, y_min = min([c[0] for c in L1] + [c[0] for c in L2]), min([c[1] for c in L1] + [c[1] for c in L2])
    L1 = [[c[0]-x_min, c[1]-y_min] for c in L1]
    L2 = [[c[0] - x_min, c[1] - y_min] for c in L2]
    return L1, L2

def degrees_generate(a, b, c):
    ba = [a[0]-b[0], a[1]-b[1]]
    bc = [c[0]-b[0], c[1]-b[1]]
    degree = 0
    if (math.sqrt(ba[0]*ba[0]+ba[1]*ba[1])*math.sqrt(bc[0]*bc[0]+bc[1]*bc[1])) > 0:
        degree = math.acos((ba[0]*bc[0]+ba[1]*bc[1])/(math.sqrt(ba[0]*ba[0]+ba[1]*ba[1])*math.sqrt(bc[0]*bc[0]+bc[1]*bc[1])))
        degree = degree*180/math.pi
    return degree

def str_Polygon_generate(L1, L2):
    res = "POLYGON (("
    for i in range(len(L1)):
        res += str(L1[i][0]) + " " + str(L1[i][1])
        if i < len(L1)-1:
            res += ", "
    res += "), ("
    for i in range(len(L2)):
        res += str(L2[i][0]) + " " + str(L2[i][1])
        if i < len(L2)-1:
            res += ", "
    res += "))"
    return res

#将多边形的内部线去掉
def PolygonCut(wkt, e=0.05):

    #数据清洗
    L1, L2 = data_clean(wkt)

    #初始标准化数据
    L1_0, L2_0 = CoordinateStandardization(L1, L2)

    L1, L2 = L1[:-1], L2[:-1]

    #循环砍点
    loop = 0
    while len(L1) >= 3:
        loop += 1
        nums0 = len(L1)

        degrees1 = [degrees_generate(L1[-1], L1[0], L1[1])] + \
                   [degrees_generate(L1[i], L1[i + 1], L1[i + 2]) for i in range(len(L1) - 2)] + \
                   [degrees_generate(L1[-2], L1[-1], L1[0])]

        print("\nloop: ", loop, "\ndegrees1: ", degrees1)

        for i in range(len(L1)):
            if degrees1[i] <= e or degrees1[i] >= 180-e:
                L1 = L1[:i]+L1[i+1:]
                break

        if len(L1) == nums0:
            break
    if len(L1) < 3:
        L1 = []
    else:
        L1.append(L1[0])

    loop = 0
    while len(L2) >= 3:
        loop += 1
        nums0 = len(L2)

        degrees2 = [degrees_generate(L2[-1], L2[0], L2[1])] + \
                   [degrees_generate(L2[i], L2[i + 1], L2[i + 2]) for i in range(len(L2) - 2)] + \
                   [degrees_generate(L2[-2], L2[-1], L2[0])]
        print("\nloop: ", loop, "\ndegrees2: ", degrees2)

        for i in range(len(L2)):
            if degrees2[i] <= e or degrees2[i] >= 180 - e:
                L2 = L2[:i] + L2[i + 1:]
                break

        if len(L2) == nums0:
            break
    if len(L2) < 3:
        L2 = []
    else:
        L2.append(L2[0])

    #最终数据
    L1_1, L2_1 = CoordinateStandardization(L1, L2)

    return str_Polygon_generate(L1, L2)

if __name__ == '__main__':
    wkt = "POLYGON (" \
          "(" \
          "11717491.32085317 4192105.64176356, 11717299.80128882 4193264.57267995, 11717006.32711409 4192915.76476822, " \
          "11717096.4439237 4193022.8728509, 11717299.8012888 4193264.5726799, 11717078.05634601 4192351.42419002, " \
          "11715386.86564284 4193357.23269148, 11715382.15507762 4193327.5158613, 11715382.1550776 4193327.5158613, " \
          "11715377.1297748 4193340.6659426, 11715390.95950833 4193383.05904149, 11715390.95950837 4193383.05904173, " \
          "11715325.6489058 4193440.2176517, 11715339.93261539 4193493.72901487, 11715339.9326154 4193493.7290149, " \
          "11715339.93261541 4193493.72901487, 11715575.92693334 4193535.72223333, 11717466.67044189 4193535.72223333, " \
          "11717515.67084537 4192091.15999688, 11717491.32085317 4192105.64176356), " \
          "(11715386.86564284 4193357.23269148, 11715388.86665743 4193369.85618956, 11715386.8656428 4193357.2326915, " \
          "11715386.86564284 4193357.23269148)" \
          ")"
    print(PolygonCut(wkt))