[从头学数学] 第241节 计算几何 导言 凸包的例子

剧情提要：
阿伟看到了一本比较有趣的书，是关于《计算几何》的，2008年由北清派出版。很好奇
它里面讲了些什么，就来看看啦。

正剧开始：
星历2016年07月21日 16:12:28, 银河系厄尔斯星球中华帝国江南行省。

[工程师阿伟]正在和[机器小伟]一起研究[计算几何]]。

<span style="font-size:18px;">#凸包问题专题
class Convex():
    #初始化
    def initVector(self, point):
        return np.array(point);

    #一般在最后输出前进行圆整
    def roundVector(self, array, precision = 3):
        #precision是小数点后的位数
        #array为[[x1, y1, (z1)], [x2, y2, (z2)], ...]这种点阵格式
        dims = np.shape(array);
        len_ = len(dims);

        if (len_ == 0):
            array = round(array, precision);            
        if (len_ == 1):
            for i in range(dims[0]):
                array[i] = round(array[i], precision);
        elif (len_ == 2):
            for i in range(dims[0]):
                for j in range(dims[1]):
                    array[i][j] = round(array[i][j], precision);
        elif (len_ == 3):
            for i in range(dims[0]):
                for j in range(dims[1]):
                    for k in range(dims[2]):
                        array[i][j][k] = round(array[i][j][k], precision);         
        return array;
    
    #点P在有向线段Vec的左侧还是右侧
    def side_P_Vec(self, Point, Vec):
        Q, R = Vec[0], Vec[1];
        P_x, P_y = Point[0], Point[1];
        deter = np.array([[1, Q[0], Q[1]],
                         [1, R[0], R[1]],
                         [1, P_x, P_y]]);

        det = np.linalg.det(deter);

        #det < 0 右侧
        #det = 0 点在直线上
        #det > 0 左侧

        return det;

    #去除重复点
    def removeDup(self, vertices):
        v_new = [];
        len_1 = len(vertices);
        len_2 = 0;

        for i in range(len_1):
            len_2 = len(v_new);
            if (len_2 < 1):
                v_new.append(vertices[i]);
            else:
                for j in range(len_2):
                    if v_new[j] == vertices[i]:
                        break;
                    if (j >= len_2-1):
                        v_new.append(vertices[i]);

        return v_new;
        
    #二维点集的凸包
    def convexHull(self, pointSet):
        pointSet = self.roundVector(pointSet);
        pointSet = self.removeDup(pointSet);
        #按照x坐标由小到大排序，x相同时考虑y坐标
        P = sorted(pointSet, key=lambda a:(a[0], a[1]));
        num = len(P);

        Lupper = [P[0], P[1]];
        
        for i in range(2, num):
            Lupper = Lupper + [P[i]];
            
            while (len(Lupper) >= 3 and self.side_P_Vec(Lupper[-3], [Lupper[-2], Lupper[-1]])>=0):
                Lupper = Lupper[:-2]+[Lupper[-1]];

        Llower = [P[-1], P[-2]];

        for i in range(num-2, -1, -1):
            Llower = Llower + [P[i]];
            while (len(Llower) >= 3 and self.side_P_Vec(Llower[-3], [Llower[-2], Llower[-1]])>=0):
                Llower = Llower[:-2]+[Llower[-1]];

        Llower = Llower[1:-1];
        
        return Lupper+Llower;

      


def tmp():        
    calc = Convex();
    shape = geo.Shape();

    #计时开始
    startTime = time.clock();    
    
    verts = shape.nStar(0, 0, 100, 7);

    verts += shape.nStar(0, 0, 102, 18);
    verts += shape.nEdge(0, 0, 150, 15);
    
    verts = calc.roundVector(verts);
    print(verts);

    result = calc.convexHull(verts);
  

    #计时结束
    endTime = time.clock();

    #打印结果
    print('result:', result);
    print('操作用时：{0:.3e} s'.format(endTime-startTime));
</span>

本节到此结束，欲知后事如何，请看下回分解。