![2a6a6c95bb47aaaba5543e0609043272.png](https://img-blog.csdnimg.cn/img_convert/2a6a6c95bb47aaaba5543e0609043272.png)
Otsu算法(双阈值)
该算法就是利用otsu算法计算出两个阈值公式
g=w0*(u0-u)^2+w1*(u1-u) ^2+ w2*(u2-u) ^2
g最大值时,就可以选出两个阈值
代码
求两个阈值
function [t1,t2]=DoubleOtsuThresh(img)
%
% Otsu 双阈值求解
% 输入 图像img,输出 最优阈值t1和t2(归一化,范围在[0,1])
%
%
BinsNum = 256;
hist = imhist(img,BinsNum);
p = hist / sum(hist); % 直方图的概率密度函数
mG= sum(p .* (1:BinsNum)'); % 全局均值
P1 = cumsum(p); % 概率分布
m1 = cumsum(p .* (1:BinsNum)')./P1; % 256*1 每个阈值的前景平均灰度
% 根据算法理论,从k2+1累加到L-1,可以先倒着累加再翻转回来
P3= cumsum(flip(p));
m3 = cumsum(flip(p) .* flip(1:BinsNum)')./P3;
P3=flip(P3);
P3=[P3(2:end) ;0]; %P3的索引用k2,则P3(1)实际上是从p(2)加到p(256),所以要移动一个
%移动后也用不到P3(1)这个值,因为k2>k1>1
m3=flip(m3);
m3=[m3(2:end) ;0];
% m2 P2为 k1*k2的索引矩阵,为 256*256 的上三角矩阵 因为k1<k2
m2=zeros(BinsNum,BinsNum);
P2=zeros(BinsNum,BinsNum);
for k1=1:BinsNum-2
for k2=k1+1:BinsNum-1
P2(k1,k2)=1-P1(k1)-P3(k2);
m2(k1,k2) = sum( (k1+1:k2)' .* p(k1+1:k2) )./ P2(k1,k2) ;
end
end
% 遍历k1,k2各种组合 求方差
variance=zeros(BinsNum,BinsNum);
for k1=1:BinsNum-2
for k2=k1+1:BinsNum-1
% variance 为256*256 的上三角矩阵 因为l1<k2
variance(k1,k2) = P1(k1)*(m1(k1)-mG)^2 + ...
P2(k1,k2)*(m2(k1,k2)-mG)^2 + ...
P3(k2)*(m3(k2)-mG)^2;
end
end
% 最大方差点即为最优阈值
[~,index] = max(variance(:));
[index_row,index_col] = ind2sub(size(variance),index);
%灰度值为[0,255] 因此需要-1
t1= (index_row-1)./(BinsNum-1);
t2= (index_col-1)./(BinsNum-1);
end
利用这两个阈值分割图像
function out=img2gray(img,t1,t2)
[M,N]=size(img);
out=zeros(M,N);
for m=1:M
for n=1:N
if img(m,n)<t1
out(m,n)=0;
else
if( img(m,n)>t1 && img(m,n)<t2)
out(m,n)=127; %灰度级可调
else
out(m,n)=255; %灰度级可调
end
end
end
end
end
主函数调用
%% 图B 彩色图
imgB =imread('window.jpg');
imgB=rgb2gray(imgB);
imgB=im2double(imgB);
%% 先高斯滤波平滑图B,去高频,再用阈值处理
gaussH=fspecial('gaussian',[5 5],10);
smoothB=imfilter(imgB,gaussH);
[t1B,t2B] =DoubleOtsuThresh(smoothB);
% 如果imgB是unit8类型,则灰度范围在[0,255] 要把算法得到的两个阈值乘255
% 如果imgB是[0,1]范围的double类型则不用乘
% T1B=t1B*255;
% T2B=t2B*255;
T1B=t1B;
T2B=t2B;
J9=img2gray(smoothB,T1B,T2B); %转为灰度图
figure
imshow(J9,[]);title('平滑+双阈值分割'); % imshow里面J9后面的[]相当于直方图均衡效果,让灰度值均匀分布,可以去掉对比一下
![8e3d13b0fdb1f59f38beb4d3b10188c0.png](https://img-blog.csdnimg.cn/img_convert/8e3d13b0fdb1f59f38beb4d3b10188c0.png)