MATLAB图像分割实战：从基础阈值到分水岭算法的全流程解析

1. 图像分割入门：为什么从阈值开始？

大家好，我是老张，在图像处理这个行当里摸爬滚打了十几年，用MATLAB处理过的图像少说也有几十万张了。今天想和大家聊聊一个非常经典，但又常被新手朋友觉得“太简单”而忽略的话题：图像分割。尤其是从最基础的阈值分割，一路讲到听起来有点“高大上”的分水岭算法。

很多朋友刚接触图像分割，可能就直接奔着深度学习、语义分割去了，结果发现效果时好时坏，调参调到怀疑人生。其实，很多工业质检、医学影像分析、甚至是我们手机里的简单抠图功能，底层逻辑都离不开这些传统但极其有效的算法。它们就像盖房子的地基，地基打牢了，上层建筑才稳。阈值分割，就是这个地基里最核心的一块砖。

那么，图像分割到底是什么呢？用大白话讲，就是把一张图片里我们感兴趣的部分（比如照片中的人、CT片里的病灶、生产线上的瑕疵）给“抠”出来。而阈值分割，就是其中最直观的一种“抠图”方法：设定一个门槛值（阈值），比这个门槛亮的像素算一类，比它暗的算另一类。听起来是不是简单得有点过分？但恰恰是这种简单，让它成为了解决大量实际问题的第一选择。我处理过很多产线上的零件图像，光照均匀、背景干净，一个简单的阈值分割，效果又快又好，根本没必要上复杂的模型。

接下来，我会手把手带你走完这个完整的流程：从最原始、需要你亲自“拍板”的人工选择阈值，到让计算机自己找最佳阈值的OTSU和自适应阈值，最后再到能处理物体“粘在一起”这种复杂情况的分水岭算法。我会分享我踩过的坑、调试的小技巧，以及如何根据你的图像特点选择最合适的方法。保证你跟着做一遍，就能把这些方法应用到自己的项目里。

2. 人工选择阈值：你的第一把分割“手术刀”

2.1 原理与直觉：用眼睛“看”出分界线

人工选择阈值，这可以说是图像分割的“原点”。它的核心思想完全依赖于我们人的视觉判断。我们拿到一张图，比如一张在白色背景下拍摄的黑色零件图，一眼就能看出零件和背景的界限在哪里。这个界限，反映在灰度值上，就是一个具体的数字。

这个过程本质上是一种二值化：你把图像中每一个像素点的灰度值，跟你自己设定的那个阈值（比如125）进行比较。大于125的，我全给你变成白色（255，代表前景或背景）；小于等于125的，我全给你变成黑色（0）。这样，一张彩色或灰度的图片，就变成了一幅只有黑和白的“剪影”。这个剪影，就是我们分割出来的目标。

我刚开始做项目的时候，特别喜欢用这个方法，因为它快，而且可控。当图像质量非常理想，前景和背景的灰度对比非常鲜明时，你几乎可以瞬间找到一个完美的阈值。但它的缺点也同样明显：主观且脆弱。说它主观，是因为你觉得合适的阈值，换一个人来看可能就觉得差一点；说它脆弱，是因为一旦光照条件稍有变化，或者背景变得复杂，你之前精心调好的那个“神奇数字”就立刻失效了。我记得有一次做一个瓶盖检测项目，上午调好的阈值，下午因为窗户阳光角度的变化，分割结果就一塌糊涂，让我吃了不少苦头。

2.2 MATLAB实战：代码与直方图的双重辅助

光说不练假把式，我们直接上MATLAB代码。我会用经典的 rice.png（一堆米粒的图像）来演示。这张图背景较暗，米粒较亮，非常适合做阈值分割。

clear all; close all; clc;
% 读取并转换为灰度图像
I = imread('rice.png');
I_gray = rgb2gray(I); % 如果原图是RGB，先转灰度

% 显示原图和它的灰度直方图
figure;
subplot(1,3,1);
imshow(I_gray);
title('原始灰度图像');

subplot(1,3,2);
imhist(I_gray); % 这是关键！直方图是我们选择阈值的“地图”
title('图像灰度直方图');
xlabel('灰度值'); ylabel('像素数量');

% 假设我们通过观察直方图，发现波谷大概在125附近
% 直方图通常有两个“峰”，一个对应背景（较暗），一个对应前景（较亮），中间的低谷就是潜在的阈值点
threshold_manual = 125;

% 进行阈值分割：大于阈值的为1（白色，前景），小于等于的为0（黑色，背景）
I_bw_manual = I_gray > threshold_manual;

subplot(1,3,3);
imshow(I_bw_manual);
title(['人工阈值分割结果: T=', num2str(threshold_manual)]);

运行这段代码，你会看到三张图。中间那张直方图是你的“决策依据”。你会看到两个明显的鼓包，左边鼓包代表暗色的背景，右边鼓包代表亮色的米粒。两个鼓包中间凹陷下去的地方，大概在120-130之间，这就是我们要找的“门槛”。选择125，分割效果通常就不错。

这里有个非常重要的技巧：imhist 函数是你的最佳战友。不要凭感觉瞎猜阈值，一定要先看直方图。一个理想的、适合人工阈值分割的图像，其直方图应该是双峰的，并且两峰之间有一个明显的低谷。这个低谷的位置，就是最佳阈值的强力候选。如果直方图是单峰的，或者跟平原一样平坦，那说明前景背景对比度太差，人工阈值就很难奏效，你得考虑后面要讲的自动方法了。

3. 自动阈值法：让计算机自己找到“黄金分割点”

手动调阈值毕竟不是长久之计，尤其是面对批量处理或者图像条件多变的情况。这时候，我们就需要把选择权交给算法，让它们根据图像自身的统计特性，自动找出那个最优的阈值。这里我重点讲两个最常用、也最实用的方法：OTSU（大津法）和自适应阈值。

3.1 OTSU算法：寻找最“分裂”的阈值

OTSU算法的思想非常巧妙，它不关心具体的灰度值是多少，而是像一个精明的裁判，寻找一个能把图像像素“分得最开”的阈值。它的评判标准是类间方差。简单类比一下：假设你要把一群人分成两组，OTSU的目标是找到一种分法，使得两组人平均身高差距最大，同时每组内部的身高又尽可能集中。这个“平均身高的差距”，就是类间方差。类间方差越大，说明两组差异越明显，这个分组方法（阈值）就越好。

算法步骤可以这样理解：

1. 遍历所有可能的灰度值（0到255），把每个灰度值都当作候选阈值T。
2. 对于每个T，把图像像素分成两组：低于T的（背景类）和高于T的（前景类）。
3. 计算这两组像素各自占的比例、平均灰度。
4. 根据这些值，算出用这个T分割时的“类间方差”。
5. 最后，选择那个能让类间方差最大的T，它就是OTSU算法找到的最佳全局阈值。

它的优点是完全自动、全局最优，对于前景和背景对比度较好的图像，效果非常稳定。在MATLAB里，调用它只需要一行代码：

% 接续之前的代码，I_gray是我们的灰度图
level_otsu = graythresh(I_gray); % 这个函数返回的是归一化的阈值，范围[0,1]
I_bw_otsu = imbinarize(I_gray, level_otsu); % 新版本推荐用imbinarize
% 或者更简洁的写法：
% I_bw_otsu = imbinarize(I_gray, ‘global’); % ‘global’ 默认就是用OTSU

figure;
subplot(1,2,1);
imshow(I_gray);
title('原始图像');
subplot(1,2,2);
imshow(I_bw_otsu);
title(['OTSU自动阈值分割，阈值=', num2str(level_otsu*255)]); % 转换回0-255范围显示

你可以对比一下之前人工选择125的结果和OTSU自动计算的结果，通常会发现OTSU选得更准，特别是当直方图的双峰不那么对称时。

3.2 自适应阈值：应对光照不均的“局部专家”

OTSU虽然好，但它是个“全局主义者”，认为整张图用一个阈值就够了。这在实际中经常会碰壁。比如，你拍一张文档，由于灯光角度，左边亮右边暗。用全局阈值，要么左边字迹消失，要么右边背景变成黑色。

这时候，就需要自适应阈值出场了。它的核心思想是“具体问题具体分析”：图像中每个像素点的阈值，不是由全局决定的，而是由它周围一小片邻居区域的灰度特性动态决定的。常用的方法是计算这个局部区域（比如一个15x15的窗口）内像素的灰度平均值，或者高斯加权平均值，然后减去一个常数C，作为该中心像素的阈值。

这样做的好处是，在亮的区域，局部平均灰度值高，阈值也高；在暗的区域，局部平均灰度值低，阈值也低。从而在整个图像上实现动态的、自适应的分割。

% 读取一张可能光照不均的图像，例如‘text.png’（MATLAB自带）
I_text = imread(‘text.png’);
if size(I_text,3)==3
    I_text_gray = rgb2gray(I_text);
else
    I_text_gray = I_text;
end

% 使用全局OTSU阈值
bw_global = imbinarize(I_text_gray, ‘global’);

% 使用自适应阈值
% ‘adaptive’ 表示自适应方法， ‘Mean’ 表示使用局部均值， 第二个参数是邻域大小（奇数）， 最后是常数C
bw_adaptive = imbinarize(I_text_gray, ‘adaptive’, ‘Sensitivity’, 0.4); 
% 注意：新版本MATLAB中，imbinarize的‘adaptive’选项通过‘Sensitivity’参数调整，值越小越敏感（背景更可能被划为前景）

% 为了更精细控制，可以使用底层函数adaptthresh先获取阈值曲面，再二值化
T = adaptthresh(I_text_gray, 0.4, ‘NeighborhoodSize’, 25, ‘Statistic’, ‘mean’);
bw_adaptive2 = imbinarize(I_text_gray, T);

figure;
subplot(2,2,1); imshow(I_text_gray); title(‘原始灰度图像（光照不均）’);
subplot(2,2,2); imshow(bw_global); title(‘全局OTSU阈值 - 效果不佳’);
subplot(2,2,3); imshow(bw_adaptive); title(‘imbinarize自适应阈值’);
subplot(2,2,4); imshow(bw_adaptive2); title(‘adaptthresh+imbinarize自适应阈值’);

关键参数解读：

• 邻域大小（NeighborhoodSize）：决定了“局部”的范围。太小了会对噪声敏感，太大了会模糊掉真正的局部光照变化。一般取图像尺寸的1/8到1/20之间的奇数。
• 常数C（或Sensitivity）：这是一个微调参数。从局部均值里减去这个值得到最终阈值。C为正数，会使得阈值更低，更多的像素被归为前景（更“敏感”）。需要根据图像内容微调。

自适应阈值的缺点是计算量比全局阈值大，因为它要对每个像素（或其邻域）进行计算。但对于光照不均、背景渐变的情况，它是无可替代的利器。

4. 分水岭算法：破解“粘连物体”分割难题

前面讲的阈值方法，无论人工还是自动，都有一个根本前提：我们能通过灰度差异把目标和背景分开。但如果图像中，我们想要的目标彼此挨在一起，甚至粘连重叠呢？比如一堆重叠的细胞、碰在一起的硬币、堆积的颗粒。这时候，阈值分割只能得到一个连成一片的白色区域，无法区分出单个物体。这就需要请出我们今天的“重磅嘉宾”——分水岭算法。

4.1 算法思想：一幅地形图的“洪水模拟”

分水岭算法的灵感来源于地理学。它把一张灰度图像想象成一幅地形图：图像中每个像素的灰度值，代表该点的海拔高度。亮度高的地方是山峰，亮度低的地方是山谷。

算法的过程，就像在这幅地形图上模拟一场洪水：

1. 预处理与求梯度：首先，对原图进行平滑去噪（比如用高斯滤波），然后计算图像的梯度幅值。梯度大的地方，意味着灰度变化剧烈，对应地形图中的“山脊”（物体的边缘）。所以，梯度图像才是我们真正要进行分水岭分割的地形图，山峰是边缘，山谷是物体内部和平坦区域。
2. 确定“注水点”：我们不可能让洪水漫无目的地流。需要先找到一些“注水点”，也就是每个物体的标记点。这些点通常位于梯度图像的局部最小值区域（即每个物体的内部深处）。寻找这些标记点是分水岭成功与否的最关键一步，做不好就会导致严重的“过分割”（一个物体被分成很多小块）或“欠分割”（多个物体没分开）。
3. 模拟洪水蔓延：从每一个标记点（山谷最低处）开始注水。水位逐渐上升，从各个标记点蔓延出的“水域”会不断扩大。
4. 形成“分水岭”：当从不同标记点蔓延出的洪水即将汇合时，我们就在它们之间筑起一道大坝，阻止它们合并。这道大坝，就是分水岭，也就是我们最终想要的分割边界。
5. 得到分割结果：当洪水淹没整个地形图，所有大坝连接起来，就形成了完整的、分割各个物体的边界线。

4.2 MATLAB实战：从距离变换到标记分水岭

理论听起来有点绕，我们直接看一个处理粘连硬币的经典例子。分水岭的难点在于标记，这里我介绍一种基于距离变换的经典标记方法。

clear all; close all; clc;
% 读取粘连的硬币图像
I = imread(‘coins.png’);
if size(I,3)==3
    I_gray = rgb2gray(I);
else
    I_gray = I;
end

% 第一步：通过阈值分割得到硬币的粗略区域（粘连在一起的）
bw = imbinarize(I_gray, ‘global’);
bw = imfill(bw, ‘holes’); % 填充孔洞
bw = bwareaopen(bw, 50); % 去除小面积噪声
subplot(2,3,1); imshow(I_gray); title(‘1. 原始灰度图像’);
subplot(2,3,2); imshow(bw); title(‘2. 全局阈值二值化结果（粘连）’);

% 第二步：计算距离变换
% 距离变换：计算二值图像中每个前景像素到最近背景像素的距离。
% 对于粘连的物体，其中心区域距离背景最远，会形成一个个“山峰”。
D = -bwdist(~bw); % 对背景求距离并取负，使背景为“高峰”，前景为“低谷”，方便后续找极小值
subplot(2,3,3); imshow(D, []); title(‘3. 距离变换（取负）’);
colorbar;

% 第三步：寻找局部极小值作为内部标记
% 直接找极小值会找到太多（受噪声和边界影响），需要先进行平滑和极值抑制
D_mod = imhmin(D, 2); % 抑制深度小于2的极小值，减少过分割
Ld = watershed(D_mod); % 对修改后的距离图进行分水岭，得到标签
bw_marker = Ld == 0; % 分水岭边界为0
% 内部标记是分水岭区域，但不是边界
internal_marker = ~bw_marker;
internal_marker = bwareaopen(internal_marker, 20); % 去除太小的标记
subplot(2,3,4); imshow(label2rgb(Ld)); title(‘4. 对距离变换分水岭（用于找标记）’);

% 第四步：将内部标记作为约束，应用到原始图像的梯度图上
% 计算梯度图（地形图）
hy = fspecial(‘sobel’); % 使用Sobel算子求梯度
hx = hy’;
grad_y = imfilter(double(I_gray), hy, ‘replicate’);
grad_x = imfilter(double(I_gray), hx, ‘replicate’);
grad_mag = sqrt(grad_x.^2 + grad_y.^2);
subplot(2,3,5); imshow(grad_mag, []); title(‘5. 梯度幅值图（地形图）’);

% 用内部标记“强制修改”梯度图：在标记点处将梯度设为最小值
grad_mod = imimposemin(grad_mag, internal_marker);

% 第五步：对修改后的梯度图进行分水岭分割
L = watershed(grad_mod);
boundary = L == 0; % 分水岭边界

% 可视化最终结果
result = I_gray;
result(boundary) = 255; % 将边界标为白色，叠加到原图

subplot(2,3,6); imshow(result); title(‘6. 最终分水岭分割结果’);
hold on;
% 可选：显示找到的内部标记点
[marker_y, marker_x] = find(internal_marker);
plot(marker_x, marker_y, ‘r.’, ‘MarkerSize’, 10);
hold off;

这段代码是分水岭算法的一个经典流程。核心技巧在于利用距离变换找到物体内部的“深谷”作为标记点，然后用这些标记点去约束梯度图上的分水岭过程，从而避免无目的的过分割。imimposemin 函数在这里起到了关键作用，它确保了洪水只从我们指定的标记点（物体内部）开始蔓延。

分水岭算法功能强大，但对噪声和初始标记非常敏感。在实际应用中，你可能需要花很多精力在预处理（去噪）和标记提取上。有时，结合边缘检测、形态学操作（开闭运算）来优化初始二值图和标记，是获得好结果的必经之路。

5. 方法对比与实战选型指南

走完了从人工阈值到分水岭的完整流程，你手里现在有了好几把“手术刀”。最后一节，我们来聊聊实战中怎么选刀，以及怎么磨刀。这比单纯学会算法更重要。

首先，我根据自己的经验，把这几种方法的特点总结成了下面这个表格，你可以快速查阅：

那么，在实际项目中该如何选择呢？我给你一个清晰的决策路径：

1. 第一步：看直方图。 打开你的图像，第一件事就是 imhist。如果直方图是漂亮的双峰形态，恭喜你，问题很简单。直接尝试 OTSU算法，它十有八九能给你一个很好的结果。如果OTSU效果不理想，但双峰依然可见，可以微调一下阈值，或者试试迭代法（代码中演示过的，通过不断逼近前景背景均值来求阈值）。

2. 第二步：检查光照均匀性。 如果OTSU分割后，图像某些部分前景没了，或者背景被误判成前景，很可能是光照不均。这时候，果断切换到自适应阈值。先试试MATLAB内置的 imbinarize(I, ‘adaptive’)，调整 ‘Sensitivity’ 参数。如果效果不够精细，就用 adaptthresh 函数生成阈值曲面，再二值化，这样可以更灵活地控制邻域大小和统计方法。

3. 第三步：判断物体是否粘连。 如果经过阈值分割（无论是全局还是自适应）后，你得到的是一大坨连在一起的白色区域，而你需要的是其中一个个独立的对象，那么分水岭算法就是你的菜。记住，分水岭的输入通常是梯度图，而成功的关键在于找到好的“标记”。多花时间在预处理上，比如用距离变换、形态学开闭运算来获取更干净的前景对象和内部标记。

4. 永远不要忘记预处理和后处理。 图像分割从来不是靠一个算法单打独斗。在分割之前，高斯滤波去噪、对比度拉伸增强，能让阈值方法更稳定。在分割之后，形态学操作（如 imopen, imclose 去除小噪点或填充空洞）、bwareaopen 过滤小面积区域、imfill 填充孔洞，这些操作能极大地提升最终结果的质量。我习惯把分割流程看作一个流水线：输入图像 -> 预处理 -> 核心分割 -> 后处理 -> 输出结果。每个环节都可能需要调试。

最后，分享一个我调试分水岭算法时的小窍门：把中间每一步的结果都 imshow 出来看看。特别是距离变换图、找到的标记点、修改前后的梯度图。可视化能让你直观地理解算法在哪里出了问题——是标记点找多了导致过分割？还是标记点没找到导致欠分割？然后有针对性地去调整前期的滤波强度、距离变换的极值抑制参数 (imhmin 的h值) 或形态学操作的核大小。图像处理很多时候就是一场和参数的“对话”，耐心点，多看中间结果，你就能摸清它的脾气。