计算一个多边形的重心点坐标(准确版)

2019-11-29 13:04:47

主要用途,就是在百度地图上画多边形,然后我们要得到这个多边形的重心。

后来我发现了,百度地图,map对象直接有这个接口。大家可以仔细看下api。

接下来,我们还是要看下这个算法问题。

参考地址 计算一个多边形的重心点坐标



背景介绍与问题分析

在之前的 《如何判断一个多边形是否合法》  一文中有提到,用无人机规划飞行路线前,往往需要框选一个多边形的区域。

而在地图控件上显示这个多边形区域时,往往会遇到这样一个需求:需要把所要测绘的多边形区域移动到地图中心。

实现这个需求的基本思路就是:获取到多边形区域的重心点坐标,然后利用地图控件的 setCenter方法,就可以把地图的显示中心移动到多边形区域重心了。那么问题来了,如何求出一个多边形的重心点坐标呢?

这里所说的重心,也常常叫几何中心

这里首先给出一个公式:

平面多边形 X 可以被剖分为 n个有限的简单图形 X_1,X_2,....X_n,这些简单图形的重心点为 C_1,面积为 A_1,那么这个平面多边形的重心点坐标为 (C_x,C_y)
C_x = \frac{\sum C_{i_x} A_i}{\sum A_i}, C_y = \frac{\sum C_{i_y} A_i}{\sum A_i}

公式参考: 维基百科

一般来说我们可以给多边形进行三角剖分,而 \sum{A_i} 即为多边形的总面积,那么这个公式可以理解为:

多边形重心横坐标 = 多边形剖分的每一个三角形重心的横坐标 * 该三角形的面积之和 / 多边形总面积

多边形重心纵坐标 = 多边形剖分的每一个三角形重心的纵坐标 * 该三角形的面积之和 / 多边形总面积

所以这里就把问题拆分成了三个小问题:

  • 求每个剖分出来的三角形的重心。

  • 求每个剖分出来的三角形的面积。

  • 求多边形的面积。

算法解析

1. 求三角形的重心


三角形重心


三角形的重心:三条中线的交点。其中重心到其中一个顶点的距离是重心到该顶点对边中点的距离的2倍。
即:GC = 2 * GP,也就是说重心坐标在 CP 线段上距离 AB 的中点 P 的 1/3 处。
假设 A,B,C 三点的坐标为:


那么通过简单坐标计算,可以得出其重心坐标为 (x,y)
x = \frac{(x_1+x_2+x_3)}{3} , y = \frac{(y_1+y_2+y_3)}{3}

2. 求三角形面积

计算三角形的面积,我们这里利用 向量积来计算,我们知道平面中的两个向量的叉乘的模等于以这两个向量为边的平行四边形的面积,那么以这个两个向量为边的三角形,则是这个平行四边形的面积的一半。

参考:向量叉积

如上图,已知平面上两点 A:(x_1,y_1),B(x_2,y_2) ,以 A,B和坐标原点 P(0,0) 构成的三角形的面积 S 为:
S=\frac{\vec{PB}\times\vec{PA}}{2} = \frac{x_2y_1 - x_1y_2 }{2}

这里给出运算草稿:


为什么这里我们会以原点作为第三个点构成三角形呢?其实是跟接下来求多边形面积是有关联的。

3. 求多边形的面积

我们在上面给出的求平面多边形重心的公式中有说到,一般我们会把多边形剖分为多个三角形。
那么这个剖分点 P 我们可以设在哪里呢?这里先给出结论:这个剖分点可以设置在多边形的内部,也可以设置到外部。

为什么这个剖分点可以设置到外部呢?我们可以通过简单的三角形情况来推广到多边形的情况。
对于三角形ABC,我们把剖分点设置在其外部 P 的一点上


如果大家还记得 《如何判断一个多边形是否合法》 一文中有讲过向量叉积是有正负之分的,并且根据上面所说的计算三角形面积,那么以 P 为剖分点,通过向量积可以得出这个三角形的面积 A 为:
A = \frac{1}{2}(\vec{PB} \times \vec{PC} + \vec{PC} \times \vec{PA} + \vec{PA} \times \vec{PB})

因为 向量PB 在 向量PA 的顺时针方向,所以 \vec{PA} \times \vec{PB} 的结果是负数的。那么上面的面积计算公式其实就可以理解为:

三角形ABC的面积 = 三角形PBC面积 + 三角形PCA面积 - 三角形PAB面积

假设这四个点的坐标为:P(x_0,y_0), A(x_1,y_1), B(x_2,y_2), C(x_3,y_3),通过上面的公式进行计算,具体的演算过程我就不给出了,这里直接给出计算结果:
A = x_1y_2-x_2y_1+x_2y_3-x_3y_2+x_3y_1-x_1y_3

我们可以发现,计算结果中没有 x_0,y_0 的项,因为它们在计算过程中给消去了,数学就是这么奇妙!所以我们可以得出一个结论,多边形的面积结果与这个剖分点的位置是无关的。那么为了计算方便,我们当然选择把这个 P 点设置到原点上啦。

那么只要我们知道多边形的每一个顶点,通过原点进行剖分成多个三角形,然后通过向量的叉乘求出每个三角的面积,最后相加,就可以求出多边形的面积了。

示例代码及解析

好了,说到这里,我们已经找到所有满足最开始的计算多边形重心点坐标的所有计算元素了。是时候上代码了,这里构建一个函数calculatePolygonGravityCenter(coordinates: [CLLocationCoordinate2D]),这个函数传入的参数是多边形在地图上的坐标点数组。

func calculatePolygonGravityCenter(coordinates: [CLLocationCoordinate2D]) -> CLLocationCoordinate2D {
    var area = 0.0 // 多边形面积
    var gravityLat = 0.0 // 重心点 latitude
    var gravityLng = 0.0 // 重心点 longitude
    for (index, coordinate) in coordinates.enumerated() {
        // 1
        let lat = coordinate.latitude        let lng = coordinate.longitude        let nextLat = coordinates[(index + 1) % coordinates.count].latitude        let nextLng = coordinates[(index + 1) % coordinates.count].longitude        // 2
        let tempArea = (nextLat * lng - nextLng * lat) / 2.0
        // 3
        area += tempArea        // 4
        gravityLat += tempArea * (lat + nextLat) / 3
        gravityLng += tempArea * (lng + nextLng) / 3
    }
    // 5
    gravityLat = gravityLat / area
    gravityLng = gravityLng / area    
    return CLLocationCoordinate2D(latitude: gravityLat, longitude: gravityLng)}

对应上面代码的注释:

  1. 拿到多边形上连续两个点的坐标,我们可以把 latitude 看做横坐标,longitude 是纵坐标。

  2. 利用向量叉乘计算这两个点与原点组成的三角形的面积。

  3. 所有面积之和得出多边形的面积,就是求公式 C_x = \frac{\sum C_{i_x} A_i}{\sum A_i} 中的 \sum A_i

  4. (lat + nextLat) / 3 是以这两个点和原点组成的三角形的重心横坐标,这样的累加gravityLat += tempArea * (lat + nextLat) / 3 其实是求公式 C_x = \frac{\sum C_{i_x} A_i}{\sum A_i} 中的 \sum C_{i_x} A_i 的值。

  5. 到这一步就简单了,直接套用公式 C_x = \frac{\sum C_{i_x} A_i}{\sum A_i}

参考资料

  1. 维基百科-Centroid

  2. 维基百科-叉积

  3. cnblogs-用向量积(叉积)计算三角形的面积和多边形面积

  4. 知乎-两个向量的叉乘为什么是面积

  5. 中国知网-任意多边形匀面重心的计算方法



  • 2019-07-26 19:31:03

    Guacamole搭建

    因项目需要,经历多天查阅各种文档,几经波折终于功德圆满,写下此篇文章供大家分享。Guacamole就个人理解而言是一个可以通过web浏览器访问远程服务器终端进行操作的可视化工具。主要由web(浏览器)、Guacamole Server(核心)、Remote Desktops(远程桌面)三大模块组成。

  • 2019-07-30 22:36:10

    使用 Spring Initializr 初始化 Spring Boot 项目

    万事开头难,你需要设置一个目录结构存放各种项目内容,创建构建文件,并在其中加入各 种依赖。Spring Boot CLI消除了不少设置工作,但如果你更倾向于传统Java项目结构,那你应该 看看Spring Initializr。

  • 2019-08-06 15:30:08

    小程序展示富文本

    然而rich-text有个问题,它不能够给内联img设置宽度100%,这样图片就容易溢出屏幕。我们只能在后台返回富文本的时候对图片的img标签进行格式化。或者前端进行格式化也行,我赞成使用后端,因为很多种情况我们不一定都能想得到。

  • 2019-08-07 09:07:32

    最全的Service Worker讲解

    Service Worker 最主要的特点是:在页面中注册并安装成功后,运行于浏览器后台,不受页面刷新的影响,可以监听和截拦作用域范围内所有页面的 HTTP 请求。 基于 Service Worker API 的特性,结合 Fetch API、Cache API、Push API、postMessage API 和 Notification API,可以在基于浏览器的 web 应用中实现如离线缓存、消息推送、静默更新等 native 应用常见的功能,以给 web 应用提供更好更丰富的使用体验。

  • 2019-08-07 09:09:19

    windows系统下定时关闭程序

    其中xxx.exe是你要关闭的进程中运行的exe,可以ctrl+alt+del打开任务管理器,进到详细信息查看 然后把.txt文件后缀改成.bat(此时要在查看一栏勾上文件拓展名,要不还是txt文档)

  • 2019-08-07 09:16:43

    一个比较完美的PWA例子

    但就目前来讲,PWA是Google主推的一项技术标准,FireFox,Chrome以及一些基于Blink的浏览器已经支持渐进式Web应用了,Edge上对渐进式Web应用的支持还在开发。Apple公司也表示会考虑在自己Safari支持PWA。然而这项功能已经进入了WebKit内核的五年计划中。长期来看,对浏览器兼容性的支持方面应该已经不算太大问题了。况且在现阶段,在不支持渐进式Web应用的浏览器中,你的应用也只是无法使用渐进式Web应用的离线功能而已,除此之外的功能均可以正常使用。

  • 2019-08-07 09:57:48

    spring data jpa 实体类中字段不与数据库表映射

    当我们使用spring data jpa开发的时候,会将实体类中的成员变量与表中的字段一一对应,当我们在实体类中加上一个不与数据库表一一对应的成员变量的时候,此时我们只要在这个成员变量上加上注解@Transient @