Matlab图像处理实战：轮廓提取与凹多边形内角计算

简介面向Matlab图像处理学习者与研究人员这份压缩包提供了一套完整的轮廓几何分析工具可获取图形边缘轮廓并将轮廓n等分进而求出n边凹多边形各内角同时结合骨料计算公式得出棱角指数AI适用于骨料形态定量评价、颗粒形状分析等场景。包内共11个文件以6个m函数脚本为核心覆盖轮廓提取、等分插值、内角计算、指数求算等关键步骤另有3张png示例图用于效果对照2个txt说明文件辅助解读代码逻辑与运行方式整包仅70KB轻量易用。已有649人学习下载适合需要快速复现骨料棱角指数算法或研究凹多边形几何特征的读者。通过参考配套博客可清晰理解各函数调用关系与参数含义借助示例图直接运行验证各脚本拆分细致便于在此基础上扩展或迁移到其他图像处理任务。1. 从图像到凹多边形内角这事比想象中麻烦但值得做透做 matlab 图像处理的人多半都遇到过这个需求给一张零件的照片要自动把外形轮廓提出来然后按弧长等分成若干段最后把每个折角的内角量出来。听起来像常规操作但一旦轮廓里有凹部事情就变味了——用常见的多边角测量方法凹进去那个角算出来往往只有几十度和实际完全对不上。这份资源的核心就是把「轮廓提取 → 弧长 n 等分 → 凹多边形内角求解」整条链路打通连凹凸点判断和角度修正这种细节都给了可跑通的代码。适合做视觉测量、形状分析或 CAD 逆向的工程师也适合被凹角问题卡住的学生。下面按我自己的拆解顺序来复现。2. 轮廓提取与坐标预处理从像素矩阵到有序边界点2.1 提取轮廓为什么用 bwboundaries 而不是找边缘再加连通域拿到一张图形第一步是把边缘轮廓变成一组有序的坐标点。常见做法是边缘检测后找连通域但这一步有个隐患边缘检测出来的点是二值化的细线可能断裂或不连续后续做等分时会出现莫名其妙的跳变。我一般直接对二值图用 bwboundaries它返回的是按连通域边界顺序排列的像素坐标天然有序省掉自己排序的麻烦。% 读取图像并转为二值图 img imread(part.png); gray rgb2gray(img); % 如果本身就是灰度图可跳过 bw imbinarize(gray, adaptive); % 自适应阈值应对光照不均 % 提取边界坐标noholes 表示只取外轮廓不要内部孔洞 B bwboundaries(bw, noholes); boundary B{1}; % 取第一个连通域的边界多数场景下直接读图、二值化、取边界就够用了。注意bwboundaries返回的坐标是[row, col]也就是[y, x]这和常规的笛卡尔坐标系相反如果后面要做向量计算第一步就得先统一为[x, y]。另外boundary是闭合的首尾是同一个点等分前要把尾点去掉。2.2 坐标预处理去重复、统一方向、转成双精度边界坐标直接来自像素网格采样步长不均匀——斜线处相邻点间距是 1.414横平竖直处是 1。这本身不影响角度计算但会影响弧长等分的精度所以要先做一次轻度平滑和去重。像素边界还有一个特点可能存在相邻点完全重复某些边界提取算法会产生不处理的话累计弧长会出现零长度段插值的时候会出 NaN。% 转换为双精度并统一为 [x, y] 顺序 boundary double(boundary(:, [2 1])); % 去重删除连续重复点 diff_idx find(any(diff(boundary), 2)); boundary boundary([1; diff_idx1], :); % 闭合处理如果首尾相同则去掉末尾点 if isequal(boundary(1,:), boundary(end,:)) boundary(end, :) []; end % 轻度平滑用滑动平均降低像素锯齿窗口取 5 k 5; boundary_smooth smoothdata(boundary, movmean, k);这里smoothdata的窗口参数值得说一句窗口太小锯齿还在窗口太大拐角会被削圆凹进去的角会被磨平导致后续算出来的内角偏大。对于一般图片轮廓我试下来窗口 5 到 7 比较稳妥具体可以在调试时画出来对比一下。等分前把坐标统一成双精度能避免后面插值时整数矩阵带来的精度丢失。3. 按弧长把轮廓 n 等分核心算法与实现细节3.1 为什么不能直接按像素索引等分很多人拿到轮廓点序列后直接boundary(1:k:end)来取等分点这个做法在轮廓接近圆形、点密度均匀时勉强能用但对真实图形基本翻车。原因在于边界点的密度不均匀斜率越大的位置相邻像素距离越远斜着走一步是 1.414平直区域距离近一步是 1凹进去的曲线部分密度更是忽高忽低。按索引等分等于把物理长度不等的段当成了等长段后续内角会系统性偏差。正确的做法是先把轮廓看成一条参数曲线计算每个点处的累计弧长再在累计弧长上做等间距采样。这个思路和 CAD 里曲线的等分完全一致只是落到 MATLAB 上需要手写一小段逻辑。% 计算相邻点之间的欧氏距离 seg_len sqrt(sum(diff(boundary_smooth).^2, 2)); % 累计弧长 cum_len [0; cumsum(seg_len)]; % 总弧长 total_len cum_len(end); % 生成 n 个等距弧长位置含首尾n 段对应 n1 个点 n 8; target_pos linspace(0, total_len, n1); % 在累计弧长上做插值重新采样得到等分点坐标 point_x interp1(cum_len, boundary_smooth(:,1), target_pos, spline); point_y interp1(cum_len, boundary_smooth(:,2), target_pos, spline); equipoints [point_x, point_y];这段代码的逻辑一句话就能讲清先算每段微元长度累加得到每个原始点的弧长位置然后构造等间距的目标弧长序列用插值函数反查坐标。interp1前两个参数分别是「原采样位置」和「原采样值」第三个参数是「目标位置」。这里务必使用spline或pchip不要用默认的linear。线性插值在轮廓弯曲处会把角切成直线段凹角直接变没。spline光滑但容易在突变处出现过冲pchip保形更好我通常优先pchip视觉上更接近原始轮廓。3.2 n 取多少合适分辨率与稳定性的平衡等分数 n 是这个流程里最敏感的参数。n 太小轮廓的弯曲特征被强行拉直凹角位置偏移n 太大等分点密集后对噪声和锯齿敏感个别点落在轮廓的微小波动上角度计算结果抖动剧烈。经验法则是先看轮廓的总像素周长n 取周长的 1/50 到 1/80比如周长 400 像素n 取 6 到 8。这里给的n 8适合中等复杂度轮廓。还有一点容易被忽略插值后等分点的首尾要闭合。linspace生成的序列天然包含起点和终点而终点恰好是起点本身因为轮廓闭合所以计算内角时最后一个角由equipoints(end)和equipoints(1)共同构成不要丢掉这个角。很多新手算完 n 个点只算了 n-1 个角就是没意识到闭合链要多算一次。4. 凹多边形内角计算向量夹角法与凹凸点判定4.1 基础内角三个相邻等分点的向量夹角拿到 n1 个等分点后内角计算就可以开始了。基本思路是对每个顶点P_i用它的前一个点P_{i-1}和后一个点P_{i1}构造两个向量然后求夹角。夹角公式是atan2加余弦定理不能用简单的acos(dot(...))后者数值不稳定角度接近 0 度或 180 度时误差会放大。% 计算多边形各顶点的内角返回角度制0~360 function angles polygon_interior_angles(pts) m size(pts, 1); angles zeros(m, 1); for i 1:m i_prev mod(i-2, m) 1; % 前一个点环形索引 i_next mod(i, m) 1; % 后一个点 v1 pts(i_prev, :) - pts(i, :); v2 pts(i_next, :) - pts(i, :); % 用 atan2 求 v1 到 v2 的有向夹角范围 (-pi, pi] ang atan2(det([v1; v2]), dot(v1, v2)); % 转换为 0 ~ 360 度 if ang 0 ang ang 2*pi; end angles(i) ang * 180 / pi; end end这里关键是理解有向夹角。atan2(det, dot)求得的是 v1 绕顶点旋转到 v2 的有向角度范围是 -180 到 180 度。对于凸多边形每个顶点的有向角方向一致全正或全负修正后直接就是内角而凹多边形内部会出现反方向顶点。上面的代码把所有负角都加 360 度得到的是「转向角」在凹点处会给出小于 180 度的角——这本身没错但它不是内角是外角的补关系需要结合凹凸性做修正。4.2 凹凸点判定与内角修正凹多边形和凸多边形的核心区别凸多边形所有顶点转向方向一致凹多边形存在个别顶点转向相反。上一节的代码如果直接用凹点的内角会被算成小角而实际凹进去的内角是大于 180 度的。判断凹凸点的方法很多常用的是多边形有向面积符号法或局部叉积法。局部叉积的做法是对顶点P_i计算边向量P_i-P_{i-1}和P_{i1}-P_i的叉积符号正负代表在闭合路径上该点是顺时针方向还是逆时针方向的转向。如果整个轮廓统一按逆时针或顺时针排列那么凸点的叉积符号相同凹点的符号相反。% 统一轮廓方向为逆时针通过有向面积判断并翻转 area sum(boundary_smooth(:,1) .* ... boundary_smooth([2:end 1],2) - ... boundary_smooth([2:end 1],1) .* boundary_smooth(:,2)) / 2; if area 0 boundary_smooth flipud(boundary_smooth); end % 等分点同样按逆时针顺序排序 equipoints sort_points_ccw(equipoints); % 计算内角修正凹点 angles zeros(n1, 1); for i 1:n1 i_prev mod(i-2, n1) 1; i_next mod(i, n1) 1; v1 equipoints(i_prev, :) - equipoints(i, :); v2 equipoints(i_next, :) - equipoints(i, :); ang atan2(det([v1; v2]), dot(v1, v2)) * 180 / pi; % 逆时针轮廓下叉积符号为负说明是凹点 cross_val v1(1)*v2(2) - v1(2)*v2(1); if cross_val 0 ang 360 - abs(ang); % 凹点的内角修正 else ang abs(ang); end angles(i) ang; end修正逻辑背后是几何直觉逆时针方向走轮廓时凸点的转向角是正的右转小于 180 度凹点的转向角是负的左转其内角应为 360 度减去转向角的绝对值。sort_points_ccw是对等分点按极角排序的辅助函数确保计算顺序和轮廓走向一致。这个修正策略在一个前提下才能用轮廓方向统一。如果原图轮廓是顺时针所有符号都会反过来凹点会被误判成凸点内角整体变成「镜像」。所以前面用有向面积做方向统一的步骤不能省。5. 避坑与常见问题这条链路里最容易翻车的四个细节5.1 现象计算出来的凹角全部小于 180 度凸角反而大于 180 度原因轮廓方向没有统一叉积符号判断与预期相反凹凸点完全对调。解决在等分前检查有向面积强制把轮廓和等分点都调整为逆时针如果图像坐标系 Y 轴向下图像常见坐标系做笛卡尔转换时记得翻转 Y 轴否则方向判断也会反。5.2 现象n 等分点落在轮廓附近但不在轮廓上画出来是锯齿折线原因插值方式用了默认的线性或者等分数太大spline 过冲把点推到轮廓外。解决先画出等分点和原始边界对比检查是否贴合优先改用pchip插值它对尖锐转角的过冲控制更好同时把等分点数量降到周长 1/60 左右再试。5.3 现象角度结果与预期相差 90 度左右而且方向性混乱原因图片坐标系下行列坐标未转换为 XY 笛卡尔坐标导致向量运算时把行当成了 y 但方向向下未翻转。解决bwboundaries返回的[row, col]在转[x, y]时要把 row 取负后才符合数学坐标系或者全程都在图像坐标系里计算只要向量运算前后统一角度结果是可以保持不变的但方向判断必须始终一致。5.4 现象同一个轮廓换一张图片角度结果波动很大原因分割阈值或图像尺寸不同导致边界位置偏移等分点对轻微变化过度敏感。解决先对二值图做形态学闭运算填补边界缺口再做边界平滑最后输出角度前打印等分点坐标肉眼核对凹点位置是否命中真实凹部定位到个别等分点漂移时可以微调插值方法而不是整体改参数。6. 进一步验证用合成凹多边形校验整套算法实际项目里最怕的是「算完了不知道对不对」。我的习惯是用合成多边形做 ground truth 验证人为构造一个顶点坐标已知的凹多边形跑完整流程看程序输出的角度能否回推出原有的顶点角度。如果吻合说明整条链路可靠再切换到真实图像。% 构造一个已知的凹五边形含一个凹点 theta linspace(0, 2*pi, 6); r 100 * ones(6, 1); r(4) 40; % 第四个点径向缩小形成一个凹部 xs r .* cos(theta); ys r .* sin(theta); % 在轮廓上密集采样模拟图像的边缘点 t_interp linspace(0, 1, 1000); cx interp1(linspace(0,1,6), xs, t_interp, spline); cy interp1(linspace(0,1,6), ys, t_interp, spline); boundary_smooth [cx cy]; % 然后走第 3、4 章相同的等分与角度计算流程r(4) 40就是在构造一个可预期的凹点。跑完对比输出的角度序列凹点那个角应该接近 240 度多一点其余凸点接近标准五边形的内角。这里提一个我个人踩过坑后的固定动作每次跑完等分我都会额外输出一张等分点和原轮廓叠加的图等分点有「悬空」立刻能发现从那以后所有角度结果我都会先做一次方向性检查再进数据表这已经成为惯性。希望你也能把这套验证流程写进自己的代码里能省下不少调试的时间——希望帮到你。本文还有配套的精品资源点击获取