简介面向图像处理初学者与需要快速实现阈值分割的开发者这份大津法二值化程序以MATLAB脚本形式给出聚焦灰度图像中前景与背景的分离问题可用于文本识别、医学图像分析、工业检测等实验场景也适合课程练习和算法入门对照。压缩包为rar格式共1个文件主体是1个m脚本包体约1KB体量轻便于直接阅读、运行和改造。已有396人学习说明其在基础图像处理练习中具备一定参考价值。代码围绕灰度直方图统计、阈值遍历、前景与背景均值及类间方差计算展开读者可据此理解最佳全局阈值的选取过程并复现二值化效果也可在脚本基础上替换输入图像、调整参数观察噪声、光照和图像质量对分割结果的影响积累算法调试与效果对比经验。1. 一批发票扫描件让我重新捡起大津法前段时间接手一批票据扫描件的预处理开头想省事直接gray 127一刀切。结果一半的图翻了车左侧被台灯照得发白数字糊成一片右边背光的角落又黑得连边框都吃掉。同一批图换个角度看阈值就得跟着改改到最后我自己都记不清哪张用的哪个数。后来把大津法二值化程序翻出来重跑一遍让它自己从灰度直方图里算阈值那批图基本一次过。这份程序做的事很朴素读图、转灰度、统计直方图、遍历 256 个候选阈值找类间方差最大的那个、回写二值图再配一个实验脚本把 Otsu 和固定阈值的结果摆在一起对比。核心算法部分只用 NumPy不装 OpenCV 也能跑完整流程OpenCV 只在对照验证时出场。适合两类人做文档扫描、OCR 前处理、工业检测预处理的开发者拿来就能接进流水线以及想把 Otsu 从公式变成能调试的代码、而不是永远调cv2.threshold的人。2. 项目定位与技术栈核心算法只留 NumPy 的取舍2.1 一份目录看全项目结构拿到资源先别急着跑先看目录。这个项目的结构非常扁平没有花哨的包管理属于「拷到哪都能跑」的类型目录长这样otsu-binarize/ ├── otsu.py # 核心算法直方图、阈值搜索、二值化 ├── main.py # 单图命令行入口 ├── experiment.py # 批量对比实验Otsu vs 固定阈值 ├── io_utils.py # 读写图、灰度化、结果落盘 ├── samples/ # 实验用图含光照不均、低对比度几类 ├── out/ # 运行产物跑完自动创建 └── requirements.txtotsu.py是全部价值所在一百行出头没有类、没有继承、没有配置对象三个纯函数gray_hist()出直方图otsu_threshold()出阈值binarize()出结果图。main.py只负责把命令行参数翻译成这三个函数的调用experiment.py负责批量跑并打印指标。这种拆法有两个好处。一是算法可被单独导入你把它塞进 Flask 接口、塞进批处理脚本、塞进 PyQt 界面都不用改一行二是出问题时定位极快阈值不对就只可能是otsu_threshold()的输入不对跟读图、跟保存格式无关。我一般接手别人的图像项目第一件事就是把「算法」和「I/O」这两层剥开剥不开的项目后面改起来基本都是血泪。2.2 技术栈与依赖为什么核心不引 OpenCV依赖清单刻意做得很轻主链路只依赖 NumPy 和 Pillow依赖用途是否必需numpy直方图统计、累积和、方差计算、数组运算必需Pillow读图、转灰度、保存结果必需opencv-python仅用于experiment.py里做结果对齐校验可选matplotlib画直方图和阈值位置调试用可选不把 OpenCV 放进主链路不是为了标新立异而是为了让阈值怎么来的这件事完全暴露在你眼前。cv2.threshold(img, 0, 255, cv2.THRESH_OTSU)一行出结果方便但你不知道它内部是用 256 桶直方图还是 128 桶不知道它遇到单峰图怎么处理出了偏差只能靠猜。自己实现一遍之后ret那个返回值是怎么算出来的、为什么有时候是 117 有时候是 118你心里就有数了。常见做法是主链路保持纯 NumPy把 OpenCV 放到验证脚本里当参照物。这样跑实验时能立刻发现「我的实现和成熟库差多少」而不是等上线后被业务方发现。2.3 大津法在二值化方案里的位置很多人一上来就问「Otsu 是不是最好的」这个问题没有统一答案得看你图的性质。把这个项目和几种常见方案放一起对比选择标准就清楚了方案阈值来源适用场景主要短板固定阈值人工指定光照稳定、批量一致的工业相机图换批次就得重调玄学感强Otsu本项目全图直方图最大化类间方差双峰明显、前景背景面积比不太悬殊光照不均、单峰图会失效自适应均值局部窗口均值减常数 C文档扫描、背景渐变空白区域容易噪点化自适应高斯局部窗口高斯加权均值减 C同上边缘更平滑计算量大块大小敏感Sauvola 等局部方法局部均值与标准差低对比度文档、退化扫描件参数多调起来费时间Otsu 的定位是「全局、自动、够快」。它的计算复杂度只跟灰度级数有关256 次遍历跟图像多大基本无关一张 4000×3000 的图跑一遍阈值搜索也就毫秒级。真正拖慢流程的是读图和写图。所以它特别适合放在流水线的第一道先把大部分质量正常的图处理掉剩下光照不均的那一小撮再交给分块 Otsu 或局部方法兜底——这个思路在后面二次开发那章会展开。3. 核心模块拆解从灰度直方图到类间方差最大3.1 灰度统计把二维图像压成一维直方图Otsu 的一切推导都建立在直方图上图像的空间信息在算阈值这一步被完全丢弃。这不是缺陷是设计它假设前景和背景的灰度分布各自成峰只要找到把两个峰分得最开的那个点就行。import numpy as np def gray_hist(image, levels256): 把图像压成灰度直方图。 image: PIL Image 或 numpy 数组任意模式 levels: 灰度级数uint8 图固定 256 返回 (灰度数组, 长度 256 的计数数组) arr np.asarray(image) # 彩图按 ITU-R BT.601 权重转灰度避免直接取平均导致对比度被压平 if arr.ndim 3: arr arr[..., :3].astype(np.float64) arr arr np.array([0.299, 0.587, 0.114]) arr np.clip(arr, 0, 255).astype(np.uint8) # bincount 比 np.histogram 快一个量级且天然给出整数计数 hist np.bincount(arr.ravel(), minlengthlevels) return arr, hist[:levels]这段里有两个容易被忽略的点。灰度化权重用的是 0.299/0.587/0.114而不是(RGB)/3后者对绿色的敏感度不足扫描件里的印章、荧光笔痕迹容易被抹平成背景色。np.bincount是给非负整数数组做计数的专用函数比np.histogram少一层区间边界计算大图上差距明显。参数levels留了口子如果你为了抗噪把图降采样到 128 级甚至 64 级这里传进去就行后续阈值搜索会自动跟着变。我一般不改256 级的信息量对这种统计型算法够用了。3.2 阈值搜索一个循环扫完 256 个候选点这是整个项目最核心的二十行。Otsu 的目标函数是类间方差把它写成关于阈值 k 的函数然后对 k 从 0 到 255 全扫一遍取最大。def otsu_threshold(hist, totalNone): 经典大津法最大化类间方差 sigma_b^2(k)。 hist: 长度 L 的灰度计数数组 total: 像素总数不传则用 hist 求和 返回最佳阈值整数0..L-1 hist np.asarray(hist, dtypenp.float64) if total is None: total hist.sum() L hist.size idx np.arange(L) p hist / total # 归一化直方图p[i] 灰度 i 的概率 omega np.cumsum(p) # omega[k] P(灰度 k)前景累计概率 mu np.cumsum(p * idx) # mu[k] 前 k 级的灰度加权和 mu_t mu[-1] # 全图平均灰度 denom omega * (1.0 - omega) # 分母 omega0 * omega1 with np.errstate(divideignore, invalidignore): sigma_b (mu_t * omega - mu) ** 2 / denom # k0 或 kL-1 时分母为 0属于无效分割直接淘汰 sigma_b[~np.isfinite(sigma_b)] -1.0 return int(np.argmax(sigma_b))逻辑拆开看是这样omega[k]是前景灰度 ≤ k的像素占比1 - omega[k]就是背景占比mu[k]是前景的灰度均值乘以占比所以前景真实均值是mu[k]/omega[k]背景真实均值是(mu_t - mu[k])/(1 - omega[k])。把类间方差按「两类均值差的平方乘以两类占比之积」展开并约分就得到上面那行(mu_t * omega - mu) ** 2 / denom。相比老老实实写双重循环分别算前景均值、背景均值这个形式是纯向量化的cumsum一次遍历就出全部前缀和。sigma_b[~np.isfinite(...)] -1.0这句是必须的。k0意味着全图都判为背景k255意味着全图都判为前景这两种情况分母为 0会产出nan或inf。np.argmax遇到nan的行为在不同 NumPy 版本里不完全一致直接置为 -1 把它踢出竞争最稳妥。这是我踩过一次坑才加上的某批全白底的图跑出来阈值是 255整张图黑成一块查了两小时才发现是这里放了nan进去。total参数暴露出来是有意的。做分块 Otsu 时每块的像素数不一样传进去能省一次求和做增量更新时甚至可以只更新直方图再重算阈值图像本身都不用重新遍历。3.3 阈值回写一个比较符号决定成败def binarize(gray, thresh, high255, low0): 按阈值把灰度图切成二值图。 gray: uint8 灰度数组 thresh: 阈值灰度 thresh 记为 high否则 low g np.asarray(gray) if g.dtype ! np.uint8: g np.clip(g, 0, 255).astype(np.uint8) # 用 而不是 与主流库的默认行为保持一致避免差 1 的对齐问题 return np.where(g thresh, high, low).astype(np.uint8)np.where会生成一个新数组内存占用是原图的两倍左右超大图要注意。真要省内存可以原地写g[g thresh] 0; g[g thresh] 255但这样会破坏输入调试阶段不建议。比较符号这里用的是也就是「大于阈值算前景」。换成结果会整体偏移一个灰度级视觉上看不出来但拿去做像素级对比校验就会报「不一致」。这个细节在 3.4 节会专门说。high/low两个参数是为了兼容反色场景。有些票据是深底浅字前景应该是黑色把high0, low255传进去即可不用改算法。3.4 与成熟库对齐差 1 到底怪谁自己实现完第一件事是拿cv2.threshold做对齐。写个校验脚本把两张图逐像素比import cv2 import numpy as np from otsu import gray_hist, otsu_threshold, binarize def check(path): img cv2.imread(path, cv2.IMREAD_GRAYSCALE) _, hist gray_hist(img) my_k otsu_threshold(hist) cv_k, cv_out cv2.threshold(img, 0, 255, cv2.THRESH_BINARY | cv2.THRESH_OTSU) my_out binarize(img, my_k) diff int(np.count_nonzero(my_out ! cv_out)) print(f{path}: mine{my_k} cv{int(cv_k)} pixel_diff{diff} f({diff / img.size:.4%})) return diff跑下来常见的三种结果对应三种不同的原因第一种阈值完全一致、像素差 0。说明图是标准双峰两边的实现路径虽然不同但收敛到同一个点。第二种阈值差 1、像素差极小万分之一量级。这是最正常的原因往往在直方图桶的归属上有的实现把灰度 i 归到「背景侧」有的归到「前景侧」等价于阈值整体平移一格。只要差异集中在这个量级不用管。第三种阈值差十几甚至几十、像素差百分之几。这就不是实现细节了得查输入cv2.imread读出来是 BGR 顺序如果你的代码里又手动做了一次 RGB 转灰度通道就串了或者一边用了INTER_AREA缩放、另一边用了默认插值直方图分布已经变了。这类问题在缩放后的图上特别容易出现我一般的做法是校验脚本里绝不缩放、绝不裁切原图进原图出。4. 跑起来环境、命令行参数与批量实验4.1 安装与最小可运行环境没什么讲究Python 3.8 以上都行虚拟环境建好之后两条命令python -m venv .venv source .venv/bin/activate # Windows 用 .venv\Scripts\activate pip install numpy pillow # 主链路 pip install opencv-python matplotlib # 跑对照实验和画直方图时才需要最小可运行只要一行输出图和阈值会一起打印出来python main.py --input samples/doc_lowlight.png --output out/跑完终端会看到类似threshold113 foreground0.284 savedout/doc_lowlight_bin.png的东西。foreground是前景像素占比这个数字比阈值本身更有诊断价值文档类图像一般在 0.05 到 0.35 之间如果跑出来是 0.9 或者 0.01几乎可以断定这张图不适合用全局 Otsu别在算法上死磕。提示第一次跑建议先拿samples/里的图别直接上手自己的业务图。样例里刻意放了光照不均、低对比度、纯色背景三类的图跑一遍就知道这套东西的边界在哪。4.2 命令行参数与入口代码main.py用 argparse参数不多但每一个都有明确用途参数类型默认值说明--inputstr必填单图路径或配合--batch传目录--outputstrout/结果目录不存在会自动创建--batchflag关开启后--input当目录处理--invertflag关反色输出深底浅字的图用--save-histflag关把直方图和阈值位置画成 png调试用--quietflag关只输出阈值不打印每张的明细入口逻辑很短核心是把参数翻译成函数调用import argparse from pathlib import Path from io_utils import load_gray, save_image from otsu import gray_hist, otsu_threshold, binarize def main(): ap argparse.ArgumentParser() ap.add_argument(--input, requiredTrue) ap.add_argument(--output, defaultout/) ap.add_argument(--batch, actionstore_true) ap.add_argument(--invert, actionstore_true) ap.add_argument(--quiet, actionstore_true) args ap.parse_args() out_dir Path(args.output) out_dir.mkdir(parentsTrue, exist_okTrue) # 单图和批处理共用一套流程只是文件列表来源不同 files sorted(Path(args.input).glob(*)) if args.batch else [Path(args.input)] for f in files: gray load_gray(f) # 统一转 uint8 灰度内部处理 EXIF 方向 _, hist gray_hist(gray) k otsu_threshold(hist) high, low (0, 255) if args.invert else (255, 0) out binarize(gray, k, highhigh, lowlow) save_image(out, out_dir / f{f.stem}_bin.png) if not args.quiet: fg float((out low).mean()) print(f{f.name}: threshold{k} foreground{fg:.3f})几个设计上的考虑。load_gray()把读图和灰度化收在一处EXIF 旋转、16 位图转 8 位、调色板图转灰度这些杂活都在里面处理掉主流程就干净了。批处理模式下用glob(*)而不是glob(*.png)是因为实际业务里目录经常混着 jpg 和 tif靠load_gray内部判断格式更省事。foreground的算法是(out low).mean()注意这里数的是前景被置为low的那部分。反色模式下low255这个统计依然是对的因为「前景」的定义跟着low走了。这种小地方不统一后面加日志分析时会给自己挖坑。4.3 实验脚本把 Otsu 和固定阈值摆在一起比项目标题里那句「实验证明很好」落到代码上就是这个脚本。它不追求严谨的学术评测只做一件朴素的事同一张图上跑两种阈值把阈值和前景占比打出来让人一眼看出固定阈值偏了多少。import numpy as np from pathlib import Path from otsu import gray_hist, otsu_threshold, binarize from io_utils import load_gray def compare(path, fixed127): gray load_gray(path) _, hist gray_hist(gray) k otsu_threshold(hist) out_otsu binarize(gray, k) out_fix binarize(gray, fixed) # 前景占比固定阈值偏离真实分布时这个数字会明显异常 fg_otsu float((out_otsu 0).mean()) fg_fix float((out_fix 0).mean()) # 两图不一致的像素比例衡量固定阈值错判了多少 mismatch float((out_otsu ! out_fix).mean()) return dict(namepath.name, kk, fg_otsufg_otsu, fg_fixfg_fix, mismatchmismatch) if __name__ __main__: rows [compare(p) for p in sorted(Path(samples).glob(*.png))] print(f{file:28}{otsu:6}{fg_otsu:9}{fg_fix:9}{mismatch:10}) for r in rows: print(f{r[name]:28}{r[k]:6}{r[fg_otsu]:9.3f} f{r[fg_fix]:9.3f}{r[mismatch]:10.3f})mismatch是这个脚本里最有价值的指标。它算的是「如果当初图省事用固定 127会有多大比例的像素被切错」。样例里的光照不均图跑出来这个值经常在 0.15 以上也就是 15% 的像素判决不同这个量级足以让 OCR 直接读不出字。脚本故意不引入 PSNR、IoU 这些需要标注真值的指标。二值化没有标准答案你说不清哪个像素「应该」是黑是白硬套指标反而误导。真正有意义的对比永远是下游任务的表现同一批图Otsu 出的二值图喂给 OCR识别率涨了几个点那才是结论。跑完这个脚本你对参数就有手感了。比如fixed从 127 调到 100、160 各跑一遍看mismatch怎么变——如果某张图在 100 到 160 之间mismatch都很小说明它本身就双峰明显用什么阈值都差不多如果调一格就剧烈变化那张图就是典型的「临界态」必须上自动阈值。5. 避坑与排查那些让阈值飘掉的细节5.1 现象阈值跑到 0 或 255整张图全黑或全白原因基本只有一个直方图是单峰的全图灰度集中在很窄的一段。纯色背景、过曝的图、全黑的图都会这样。这时候类间方差的最大值出现在端点而端点又被 3.2 节那句sigma_b[~np.isfinite(...)] -1.0排除了argmax只能在次优的无效点上随便挑一个。解决分两层。第一层是入口拦截算完foreground之后加个判断占比低于 0.01 或高于 0.99 就直接标记为「不适合全局 Otsu」走回退分支比如保留固定阈值 127或者直接判为异常图打回重拍。第二层是给otsu_threshold()加一个返回值把最大类间方差和次大值一起吐出来两者比值小于某个经验值比如 1.05就认为峰不够分明。这个判据比单纯看前景占比可靠。5.2 现象彩色图跑出来一片糊边缘全是噪点原因是根本没转灰度但代码「没报错」。如果你的binarize()收到的是三维数组np.where(g thresh, 255, 0)会逐通道独立比较最后得到一个三通道的「伪二值图」肉眼看是彩色的黑白片存成 jpg 之后又会被压缩算法掺进中间色。解决办法是在load_gray()里强校验assert arr.ndim 2不是二维就直接抛异常别给下游留幻想。另外 OpenCV 用户要特别注意cv2.imread默认是 BGR 顺序如果你在转灰度前先做了颜色判断比如「红色印章区域跳过」通道顺序错了会把蓝色当成红色这个坑我在别的项目里踩过排查了半天。5.3 现象和 OpenCV 结果差 1逐像素对比报不一致前面 3.4 节提过根因是比较符号。cv2.threshold的THRESH_BINARY判定是src threshTHRESH_BINARY_INV是src thresh不同版本、不同 flag 组合下边界归属会有差异。你自己实现时用还是直接决定差 1 还是差 0。解决方式不是去猜哪个「对」而是把比较符号提成参数校验脚本里两个都试一遍取像素差小的那个。另外cv2.threshold返回的ret是float打印时int()截断会掩盖小数部分——如果你的实现返回整数而库返回 117.0那多半是库内部用了浮点直方图不影响结果别在这上面浪费时间。5.4 现象光照不均的图左半边全黑、右半边全白这是全局阈值的固有限制不是 bug。台灯照亮的区域整体灰度抬高 40 到 60暗部整体压低两个区域各自的双峰被叠在一起直方图变成一个宽包Otsu 找不到中间点切出来的结果就是「一半一半」。解决靠预处理或者分块。预处理方向是先把不均匀的光照压平形态学顶帽cv2.morphologyEx(gray, cv2.MORPH_TOPHAT, kernel)能扣掉大面积的亮背景CLAHE 做局部对比度均衡也有效两者叠加对扫描件效果很明显。分块方向就是 6.1 节的 tile Otsu块大小取 64 到 128重叠 8 到 16 像素。两条路我都试过实测下来「先 CLAHE 再全局 Otsu」的代码量最小、参数最少优先试这个。5.5 现象大图算出来的阈值和别人不一样小图却对得上原因是累加溢出或者精度丢失。np.bincount出来的计数是int64安全但如果你自己写了hist np.zeros(256, dtypenp.uint8)然后循环累加像素数一超过 255 就回绕直方图整个是错的。另一种情况是用float32做cumsum几百万像素累加下来末位精度会飘argmax在平坦区段挑到的点可能差好几格。解决办法很简单otsu_threshold()开头那句hist.astype(np.float64)就是干这个的别省。自己写统计代码时也统一用np.int64计数、np.float64算概率。这类问题最恶心的地方在于小图上完全复现不了等上线跑了批大图才发现那时候日志里只有阈值根本看不出是直方图错了。6. 二次开发把全局 Otsu 改成能落地分块版本全局 Otsu 的短板很明确改造方向也就那么几条。最见效的是分块其次是预处理再往上才是多阈值扩展。分块 Otsu 的核心是「每块自己算阈值重叠区做加权平均」。直接拼块会看到明显的方格接缝加重叠就是为了把接缝抹平import numpy as np from otsu import otsu_threshold def tile_otsu(gray, tile96, overlap16): 分块大津法逐块求阈值并二值化重叠区累加后取平均。 tile: 块边长建议 64~128太小会放大小噪声 overlap: 相邻块重叠像素建议取 tile 的 1/8 左右 h, w gray.shape step max(tile - overlap, 1) acc np.zeros((h, w), dtypenp.float32) # 二值结果累加 cnt np.zeros((h, w), dtypenp.float32) # 每像素被覆盖次数 for y in range(0, h, step): for x in range(0, w, step): y2, x2 min(y tile, h), min(x tile, w) patch gray[y:y2, x:x2] if patch.size tile * tile // 4: # 边角碎块跳过避免阈值不稳 continue hist np.bincount(patch.ravel(), minlength256) k otsu_threshold(hist) acc[y:y2, x:x2] np.where(patch k, 255.0, 0.0) cnt[y:y2, x:x2] 1.0 # 没被任何块覆盖的像素理论上不该有按背景处理 return (acc / np.maximum(cnt, 1.0)).astype(np.uint8)tile和overlap这两个参数是调出来的不是想出来的。tile太小比如 32块内像素太少直方图统计不稳空白块容易整块判成前景出来的图一片麻点tile太大比如 512又退化回全局阈值光照不均的图照样救不回来。overlap决定接缝的平滑程度取tile的八分之一到四分之一之间比较合适太小接缝明显太大计算量翻倍。性能上这个双层循环是纯 Python 的一块一块跑2000×3000 的图用 96 的块大概几百次迭代单张耗时在几十毫秒量级批量处理够用。真要压榨可以用积分图把每块的直方图做成增量更新或者直接用cv2.integral加速块内求和。不过说实话这种场景下读图才是瓶颈算法优化到毫秒级意义不大我一般先保证结果对性能等真成了瓶颈再说。另外两条改造路线可以顺带记一下。预处理路线cv2.createCLAHE(clipLimit2.0, tileGridSize(8,8))先把光照压平再走全局 Otsu代码三行效果对文档类图像很稳。多阈值路线把otsu_threshold()里的单阈值搜索扩展成两类以上复杂度从 O(L) 涨到 O(L²) 甚至更高L256 时还扛得住但收益主要体现在有三类以上目标的场景比如票据上的表格线、印刷字、印章普通的黑白分离用不上。判断标准很简单先看全局 Otsu 的foreground是否落在合理区间落进去了就别折腾落不进去再上分块。这套东西我现在的习惯是任何图像预处理脚本跑完必须打印三样东西——阈值、前景占比、处理耗时。阈值告诉你算法在想什么前景占比告诉你这张图适不适合这个算法耗时告诉你批处理要跑多久。三样少一样等出问题回头查日志就是抓瞎。吃过一次亏之后我每次给新图源接管道都强制先跑一遍experiment.py看mismatch分布确认那批图不是「临界态」再往下走。希望帮到你。本文还有配套的精品资源点击获取 用NumPy实现大津法二值化：从灰度直方图到分块Otsu实践
简介面向图像处理初学者与需要快速实现阈值分割的开发者这份大津法二值化程序以MATLAB脚本形式给出聚焦灰度图像中前景与背景的分离问题可用于文本识别、医学图像分析、工业检测等实验场景也适合课程练习和算法入门对照。压缩包为rar格式共1个文件主体是1个m脚本包体约1KB体量轻便于直接阅读、运行和改造。已有396人学习说明其在基础图像处理练习中具备一定参考价值。代码围绕灰度直方图统计、阈值遍历、前景与背景均值及类间方差计算展开读者可据此理解最佳全局阈值的选取过程并复现二值化效果也可在脚本基础上替换输入图像、调整参数观察噪声、光照和图像质量对分割结果的影响积累算法调试与效果对比经验。1. 一批发票扫描件让我重新捡起大津法前段时间接手一批票据扫描件的预处理开头想省事直接gray 127一刀切。结果一半的图翻了车左侧被台灯照得发白数字糊成一片右边背光的角落又黑得连边框都吃掉。同一批图换个角度看阈值就得跟着改改到最后我自己都记不清哪张用的哪个数。后来把大津法二值化程序翻出来重跑一遍让它自己从灰度直方图里算阈值那批图基本一次过。这份程序做的事很朴素读图、转灰度、统计直方图、遍历 256 个候选阈值找类间方差最大的那个、回写二值图再配一个实验脚本把 Otsu 和固定阈值的结果摆在一起对比。核心算法部分只用 NumPy不装 OpenCV 也能跑完整流程OpenCV 只在对照验证时出场。适合两类人做文档扫描、OCR 前处理、工业检测预处理的开发者拿来就能接进流水线以及想把 Otsu 从公式变成能调试的代码、而不是永远调cv2.threshold的人。2. 项目定位与技术栈核心算法只留 NumPy 的取舍2.1 一份目录看全项目结构拿到资源先别急着跑先看目录。这个项目的结构非常扁平没有花哨的包管理属于「拷到哪都能跑」的类型目录长这样otsu-binarize/ ├── otsu.py # 核心算法直方图、阈值搜索、二值化 ├── main.py # 单图命令行入口 ├── experiment.py # 批量对比实验Otsu vs 固定阈值 ├── io_utils.py # 读写图、灰度化、结果落盘 ├── samples/ # 实验用图含光照不均、低对比度几类 ├── out/ # 运行产物跑完自动创建 └── requirements.txtotsu.py是全部价值所在一百行出头没有类、没有继承、没有配置对象三个纯函数gray_hist()出直方图otsu_threshold()出阈值binarize()出结果图。main.py只负责把命令行参数翻译成这三个函数的调用experiment.py负责批量跑并打印指标。这种拆法有两个好处。一是算法可被单独导入你把它塞进 Flask 接口、塞进批处理脚本、塞进 PyQt 界面都不用改一行二是出问题时定位极快阈值不对就只可能是otsu_threshold()的输入不对跟读图、跟保存格式无关。我一般接手别人的图像项目第一件事就是把「算法」和「I/O」这两层剥开剥不开的项目后面改起来基本都是血泪。2.2 技术栈与依赖为什么核心不引 OpenCV依赖清单刻意做得很轻主链路只依赖 NumPy 和 Pillow依赖用途是否必需numpy直方图统计、累积和、方差计算、数组运算必需Pillow读图、转灰度、保存结果必需opencv-python仅用于experiment.py里做结果对齐校验可选matplotlib画直方图和阈值位置调试用可选不把 OpenCV 放进主链路不是为了标新立异而是为了让阈值怎么来的这件事完全暴露在你眼前。cv2.threshold(img, 0, 255, cv2.THRESH_OTSU)一行出结果方便但你不知道它内部是用 256 桶直方图还是 128 桶不知道它遇到单峰图怎么处理出了偏差只能靠猜。自己实现一遍之后ret那个返回值是怎么算出来的、为什么有时候是 117 有时候是 118你心里就有数了。常见做法是主链路保持纯 NumPy把 OpenCV 放到验证脚本里当参照物。这样跑实验时能立刻发现「我的实现和成熟库差多少」而不是等上线后被业务方发现。2.3 大津法在二值化方案里的位置很多人一上来就问「Otsu 是不是最好的」这个问题没有统一答案得看你图的性质。把这个项目和几种常见方案放一起对比选择标准就清楚了方案阈值来源适用场景主要短板固定阈值人工指定光照稳定、批量一致的工业相机图换批次就得重调玄学感强Otsu本项目全图直方图最大化类间方差双峰明显、前景背景面积比不太悬殊光照不均、单峰图会失效自适应均值局部窗口均值减常数 C文档扫描、背景渐变空白区域容易噪点化自适应高斯局部窗口高斯加权均值减 C同上边缘更平滑计算量大块大小敏感Sauvola 等局部方法局部均值与标准差低对比度文档、退化扫描件参数多调起来费时间Otsu 的定位是「全局、自动、够快」。它的计算复杂度只跟灰度级数有关256 次遍历跟图像多大基本无关一张 4000×3000 的图跑一遍阈值搜索也就毫秒级。真正拖慢流程的是读图和写图。所以它特别适合放在流水线的第一道先把大部分质量正常的图处理掉剩下光照不均的那一小撮再交给分块 Otsu 或局部方法兜底——这个思路在后面二次开发那章会展开。3. 核心模块拆解从灰度直方图到类间方差最大3.1 灰度统计把二维图像压成一维直方图Otsu 的一切推导都建立在直方图上图像的空间信息在算阈值这一步被完全丢弃。这不是缺陷是设计它假设前景和背景的灰度分布各自成峰只要找到把两个峰分得最开的那个点就行。import numpy as np def gray_hist(image, levels256): 把图像压成灰度直方图。 image: PIL Image 或 numpy 数组任意模式 levels: 灰度级数uint8 图固定 256 返回 (灰度数组, 长度 256 的计数数组) arr np.asarray(image) # 彩图按 ITU-R BT.601 权重转灰度避免直接取平均导致对比度被压平 if arr.ndim 3: arr arr[..., :3].astype(np.float64) arr arr np.array([0.299, 0.587, 0.114]) arr np.clip(arr, 0, 255).astype(np.uint8) # bincount 比 np.histogram 快一个量级且天然给出整数计数 hist np.bincount(arr.ravel(), minlengthlevels) return arr, hist[:levels]这段里有两个容易被忽略的点。灰度化权重用的是 0.299/0.587/0.114而不是(RGB)/3后者对绿色的敏感度不足扫描件里的印章、荧光笔痕迹容易被抹平成背景色。np.bincount是给非负整数数组做计数的专用函数比np.histogram少一层区间边界计算大图上差距明显。参数levels留了口子如果你为了抗噪把图降采样到 128 级甚至 64 级这里传进去就行后续阈值搜索会自动跟着变。我一般不改256 级的信息量对这种统计型算法够用了。3.2 阈值搜索一个循环扫完 256 个候选点这是整个项目最核心的二十行。Otsu 的目标函数是类间方差把它写成关于阈值 k 的函数然后对 k 从 0 到 255 全扫一遍取最大。def otsu_threshold(hist, totalNone): 经典大津法最大化类间方差 sigma_b^2(k)。 hist: 长度 L 的灰度计数数组 total: 像素总数不传则用 hist 求和 返回最佳阈值整数0..L-1 hist np.asarray(hist, dtypenp.float64) if total is None: total hist.sum() L hist.size idx np.arange(L) p hist / total # 归一化直方图p[i] 灰度 i 的概率 omega np.cumsum(p) # omega[k] P(灰度 k)前景累计概率 mu np.cumsum(p * idx) # mu[k] 前 k 级的灰度加权和 mu_t mu[-1] # 全图平均灰度 denom omega * (1.0 - omega) # 分母 omega0 * omega1 with np.errstate(divideignore, invalidignore): sigma_b (mu_t * omega - mu) ** 2 / denom # k0 或 kL-1 时分母为 0属于无效分割直接淘汰 sigma_b[~np.isfinite(sigma_b)] -1.0 return int(np.argmax(sigma_b))逻辑拆开看是这样omega[k]是前景灰度 ≤ k的像素占比1 - omega[k]就是背景占比mu[k]是前景的灰度均值乘以占比所以前景真实均值是mu[k]/omega[k]背景真实均值是(mu_t - mu[k])/(1 - omega[k])。把类间方差按「两类均值差的平方乘以两类占比之积」展开并约分就得到上面那行(mu_t * omega - mu) ** 2 / denom。相比老老实实写双重循环分别算前景均值、背景均值这个形式是纯向量化的cumsum一次遍历就出全部前缀和。sigma_b[~np.isfinite(...)] -1.0这句是必须的。k0意味着全图都判为背景k255意味着全图都判为前景这两种情况分母为 0会产出nan或inf。np.argmax遇到nan的行为在不同 NumPy 版本里不完全一致直接置为 -1 把它踢出竞争最稳妥。这是我踩过一次坑才加上的某批全白底的图跑出来阈值是 255整张图黑成一块查了两小时才发现是这里放了nan进去。total参数暴露出来是有意的。做分块 Otsu 时每块的像素数不一样传进去能省一次求和做增量更新时甚至可以只更新直方图再重算阈值图像本身都不用重新遍历。3.3 阈值回写一个比较符号决定成败def binarize(gray, thresh, high255, low0): 按阈值把灰度图切成二值图。 gray: uint8 灰度数组 thresh: 阈值灰度 thresh 记为 high否则 low g np.asarray(gray) if g.dtype ! np.uint8: g np.clip(g, 0, 255).astype(np.uint8) # 用 而不是 与主流库的默认行为保持一致避免差 1 的对齐问题 return np.where(g thresh, high, low).astype(np.uint8)np.where会生成一个新数组内存占用是原图的两倍左右超大图要注意。真要省内存可以原地写g[g thresh] 0; g[g thresh] 255但这样会破坏输入调试阶段不建议。比较符号这里用的是也就是「大于阈值算前景」。换成结果会整体偏移一个灰度级视觉上看不出来但拿去做像素级对比校验就会报「不一致」。这个细节在 3.4 节会专门说。high/low两个参数是为了兼容反色场景。有些票据是深底浅字前景应该是黑色把high0, low255传进去即可不用改算法。3.4 与成熟库对齐差 1 到底怪谁自己实现完第一件事是拿cv2.threshold做对齐。写个校验脚本把两张图逐像素比import cv2 import numpy as np from otsu import gray_hist, otsu_threshold, binarize def check(path): img cv2.imread(path, cv2.IMREAD_GRAYSCALE) _, hist gray_hist(img) my_k otsu_threshold(hist) cv_k, cv_out cv2.threshold(img, 0, 255, cv2.THRESH_BINARY | cv2.THRESH_OTSU) my_out binarize(img, my_k) diff int(np.count_nonzero(my_out ! cv_out)) print(f{path}: mine{my_k} cv{int(cv_k)} pixel_diff{diff} f({diff / img.size:.4%})) return diff跑下来常见的三种结果对应三种不同的原因第一种阈值完全一致、像素差 0。说明图是标准双峰两边的实现路径虽然不同但收敛到同一个点。第二种阈值差 1、像素差极小万分之一量级。这是最正常的原因往往在直方图桶的归属上有的实现把灰度 i 归到「背景侧」有的归到「前景侧」等价于阈值整体平移一格。只要差异集中在这个量级不用管。第三种阈值差十几甚至几十、像素差百分之几。这就不是实现细节了得查输入cv2.imread读出来是 BGR 顺序如果你的代码里又手动做了一次 RGB 转灰度通道就串了或者一边用了INTER_AREA缩放、另一边用了默认插值直方图分布已经变了。这类问题在缩放后的图上特别容易出现我一般的做法是校验脚本里绝不缩放、绝不裁切原图进原图出。4. 跑起来环境、命令行参数与批量实验4.1 安装与最小可运行环境没什么讲究Python 3.8 以上都行虚拟环境建好之后两条命令python -m venv .venv source .venv/bin/activate # Windows 用 .venv\Scripts\activate pip install numpy pillow # 主链路 pip install opencv-python matplotlib # 跑对照实验和画直方图时才需要最小可运行只要一行输出图和阈值会一起打印出来python main.py --input samples/doc_lowlight.png --output out/跑完终端会看到类似threshold113 foreground0.284 savedout/doc_lowlight_bin.png的东西。foreground是前景像素占比这个数字比阈值本身更有诊断价值文档类图像一般在 0.05 到 0.35 之间如果跑出来是 0.9 或者 0.01几乎可以断定这张图不适合用全局 Otsu别在算法上死磕。提示第一次跑建议先拿samples/里的图别直接上手自己的业务图。样例里刻意放了光照不均、低对比度、纯色背景三类的图跑一遍就知道这套东西的边界在哪。4.2 命令行参数与入口代码main.py用 argparse参数不多但每一个都有明确用途参数类型默认值说明--inputstr必填单图路径或配合--batch传目录--outputstrout/结果目录不存在会自动创建--batchflag关开启后--input当目录处理--invertflag关反色输出深底浅字的图用--save-histflag关把直方图和阈值位置画成 png调试用--quietflag关只输出阈值不打印每张的明细入口逻辑很短核心是把参数翻译成函数调用import argparse from pathlib import Path from io_utils import load_gray, save_image from otsu import gray_hist, otsu_threshold, binarize def main(): ap argparse.ArgumentParser() ap.add_argument(--input, requiredTrue) ap.add_argument(--output, defaultout/) ap.add_argument(--batch, actionstore_true) ap.add_argument(--invert, actionstore_true) ap.add_argument(--quiet, actionstore_true) args ap.parse_args() out_dir Path(args.output) out_dir.mkdir(parentsTrue, exist_okTrue) # 单图和批处理共用一套流程只是文件列表来源不同 files sorted(Path(args.input).glob(*)) if args.batch else [Path(args.input)] for f in files: gray load_gray(f) # 统一转 uint8 灰度内部处理 EXIF 方向 _, hist gray_hist(gray) k otsu_threshold(hist) high, low (0, 255) if args.invert else (255, 0) out binarize(gray, k, highhigh, lowlow) save_image(out, out_dir / f{f.stem}_bin.png) if not args.quiet: fg float((out low).mean()) print(f{f.name}: threshold{k} foreground{fg:.3f})几个设计上的考虑。load_gray()把读图和灰度化收在一处EXIF 旋转、16 位图转 8 位、调色板图转灰度这些杂活都在里面处理掉主流程就干净了。批处理模式下用glob(*)而不是glob(*.png)是因为实际业务里目录经常混着 jpg 和 tif靠load_gray内部判断格式更省事。foreground的算法是(out low).mean()注意这里数的是前景被置为low的那部分。反色模式下low255这个统计依然是对的因为「前景」的定义跟着low走了。这种小地方不统一后面加日志分析时会给自己挖坑。4.3 实验脚本把 Otsu 和固定阈值摆在一起比项目标题里那句「实验证明很好」落到代码上就是这个脚本。它不追求严谨的学术评测只做一件朴素的事同一张图上跑两种阈值把阈值和前景占比打出来让人一眼看出固定阈值偏了多少。import numpy as np from pathlib import Path from otsu import gray_hist, otsu_threshold, binarize from io_utils import load_gray def compare(path, fixed127): gray load_gray(path) _, hist gray_hist(gray) k otsu_threshold(hist) out_otsu binarize(gray, k) out_fix binarize(gray, fixed) # 前景占比固定阈值偏离真实分布时这个数字会明显异常 fg_otsu float((out_otsu 0).mean()) fg_fix float((out_fix 0).mean()) # 两图不一致的像素比例衡量固定阈值错判了多少 mismatch float((out_otsu ! out_fix).mean()) return dict(namepath.name, kk, fg_otsufg_otsu, fg_fixfg_fix, mismatchmismatch) if __name__ __main__: rows [compare(p) for p in sorted(Path(samples).glob(*.png))] print(f{file:28}{otsu:6}{fg_otsu:9}{fg_fix:9}{mismatch:10}) for r in rows: print(f{r[name]:28}{r[k]:6}{r[fg_otsu]:9.3f} f{r[fg_fix]:9.3f}{r[mismatch]:10.3f})mismatch是这个脚本里最有价值的指标。它算的是「如果当初图省事用固定 127会有多大比例的像素被切错」。样例里的光照不均图跑出来这个值经常在 0.15 以上也就是 15% 的像素判决不同这个量级足以让 OCR 直接读不出字。脚本故意不引入 PSNR、IoU 这些需要标注真值的指标。二值化没有标准答案你说不清哪个像素「应该」是黑是白硬套指标反而误导。真正有意义的对比永远是下游任务的表现同一批图Otsu 出的二值图喂给 OCR识别率涨了几个点那才是结论。跑完这个脚本你对参数就有手感了。比如fixed从 127 调到 100、160 各跑一遍看mismatch怎么变——如果某张图在 100 到 160 之间mismatch都很小说明它本身就双峰明显用什么阈值都差不多如果调一格就剧烈变化那张图就是典型的「临界态」必须上自动阈值。5. 避坑与排查那些让阈值飘掉的细节5.1 现象阈值跑到 0 或 255整张图全黑或全白原因基本只有一个直方图是单峰的全图灰度集中在很窄的一段。纯色背景、过曝的图、全黑的图都会这样。这时候类间方差的最大值出现在端点而端点又被 3.2 节那句sigma_b[~np.isfinite(...)] -1.0排除了argmax只能在次优的无效点上随便挑一个。解决分两层。第一层是入口拦截算完foreground之后加个判断占比低于 0.01 或高于 0.99 就直接标记为「不适合全局 Otsu」走回退分支比如保留固定阈值 127或者直接判为异常图打回重拍。第二层是给otsu_threshold()加一个返回值把最大类间方差和次大值一起吐出来两者比值小于某个经验值比如 1.05就认为峰不够分明。这个判据比单纯看前景占比可靠。5.2 现象彩色图跑出来一片糊边缘全是噪点原因是根本没转灰度但代码「没报错」。如果你的binarize()收到的是三维数组np.where(g thresh, 255, 0)会逐通道独立比较最后得到一个三通道的「伪二值图」肉眼看是彩色的黑白片存成 jpg 之后又会被压缩算法掺进中间色。解决办法是在load_gray()里强校验assert arr.ndim 2不是二维就直接抛异常别给下游留幻想。另外 OpenCV 用户要特别注意cv2.imread默认是 BGR 顺序如果你在转灰度前先做了颜色判断比如「红色印章区域跳过」通道顺序错了会把蓝色当成红色这个坑我在别的项目里踩过排查了半天。5.3 现象和 OpenCV 结果差 1逐像素对比报不一致前面 3.4 节提过根因是比较符号。cv2.threshold的THRESH_BINARY判定是src threshTHRESH_BINARY_INV是src thresh不同版本、不同 flag 组合下边界归属会有差异。你自己实现时用还是直接决定差 1 还是差 0。解决方式不是去猜哪个「对」而是把比较符号提成参数校验脚本里两个都试一遍取像素差小的那个。另外cv2.threshold返回的ret是float打印时int()截断会掩盖小数部分——如果你的实现返回整数而库返回 117.0那多半是库内部用了浮点直方图不影响结果别在这上面浪费时间。5.4 现象光照不均的图左半边全黑、右半边全白这是全局阈值的固有限制不是 bug。台灯照亮的区域整体灰度抬高 40 到 60暗部整体压低两个区域各自的双峰被叠在一起直方图变成一个宽包Otsu 找不到中间点切出来的结果就是「一半一半」。解决靠预处理或者分块。预处理方向是先把不均匀的光照压平形态学顶帽cv2.morphologyEx(gray, cv2.MORPH_TOPHAT, kernel)能扣掉大面积的亮背景CLAHE 做局部对比度均衡也有效两者叠加对扫描件效果很明显。分块方向就是 6.1 节的 tile Otsu块大小取 64 到 128重叠 8 到 16 像素。两条路我都试过实测下来「先 CLAHE 再全局 Otsu」的代码量最小、参数最少优先试这个。5.5 现象大图算出来的阈值和别人不一样小图却对得上原因是累加溢出或者精度丢失。np.bincount出来的计数是int64安全但如果你自己写了hist np.zeros(256, dtypenp.uint8)然后循环累加像素数一超过 255 就回绕直方图整个是错的。另一种情况是用float32做cumsum几百万像素累加下来末位精度会飘argmax在平坦区段挑到的点可能差好几格。解决办法很简单otsu_threshold()开头那句hist.astype(np.float64)就是干这个的别省。自己写统计代码时也统一用np.int64计数、np.float64算概率。这类问题最恶心的地方在于小图上完全复现不了等上线跑了批大图才发现那时候日志里只有阈值根本看不出是直方图错了。6. 二次开发把全局 Otsu 改成能落地分块版本全局 Otsu 的短板很明确改造方向也就那么几条。最见效的是分块其次是预处理再往上才是多阈值扩展。分块 Otsu 的核心是「每块自己算阈值重叠区做加权平均」。直接拼块会看到明显的方格接缝加重叠就是为了把接缝抹平import numpy as np from otsu import otsu_threshold def tile_otsu(gray, tile96, overlap16): 分块大津法逐块求阈值并二值化重叠区累加后取平均。 tile: 块边长建议 64~128太小会放大小噪声 overlap: 相邻块重叠像素建议取 tile 的 1/8 左右 h, w gray.shape step max(tile - overlap, 1) acc np.zeros((h, w), dtypenp.float32) # 二值结果累加 cnt np.zeros((h, w), dtypenp.float32) # 每像素被覆盖次数 for y in range(0, h, step): for x in range(0, w, step): y2, x2 min(y tile, h), min(x tile, w) patch gray[y:y2, x:x2] if patch.size tile * tile // 4: # 边角碎块跳过避免阈值不稳 continue hist np.bincount(patch.ravel(), minlength256) k otsu_threshold(hist) acc[y:y2, x:x2] np.where(patch k, 255.0, 0.0) cnt[y:y2, x:x2] 1.0 # 没被任何块覆盖的像素理论上不该有按背景处理 return (acc / np.maximum(cnt, 1.0)).astype(np.uint8)tile和overlap这两个参数是调出来的不是想出来的。tile太小比如 32块内像素太少直方图统计不稳空白块容易整块判成前景出来的图一片麻点tile太大比如 512又退化回全局阈值光照不均的图照样救不回来。overlap决定接缝的平滑程度取tile的八分之一到四分之一之间比较合适太小接缝明显太大计算量翻倍。性能上这个双层循环是纯 Python 的一块一块跑2000×3000 的图用 96 的块大概几百次迭代单张耗时在几十毫秒量级批量处理够用。真要压榨可以用积分图把每块的直方图做成增量更新或者直接用cv2.integral加速块内求和。不过说实话这种场景下读图才是瓶颈算法优化到毫秒级意义不大我一般先保证结果对性能等真成了瓶颈再说。另外两条改造路线可以顺带记一下。预处理路线cv2.createCLAHE(clipLimit2.0, tileGridSize(8,8))先把光照压平再走全局 Otsu代码三行效果对文档类图像很稳。多阈值路线把otsu_threshold()里的单阈值搜索扩展成两类以上复杂度从 O(L) 涨到 O(L²) 甚至更高L256 时还扛得住但收益主要体现在有三类以上目标的场景比如票据上的表格线、印刷字、印章普通的黑白分离用不上。判断标准很简单先看全局 Otsu 的foreground是否落在合理区间落进去了就别折腾落不进去再上分块。这套东西我现在的习惯是任何图像预处理脚本跑完必须打印三样东西——阈值、前景占比、处理耗时。阈值告诉你算法在想什么前景占比告诉你这张图适不适合这个算法耗时告诉你批处理要跑多久。三样少一样等出问题回头查日志就是抓瞎。吃过一次亏之后我每次给新图源接管道都强制先跑一遍experiment.py看mismatch分布确认那批图不是「临界态」再往下走。希望帮到你。本文还有配套的精品资源点击获取