简介这套 MATLAB 源码包面向初学者专注于多边形面积与形心的快速求解适用于数据分析、图像处理、物理建模等需要几何度量的实际场景。压缩包仅含 2 个 .m 脚本文件整体大小约 1KB代码轻量精简便于直接阅读、调试和复用已有 1055 人学习下载。两个源码函数分别封装面积计算与形心计算逻辑均以多边形顶点坐标矩阵作为入口参数内部通过遍历顶点并结合面积公式与形心公式完成求解最终返回面积数值和形心坐标使用方式非常直观。对于需要批量处理多个多边形的用户只需维护好顶点矩阵并按对应格式传入即可快速获得结果省去手动推导公式的麻烦。资源适合正在学习 MATLAB 几何计算或需要在图像分割、区域测量中集成面积与形心工具的开发者是一份务实、可直接套用的参考实现。1. “MATLAB求形心源码”不是下载问题是计算口径问题看起来像是一个压缩包下载需求拿一份源码输入坐标输出面积和形心。真正动过手的人知道问题从来不在“调用哪个函数”而在三个绕不开的口径顶点顺序是否闭合、自相交区域要不要排除、坐标数值很大时鞋带公式还稳不稳。形心又叫质心在均匀密度下等同于重心配重设计、截面抗弯、图像目标定位里每天都在用而面积计算和形心求解又是一对孪生操作先求面积再对面积加权。这篇按“数学公式 → 可运行的 MATLAB 函数 → 图像与网格场景 → 结果自检”推进适合要立刻算数、也想知道边界条件的工程师。新手照着抄能出数老手看参数和坑。2. 理论先落地鞋带公式与面积加权形心别用均值当形心2.1 为什么平均顶点坐标不是形心第一反应多半是把多边形顶点的 x 和 y 各取平均当作形心。对一个三角形、一个轴对齐的矩形这个结果恰好是对的对稍微复杂一点的 L 形截面、凹多边形均值就会偏离实际形心。原因要回到定义形心是对面积元素做加权平均的结果权重是每个微小面积块的大小而不是轮廓上几个控制点。L 形截面左下方有一大块面积右下是镂空/非图形区域把全部顶点坐标一平均均值点会明显向没有面积的方向偏移这就是静矩与几何均值在本质上的差别。在材料力学里截面形心的定义是静矩除以面积。对匀质薄板静矩 Sx ∫ y dASy ∫ x dA形心坐标 Cx Sy/ACy Sx/A。连续的积分在计算机里没法精确计算于是把边界离散成多边形用“面积加权”逼近这个积分。理解到这一层才能接受后面的公式不是“另一种算法”而是同一个积分在离散轮廓下的精确闭合解。2.2 鞋带公式多边形面积的标准算法多边形面积最常用的算法是鞋带公式也叫高斯面积公式。设顶点按逆时针顺序排列为 (x1,y1), (x2,y2), …, (xn,yn)最后一点不与第一个点重复则有向面积为A 1/2 × Σ(xi × y(i1) − x(i1) × yi)其中 i1 超过 n 时回到 1。每一项 (xi×y(i1) − x(i1)×yi) 的几何意义是两个相邻顶点与原点构成的平行四边形的有向面积。把多边形三角剖分后这些项会在内部相互抵消只剩下多边形本身的带符号面积。顶点按逆时针排列时 A 为正顺时针时为负MATLAB 内置的 polyarea 会取绝对值而后面计算形心时分母上的面积需要符号一致所以自定义函数要保留这个符号。方法与工具复杂度适用场景是否返回形心polyareaO(n)内置函数只算面积不返回形心否鞋带公式自写循环/向量化O(n)同时求面积与形心适合闭合多边形是三角剖分加每个三角形加权O(n log n)凹多边形、带孔区域拆分为子区域是蒙特卡洛随机点法O(N)任意闭合曲线、无法写成多边形近似鞋带公式对简单多边形边与边不相交是精确解比蒙特卡洛算法更可靠。对自相交的多边形公式算出的“有向面积”会部分抵消数值上依然有定义但不再是几何意义上的封闭区域面积。因此在把数据交给函数之前先确认轮廓是简单多边形。2.3 面积加权形心公式与最小可运行脚本在鞋带公式的基础上面积加权形心公式为Cx Σ( (xi x(i1)) × (xi×y(i1) − x(i1)×yi) ) / (6 × A) Cy Σ( (yi y(i1)) × (xi×y(i1) − x(i1)×yi) ) / (6 × A)这个式子可以直接写成 MATLAB 循环方便与公式字面逐项对照x [0 1 1 0]; % 单位正方形的四个顶点逆时针 y [0 0 1 1]; n length(x); A 0; cx 0; cy 0; for i 1:n j mod(i, n) 1; % 下一个顶点末点回到首点 crossTerm x(i)*y(j) - x(j)*y(i); % 相邻顶点叉积即鞋带项 A A crossTerm; cx cx (x(i) x(j)) * crossTerm; cy cy (y(i) y(j)) * crossTerm; end A A / 2; cx cx / (6 * A); cy cy / (6 * A); fprintf(A%.6f, C(%.6f, %.6f)\n, A, cx, cy);这段脚本里 crossTerm 是核心它两用既参与面积累加又是形心分子上的权重。对单位正方形跑一遍输出 A1C(0.5, 0.5)和直觉一致。把顶点改成三角形 (0,0)、(1,0)、(0,1)输出 A0.5C(1/3, 1/3)与三角形形心是顶点均值的结论吻合。提示循环写法只适合理解公式。顶点数到了百万级比如有限元网格的外轮廓循环会有明显耗时下一章给出向量化版本。3. 工程化实现从 polyarea 到自定义形心函数 centroidPolygon3.1 最小实现circshift 向量化版本的函数标题里的“源码”核心资产通常就是一个能直接调用的函数。把上面的循环压成向量运算用 circshift 一次取得“下一个点”避免 for 循环和 mod 索引function [A, cx, cy] centroidPolygon(x, y) % centroidPolygon 计算闭合多边形的有向面积与形心 % 输入: x, y 为顶点坐标行向量或列向量均可 % 输出: A 带符号面积cx, cy 形心坐标 x x(:); y y(:); % 统一转成列向量 if x(1) x(end) y(1) y(end) % 如果首尾重复则去掉尾点 x(end) []; y(end) []; end xp circshift(x, -1); % 后移一位得到下一个顶点 yp circshift(y, -1); cross x .* yp - xp .* y; % 鞋带项逐点叉积 A 0.5 * sum(cross); cx sum((x xp) .* cross) / (6 * A); cy sum((y yp) .* cross) / (6 * A); end函数粗看和公式一一对应但有三处工程化处理第一输入统一转成列向量行向量列向量混用不会出错第二末尾自动容忍首尾重复的写法很多人从 CAD 里导出坐标时最后一个点和第一个点相同不处理会多算一项第三cross 是向量参与运算circshift 也在编译优化快的路径上百万顶点毫秒级算完。3.2 输入约定与三个常见坑这个函数对输入有一个隐含约定顶点必须按边界顺序单向排列不能跳点。下面是常见的三个坑及规避方式。坑一是顶点顺序乱。如果坐标列表是散点而不是有序轮廓鞋带公式计算面积时连线会来回穿过图形内部结果完全无意义。判定顺序是否正确的办法是检查有向面积符号逆时针为正时可用if A 0, x x(end:-1:1); y y(end:-1:1); end强制翻转正向的判定标准不一致时符号只影响面积正负不影响形心位置所以计算前统一逆时针更利于调试。坑二是自相交。图形像“8”字形交叉时鞋带公式会把内部交叠区域抵消形心也不知道偏向哪边。此时要先把轮廓拆成若干简单多边形分别调用函数再用面积加权组合形心。拆分的判断依据通常是交点检测MATLAB 里可用 polyshape 对象配合 intersect 检查常见做法是“如果对象是 polyshape直接使用其内置的 centroid 属性”但要确认不涉及需要排除的自相交环。坑三是输入里混入重复相邻点。两个完全相同的连续坐标不会影响面积却会多产生一个零长度的叉积项极端情况下放大浮点误差。简单去重可用idx [true; diff(x) ~ 0 | diff(y) ~ 0]; x x(idx); y y(idx);完成。3.3 带孔洞的组合区域多层多边形环单一闭合轮廓是简单情形不少工程问题是“外轮廓挖去几个孔洞”比如截面要算净面积和净形心。处理这类图形工程上用的是环分解外环按逆时针给出孔洞按顺时针给出所有环的有向面积均参与求和总面积 Σ有向面积形心 Σ(各环面积 × 各环形心) / 总面积。outerX [0 5 5 0]; % 外矩形逆时针 outerY [0 0 4 4]; innerX [1 3 3 1]; % 内孔顺时针 innerY [1 1 3 3]; [Ao, cxo, cyo] centroidPolygon(outerX, outerY); [Ai, cxi, cyi] centroidPolygon(innerX, innerY); A Ao Ai; % 内孔顺时针Ai 为负数 cx (cxo * Ao cxi * Ai) / A; % 面积加权组合 cy (cyo * Ao cyi * Ai) / A;这段逻辑的关键在符号内环按顺时针输入后centroidPolygon 返回的 Ai 是负值所以总面积用Ao Ai就是净面积加权时也自然把内孔的面积权重视为负。上例中外矩形面积 20内孔面积 4净面积 16输出形心 (2.625, 2)比外矩形几何中心略靠近孔洞反方向符合直觉。组合法对任意多个孔洞同样适用只要每个环都是简单多边形且嵌套方向满足“外逆内顺”。注意如果把内环也按逆时针输入Ai 变正值总面积会变成 24形心也会被错误地往孔洞中心拉偏。多环数据务必提前规范方向再进函数。4. 进阶图像掩膜与网格单元的形心计算4.1 像素级二值图像用 regionprops 求形心多边形能精确描述 CAD 轮廓但很多场景拿到的是图像掩膜比如摄像头拍到工件的二值图或医学图像分割结果。此时最优方案不是把边界点抠出来再调用鞋带公式而是直接使用图像处理工具箱里的 regionpropsbw imread(part_mask.png) 0; % 读入灰度图并二值化 stats regionprops(bw, Area, Centroid); pixelArea stats.Area; cxPixel stats.Centroid(1); % 像素坐标x 向右 cyPixel stats.Centroid(2); % 像素坐标y 向下regionprops 的原理是统计所有前景像素的坐标均值等价于对每个像素做单位权重累加像素级面积加权形心天然包含了内部填充区域不用先提轮廓。它有三个容易被忽略的细节坐标原点是左上角y 轴向下若要映射到世界坐标需要做y_world (cyPixel - 0.5) * scale之类的变换统计的是连通域而不是整幅图多目标时每个连通域一个结果输出面积单位是像素个数换算物理面积要乘像素尺度。和鞋带公式做个对比同样一个旋转矩形多边形的鞋带公式得到的是精确形状面积regionprops 得到的是栅格化后的像素面积两者差异就是像素分辨率引入的离散误差。分辨率越高越接近但计算量也越大。4.2 网格单元三角剖分与四边形细分的形心有限元、图像剖分、点云重建里形状是大量网格单元的组合。此时用鞋带公式对每个单元循环一圈也行但三角形和四边形有更轻量的写法。三角形的形心就是三个顶点坐标的平均值这是线性单元的特殊性质面积用叉积的一半求。四边形的形心不能简单取四个顶点均值需要沿一条对角线切成两个三角形再用三角形面积加权。tri [0 0; 2 0; 1 1.5]; % 三个顶点 ctrTri mean(tri, 1); % 三节点均值 areaTri 0.5 * abs(det([tri(2,:)-tri(1,:); tri(3,:)-tri(1,:)])); quad [0 0; 2 0; 2 1; 0 1]; % 四边形 [At1, ct1] triangleCentroid(quad(1,:), quad(2,:), quad(3,:)); [At2, ct2] triangleCentroid(quad(1,:), quad(3,:), quad(4,:)); ctrQuad (ct1 * At1 ct2 * At2) / (At1 At2);四边形均值形心和三角形加权结果不一致的典型场景是严重歪斜的四边形均值点偏向两个长边夹角较小的顶点而三角形剖分后的面积加权形心更贴近真实几何。如果你的网格全是质量较好的矩形或平行四边形均值法和三角形剖分结果一样网格畸变严重时必须用剖分法。对大规模网格这个操作可以批量并行化把所有四边形顶点组织成 n×4×2 的数组逐层调用向量化函数避免 for 循环。4.3 数值稳定性与量纲大坐标平移后计算鞋带公式的每一项包含 xi×y(i1) 和 x(i1)×yi 两个乘积。坐标值在几千以内时双精度毫无压力如果坐标来自 UTM 投影或工程控制网x、y 数值达到 1e5 甚至 1e6 量级乘积会到 1e10两个大数相减后有效数字被快速消耗面积和形心的相对误差可能达到不可接受的程度。这是鞋带公式的典型数值弱点。标准的缓解手法是坐标平移先取一个参考点通常取顶点均值或第一个顶点把坐标整体平移到原点附近计算完毕再把参考点加回形心。平移不改变面积形心平移量等于参考点坐标。x0 mean(x); y0 mean(y); [A, cxLocal, cyLocal] centroidPolygon(x - x0, y - y0); cx cxLocal x0; cy cyLocal y0;另一个量纲相关的问题是这套几何算法只适用于平面直角坐标系。如果坐标是经纬度直接代入会得到“经纬度单位下的面积”数值没有任何物理意义。遇到地理坐标数据常见做法是先投影到合适的平面坐标系再算面积和形心这和 GIS 里的椭球面积计算是两套口径不能混用。同样的单位不一致时一部分顶点是毫米、另一部分是米得到的结果是错乱的进入函数前先统一单位。5. 三种自检方法组合法、蒙特卡洛与栅格化对照5.1 组合法自检子区域面积加权应与整体一致这是最符合工程师习惯的一种验证。把多边形沿一条对角线拆成两个或多个子多边形各子多边形分别调用 centroidPolygon再用面积加权合成结果应等于不拆分时的整体结果。误差来源只有浮点累积通常不超过 1e-12 量级。以矩形 (0,0)、(2,0)、(2,1)、(0,1) 为例x [0 2 2 0]; y [0 0 1 1]; [A, cx, cy] centroidPolygon(x, y); % 整体结果 % 拆成两个三角形 tri1x [0 2 0]; tri1y [0 0 1]; tri2x [2 2 0]; tri2y [0 1 1]; [A1, cx1, cy1] centroidPolygon(tri1x, tri1y); [A2, cx2, cy2] centroidPolygon(tri2x, tri2y); % 面积加权合成 cxCheck (cx1 * A1 cx2 * A2) / (A1 A2); cyCheck (cy1 * A1 cy2 * A2) / (A1 A2);配合第三节里带孔洞的多环组合这个自检可以延展到任意复杂截面。凡是计算形心的代码第一次写完后都应该先做一次这样的拆合验证。5.2 蒙特卡洛自检不用公式直接验证结果蒙特卡洛是另一种独立的验证思路在包围盒内随机撒点统计落在多边形内部的比例来估计面积和质量分布。它不依赖鞋带公式能独立发现公式书写错误或顶点顺序问题。注意它的误差是统计性的验证对比时设定容差不要要求完全相等。rng(42); x [0 4 4 0 0]; y [0 0 3 3 0]; % 4x3 矩形 N 5e6; xmin min(x) - 0.5; xmax max(x) 0.5; ymin min(y) - 0.5; ymax max(y) 0.5; u xmin (xmax - xmin) * rand(N, 1); v ymin (ymax - ymin) * rand(N, 1); inside inpolygon(u, v, x, y); areaMC mean(inside) * (xmax - xmin) * (ymax - ymin); cxMC mean(u(inside)); cyMC mean(v(inside));inpolygon 是 MATLAB 自带的多边形包含判断内部点集合的均值就是形心的无偏估计。对矩形结果areaMC 与 12 的偏差通常在 0.1% 以内cxMC、cyMC 与 (2,1.5) 的偏差类似。这套方法还可以直接搬到任意封闭曲线只要离散成多边形调用 inpolygon 即可不需要解析表达式。5.3 栅格化对照回到像素口径验证最后一招是把多边形栅格化成二值图像再交给 regionprops 反算一次。逻辑上它是从解析解到离散解的过渡验证鞋带公式的宏观正确性[xg, yg] meshgrid(0:0.01:4, 0:0.01:3); bw inpolygon(xg, yg, x, y); stats regionprops(bw, Area, Centroid); areaGrid stats.Area * 0.01^2; % 每个像素对应的物理面积 cxGrid stats.Centroid(1) * 0.01 - 0.005; cyGrid stats.Centroid(2) * 0.01 - 0.005;栅格分辨率为 0.01 时面积误差来自边界像素的锯齿通常能到 0.5% 以内。把三种自检放到同一个脚本里打印出解析值、蒙特卡洛值、栅格化值三者互相印证需要批量计算时也可以把第 3 章的 centroidPolygon 存成独立 m 文件放进工程路径后续每次调用就退化成一行函数调用。本文还有配套的精品资源点击获取