数字图像处理11：阈值分割（基本全局阈值处理，Otsu 的最佳全局阈值，图像平滑改善全局阈值处理，图像分块的可变阈值）
阈值分割
从背景中提取物体的一种明显方法是，选择一个将这些模式分开的阈值 T。然后， f ( x , y ) > T f(x,y)>T f(x,y)>T的任何点 (x, y) 称为个对象点，否则该点称为背景点。
灰度阈值的成功与否直接关系到可区分直方图模式的波谷的宽度和深度。
而影响波谷特性的关键因素有以下几点：
(1)波峰间的间隔 (波峰离得越远，分离这些模式的机会越大);
(2)图像中的噪声内容 (模式随噪声的增加而展宽);
(3)物体和背景的相对尺寸;
(4)光源的均匀性;
(5)图像反射特性的均匀性。
基本全局阈值处理
当物体和背景像素的灰度分布十分明显时，可以用适用于整个图像的单个 (全局) 阈值。
实验思想
在大多数应用中，通常图像之间有较大变化，即使全局阈值是种合适的方法，也需要有能对每幅图像自动估计阈值的算法。下面的迭代算法可用于这一目的:
1.为全局阈值 T 选择一个初始估计值。
2.使用用初始的 T 分割该图像。这将产生两组像素:
G1 由灰度值大于 T 的所有像素组成，G2 由所有小于等于 T 的像素组成。
3.对 G1 和 G2 的像素分别计算平均灰度值 (均值)m1 和 m2。
4.计算一个新的阈值: T = (m1 + m2)/2
5.重复步骤 2 到步骤 4，直到连续迭代中的 T 值间的差小于个预定义的参数 dT 为止。
代码
主函数：
%% 运行 main_Global.m clc;
clear;
close all;%% 课本图 10.34
im = imread('Fig1038(a)(noisy_fingerprint).tif'); % 原始图像 uint8
% 得到图像的行数、列数、以及每像素的维数（防止出现RGB图像）
[line,row,v]=size(im);
im = im(:,:,1);
t0 = sum(im(:))/(line*row);
[im1,NK] = Global_Threshold(im,0,t0);
x = 1:256;%% 将结果保存到当前目录下的result文件夹下
imwrite(im, sprintf('result/%s.jpg','2_im1'));
imwrite(im1, sprintf('result/%s.jpg','2_im2'));%% 显示图像
figure(1);
subplot(131); imshow(im); title('原图'); axis on
subplot(133); imshow(im1); title('全局阈值'); axis on
subplot(132); plot(x,NK); title('直方图'); axis on
功能函数：
function [im1,NK] = Global_Threshold(im,dT,T0)[h,w] = size(im);%% 得到直方图和累计直方图
[NK,CH] = myHisteq(im);G1sum = 0;
G2sum = 0;
dTemp = 255;while(dTemp<dT)for i = 1:T0G1sum = G1sum + (i-1)*NK(i);endG1Ave = G1sum/CH(T0);for i = T0+1:256G2sum = G2sum + (i-1)*NK(i);endG2Ave = G2sum/(CH(256)-CH(T0));temp = round((G1Ave+G2Ave)/2);dTemp = abs(temp - T0);T0 = temp;
end
im1 = zeros(h,w);
for i =1:hfor j = 1:wif(im(i,j)>=T0)im1(i,j) = 1;endend
end
im1 = logical(im1);
实验结果
用 Otsu 的最佳全局阈值
Otsu 方法是最佳的，因为它使得类间方差最大化。其基本思想是，适当的阈值化的类就其像素灰度值而言，应当是截然不同的，相反地，就其灰度值而言，给出最佳类间分离的阈值 将是最佳的阈值。除了其最佳性之外，Otsu 方法还有一个重要的特性，即它完全以在一幅图像的直方图上执行计算为基础，直方图是很容易得到的一维阵列。
实验思想
Otsu 算法计算步骤如下:
1.计算输人图像的归一化直方图。使用 pi,i=0,1,2,.,L-1 表示该直方图的各个分量。
2. 对于 k=0, 1,2,.,L-1，计算累积和 P1(k) 和 P2(k), 计算公式为：
P 1 ( k ) = ∑ i = 0 k p i = 1 − P 2 ( k ) P1(k) = \sum_{i=0}^{k} pi = 1 - P2(k) P1(k)=i=0∑kpi=1−P2(k)
3.对于 k=0,1,2…,L-1, 计算累积均值 m1(k) 和 m2(k)。计算公式为：
m 1 ( k ) = 1 P 1 ( k ) ∗ ∑ i = 0 k ( i ∗ p i ) m1(k) = \frac{1}{P1(k)} * \sum_{i=0}^{k} (i*pi) m1(k)=P1(k)1∗i=0∑k(i∗pi) m 2 ( k ) = 1 P 2 ( k ) ∗ ∑ i = k + 1 L − 1 ( i ∗ p i ) m2(k) = \frac{1}{P2(k)} * \sum_{i=k+1}^{L-1} (i*pi) m2(k)=P2(k)1∗i=k+1∑L−1(i∗pi)
4.计算全局灰度均值 mg。计算公式为： m g = ∑ i = 0 L − 1 ( i ∗ p i ) mg = \sum_{i=0}^{L-1} (i*pi) mg=i=0∑L−1(i∗pi)
5.对于 k=0,1,2…,L-1, 计算类间方差 b 2 ( k ) b^2(k) b2(k)。计算公式：
b 2 ( k ) = P 1 ( k ) ∗ ( m 1 ( k ) − m g ) 2 + P 2 ( k ) ∗ ( m 2 ( k ) − m g ) 2 b^2(k) = P1(k)*(m1(k)-mg)^2 + P2(k)*(m2(k)-mg)^2 b2(k)=P1(k)∗(m1(k)−mg)2+P2(k)∗(m2(k)−mg)2
6.得到 Otsu 阈值 k，即使得 b2(k) 最大的 k 值。如果最大值不唯一，用相应检测到的各个值 k 的平均得到 k*。
代码
主函数：
%% 运行 main_Otsu.m clc;
clear;
close all;%% 课本图 10.39
im = imread('Fig1039(a)(polymersomes).tif'); % 原始图像 uint8
% 得到图像的行数、列数、以及每像素的维数（防止出现RGB图像）
[line,row,v]=size(im);
im = im(:,:,1);[im1,NK] = Global_Threshold(im,0,sum(im(:))/(line*row));
[im2,NK] = my_Otsu(im);
x = 1:256;%% 将结果保存到当前目录下的result文件夹下
imwrite(im, sprintf('result/%s.jpg','3_im1'));
imwrite(im1, sprintf('result/%s.jpg','3_im2'));
imwrite(im2, sprintf('result/%s.jpg','3_im3'));%% 显示图像
figure(1);
subplot(221); imshow(im); title('原图'); axis on
subplot(222); plot(x,NK); title('直方图'); axis on
subplot(223); imshow(im1); title('全局阈值'); axis on
subplot(224); imshow(im2); title('otsu阈值'); axis on功能函数：
function [im1,NK] = my_Otsu(im)[h,w] = size(im);
%% 像素值从0到255，下标从1到256
%% 得到直方图和累计直方图
[NK,CH] = myHisteq(im);
%归一化
nK = NK/(h*w);
cH = CH/(h*w);%% 计算P1k P2k
P1 = cH;
P2 = 1-cH;%% 累计均值的计算，仅有255个级别，比像素点取值总数小1
m1 = zeros(255,1);
for k = 1:255for j = 1:km1(k) = m1(k) + j*nK(j);endm1(k) = m1(k)/cH(k);
end
m2 = zeros(255,1);
for k = 1:255for j = k+1:256m2(k) = m2(k) + j*nK(j);endm2(k) = m2(k)/(cH(256)-cH(k));
end
%% 全局灰度均值
mG = 0;
for j = 1:256mG = mG + j*nK(j);
end
mG = mG/cH(256);
%% 计算类间方差
phi2 = zeros(255,1);
for k = 1: 255phi2(k) = P1(k)*(m1(k)-mG)^2 + P2(k)*(m2(k)-mG)^2;
end%% 最优phi2值
m = max(phi2);%% 防止存在相等的最优值
num = 0;
k = zeros(10,0);
for i = 1:255if(phi2(i) == m) num = num+1;k(num) = i;end
end
Bestk = sum(k)/num;
%% 得到最优二值化图像
im1 = zeros(h,w);
for i =1:hfor j = 1:wif(im(i,j)>=Bestk)im1(i,j) = 1;endend
end
im1 = logical(im1); 实验结果
图像平滑改善全局阈值处理
这道题未在本次实验要求范围内，但给了这张图片，于是进行了下面的实验。
实验思想
噪声会将简单的阈值处理问题变得不可解决，当噪声不能在源头减少，并且阈值处理又是所以选择的分割方法时，那么通常能增强性能的一种技术是，在阈值处理之前平滑图像。
部分核心代码
主函数：
%% 运行 main_SmoothOtsu.m clc;
clear;
close all;%% 课本图 10.39
im = imread('Fig1040(a)(large_septagon_gaussian_noise_mean_0_std_50_added).tif'); % 原始图像 uint8
% 得到图像的行数、列数、以及每像素的维数（防止出现RGB图像）
[line,row,v]=size(im);
im = im(:,:,1);
[im1,NK1] = my_Otsu(im);
im2 = myAverage(im);
[im3,NK2] = my_Otsu(im2);
x = 1:256;%% 将结果保存到当前目录下的result文件夹下
imwrite(im, sprintf('result/%s.jpg','4_im1'));
imwrite(im1, sprintf('result/%s.jpg','4_im2'));
imwrite(im2, sprintf('result/%s.jpg','4_im3'));
imwrite(im3, sprintf('result/%s.jpg','4_im4'));%% 显示图像
figure(1);
subplot(231); imshow(im); title('原图'); axis on
subplot(232); plot(x,NK1); title('直方图'); axis on
subplot(233); imshow(im1); title('全局阈值'); axis on
subplot(234); imshow(im2); title('原图'); axis on
subplot(235); plot(x,NK2); title('直方图'); axis on
subplot(236); imshow(im3); title('全局阈值'); axis on
功能函数前面都已经给出过
function [img_2] = myAverage(img_1)size_1 = size(img_1);
h = size_1(1);
w = size_1(2);
img_2 = zeros(h, w);%%边缘延拓四行四列
a = img_1(2:-1:1,:);
b = img_1(h:-1:h-1,:);
img_1 = [a;img_1;b];
c = img_1(:,2:-1:1);
d = img_1(:,w:-1:w-1);
img_1 = double([c,img_1,d]);%5X5均值模板
L = 1/25*[1 1,1,1 1;1,1,1 1 1;1,1,1 1 1;1,1,1 1 1;1,1,1 1 1];for i= 1:hfor j = 1:wim = [img_1(i,j) img_1(i,j+1) img_1(i,j+2) img_1(i,j+3) img_1(i,j+4);...img_1(i+1,j) img_1(i+1,j+1) img_1(i+1,j+2) img_1(i+1,j+3) img_1(i+1,j+4);...img_1(i+2,j) img_1(i+2,j+1) img_1(i+2,j+2) img_1(i+2,j+3) img_1(i+2,j+4);...img_1(i+3,j) img_1(i+3,j+1) img_1(i+3,j+2) img_1(i+3,j+3) img_1(i+3,j+4);...img_1(i+4,j) img_1(i+4,j+1) img_1(i+4,j+2) img_1(i+4,j+3) img_1(i+4,j+4);...];img_2(i,j) = round(sum(sum(L.*im)));end
end
img_2 =uint8(img_2);
end
实验结果
图像分块的可变阈值
由于噪声或者光照不均匀等等一系列问题，导致无法直接进行阈值分割，图像平滑和边缘信息有益于阈值处理，但是在更经常的情况下，上述的两种做法效果不明显，只能使用可变阈值来进行解决。
实验思想
可变阈值处理最简单的方法之一是，把一幅图像分成不重叠的矩形。这种方法用于补偿光照和或反射的不均匀性。选择的矩形要足够小，以便每个矩形的光照都近似是均匀的。
代码
主函数：
%% 运行 main_divideOtsu.m clc;
clear;
close all;%% 课本图 10.39
im = imread('Fig1046(a)(septagon_noisy_shaded).tif'); % 原始图像 uint8
% 得到图像的行数、列数、以及每像素的维数（防止出现RGB图像）
[line,row,v]=size(im);
im = im(:,:,1);[im1,NK] = Global_Threshold(im,0,sum(im(:))/(line*row));
[im2,NK] = my_Otsu(im);
[im3,im4] = my_divideOtsu(im);
x = 1:256;%% 将结果保存到当前目录下的result文件夹下
imwrite(im, sprintf('result/%s.jpg','5_im1'));
imwrite(im1, sprintf('result/%s.jpg','5_im2'));
imwrite(im2, sprintf('result/%s.jpg','5_im3'));
imwrite(im3, sprintf('result/%s.jpg','5_im4'));
imwrite(im4, sprintf('result/%s.jpg','5_im5'));%% 显示图像
figure(1);
subplot(231); imshow(im); title('原图'); axis on
subplot(232); plot(x,NK); title('直方图'); axis on
subplot(233); imshow(im1); title('全局阈值'); axis on
subplot(234); imshow(im2); title('otsu阈值'); axis on
subplot(235); imshow(im3); title('划分区域'); axis on
subplot(236); imshow(im4); title('可变阈值'); axis on
功能函数：
function [imd,imv] = my_divideOtsu(im)
[h,w] = size(im);
h1 = floor(h/2);
w1 = floor(w/3);
im1 = im(1:h1,1:w1);
im2 = im(1:h1,w1+1:2*w1);
im3 = im(1:h1,2*w1+1:w);
im4 = im(h1+1:h,1:w1);
im5 = im(h1+1:h,w1+1:2*w1);
im6 = im(h1+1:h,2*w1+1:w);
imd = im;
imd(h1,:) = 255;
imd(:,w1) = 255;
imd(:,2*w1) = 255;[im1,n] = my_Otsu(im1);
[im2,n] = my_Otsu(im2);
[im3,n] = my_Otsu(im3);
[im4,n] = my_Otsu(im4);
[im5,n] = my_Otsu(im5);
[im6,n] = my_Otsu(im6);
imv = [im1,im2,im3;im4,im5,im6];实验结果
实验小结
在上面的图像分割中，通过不同的阈值处理方法，得到良好的图像分割结果。不难发现， 对于不同的噪声环境和光照影响等等，一种单一的阈值分割方法往往无法出色的完成任务。
在今后的实验中，要合理的根据自己的图像样本来选择合理的算法。例如是否需要对图像进行平滑去噪，是否需要图像分区来使光照近似均匀。不仅要掌握上面的处理方法，也明白了为什么对上面的图像进行这样的处理和分析。
数字图像处理11：阈值分割（基本全局阈值处理，Otsu 的最佳全局阈值，图像平滑改善全局阈值处理，图像分块的可变阈值）相关推荐
- 自适应阈值分割—大津法（OTSU算法）C++实现
转自:https://blog.csdn.net/dcrmg/article/details/52216622 大津法是一种图像灰度自适应的阈值分割算法,是1979年由日本学者大津提出,并由他的名字命 ...
- 【MATLAB数字图像处理】四叉树分割图像
clc clear all dir = 'image/56.jpg'; img = imread(dir); imggray = rgb2gray(img);%灰度处理 imggray = imres ...
- 数字图像处理（第四版）胡学龙：编程实现图3.5中不同采样率图像的显示效果
编程实现图3.5中不同采样率图像的显示效果(纯手敲,代码+运行结果+1对1教你怎么运行) 链接:
- 数字图像处理第十章笔记——图像分割
目录 引言 一.基础知识 二. 点.线和边缘检测 2.1 背景知识 2.2 孤立点检测 2.3 线检测 2.4 边缘检测 2.5 基本边缘检测.更先进的边缘检测 三.阈值处理 3.1 基础知识 3.2 ...
- 数字图像处理实验（六）|图像分割{阈值分割、直方图法、OTUS最大类间方差法(edge、im2dw、imfilter、imresize)、迭代阈值法、点检测}（附matlab实验代码和截图）
文章目录 一.实验目的 二.实验原理 (一) 阈值分割 1. 直方图法 2.OTSU法(最大类间方差法)确定阈值 3. 迭代阈值法 4. 点检测 (二)边缘检测 三.实验内容 (一)阈值分割 1. 直 ...
- 图像处理合集：图像基础操作（图像翻转、图像锐化、图像平滑等）、图像阈值分割（边缘检测、迭代法、OSTU、区域增长法等）、图像特征提取（图像分割、灰度共生矩阵、PCA图像压缩）
文章目录 说明 一.图像锐化或增强相关 1. 图像点处理 1.1 图像翻转 1.2 幂运算和对数运算 2. 直方图处理 3. 图像平滑 4. 图像锐化 5. 图像增强 二.图像阈值分割 1. 边缘检测 ...
- 【医学图像处理】 2 灰度直方图、图像二值化（阈值分割）
文章目录 1 灰度直方图 1.1 直方图理解 1.2 直方图计算 1.3 直方图均衡化 1.3.1 全局均衡化 1.3.2 自适应(局部)均值化 2 图像二值化(阈值分割) 2.1 二值化理解 2.2 ...
- 基于MATLAB改进Otsu阈值分割的车道线检测
基于MATLAB改进Otsu阈值分割的车道线检测 摘要:在判断车道偏离以防止车辆碰撞等危害时,车道标线检测需要通过图像处理来进行,检测方法是否适用于各种背景环境条件以及检测的及时性至关重要传统的Ots ...
- 图像二值化方法及适用场景分析（OTSU Trangle 自适应阈值分割）
图像二值化 应用场景 二值图像定义 阈值获取的方法 手动阈值法 自动阈值法 灰度均值法 基于直方图均值法 OTSU Triangle 自适应均值阈值分割方法 总结 参考文献 应用场景 二值图像处理与分 ...
最新文章
- Java基础学习三：循环结构的使用
- TMG 2010 建立站对站***隧道
- Python中threading的join和setDaemon的区别及用法 例子
- 将Python源码编译成pyc和pyo文件
- 学会阅读硬件的原理图、数据手册大全
- 职称计算机word模拟题,2017年职称计算机考试Word2003模拟题及答案(1)
- 用户行为分析最重要的3个点 渠道转化留存
- LINUX 查看分区UUID的两种方法
- 为猿七年有余，痒否？痛否？
- 32位系统支持多大内存 Windows32位/64位系统最大支持内存详解
- CTF-网络信息安全攻防学习平台(脚本关)
- mysql jail_FreeNAS：如何在Jail里面安装软件？
- 屋里的大象：粒子物理学有自己的死神
- lvs工作在第几层_四层负载均衡——LVS
- 会话技巧---英文单词
- 007. RPX服务端和设计端说明
- 自然语言处理之hmm（隐马尔可夫模型）
- 如何通过AppStore变态审核:看同行经验
- 二叉树 二叉树遍历 通过二叉树遍历求得二叉树
- Python采集ppt素材模板 (多线程版本)，答辩、演讲再也不怕没有好用的PPT模板了（含完整源代码）