简介针对生物识别、医学诊断与人机交互中的瞳孔定位任务这份基于MATLAB的开发包将阈值分割与灰度分布特性相结合为图像处理学习者和算法研究人员提供了一套可复现的参考实现。压缩包共5个文件包含4个.m源码脚本与1张示例眼部图像整体仅563KB内容覆盖灰度直方图统计、自动阈值选取、图像二值化、形态学腐蚀膨胀及瞳孔边界轮廓提取等关键步骤。已有447人浏览学习适合需要快速完成课程设计或启动相关课题实验的读者。代码中通过定位直方图首峰之后的低谷完成初始分割并支持Otsu与Isodata等自适应方法优化阈值同时借助形态学开闭运算去除眼睫毛与光斑干扰利用连通域分析提取平滑的瞳孔边界最终输出精确的瞳孔中心坐标与半径。配合示例图像可完整观察从灰度分布分析到瞳孔轮廓提取的全流程且模块化设计便于将分割思路迁移至虹膜识别、视线估计、疲劳驾驶监测等下游应用。1. 为什么瞳孔定位要先看灰度分布而不是直接上边缘检测在虹膜识别、眼动追踪和驾驶员疲劳监测这类应用里瞳孔定位的精度直接决定后端特征提取的成败。一上来就用 Canny 边缘检测加 Hough 圆拟合的做法很常见可惜在近红外光照下角膜反光点、睫毛投影和虹膜纹理都会产生强边缘响应圆拟合的结果经常把上下眼睑当成候选圆。反过来想瞳孔是图像里灰度最低的连通区域它的灰度值范围和虹膜、皮肤有明显分离先通过灰度直方图摸清分布特性再用阈值分割把暗区切出来最后用质心法定位这条链路在光照波动、反光干扰下反而更稳。这篇文章会从灰度分布分析出发逐步拆解阈值分割和瞳孔质心定位在 MATLAB 里的完整落地流程并给出每一步的参数选型和踩坑点。2. 瞳孔图像预处理与Matlab灰度直方图分析2.1 图像灰度化与感兴趣区域裁剪原始的眼部图像可能是 RGB 彩色图也可能是近红外相机直接输出的 8 位灰度图。最稳妥的起点是统一转成灰度单通道再做中值滤波。中值滤波在去除椒盐噪声的同时能保住边缘轮廓对瞳孔这种低灰度区域的干扰比均值滤波小很多。处理流程里还有一个容易被忽略的环节——ROI 裁剪。整幅人脸图像里瞳孔只占很小面积直方图统计时皮肤和头发的高灰度像素会把暗部峰压扁导致阈值选取偏离。裁出包含眼球的矩形区域让灰度分布集中在瞳孔和虹膜附近后续所有处理都在这个小范围内完成。% 读取图像并灰度化 I_raw imread(eye_ir.jpg); if size(I_raw, 3) 3 I_gray rgb2gray(I_raw); else I_gray I_raw; end % 3x3 中值滤波去掉近红外图像中的孤立噪点 I_roi0 medfilt2(I_gray, [3 3]); % 裁剪 ROI矩形格式为 [x y width height] I_roi imcrop(I_roi0, [60 80 240 180]);中值滤波窗口选 3x3 而不是更大的 5x5原因在于瞳孔边缘的梯度信息直接关系到后续质心定位精度窗口过大会让角膜缘的过渡带被磨平瞳孔面积被低估。ROI 坐标不推荐手写固定先对全图做一个粗略的暗区检测取最大连通域的包围盒再放大 1.2 倍是标定 ROI 的通用做法。2.2 用 imhist 观察瞳孔与背景的灰度分布特性灰度分布是整个阈值分割流程的地基。把 ROI 的灰度直方图画出来后能看到两类典型形态一类是双峰分布低灰度区的瞳孔峰和高灰度区的虹膜皮肤峰之间有明显谷底说明阈值很好选另一类是单峰带长尾分布瞳孔峰和睫毛阴影叠在一起谷底模糊说明照明或反光有问题。观察时不能只看形状还需要记录峰的灰度位置和像素占比用数值说话。% 统计 ROI 的灰度直方图压缩到 255 级 [counts, grayLevels] imhist(I_roi, 255); % 定位暗部峰只观察 0~51 灰度区间 dark_bin floor(255 * 0.2); dark_counts counts(1:dark_bin); [max_dark_count, dark_peak_pos] max(dark_counts); fprintf(瞳孔峰灰度级: %d像素数: %d\n, ... grayLevels(dark_peak_pos), max_dark_count); % 计算 ROI 内暗像素占比用于判断瞳孔面积是否合理 dark_ratio sum(dark_counts) / numel(I_roi); fprintf(暗区像素占比: %.2f\n, dark_ratio);暗部峰位置反映照明强度峰越靠左说明瞳孔区域越暗红外补光充足峰靠近 60 以上说明环境光干扰明显后续阈值需要放宽。暗区占比则能辅助判断 ROI 是否选得合适正常瞳孔占比在 5%~20%如果远低于这个区间说明 ROI 里包含太多眼白和皮肤应该重新裁剪。2.3 自动确定阈值区间Otsu 与双峰法取舍确定阈值有两条路线。Otsu 法通过最大化类间方差自动算出全局阈值瞳孔与背景对比强烈时结果非常准确几乎不用人工干预。但当瞳孔内部出现大面积反光点反光块的亮度接近虹膜Otsu 会把阈值往上抬导致瞳孔边缘被吞掉。双峰法直接在直方图上找两个峰之间的谷底可解释性强但谷底位置受睫毛和噪声影响大经常需要平滑处理。判断维度Otsu 法双峰法适用光照均匀环形补光光照不均匀或存在大面积反光反光敏感度阈值被高光拉高谷底被噪声填平人工干预程度几乎免调需要根据峰位置微调计算开销单次直方图遍历平滑加谷底搜索把两种方法结合使用是工程上更常见的做法。% Otsu 自动阈值 level_otsu graythresh(I_roi); T_otsu round(level_otsu * 255); fprintf(Otsu 阈值: %d\n, T_otsu); % 双峰法先高斯平滑直方图再找谷底 smoothed_counts smoothdata(counts, gaussian, 5); [~, dark_peak] max(smoothed_counts(1:30)); [~, bright_peak] max(smoothed_counts(80:200)); bright_peak bright_peak 79; valley_segment smoothed_counts(dark_peak:bright_peak); [~, valley_offset] min(valley_segment); T_valley dark_peak valley_offset - 1; fprintf(双峰法阈值: %d\n, T_valley);平滑窗口选 5 而不是更大的 15是为了保留谷底的细节。验算逻辑是如果 T_otsu 与 T_valley 相差不超过 15 个灰度级直接用 Otsu 的结果因为两者互相印证偏差超过 20 时大概率是反光把 Otsu 拉高取两者中的较小值更利于保住瞳孔边缘。平滑函数 smoothdata 在 MATLAB R2017a 之后可用旧版本用户可以用 filter 或 sgolayfilt 代替。3. 在Matlab中用阈值分割分离瞳孔区域的完整实现3.1 全局阈值与自适应阈值选型拿到阈值后进入分割阶段。全局阈值速度快、可解释性强适合照明均匀的近红外条件自适应阈值则能处理光照渐变但计算开销大还会把虹膜纹理放大成噪声。选型时先看直方图形态双峰清晰用全局谷底模糊且暗部峰宽平优先考虑自适应。值得注意的是瞳孔定位场景里的自适应阈值必须设置暗前景选项否则默认模式下会分割出反光点而不是瞳孔本身。3.2 imbinarize 与手动阈值的参数设定MATLAB 从 R2016b 开始主推 imbinarize功能上覆盖了旧接口 im2bw 和 graythresh 的组合。全局模式直接传阈值归一化值自适应模式下需要重点关注 Sensitivity 参数。这个参数控制像素判为前景的灵敏度数值越大越容易把睫毛和阴影也判成瞳孔区域。瞳孔定位的经验区间是 0.80~0.85超过 0.9 几乎必然引入噪声。% 全局阈值分割前景设置为暗区 BW_global imbinarize(I_roi, T_otsu / 255); % 自适应阈值分割ForegroundPolarity 设为 dark BW_adapt imbinarize(I_roi, adaptive, ... ForegroundPolarity, dark, ... Sensitivity, 0.82); figure(Name, 阈值分割对比); subplot(1, 2, 1); imshow(BW_global); title(全局阈值, FontSize, 10); subplot(1, 2, 2); imshow(BW_adapt); title(自适应阈值, FontSize, 10);这里的 ForegroundPolarity 是瞳孔分割最容易踩的坑默认值是 bright 表示亮区作为前景用默认值运行会得到一团白色区域。Sensitivity 的调节参考像素占比分割结果中对暗区占比做统计如果得到的瞳孔面积明显大于暗部峰像素数说明灵敏度偏高需要往下调。3.3 形态学处理去掉睫毛和反光点二值分割结果里通常伴随两类杂质。睫毛形成细长前景区域与瞳孔相连或紧贴边缘反光点则是在瞳孔内部形成小孔洞或碎块。形态学开运算先腐蚀后膨胀能把细长结构断开去掉imfill 再补充填充孔洞。结构元素的选择决定这类操作的效果圆盘半径 5 适合瞳孔直径在 100 像素以上的图像小图需要相应缩小。% 开运算去除睫毛细线 se strel(disk, 5); BW_open imopen(BW_global, se); % 填充瞳孔内部的反射孔洞 BW_holes imfill(BW_open, holes); % 保留面积最大的连通域去掉眼睑边缘残留 BW_final bwareafilt(BW_holes, 1); figure(Name, 形态学后处理); subplot(1, 3, 1); imshow(BW_global); title(分割原始结果, FontSize, 10); subplot(1, 3, 2); imshow(BW_open); title(开运算去睫毛, FontSize, 10); subplot(1, 3, 3); imshow(BW_final); title(填孔最大连通域, FontSize, 10);用圆形结构元素对任意角度的睫毛都能有效处理。但如果睫毛方向高度统一比如垂直向下延伸改用 strel(line, 10, 90) 的效果会更好。bwareafilt 保留最大连通域的操作隐含了一个前提——假设瞳孔是图像中面积最大的暗连通域如果眼睑闭合导致大面积暗色区域出现这个假设就会失效需要在连通域筛选时加上形状校验。4. 瞳孔质心定位与灰度重心法校正4.1 连通域标记与瞳孔区域提取二值分割完成后先把所有连通域标记出来再用区域属性筛选目标。这里的筛选条件不只看面积还要看形状指标。瞳孔形状接近圆而睫毛残段通常是细长的通过偏心率属性就能有效区分。regionprops 一次调用能拿到所有候选区域的几何参数筛选过程用表格列出更清晰。% 标记连通域 labeled bwlabel(BW_final); % 提取区域属性包括面积、质心、偏心率 stats regionprops(labeled, Area, Centroid, Eccentricity); % 面积排序取最大的候选域做圆度校验 [~, idx] sort([stats.Area], descend); for k 1:min(3, numel(stats)) if stats(idx(k)).Eccentricity 0.6 pupil_centroid stats(idx(k)).Centroid; fprintf(候选瞳孔中心: (%.1f, %.1f)面积: %d\n, ... pupil_centroid(1), pupil_centroid(2), stats(idx(k)).Area); break; end end偏心率小于 0.6 意味着该连通域接近椭圆这是瞳孔的基本形状约束。循环从面积最大的候选往下查最先满足条件的即为目标。area 过小的候选区域可能是反光残留离心率过大的可能是睫毛或眼睑边缘这两个条件叠加后误检率下降很多。4.2 质心计算regionprops 与灰度加权对比regionprops 返回的质心是二值区域的几何中心不包含像素亮度信息。当瞳孔内部存在反光点时几何质心依然位于瞳孔中央但边缘如果因为阈值选取不当被不对称切掉一块几何质心就会向切掉的相反方向偏移。灰度重心法把每个像素的灰度值当作权重对瞳孔定位有两个方向的影响暗区权重小反光点权重大如果反光点位于瞳孔边缘质心会被明显拉偏。% 灰度重心法以灰度值的倒数作为权重暗像素权重更大 I_double double(I_roi); % 防止除零加一个小常量 inv_gray 1 ./ (I_double 1); % 只在分割区域内累加 mask double(BW_final); sum_w sum(sum(inv_gray .* mask)); weighted_x sum(sum(inv_gray .* mask .* (1:size(I_roi, 2)))) / sum_w; weighted_y sum(sum(inv_gray .* mask .* (1:size(I_roi, 1)))) / sum_w; fprintf(灰度重心: (%.1f, %.1f)\n, weighted_x, weighted_y);对瞳孔定位而言反光点 q 是误差源因为反光的亮度接近虹膜在二值图上被排除掉但灰度重心里反光点保留了权重把质心拉向反光一侧。更稳妥的组合是二值质心负责粗定位灰度重心作为校验两者距离超过 3 像素时回到分割阶段检查阈值和形态学参数。中心偏移的另一个原因是不对称的角膜反光这种情况下闭运算处理比调整质心算法更有效。4.3 光照不均下的局部阈值补救全局阈值跑通流程后遇到光照不均的数据会暴露问题。瞳孔一侧边缘被切开二值质心偏移。此时回到分割阶段做局部处理先用 imtophat 顶帽变换校正背景亮度再做全局阈值。顶帽变换能移除不均匀照明造成的低频背景让瞳孔区域和边缘的亮度梯度恢复均匀。% 顶帽变换校正不均匀光照结构元素大小约等于瞳孔直径 se_light strel(disk, 30); I_corrected imtophat(I_roi, se_light); % 校正后重新做全局阈值 level_corr graythresh(I_corrected); BW_corr imbinarize(I_corrected, level_corr); BW_corr imfill(BW_corr, holes); BW_corr bwareafilt(BW_corr, 1); figure(Name, 光照校正); subplot(1, 2, 1); imshow(I_roi); title(原始 ROI, FontSize, 10); subplot(1, 2, 2); imshow(BW_corr); title(顶帽变换后分割, FontSize, 10);顶帽变换结构元素的直径必须大于瞳孔直径否则会把瞳孔本身当作背景减掉。30 像素圆盘适用瞳孔直径在 20~25 像素的图像。顶帽变换后的直方图通常呈现更尖锐的双峰用它辅助设计阈值不仅解决本张图像的问题还能在批量处理时提高自动阈值的稳定性。5. 一个实战技巧用灰度分布曲线验证分割阈值是否合理分割和定位做完后需要一套验证方法确认结果可复用而不是手动调参数恰好拟合某一张图。常见的验证手段是人工检查每张图的分割掩膜但在批量数据上逐张看效率太低。我常用一个基于灰度分布的客观校验法把分割得到的瞳孔掩膜映射回原始 ROI统计掩膜内外的灰度分布看两组的分离程度。% 掩膜内外灰度分布对比 pupil_pixels I_roi(BW_final); background_pixels I_roi(~BW_final); % 计算分离度指标类间方差比 var_between (mean(pupil_pixels) - mean(background_pixels))^2; var_within var(double(pupil_pixels)) var(double(background_pixels)); separation_ratio var_between / (var_within eps); fprintf(灰度分离度: %.2f\n, separation_ratio); fprintf(瞳孔平均灰度: %.1f背景平均灰度: %.1f\n, ... mean(pupil_pixels), mean(background_pixels));分离度指标大于 2 说明两组灰度区间明显拉开阈值选得合理小于 1.5 则需要检查阈值是否偏移或 ROI 是否包含了过多干扰区域。这个指标的优点是全程量化不依赖主观视觉判断可以用在自动化的批量质量筛选里。另一个有用的校验是观察瞳孔平均灰度与暗部峰的偏差正常情况下瞳孔平均灰度应在暗部峰灰度级附近偏差超过 15 说明分割混入了虹膜像素。把这两条校验写入批处理脚本配合阈值微调能让分割模块在多组数据集上稳定复用。本文还有配套的精品资源点击获取