资讯动态

基于Matlab图像处理的不规则颗粒粒径与形态测量系统

发布时间:2026/10/15 23:27:53 来源:尧图企业网站定制
1. 不规则颗粒测量为什么不能靠“目测”在岩土工程、粉体材料、制药工艺和食品分级这些行当里待过的人应该都有同一个感受手头那些颗粒样品几乎没有一个长成标准圆的。碎石、药片、砂粒、催化剂载体形状千奇百怪有的像土豆有的带棱角有的甚至带着明显的内凹。过去大家习惯了用筛分法过一遍标准筛得到一个网目区间或者拿游标卡尺手动量几个“最宽处”然后直接当成粒径去参与后续计算。但这些做法都有一个致命问题筛分只能告诉你颗粒能不能通过某个孔径卡尺只能量跨度它们都丢失了颗粒真实的二维形态信息。而工程上很多关键参数其实恰恰取决于这些形态信息。比如混凝土骨料的表面积直接关系到水泥浆用量粉体颗粒的周长和面积影响比表面积进而影响反应速率或溶解速率制药行业里颗粒圆形度与粉体流动性、压片质量强相关。这些指标的共同特征是它们依赖的不是“直径”而是颗粒在图像平面上的投影面积和投影周长。更麻烦的是不规则颗粒的粒径根本没有统一答案——同一个颗粒用面积等效、周长等效、最大费雷特直径、最小费雷特直径来定义会得到完全不同的“粒径”。所以在实际项目里我更喜欢用Matlab搭一套基于图像处理的颗粒测量系统一次性输出颗粒数量、投影面积、周长、等效粒径、圆形度等一组指标然后把数据导出成表格再交给后续的统计分析。这不只是“省人力”的问题而是把测量标准从“人眼主观判定”变成“可复现的算法逻辑”。同一批图像不管谁来跑只要参数固定结果就是一致的。这种可复现性在项目评审、论文支撑和工艺对比里都是硬需求。这套方案适合谁用说实话门槛不高。你只需要对Matlab的基本语法有一点了解哪怕只会写简单脚本也能照着下面这套流程把系统跑起来。它特别适合做材料检测、岩土分析、化工粉体研究的学生和工程技术人员也适合那些还在用“人工数数、人工量尺”的实验室想快速升级成数字化测量的人。2. 系统整体设计与测量参数换算逻辑2.1 从像素到物理尺寸先建立标尺很多人第一次做图像测量最容易犯的错误就是直接用图像里的像素数当作“面积”。一颗颗粒在图像里占了3000个像素这只是一个无量纲的数字跟真实的平方毫米、平方微米完全不搭边。要得到物理尺寸必须先完成像素标定。最简单的标定方法是在拍摄颗粒时同场景放入一个已知尺寸的参照物。常见操作是放一段直尺、一枚标准圆片或者工业上放的标定板。拍完之后我们在Matlab里量出参照物在图像中占据的像素长度就能算出一个比例系数其中k表示每个像素对应的物理尺寸单位是毫米每像素。比如直尺上10毫米在图像中对应500像素那么k0.02mm/px。后续所有面积、周长、粒径结果都要乘上这个系数的平方或一次方才能换算成真实物理单位。这里有一个很容易被忽略的细节标尺必须和颗粒放在同一个焦平面上并且镜头光轴要尽量垂直于载物面。如果标尺放得比颗粒高一点或者低一点都会造成比例系数的轻微偏差如果镜头带透视角度图像的边缘区域会有畸变标尺区域和颗粒区域的比例系数就不再一致。我在实际项目中一般要求拍摄时把样品铺平标尺尽量贴近样品层并且避免使用广角端拍摄。对于高精度的需求建议拍摄前用标定板做整幅图像的畸变校正但不做校正至少要保证标尺和样品同层。2.2 一组核心测量指标与计算公式测量系统输出的指标不是越多越好关键是每个指标都有明确的物理含义和适用范围。我在这套系统里设计了五个核心输出每个都对应一个具体参数或计算公式。第一是颗粒数量。这个直接对二值化图像做连通域标记有多少个独立的连通区域就是多少个颗粒。第二是投影面积。二值图中每个连通区域的像素总数乘以k的平方就得到物理面积。第三是周长。Matlab的regionprops函数里自带Perimeter属性它会计算区域边界的长度单位是像素乘上k就得到物理周长。第四是等效粒径。这个参数使用的是“面积等效圆”的概念把颗粒看成与实际投影面积相同的圆此时圆的直径就是等效粒径公式为。第五是圆形度也叫形状因子公式为。圆形度越接近1颗粒越接近圆形碎石颗粒的圆形度一般只有0.7左右而球形颗粒可以到0.95以上。这些指标之间是有内在关联的。等效粒径适合用在粉体粒度分布统计里圆形度适合用在流动性评估里。而周长和面积则是计算比表面积的基础数据。在实际报告中我通常会把颗粒数量、平均面积、平均等效粒径、平均圆形度这四项作为一页摘要再附上完整的数据表用于后续处理。有一个经验值是在像素分辨率不足的情况下周长的测量误差往往比面积大。因为周长对边界细节很敏感一个颗粒边界只要出现几个像素的锯齿波动周长就会显著偏大。所以我建议拍摄时尽量提高分辨率让单个颗粒至少占几百个像素以上否则周长数据不太可信。3. 图像前处理这一步决定了测量可信度不少人拿到图片第一反应就是直接二值化然后开始数颗粒。这样做运气好能得到结果运气不好整个测量就是错的。图像前处理是整个系统里最枯燥、却最决定上限的环节。我在迭代这套系统时前处理代码占了一半以上的调试时间。3.1 灰度化与光照校正显微照片或者工业相机拍回来的图片一般是彩色RGB图像。对测量来说颜色信息并不关键真正有用的是前景和背景的灰度差异。所以第一步是转灰度。这里有一个选择直接使用rgb2gray还是取单通道我的经验是别盲用默认函数。有时候颗粒带有明显色彩比如红色或黄色颗粒它在蓝色通道的灰度和背景非常接近但在红色通道却有很好的区分度。这时候只取某个颜色通道分割效果远好于直接转灰度。建议做法是先分别显示三个通道的灰度图人工看一眼选对比度最强的通道作为后续处理的基础。光照不均是实际拍摄里最常见的问题。中心亮、四周暗或者某一个角落有阴影都会导致同一张图里颗粒和背景的灰度关系不一致全局阈值根本没法用。解决办法有两个一是使用自适应阈值方法通过imbinarize的adaptive选项按局部区域计算阈值二是先做背景校正通过imopen提取背景光照分布再与原图相减。我自己更推荐先做背景校正因为它能把后续的处理拉回到“全局阈值可用”的简单状态。背景校正的典型做法是这样用imopen对灰度图做形态学开运算结构元素选一个比颗粒大很多的圆盘比如颗粒直径大概是20像素结构元素就选50到100像素半径开运算提取出的结果就是背景光照分布然后原图减去这个背景光照不均就被压平了。这一步之后图像看起来更像一个“均匀背景上的颗粒集合”后面不管用Otsu还是固定阈值效果都会稳定很多。3.2 二值化与孔洞填充二值化是测量系统的分水岭这一步做得不好后面全是白费。Matlab里最常用的全局阈值方法是Otsu方法函数是graythresh。但我建议不要直接把它当成万能钥匙。Otsu假设图像只有两类像素且类间方差最大可如果颗粒之间灰度差异大或者背景有噪声Otsu会倾向于把多数像素归到大类导致小颗粒被吃掉。有两个更稳妥的做法。一种是在背景校正后手动选择一个小范围阈值区间通过im2bw配合自己调整的阈值实时看二值图效果。这么做看起来“不自动化”但对一批固定拍摄条件下采集的图片效果非常稳。另一种是用imbinarize配合adaptive参数适用于光照不均又不想做背景校正的场景。自适应阈值会把噪声也放大所以一般要结合后面的去噪步骤一起用。二值化之后有一个操作必须做imfill填充孔洞。很多颗粒表面有纹理、凹坑或者反光造成的暗斑二值化后会形成内部空洞。如果不填充面积会偏小周长会因为边界凹凸而严重偏大。填充孔洞本身就是一种形态学修复它把连通的背景区域从颗粒内部补上实现方式很简单bw imfill(bw, holes);fill holes的处理只针对封闭区域颗粒表面的开孔不会受到影响所以安全性很高。做完这一步颗粒的基本形态就算还原出来了。3.3 去噪与边缘保留的取舍二值图里经常有一些零星的小噪声区域通常是几像素到几十像素的杂点可能是灰尘、反光点、传感器坏点。直接去掉它们最省事的方法是bwareaopen把小于指定像素阈值的连通域删掉。阈值怎么选我一般取颗粒最小投影面积的十分之一作为初值。比如预计最小颗粒占500像素那我就删掉所有小于50像素的连通域。这个值宁可设小一点也不能设大设大了会把真实的小颗粒也删掉。删除之后如果发现还有少量杂点可以用中值滤波或者形态学开运算再做一次轻量清理但要注意开运算的窗口大小不能超过最小颗粒半径的一半否则颗粒边缘会被“削”掉一圈。这一环节有个深刻教训边缘和噪声的取舍本质上是测量的“分辨率”。你为了去噪做得越狠颗粒的周长损失就越大。因为周长的定义依赖边界细节你去掉了一个边缘突起等于少算了一段边界。所以我现在的做法是去噪只删面积噪声不做边缘平滑。真正的边缘平滑留到测量之后如果发现周长曲线太毛糙再做形态学闭运算但结果要人工复核不能闭着眼交给算法。4. 颗粒分割与核心测量实现4.1 regionprops一次拿到所有形态参数前处理做完图像已经变成清晰的二值掩膜白色是颗粒黑色是背景。接下来就是测量系统的心脏部分。Matlab里regionprops可谓是测量利器一行代码可以同时获取面积、周长、质心、长短轴、外接矩形等几十个属性。先做连通域标记[labeled, n] bwlabel(bw);labeled是标签矩阵每个连通的白色区域有了独立的整数编号n就是区域总数。然后调用regionpropsstats regionprops(labeled, Area, Perimeter, Centroid, ... MajorAxisLength, MinorAxisLength, EquivDiameter);这里我要特别强调一下Perimeter属性的计算方式。Matlab默认的周长是基于边界链码计算的通俗说它会沿着区域边界的每个像素中心走一圈把邻接边的长度累加起来。这个算法对光滑目标很准确但对带锯齿边界的颗粒会高估周长因为锯齿边界的曼哈顿距离大于光滑弧线。如果对周长精度要求高可以把边界提取出来后自己用多边形长度公式重新算一遍或者用像元比例折减的经验修正。不过大多数工程应用里默认的Perimeter配合圆形度计算已经能反映形状差异了不必过度纠结。等效粒径可以直接用EquivDiameter属性它计算的就是面积等效圆的直径。但对于不规则颗粒我建议同时输出最大和最小费雷特直径也就是MajorAxisLength和MinorAxisLength属性。这两个值描述的是颗粒外接椭圆的长短轴能反映出颗粒的拉伸程度。假如某个样本的长短轴之比普遍大于1.5说明颗粒形状偏针片状这在骨料质量评估里是一个非常重要的不合格指标。4.2 粘连颗粒的分水岭分割实际拍摄中颗粒很少能一个个完美散开。尤其是粉体样品颗粒之间相互堆叠、紧挨着二值化后就会连成一个大块。如果不做粘连分割一个大块的面积会被当成一个颗粒计数严重偏低粒径分布完全失真。解决粘连颗粒的经典算法是分水岭分割。核心思想是把二值图当作地形图离背景越远的地方海拔越高然后模拟注水过程不同水源相遇的地方形成分水岭也就是颗粒的分界线。Matlab里的基本流程是这样D -bwdist(~bw); mask imextendedmin(D, 2); D imimposemin(D, mask); L watershed(D); bw_split bw; bw_split(L 0) 0;这个过程并不复杂但有两个前置条件必须满足。第一颗粒之间必须存在哪怕一个像素的窄连接分水岭才有发挥空间如果颗粒完全重叠成一个实心团分水岭也无能为力。第二内部标记点设置很关键。默认的watershed对粘连物体经常产生过度分割把一个颗粒切开成两半。为了抑制过度分割我们要通过imextendedmin先找到颗粒内部的“种子点”然后强制这些种子点成为局部极小值。这样分水岭只会在颗粒粘连的位点分割而不会把单颗粒内部也切开。过度分割是分水岭最常踩的坑。我见过不少初学者直接调用watershed然后发现结果惨不忍睹。后来我的经验是先用imerode对二值图做腐蚀腐蚀后的连通域数量大致相当于真实颗粒数量再把腐蚀后的区域作为标记进行分水岭。这个“先腐蚀、后分水岭”的组合拳对大部分粘连颗粒效果都很好。分水岭处理之后不要忘记对分割后的结果再跑一次regionprops。因为分水岭会在连通域之间切出新的边界重新标记后统计才是最终结果。4.3 像素折算与结果导出分割完成测量指标也都拿到了现在是换算和导出的环节。把所有像素单位的Area和Perimeter乘上k的平方和k把长度指标乘上k就能得到物理单位的结果。为了方便后续分析我会把结果存成一个表格结构每一行是一个颗粒每一列是面积、周长、等效粒径、圆形度、长短轴等指标然后用writetable函数导出成Excel或CSV。这一步最容易被忽视的是把面积和周长从像素值折算成物理值之前仔细检查比例系数k有没有写反。我吃过亏标尺量出来是2毫米对应100像素k0.02mm/px结果代码里粗心写成了20导致所有面积放大了100万倍整个数据直接没法用。所以后来我在代码开头定义常量时就写清楚mmPerPixel 0.02;接着所有计算都用mmPerPixel和mmPerPixel^2并且先跑一次“冒烟测试”——手动量一个颗粒确认输出值在合理范围再批量处理所有图片。导出后的数据不要直接拿去写报告。先用直方图看一眼分布比如等效粒径分布、圆形度分布如果发现某个颗粒的面积异常大很可能是有颗粒没分割开如果某个颗粒的面积异常小可能是噪声没删干净。这一步人工复核很关键算法不会告诉你它哪里算错了。5. 计数不是“数一下”那么简单5.1 标签矩阵与连通域计数当所有颗粒都成功分割计数其实就写在标签矩阵里了。bwlabel返回的n就是连通域数量这也是系统输出的颗粒总数。但要注意这个计数依赖一个前提分割完全正确。前处理阶段如果少删了一个噪声点计数就多一如果没分割开两个粘连颗粒计数就少一。所以更稳妥的做法是不直接用n这个数字而是从regionprops的结果里统计。例如设定面积范围阈值只在[minArea, maxArea]区间内的连通域才算颗粒。这个思路可以自动排除大块未分割区域和细小噪声valid stats([stats.Area] minArea [stats.Area] maxArea); count numel(valid);这个滤波器放在批量自动处理环节里极为有用因为不同批次的图像噪声水平不同一个固定的面积范围可以帮你剔除大部分干扰不让计数随噪声波动。5.2 边界处理的坑与人工复核方案图像边缘是计数的重灾区。一颗颗粒本来在视野内但有很大一部分被图像边框切掉了二值化后它仍然是一个区域会被当成一个完整的颗粒参与统计。这时候面积和周长都严重偏低等效粒径也无意义。对于这类边界颗粒我建议在编写测量代码时就做一个“边界抑制”把与图像边缘相接的连通域单独标记出来要么剔除要么额外备注。实现方法非常简单用bwlabel后再检测每个连通域是否与图像边缘像素有交集。如果是为了统计总颗粒数边界颗粒应该单独列出来不计入面积统计的平均值如果只是估算总量且采样量足够大也可以直接剔除。设置一个“边缘颗粒报告”的输出比单纯删掉更科学因为它在报告中保留了这部分信息。总之计数这个看似最简单的功能最后往往是最容易出错的地方。宁可多花十分钟做人工复核也不要直接信任代码输出的那个n。6. 常见问题与排查技巧实录6.1 颗粒边缘发虚导致面积偏大拍照时如果对焦不够锐利或者颗粒表面反光二值化后颗粒边缘会有一圈半透明的灰度过渡带。这一圈过渡带会被二值化纳入颗粒区域导致每个颗粒的面积普遍偏大几个百分点周长边缘额外增加更多。排查方法是看面积统计分布如果所有颗粒的面积都比手工测量的预期值大而且大出的比例比较稳定基本可以判定是边缘发虚。解决办法有两个方向一是提高拍摄对焦质量这是治本二是二值化时调高阈值让过渡带中较暗的部分被划到背景减小膨胀区域。另一个技巧是用imerode对整幅二值图腐蚀一圈然后测量再根据边缘过渡宽度加回一个平均补偿值。这种方法比较“粗”适合对绝对精度要求不高的场景。6.2 粘连严重时计数严重偏低样品如果完全没有分散颗粒抱成一团分水岭的改进能力也有限。这时候算法层面能做的已经到头了需要回到制样环节。我一般建议在拍摄前做分散处理例如用超声波分散、搅拌粉体或者用洗瓶冲散样品。另一个思路是降低样品浓度让颗粒在视野内的分布密度降低重叠概率自然下降。测量系统再厉害也补不了制样的坑。如果必须在算法层面硬扛可以尝试对二值图做距离变换后的局部极大值检测把“颗粒核”找出来然后基于这些核做分水岭。但这个方法对颗粒尺寸相近的情况效果好尺寸差异大时小颗粒的核容易被大颗粒的核压制测出的小颗粒数会偏低。6.3 标尺像素换算误差放大比例系数k的误差会等比放大到所有测量指标。面积是k的平方项所以差之毫厘谬以千里。我见过有人把标尺拍“斜”了本该平行的刻度在图像里呈一定角度量出来的像素长度偏短k偏小最后面积普遍偏小误差被平方放大。正确做法是量标尺的时候要在刻度线的最两端垂直量取尽量用像素插值的方式取小数位。如果标尺上10毫米对应496和504像素之间波动说明标尺放得不够正。拍摄时还要注意把标尺尽量放在画面中央附近避开镜头畸变大的区域。条件允许的话可以在整幅图像的不同位置多次量取同一个标尺算出平均值能在一定程度上校正畸变带来的误差。6.4 批量图片的自动化测量一批样品可能拍了几十张甚至上百张图手动一张张处理肯定不现实。所以系统一定要设计成批量脚本读入文件夹内所有图片逐张处理每张图片的结果以文件名标识存入结构体最后汇总成一个大表格。批量处理时要格外注意不同图片之间光照条件的微小差异。如果第一张图的二值化阈值能分割得很好不代表第二张也能。这时候要么先把所有图片统一做背景校正要么为每张图片单独计算自适应阈值。批量处理还有一个常见的坑某一张图片里出现了意外的异物比如拍摄时手影、镊子尖端它会直接污染最终的统计结果。我建议在汇总统计前先做一遍分布检查把所有颗粒的面积、周长画成散点图一眼就能看到异常群体然后锁定是哪张图片出了问题。7. 实测效果与个人心得这套系统在我自己经手的骨料粒径分析项目里跑过一段时间。样品是机制砂和碎石颗粒数量从几十颗到几百颗不等图像分辨率在1200万像素级别颗粒最小投影面积大约500像素。实测下来计数准确率在颗粒充分分散时可以达到95%以上主要误差来源是少数极小的边缘颗粒。面积和等效粒径的重复性测试同一批图片在相同参数下跑十次结果几乎完全一致区别只在浮点噪声级别。我踩过几次坑之后对这套系统最大的体会是算法不是瓶颈制样和图像采集才是。二值化和分割再精巧也救不了一张对焦模糊、光照不均的照片。反过来如果样品铺得平整、光照均匀、标尺放置正确哪怕只用一个简单的Otsu阈值加regionprops也能得到很漂亮的结果。测量系统的本质是把好的图像采集转化成数字。而图像采集的功夫需要花在拍摄台上。另一个值得分享的经验是不要把周长测量想得太完美。对于不规则颗粒的实际需求周长这个指标更适合做相对比较而不是绝对测量。因为周长对边界细节的敏感度太高两个形态几乎相同的颗粒只要图像分辨率不同周长就可能差出几个百分点。相比之下面积和等效粒径的稳定性要好很多更适合作为跨图像对比的核心指标。最后再分享一个小技巧整套流程里把每一步的中间结果都保存下来包括灰度图、背景校正图、二值图、分水岭标记图。不要怕占用空间这些中间图是排查问题的关键。当最终数据出现异常时你能一眼定位到是前处理丢了颗粒还是分割阶段过分割了。市面上的通用图像测量软件往往会把这些中间过程封装起来出了问题你无从下手自己用Matlab搭这套流程最大的优势恰恰就是每一步都能看得见、摸得着参数调整也完全可控。如果你正好需要在一批不规则颗粒样品上做面积、周长和计数的统计分析这套思路可以直接参考着上。

读完文章,也想定制专属网站?

尧图设计师 24 小时内与您沟通定制方案

免费获取报价 →
↑