Graham-Scan算法计算凸包的Python代码实现

对于一个点集P来讲，它的凸包就是一个凸多边形Q，其中满足P中的每个点都在Q的边界上或内部。就像下图所示

凸包的计算算法有好多种，wiki和算法导论第33章中都有比较详细的介绍，比如下面是算法导论中给出的Graham-Scan算法计算凸包的伪代码。

现在网上已经有了好多计算点集凸包的优秀代码，比如这篇文章，作者在文中使用了一个动画来表示了Graham-Scan算法计算凸包的过程，并给出了python程序的实现，十分有助于学习者对算法的理解。最近有个东西需要使用到凸包，本着“不要重复造轮子”的原则，我在开始的时候直接使用了作者文中的程序。开始使用的时候没有发现什么问题，比如下图所示的效果。

但是，当我在使用其他的数据进行测试的时候，发现程序执行的效果并不太好，下图是使用著名的Iris数据集的前两个维度测试的效果，这明显不是一个凸包。

没有找到原因，于是开始自己写程序，发现结果也是这样的。后来将执行过程中栈的状态打印输出，才发现我们两个的程序中都存在一个小问题。当点集中出现两个完全相同的点的时候，算法就是失效了。上面伪代码中序号8所示的while循环的判断条件不应该是两个向量的叉乘大于0，而应该是大于等于0。否则的话，对于点集中出现连续两个完全相同的点的时候，栈只能弹出其中一个点，而另一个点会仍然在栈中。所以会出现上图所示的结果。修改之后的结果如下图所示，这才是一个正确的凸包。

具体的Python实现代码如下

import matplotlib.pyplot as plt
import mathimport sklearn.datasets as datasets"""
使用Graham扫描法计算凸包
网上的代码好多运行效果并不好
算法参见《算法导论》第三版 第605页
"""def get_bottom_point(points):"""返回points中纵坐标最小的点的索引，如果有多个纵坐标最小的点则返回其中横坐标最小的那个:param points::return:"""min_index = 0n = len(points)for i in range(0, n):if points[i][1] < points[min_index][1] or (points[i][1] == points[min_index][1] and points[i][0] < points[min_index][0]):min_index = ireturn min_indexdef sort_polar_angle_cos(points, center_point):"""按照与中心点的极角进行排序，使用的是余弦的方法:param points: 需要排序的点:param center_point: 中心点:return:"""n = len(points)cos_value = []rank = []norm_list = []for i in range(0, n):point_ = points[i]point = [point_[0]-center_point[0], point_[1]-center_point[1]]rank.append(i)norm_value = math.sqrt(point[0]*point[0] + point[1]*point[1])norm_list.append(norm_value)if norm_value == 0:cos_value.append(1)else:cos_value.append(point[0] / norm_value)for i in range(0, n-1):index = i + 1while index > 0:if cos_value[index] > cos_value[index-1] or (cos_value[index] == cos_value[index-1] and norm_list[index] > norm_list[index-1]):temp = cos_value[index]temp_rank = rank[index]temp_norm = norm_list[index]cos_value[index] = cos_value[index-1]rank[index] = rank[index-1]norm_list[index] = norm_list[index-1]cos_value[index-1] = temprank[index-1] = temp_ranknorm_list[index-1] = temp_normindex = index-1else:breaksorted_points = []for i in rank:sorted_points.append(points[i])return sorted_pointsdef vector_angle(vector):"""返回一个向量与向量 [1, 0]之间的夹角， 这个夹角是指从[1, 0]沿逆时针方向旋转多少度能到达这个向量:param vector::return:"""norm_ = math.sqrt(vector[0]*vector[0] + vector[1]*vector[1])if norm_ == 0:return 0angle = math.acos(vector[0]/norm_)if vector[1] >= 0:return angleelse:return 2*math.pi - angledef coss_multi(v1, v2):"""计算两个向量的叉乘:param v1::param v2::return:"""return v1[0]*v2[1] - v1[1]*v2[0]def graham_scan(points):# print("Graham扫描法计算凸包")bottom_index = get_bottom_point(points)bottom_point = points.pop(bottom_index)sorted_points = sort_polar_angle_cos(points, bottom_point)m = len(sorted_points)if m < 2:print("点的数量过少，无法构成凸包")returnstack = []stack.append(bottom_point)stack.append(sorted_points[0])stack.append(sorted_points[1])for i in range(2, m):length = len(stack)top = stack[length-1]next_top = stack[length-2]v1 = [sorted_points[i][0]-next_top[0], sorted_points[i][1]-next_top[1]]v2 = [top[0]-next_top[0], top[1]-next_top[1]]while coss_multi(v1, v2) >= 0:stack.pop()length = len(stack)top = stack[length-1]next_top = stack[length-2]v1 = [sorted_points[i][0] - next_top[0], sorted_points[i][1] - next_top[1]]v2 = [top[0] - next_top[0], top[1] - next_top[1]]stack.append(sorted_points[i])return stackdef test1():points = [[1.1, 3.6],[2.1, 5.4],[2.5, 1.8],[3.3, 3.98],[4.8, 6.2],[4.3, 4.1],[4.2, 2.4],[5.9, 3.5],[6.2, 5.3],[6.1, 2.56],[7.4, 3.7],[7.1, 4.3],[7, 4.1]]for point in points:plt.scatter(point[0], point[1], marker='o', c='y')result = graham_scan(points)length = len(result)for i in range(0, length-1):plt.plot([result[i][0], result[i+1][0]], [result[i][1], result[i+1][1]], c='r')plt.plot([result[0][0], result[length-1][0]], [result[0][1], result[length-1][1]], c='r')plt.show()def test2():"""使用复杂一些的数据测试程序运行效果:return:"""iris = datasets.load_iris()data = iris.datapoints_ = data[:, 0:2]points__ = points_[0:50, :]points = points__.tolist()temp_index = 0for point in points:plt.scatter(point[0], point[1], marker='o', c='y')index_str = str(temp_index)plt.annotate(index_str, (point[0], point[1]))temp_index = temp_index + 1result = graham_scan(points)print(result)length = len(result)for i in range(0, length-1):plt.plot([result[i][0], result[i+1][0]], [result[i][1], result[i+1][1]], c='r')plt.plot([result[0][0], result[length-1][0]], [result[0][1], result[length-1][1]], c='r')# for i in range(0, len(rank)):plt.show()if __name__ == "__main__":test2()

Graham-Scan算法计算凸包的Python代码实现相关推荐

凸包Graham Scan算法实现
凸包Graham Scan算法实现凸包算法实现点集合中搜索凸包顶点的功能,可以处理共线情况,可以输出共线点也可以不输出而只输出凸包顶点.经典的Graham Scan算法,点排序使用极角排序方式,并对 ...
计算几何--凸包之graham scan算法
Graham scan算法主要步骤: 找出所有已知点的y值最小,如果相同,取x值最小的点,作为基准点s. 以s为基准,所有的点按照与X轴夹角从小到大排序. 使用两个栈,一个记录已访问的点,一个记录未访 ...
计算几何入门 1.6：凸包的构造——Graham Scan算法
上文简要分析出了凸包构造问题算法的下界:O(nlogn),在此就引入一种下界意义上最优的算法:Graham Scan算法.这种算法可以保证在最坏情况下时间复杂度也不超过nlogn.我们先大致了解一下算 ...
graham算法 java_凸包Graham Scan算法实现
凸包算法实现点集合中搜索凸包顶点的功能,可以处理共线情况,可以输出共线点也可以不输出而只输出凸包顶点.经典的Graham Scan算法,点排序使用极角排序方式,并对共线情况做特殊处理.一般算法是将共线 ...
Graham Scan算法
预备知识:凸集:集合S中任意两点的连线都在集合S中,如果简单的理解,可以理解为凸边形凸包:对于给定集合X,所有包含X的凸集的交集,简单的理解,就是包含X的最小凸集,或者就是最外圈的点连起来极角排序 ...
用通俗易懂的方式讲解：主成分分析(PCA)算法及案例（Python 代码）
文章目录知识汇总加入方式一.引入问题二.数据降维三.PCA基本数学原理 3.1 内积与投影 3.2 基 3.3 基变换的矩阵表示 3.4 协方差矩阵及优化目标 3.5 方差 3.6 协方差 ...
算法学习之模拟退火算法路径规划（python代码实现）
模拟退火算法路径规划(python代码实现) 一.引言二.算法介绍以及伪代码 1.算法通俗介绍 2.路径规划算法伪代码三.算法流程及代码实现 1.地图创建 2.初始化路径小结 3.计算适应度值 ...
kmeans算法详解和python代码实现
kmeans算法详解和python代码实现 kmeans算法无监督学习和监督学习监督学习: 是通过已知类别的样本分类器的参数,来达到所要求性能的过程简单来说,就是让计算机去学习我们已经创建好了的 ...
编辑距离算法详解和python代码
编辑距离(Levenshtein Distance)算法详解和python代码最近做NLP用到了编辑距离,网上学习了很多,看到很多博客写的有问题,这里做一个编辑距离的算法介绍,步骤和多种python ...

Graham-Scan算法计算凸包的Python代码实现

Graham-Scan算法计算凸包的Python代码实现相关推荐

最新文章

热门文章