1 简介
利用Tsallis熵的非广延性,提出了三维Tsallis熵多阈值分割方法。并以最大三维Tsallis熵为准则,利用改进粒子群优化算法进行多阈值优化搜索。该方法不仅考虑了图像像元点的灰度分布信息和像元点之间的灰度相关信息,而且考虑了目标和背景之间的相互关系,在灰度级上有不同的反应,具有很强的抗噪声能力。试验结果表明,该算法不仅能对图像进行正确的分割,而且还具有稳定性高,易于实现等优点。
2 部分代码
%*********************************************
%Image thresholding using Tsallis entropy
%2021-3-26
%*********************************************
clear all;
close all;
clc
format short g;
% disp(strcat('***************apply on whiteRose.jpg'));
% I=imread('whiteRose.jpg');
I=imread('C:\Users\lenovo\Desktop\113172244Tsallis-entropy\Tsallis entropy/cell.jpg');
% I=imread('ndmz.png');
% I=imread('star.jpg');
% I=imread('cell2.jpg');
imgray = rgb2gray(I);%灰度转化
thresh=graythresh(imgray);%求阈值
I2=im2bw(imgray,thresh);%二值图转化
subplot(2,3,1)
imshow(I);title('原图');
subplot(2,3,2)
imshow(imgray);title('灰度图');
subplot(2,3,3)
imhist(imgray);title('灰度直方图');
subplot(2,3,4)
omega=num2str(uint8(round(thresh*255)));
imshow(I2);title(['OTSU 阈值：',omega]);
disp(strcat('OTSU 计算灰度阈值：',num2str(uint8(round(thresh*255)))));
thresh1=thresh
[m,n]=size(imgray);
N=m*n ; %灰度值的总数
%****************************************
%KAPUR_entroy最大熵值法求阈值 figure(1)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[count,x]=imhist(imgray,256);
A=[x,count];
p=count/N;%求每个灰度值出现的概率
%用'1985-histo'文献里提到的新算法计算阈值
a=min(find(p~=0.0));
P=p(p~=0);%提取p不为零的项，剔除奇点
p=p(a:256);%把前几项为零的奇点剔除
b=length(p);
ps=cumsum(p);%概率1-T求和
PS=cumsum(P);
Hs=-cumsum(P.*log(P));
HA=log(PS)+Hs./PS;
HB=log(1-PS)+(Hs(end)-Hs)./(1-PS);
L2=length(HA);%剔除1-PS为零的点
HA=HA(1:(L2)-1);
HB=HB(1:(L2)-1);
Ts=HA+HB;
s=mean(find(Ts==max(Ts))+a-2);
disp(strcat('KAPUR_entroy计算灰度阈值：',num2str(s)));
thresh2=double(s)/255
I3=im2bw(imgray,double(s)/255);
subplot(2,3,5);imshow(I3);title(['KAPUR 阈值：',num2str(s)]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Image thresholding using Tsallis entropy
q=0.5;
for T=1:1:b
sumA=0.0;sumB=0.0;
for i=1:1:T
sumA=sumA+(p(i)/ps(T))^q;
end
for j=T+1:1:b
sumB=sumB+(p(j)/(1-ps(T)))^q;
end
c(T)=sumA;
d(T)=sumB;
end
L=min(find(d==0));%剔除d=0的奇点
c=c(1:L-1);
d=d(1:L-1);
g=[c',d'];
Sqa=(1-c)./(q-1);
Sqb=(1-d)./(q-1);
g=[Sqa',Sqb'];
ST=Sqa+Sqb+(1-q).*Sqa.*Sqb;%Tsallis entropy定义式
t=mean(find(ST==max(ST))+a-2);
disp(strcat('Tsallis entropy计算灰度阈值：',num2str(t)));
thresh3=double(t)/255
I3=im2bw(imgray,thresh3);
subplot(2,3,6);imshow(I3);title(['Tsallis 阈值：',num2str(t)]);
%*****************************************
figure(2) %%%*************不同q值下的情况*******
%*****************************************
% disp(strcat('***************apply on star.jpg'));
% I=imread('whiteRose.jpg');
% I=imread('ROSE2.jpg');
% I=imread('ndmz.png');
% I=imread('star.jpg');
% I=imread('inftest.jpg');
% I=imread('kejian.jpg');
I=imread('fengz.jpg');
% I=imread('cell.jpg');
% I=imread('xingzuo.jpg');
imgray=rgb2gray(I);%灰度转化
thresh=graythresh(imgray);%求阈值
I2=im2bw(imgray,thresh);%二值图转化
subplot(2,3,1)
imshow(I);title('原图');
subplot(2,3,2)
imshow(imgray);title('灰度图');
subplot(2,3,3)
imhist(imgray);title('灰度直方图');
subplot(2,3,4)
omega=num2str(uint8(round(thresh*255)));
imshow(I2);title(['OTSU 阈值：',omega]);
disp(strcat('OTSU 计算灰度阈值：',num2str(uint8(round(thresh*255)))));
thresh1=thresh
[m,n]=size(imgray);
N=m*n ; %灰度值的总数
%****************************************
%KAPUR_entroy最大熵值法求阈值 figure(2)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[count,x]=imhist(imgray,256);
A=[x,count];
p=count/N;%求每个灰度值出现的概率
%用'1985-histo'文献里提到的新算法计算阈值
a=min(find(p~=0.0));
P=p(p~=0);%提取不为零的项
p=p(a:256);
b=length(p);
ps=cumsum(p);%概率1-T求和
PS=cumsum(P);
Hs=-cumsum(P.*log(P));
HA=log(PS)+Hs./PS;
HB=log(1-PS)+(Hs(end)-Hs)./(1-PS);
L2=length(HA);%剔除1-PS为零的奇点
HA=HA(1:(L2)-1);
HB=HB(1:(L2)-1);
Ts=HA+HB;
s=mean(find(Ts==max(Ts))+a-2);
disp(strcat('KAPUR_entroy计算灰度阈值：',num2str(s)));
thresh2=double(s)/255
I3=im2bw(imgray,double(s)/255);
subplot(2,3,5);imshow(I3);title(['KAPUR 阈值：',num2str(s)]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Image thresholding using Tsallis entropy ****figure(2)
q=[0.1,0.5,1,1.5,2];
for n=1:1:5
if q(n)~=1%计算q/=1时候的情况Tsallis entropy
for T=1:1:b
sumA=0.0;sumB=0.0;
for i=1:1:T
sumA=sumA+(p(i)/ps(T))^q(n);
end
for j=T+1:1:b
sumB=sumB+(p(j)/(1-ps(T)))^q(n);
end
c(T)=sumA;
d(T)=sumB;
end
L=min(find(d==0));%处理d=0的奇点
c=c(1:L-1);
d=d(1:L-1);
f=[c',d'];
Sqa=(1-c)./(q(n)-1);
Sqb=(1-d)./(q(n)-1);
h=[Sqa',Sqb'];
ST=Sqa+Sqb+(1-q(n)).*Sqa.*Sqb;
ss=[f,h,ST'];
find(ST==max(ST));
t(n)=mean(find(ST==max(ST))+a-2);
else
%****************************
%处理q=1的情况
for T=1:1:b
sumA=0.0;sumB=0.0;
for i=1:1:T
ft=p(i)/ps(T);
if ft==0
sumA=sumA;%处理ft=0的奇点
else
sumA=sumA+(p(i)/ps(T))*log(p(i)/ps(T));%shannon entropy sum(pi*ln(pA)),pA=pi/ps(T)
end
end
for j=T+1:1:b
segam2=(p(j)/(1-ps(T)));
if segam2==0
sumB=sumB;
else
sumB=sumB+(p(j)/(1-ps(T)))*log((p(j)/(1-ps(T))));
end
end
c(T)=sumA;
d(T)=sumB;
S=-(c+d); %shannon entropy
end
% f=[p,c',d'];
t(n)=mean(find(S==max(S))+a-2);
end
end
t;
disp(strcat('Tsallis entropy计算灰度阈值：',num2str(t(1))));
thresh3=double(t(1))/255
I3=im2bw(imgray,thresh3);
subplot(2,3,6);imshow(I3);title(['Tsallis 阈值：',num2str(t(1)),' q=',num2str(q(1))]);
figure(3)
subplot(2,3,1)
imshow(I);title('原图');
thresha=double(t(1))/255;
Ia=im2bw(imgray,thresha);
subplot(2,3,2);imshow(Ia);title(['Tsallis 阈值：',num2str(t(1)),' q=',num2str(q(1))]);
threshb=double(t(2))/255;
Ib=im2bw(imgray,threshb);
subplot(2,3,3);imshow(Ib);title(['Tsallis 阈值：',num2str(t(2)),' q=',num2str(q(2))]);
threshc=double(t(3))/255;
Ic=im2bw(imgray,threshc);
subplot(2,3,4);imshow(Ic);title(['Tsallis 阈值：',num2str(t(3)),' q=',num2str(q(3))]);
threshd=double(t(4))/255;
Id=im2bw(imgray,threshd);
subplot(2,3,5);imshow(Id);title(['Tsallis 阈值：',num2str(t(4)),' q=',num2str(q(4))]);
threshe=double(t(5))/255;
Ie=im2bw(imgray,threshe);
subplot(2,3,6);imshow(Ie);title(['Tsallis 阈值：',num2str(t(5)),' q=',num2str(q(5))]);
3 仿真结果
4 参考文献
[1]邸秋艳. 基于Tsallis熵的阈值图像分割方法研究[D]. 燕山大学.
博主简介：擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机等多种领域的Matlab仿真，相关matlab代码问题可私信交流。
部分理论引用网络文献，若有侵权联系博主删除。
