简介：面向遥感影像处理与地理信息系统应用，这份专题资料聚焦全色影像中的水体提取，以阈值分割为核心手段，适合图像处理初学者、测绘或环境监测从业者快速上手。压缩包内共2个文件，包括1个Matlab脚本与1个BMP结果示例图，整体大小2.18MB，脚本简洁可直接运行，输出图像便于直观验证分割效果。目前已有1270人学习，属于轻量实用的入门级资源。通过该资源，既能理解阈值分割的基本原理，也能看到从影像预处理、灰度特征选择、阈值设定到二值化分割、后处理去噪的完整流程；同时对Otsu自动阈值选取、连通区域分析、边缘修补等常见操作会有具体感知。即使仅对全色影像做快速初步的水体识别，这套方法也足够简单高效，可迁移至其他单波段影像分析任务。

1. 阈值分割做水体提取：全色影像里把河流与湖泊切出来的最小可行方案

做地理信息相关项目时，经常会遇到一个挺尴尬的需求：手里只有一张高分辨率全色影像，没有多光谱波段，却要快速把水体范围圈出来。全色影像只有一个波段，看不到近红外，教科书里的 NDWI、MNDWI 这类多光谱指数全都用不上。这个专题包里的 water.m 和 out.bmp，解决的就是这个单波段场景下的水体提取问题——用最基础的阈值分割，把影像里灰度值偏低的像素判定为水体，直接输出一张黑白分明的水体分布图。整个方案不依赖任何第三方库，一个 Matlab 脚本几步就走完全流程，特别适合做前期的快速摸底，或者给后续更精细的监督分类提供一个干净的初始掩膜。

2. 全色影像与阈值分割的原理：为什么单波段也能把水提出来

2.1 全色影像的水体灰度特征

全色影像本质上就是一张单波段的灰度图，每个像素只记录一个辐射亮度值。但恰恰是这一个值，在水体识别里反而给了我们一个稳定的物理依据：水体在可见光波段对电磁波的吸收率远高于土壤、植被和建筑物，所以在全色影像上，水体区域通常表现为大面积的暗色块，灰度值明显低于周边地物。也就是说，水在影像里天生是“偏黑”的。

这个特性让阈值分割有了用武之地。阈值分割的思路非常直接：选定一个灰度阈值 T，遍历影像每个像素，灰度值大于等于 T 的归为一类，小于 T 的归为另一类。在水体提取的语境下，低于阈值的判为水体，高于阈值的判为背景。这样一个线性决策规则虽然简单，但它在水体这种“低灰度同质区域”上往往表现得相当稳健。河流、湖泊内部的纹理差异很小，不像城市区域那样有密集的边缘，所以单靠灰度阈值就能把大片水域从地面背景里剥离开。

当然，单波段阈值分割也不是万能药。它最怕两种干扰：一是水体表面有强反射的漂浮物（比如船只、白色波浪），局部灰度被拉高，导致漏提；二是山体阴影和建筑物阴影的灰度值与水体重叠，它们同样很暗，会被错误地划进水体。后处理阶段要专门处理这些噪声，这部分在后面的代码里会具体讲到。

2.2 固定阈值与 Otsu 自适应阈值的选择

阈值分割首先要回答一个问题：T 取多少。最省事的做法是人工盯着一张测试影像，反复试出来一个固定阈值，比如 40。这样做的优点是快、逻辑透明，缺点是一批影像的光照条件变了，这个 40 就成了拍脑袋值。不同时间、不同季节拍到的同一片水域，灰度分布可能差出好几个数量级，固定阈值很容易整片失效。

更稳妥的做法是让算法自己算阈值，其中最经典、在 Matlab 里一行就能调用的就是 Otsu 方法。Otsu 的核心思想是遍历 0 到 255 之间的每个灰度级，找到那个能让分割出来的两类像素“类内方差最小、类间方差最大”的阈值。用人话说，就是让水体类尽量黑得整齐、背景类尽量亮得整齐，两类之间的差距尽量拉大。这样即使整幅影像偏暗或偏亮，Otsu 都会根据直方图的实际分布自适应地调整阈值，保证分割不会被整体灰度偏移带跑。

在实际开发里，我一般会把 Otsu 当作默认值，再用一个开关保留手工指定阈值的入口。原因很朴素：Otsu 在大多数情况下已经够好，但如果数据里有大面积的高亮噪声，它的阈值会被噪声拉偏，这时手工阈值反而是最快的兜底手段。代码里同时留了 autoMode 和 manualThresh 两个参数，就是为了让使用者不用改算法，只改一个数字就能在两种模式间切换。

3. 核心代码拆解：water.m 里一个像素从灰度到水与非水的完整决策路径

3.1 主脚本的整体流程与三个核心函数

water.m 这个脚本的完整逻辑可以拆成五步：读取图像 → 预处理去噪 → 计算阈值 → 二值化分割 → 后处理优化。在 Matlab 里，这个过程非常紧凑，核心代码可以按顺序写成下面这样一个完整的脚本。为了便于理解，我把每一步用注释标了出来。

%% water.m - 全色影像水体提取主脚本 clear; clc; close all; % ========== 1. 读取影像 ========== srcFile = 'demo.tif'; % 输入的全色影像路径 img = imread(srcFile); if size(img, 3) == 3 % 偶尔会有三通道伪装的全色图 img = rgb2gray(img); end img = double(img); % 转double，避免后续计算溢出 % ========== 2. 预处理：中值滤波去噪 ========== medWinSize = 3; % 滤波窗口，3x3即可 imgF = medfilt2(img, [medWinSize medWinSize]); % ========== 3. 计算分割阈值 ========== autoMode = true; % true用Otsu，false用手动阈值 manualThresh = 50; % 手动阈值，仅在autoMode=false时生效 if autoMode level = graythresh(uint8(imgF)); % Otsu，返回归一化阈值 threshVal = level * 255; % 换算到0-255灰度域 else threshVal = manualThresh; end % ========== 4. 阈值分割 ========== bw = imgF > threshVal; % 大于阈值是背景，小于等于阈值是水体 bw = ~bw; % 取反，让水体变成白色(1) % ========== 5. 后处理 ========== bw = bwareaopen(bw, 80); % 删除小于80像素的孤立碎块 se = strel('disk', 3); % 半径为3的圆盘结构元素 bw = imclose(bw, se); % 闭合操作，填补内部小孔洞 bw = imfill(bw, 'holes'); % 继续补洞 % ========== 6. 保存与显示 ========== imwrite(uint8(bw) * 255, 'out.bmp'); % 输出8位二值bmp figure; imshow(bw); title('水体提取结果');

这段脚本的核心逻辑并不复杂，重点在于几个参数的取值和它对应的效果。medWinSize 取 3 是因为水体轮廓通常比较平滑，3x3 的中值滤波已经能干掉大部分独立噪点，同时不会把细小的支流水道抹掉；如果影像里有明显的椒盐噪声，可以把它调到 5，但要观察支流会不会因此断线。graythresh 是 Matlab 自带的 Otsu 实现，返回的 level 是一个 0 到 1 之间的归一化值，所以要乘回 255 才能和原始灰度值比较。

3.2 二值化方向与后处理的设计逻辑

很多人第一次写阈值分割时会卡在一个问题上：到底大于阈值是水体还是小于阈值是水体。这个其实取决于影像自身的特性。在全色影像里水体是暗的，所以小于阈值为水体。但代码里我先写了imgF > threshVal得到背景，再用bw = ~bw反转，这样做的目的是让读者一眼就能理解“大于阈值的是背景这一语义”，避免直写imgF < threshVal带来的理解歧义。实际跑代码时，这两种写法得到的结果完全一样。

接下来是后处理的三连招。bwareaopen(bw, 80)的作用是删除面积小于 80 个像素的连通区域。这一步处理的噪声主要是两类：水体外面零星的暗色噪点被判成的小块水体，以及水体边缘因为灰度过渡被切断产生的碎屑。strel('disk', 3)配合imclose做的是形态学闭运算，先膨胀再腐蚀，用于填充水体内部的小空洞。水体表面如果有高亮的船只或浮萍，在二值图像里会形成黑洞，闭运算可以把它周边的水体区域连接起来，把洞补上。最后imfill(bw, 'holes')再兜底一次，把任何剩余的封闭孔洞填掉。

实际效果可以这样预估：主要河道的轮廓会非常完整，宽度在 10 个像素以上的水道基本不会被中值滤波破坏；零星的小水洼如果面积小于 80 像素，会被 bwareaopen 误删。这是后处理带来的必然取舍——为了去掉噪声，必然要牺牲一些极小的真实水体。

3.3 输出结果的规格与后续使用方式

脚本最后直接将二值矩阵写成了 8 位 BMP 文件。这一步有一个容易被忽略的细节：直接用imwrite(bw, 'out.bmp')会把逻辑真(1)和假(0)写入文件，在很多看图软件里会显示成一张近乎全黑的图，因为 1 的灰度太暗。代码里用uint8(bw) * 255把 1 映射到 255，这样水体就是纯白色，背景是纯黑色，看起来直观，也方便直接送入 GIS 软件做矢量化。

输出的 out.bmp 虽然是图像格式，但它的本质是一个栅格掩膜。在 ArcGIS 或 QGIS 里，可以直接用「栅格计算器」或者「栅格转面」工具，把白色区域转成矢量多边形，再统计面积、计算与其它图层的叠置关系。这也是这个专题最常用的落地链路：从像素级的二值判定，到图斑级的业务数据。

4. 跑起来再说：环境准备、数据输入与第一次运行全流程

4.1 运行环境和必要假设

这个项目对运行环境的要求很低。只要是近十年的 Matlab 版本都可以直接打开 water.m 运行，依赖的全部函数 medfilt2、graythresh、bwareaopen、strel、imclose、imfill 都是 Image Processing Toolbox 里的标准函数，不需要额外安装任何第三方工具包。如果你用的是 Octave，大部分函数也是兼容的，但 imclose 和 bwareaopen 的行为在个别版本上会有细微差异，建议优先用 Matlab。

输入数据方面，water.m 里把影像路径写死成了demo.tif，所以第一次运行前需要把你自己的全色影像放到脚本同一目录下，并改名为 demo.tif，或者把脚本里的srcFile改成你的实际文件名。这里要特别注意：如果输入影像的位深不是 8 位而是 16 位，graythresh 会返回一个基于 0 到 65535 动态范围的阈值，后面乘 255 这一步就会出错，需要改成level * 65535。

4.2 从零跑通全流程的逐步操作

第一次运行建议按下面这个顺序操作，每跑完一步都停下来看一眼中间结果，不要直接跳到最终输出。这样做的好处是，如果最终的水体图有问题，你能立刻定位到是预处理、阈值还是后处理哪个环节出了问题。

第一步，准备数据。在 Matlab 当前目录建一个专用的文件夹，比如 waterExtract，把全色影像放进去，重命名为 demo.tif。第二步，打开脚本。双击 water.m，在编辑器里确认当前目录和文件路径一致，然后直接点运行。第三步，看命令行窗口的输出。如果脚本执行成功，工作区里会多出bw、imgF、threshVal这几个变量，其中threshVal的值就是这次分割实际用到的阈值。

第四步，查看图窗。脚本运行后弹出来的图窗里，白色区域代表被识别为水体的部分。把它和原始影像对照看：如果水体位置被正确圈出且边界没有明显的杂乱碎块，说明整套参数是有效的。第五步，检查输出的 out.bmp 文件，确认它在文件管理器里存在，且文件大小不是一个完全空白的尺寸。

# 如果是在服务器上没有图形界面的环境，可以用Matlab命令行批处理运行： matlab -batch "run('water.m')"

上面这行命令适合把提取流程接入自动化流水线的场景。-batch模式会在不打开桌面界面的情况下执行脚本，执行完自动退出。此时脚本里的 figure 命令会因为没有图形窗口而被跳过，但 out.bmp 依然会正常写出。

4.3 不同场景下的参数调整对照

跑通了默认参数后，最常做的调整集中在三个地方：中值滤波窗口、阈值模式、后处理的面积阈值。下面这张表总结了常见情况下的取值建议和预期效果，遇到提取结果不理想时，可以对照排查。

调参时有个经验可以参考：如果第二次跑出来的结果比第一次还差，不要反复微调单一参数，而是把滤波窗口和面积阈值同时恢复到默认值，只动阈值这一个变量，这样能更快定位瓶颈。

5. 避坑记录：水体提取最容易翻车的四个操作细节

5.1 out.bmp 存完打开是黑的或全白的

现象：脚本提示执行成功，out.bmp 也生成了，但打开看要么整张都是纯黑，要么整张都是纯白。这种情况在第一次跑这个项目的开发里非常常见，不过原因通常不是算法错了。

原因：二值矩阵里的值为逻辑 1 和 0，直接用 imwrite 存储时，很多看图软件会把 1 显示成接近黑色的灰度，导致视觉上“消失”了。另一种可能是代码里漏了取反操作，把水体设成了 0、背景设成了 1，结果整幅图只有零星的白点，看起来像是全黑。

解决：参考 3.3 节的写法，输出前先转 uint8 并乘以 255，强制把逻辑值映射到 0 和 255 两个灰度级。取反操作放在后处理之前完成，确保水体是前景白色。判断是否修正成功的标准，是用 Matlab 读回 out.bmp，检查像素值只有 0 和 255 两种取值。

5.2 山体阴影和深色土壤被错判成水体

现象：提取结果里，山坡背阴面和低洼湿地在图上被标成了白色，实际的水体区域反而没有完整连成片。

原因：单波段影像里，山体阴影的灰度值和水体非常接近，Otsu 方法只看灰度直方图，它没有地形知识，无法区分“暗色的水体”和“暗色的阴影”。这是全色影像水体提取的固有缺陷，不是代码 bug。

解决：可以在后处理链路里增加一个形状约束。水体的连通域通常形状比较规整，边缘平滑，而山体阴影往往呈现出沿坡向伸展的细长形态。用regionprops计算每个连通域的短长轴比，把离心率过大的区域丢弃。具体思路是给bwareaopen之后再加一步筛选，只保留面积大于某阈值且离心率低于某阈值的连通域。

5.3 提取结果里水体内部全是碎洞

现象：水体部分整体轮廓能看出来了，但内部星星点点全是黑色的小洞，整片湖泊像被虫蛀过。

原因：水体表面存在高亮物体，比如白色波浪纹理、船只、水面的漂浮垃圾。这些物体反射率高，灰度值接近陆地背景，二值化后它们在白色水体区域内部形成了小块的黑色孤岛。

解决：形态学闭合操作是标准解法。把 3.1 节代码里的strel('disk', 3)改成strel('disk', 5)往往就能明显改善，因为更大的结构元素能把更宽范围的白色水体连接起来，覆盖掉更大尺寸的黑色孤岛。如果改大半径后支流被磨掉了，可以改用imfill加连通域面积限制的组合，而不是一味加大结构元素。

5.4 同一幅影像左侧能提出水体右侧全是噪声

现象：把影像从中间切开，左半部分的分割结果相当干净，右半部分却出现了一大片连续的高亮噪声，阈值好像只在局部有效。

原因：这是光照不均造成的。全色影像拍摄时如果云层遮挡不均匀，或者镜头存在渐晕效应，会导致同一影像不同区域的灰度基准不一致。左边水体灰度均值 30，右边背景灰度均值 45，那全局一个阈值必然顾此失彼。

解决：分块阈值分割是应对这个问题的一个不复杂且有效的方案。把影像按网格切分成若干小块，在每个小块内独立运行 Otsu 计算阈值，再把所有小块的分割结果拼接回来。需要注意两点：一是块的大小要适中，至少要能包含一定比例的水体和背景，两个 128x128 的块通常够用；二是块与块之间存在阈值跳变，在拼接后再做一次中值滤波或者形态学开闭操作，消除块边界的接缝痕迹。

6. 二次开发指南：从单阈值到多策略水体识别的改造路线

6.1 把脚本改造成可配置的函数式模块

water.m 目前是脚本形式，所有参数写死在代码里。如果要接入到更大的业务链路里，第一步建议把它改造成函数。我一般会这样改：把影像路径、自动/手动模式、阈值、滤波窗口大小、面积阈值、闭合结构元素半径这六个参数全部提升为函数入参，输出除了 bmp 文件之外，还把 bw 矩阵返回出来，方便后续计算面积和叠加分析。

function [bw, threshVal] = waterExtract(srcFile, varargin) % 输入: srcFile 影像路径 % 'AutoMode', true/false % 'ManualThresh', 数值 % 'MedWin', 3/5 % 'MinArea', 像素个数 % 'CloseRadius', 结构元素半径 % 输出: bw 二值掩膜, threshVal 实际使用的阈值 p = inputParser; addParameter(p, 'AutoMode', true); addParameter(p, 'ManualThresh', 50); addParameter(p, 'MedWin', 3); addParameter(p, 'MinArea', 80); addParameter(p, 'CloseRadius', 3); parse(p, varargin{:}); % ... 主逻辑与3.1节相同 end

这样改造之后，调用方只需要一行[bw, th] = waterExtract('demo.tif', 'AutoMode', false, 'ManualThresh', 40);就能完成提取，并拿到实际使用的阈值。调研阶段可以写一个 for 循环，遍历多组参数组合，把每组输出的 out.bmp 拼成一张对比图，快速筛选出最适合当前数据源的参数组合。

6.2 分块 Otsu 的具体改造实现

光照不均问题的另一个彻底解法，是把 5.4 节提到的分块思路落地成代码。核心步骤如下：先把图像切成 n x n 的块，每个块独立跑 graythresh，得到 n 个阈值组成一个矩阵，再用插值把矩阵插回原图大小，生成一个每个像素都有独立阈值的阈值网格。最后直接用imgF > threshGrid把原图和新阈值网格逐像素比较，这样既保留了局部的光照适应性，也避免了块边界跳变。

在实现上要注意两个细节。块大小建议取 256x256，太小会让背景类数量不足，Otsu 就退化成找单峰噪声；太大又失去了分块的意义。插值建议用线性插值，不要用最近邻，否则阈值网格的接缝会直接映射到输出结果上成为明显的块状痕迹。

6.3 改造完成后的验证方法

每次改完代码，都要用一套固定的验证流程来防止“新改没把旧功能改坏”。我自己的习惯是准备三张测试影像：一张是宽而平静的干流河道，用来验证大水面完整性；一张是细小支流密集的区域，用来验证窄水体的保真度；一张是包含大片山体阴影的丘陵区，用来验证抗干扰能力。每跑一次改造版本，就把三张影像的输出结果和上一版并列对比，重点看水体边界偏移是否超过 2 个像素、支流断线是否增加、阴影误判区是否扩大。这组数据固定下来以后，每次改动都有了量化的评判标准，而不是凭感觉说“好像更好了一点”。

从那时起，我每次改分割类的算法，都会强制自己先跑一遍三张验证图的基线结果，再做改动。这个习惯帮我挡掉了不少自我感觉良好但实际在退化的情况。这个专题虽然只用了阈值分割这一招，但它的工程链路——预处理、自适应阈值、形态学后处理、评估验证——几乎可以原样复用到更多复杂的提取任务里。希望帮到你。

本文还有配套的精品资源，点击获取