32、找凸包

import random as rd
import matplotlib.pyplot as plt
import functools
import math
rd.seed(10)

import time
#记录每个函数执行的时间
loop = 1
def timeLog(f):
    def wrapper(*args, **kw):
        global loop
        now = time.perf_counter()
        res = f(*args, **kw)
        print("%s-%s:"%(loop, f.__name__), time.perf_counter()-now)
        loop += 1
        return res
    return wrapper

#先看纵坐标，再看横坐标
def cmp_2d(x,y):
    if x[1] < y[1]:
        return -1
    elif x[1] > y[1]:
        return 1
    elif x[0] < y[0]:
        return -1
    elif x[0] > y[0]:
        return 1
    else:
        return 0

def points_generate(x_min, y_min, x_max, y_max, nums):

    points = []
    for i in range(nums):
        points.append([int((x_max-x_min) * rd.random() + x_min),
                     int((y_max-y_min) * rd.random() + y_min)])
    return points

def clusters_plot(points, boundary):
    x = []
    y = []
    x_b = []
    y_b = []
    for i in range(len(points)):
        x.append(points[i][0])
        y.append(points[i][1])
    for i in range(len(boundary)):
        x_b.append(boundary[i][0])
        y_b.append(boundary[i][1])
    #点、线展示
    plt.scatter(x, y)
    plt.plot(x_b, y_b, color='r')

    plt.show()
    # 保存
    # plt.scatter(x, y, c=col)
    # for i in range(len(boundary)-1):
    #     plt.plot(x[i:i+2],y[i:i+2], color = 'r')
    # plt.savefig("result-1-clusters_DelvFront.jpg")

def L_Region_plot(points, L_boundary):

    x = []
    y = []
    for i in range(len(points)):
        x.append(points[i][0])
        y.append(points[i][1])
    plt.scatter(x, y)

    for boundary in L_boundary:
        x_b = []
        y_b = []
        for i in range(len(boundary)):
            x_b.append(boundary[i][0])
            y_b.append(boundary[i][1])
        plt.plot(x_b, y_b, color='r')

    plt.show()

@timeLog
def boundary_find(points):

    if len(points) == 0:
        return []

    points = sorted(points, key=functools.cmp_to_key(cmp_2d))

    points1 = [points[0]]
    for i in range(1,len(points)):
        if points[i][0] != points1[-1][0] or points[i][1] != points1[-1][1]:
            points1.append(points[i])

    points = points1

    #计算每两个点之间的角度
    LL_rad = []
    for i in range(len(points)):
        L_rad = []
        for j in range(len(points)):
            x1 = points[i][0]
            y1 = points[i][1]
            x2 = points[j][0]
            y2 = points[j][1]
            if x1 == x2 and y1 == y2:
                L_rad.append(0)
                continue
            if x2 > x1:
                if y2 >= y1:
                    L_rad.append(math.atan((y2 - y1) / (x2 - x1)))
                elif y2 < y1:
                    L_rad.append(2*math.pi + math.atan((y2 - y1) / (x2 - x1)))
            elif x2 < x1:
                if y2 >= y1:
                    L_rad.append(math.pi + math.atan((y2 - y1) / (x2 - x1)))
                elif y2 < y1:
                    L_rad.append(math.pi + math.atan((y2 - y1) / (x2 - x1)))
            else:
                if y2 > y1:
                    L_rad.append(math.pi/2)
                else:
                    L_rad.append(3*math.pi/2)
        LL_rad.append(L_rad)

    LL_rad = [[round(1000*LL_rad[i][j])/1000 for j in range(len(LL_rad[i]))] for i in range(len(LL_rad))]

    # 先找到纵坐标最小的点
    points_used = [0 for i in range(len(points))]
    p0 = points[0]
    idex0 = 0
    rad0 = 0
    for i in range(len(points)):
        if p0[1] > points[i][1] or (p0[1] == points[i][1] and p0[0] > points[i][0]):
            p0 = points[i]
            idex0 = i

    boundary = []
    boundary.append(p0)
    points_used[idex0] = 1

    loop = 0
    while loop < len(points):
        loop += 1
        idexs = [i for i in range(len(points)) if LL_rad[idex0][i] >= rad0 and i != idex0]
        if len(idexs) == 0:
            break
        idex1 = idexs[0]
        for i in idexs:
            if LL_rad[idex0][idex1] > LL_rad[idex0][i]:
                idex1 = i
        boundary.append(points[idex1])
        if points_used[idex1] == 1:
            break
        points_used[idex1] = 1
        rad0 = LL_rad[idex0][idex1]
        idex0 = idex1

    return boundary

#余弦值计算
def cosSinCount(x1, y1, x2, y2):
    dx, dy = x2-x1, y2-y1
    r = math.sqrt(dx*dx+dy*dy)
    cos = dx/r if r > 0 else dx
    sin = dy/r if r > 0 else dy
    return cos, sin

#夹角统计
def radCount(x1, y1, x2, y2):
    cos, sin = cosSinCount(x1, y1, x2, y2)
    if cos == 0 and sin == 0:
        return 0
    if cos > 0:
        if sin >= 0:
            return math.asin(sin)
        else:
            return math.asin(sin) + 2 * math.pi
    elif cos < 0:
        return math.pi - math.asin(sin)
    else:
        if sin > 0:
            return math.pi / 2
        else:
            return 3 * math.pi / 2

@timeLog
def boundary_find2(points):

    if len(points) <= 3:
        return points

    #去除重复的点
    points1 = set([(p[0], p[1]) for p in points])
    points1 = [p for p in points1]

    #先找y值最小的点（多个取x最小）P0
    p0 = points1[0]
    for p in points1[1:]:
        if p[1] < p0[1] or (p[1] == p0[1] and p[0] < p0[0]):
            p0 = p
    points1.remove(p0)

    #其他点按照与P0形成夹脚的余弦值进行从大到小排序，设为P1,P2,...,Pn
    points2 = [[p[0], p[1], cosSinCount(p0[0], p0[1], p[0], p[1])[0]] for p in points1]
    points2.sort(key= lambda p:p[2], reverse=True)

    #P1一定在凸包上，之后逐个加进来，如果往外拐则进行更新，直至绕回来，复杂度为O(nlogn)+O(n)
    boundary = [p0, points2[0]]
    rads = [radCount(boundary[0][0], boundary[0][1], boundary[1][0], boundary[1][1])]
    for p in points2[1:]:
        radi = radCount(boundary[-1][0], boundary[-1][1], p[0], p[1])
        while radi < rads[-1]:
            boundary.pop()
            rads.pop()
            radi = radCount(boundary[-1][0], boundary[-1][1], p[0], p[1])
        boundary.append(p)
        rads.append(radi)
    boundary.append(p0)
    boundary = [[p[0], p[1]] for p in boundary]
    return boundary


if __name__ == "__main__":
    points = points_generate(0,0,100,100,10)
    points = points_generate(0, 0, 1000, 1000, 1000)
    # points = [[1, 1], [2, 2], [1, 1], [3, 3], [2, 2]]
    # points = [[32.82964904320003, 23.95932630554721], [42.147973315233166, 36.412172357406526], [52.74098811758407, 14.393079909900852], [4.642057278225897, 28.317158594325143]]
    # points = [[113.166123,41.643543],[113.166124,41.643544],[113.155905,41.715037],[113.154363,41.71566],[113.118134,41.711303],[113.117487,41.711158],
    #           [113.117198,41.710728],[113.115037,41.710429],[113.105046,41.707767],[113.102645,41.707553],[113.105096,41.705223],[113.105209,41.705104],
    #           [113.097391,41.702161],[112.981497,41.67395],[112.943379,41.668826],[112.920114,41.672687],[112.920305,41.673401],[112.920053,41.673463],
    #           [112.920051,41.673461],[112.919231,41.673812],[112.918151,41.674616],[112.918151,41.674611],[112.917452,41.673712],[112.915135,41.669412]]
    print(points)
    boundary = boundary_find(points)
    print(boundary)
    # clusters_plot(points, boundary)
    boundary2 = boundary_find2(points)
    print(boundary2)
    # clusters_plot(points, boundary2)