MATLAB 粒子群图像分割算法

%对累积分布概率 平均值，减少其运算次数
clc
clear all
tic
a=imread（‘C:\Users\Administrator\Desktop\2018-3-8-FC-分配给学生的部分文章\魏晶茹/lena.jpg‘）;     %Read a image    
[mn]=size（a）; 

imgmax=max（max（a））;
 
a0=double（a）; 
 h=1;                           
 a1=zeros（mn）; 
% 计算平均领域灰度的一维灰度直方图 
for i=1:m 
   for j=1:n 
        for k=-h:h 
            for w=-h:h; 
                 p=i+k; 
               q=j+w; 
                if （p<=0）||（ p>m） 
                    p=i; 
                end 
                if （q<=0）||（q>n） 
                     q=j; 
                end 
                 a1（ij）=a0（pq）+a1（ij）; 
             end 
        end 
        a2（ij）=uint8（1/9*a1（ij））; 
       % a2（ij）=abs（a0（ij）-a2（ij））;
    end 
 end 


gradmax=max（max（a2））;
fxy=zeros（mn）; 
% 计算二维直方图 
for i=1:m 
    for j=1:n 
        c=a0（ij）; 
        d=double（a2（ij））; 
        fxy（c+1d+1）=fxy（c+1d+1）+1; 
    end 
end 
Pxy=fxy/m/n; %二维直方图每个坐标点的概率
P0=zeros（mn）; 
Ui=zeros（mn）;  
Uj=zeros（mn）;  
P0（11）=Pxy（11）; 
for i=2:m
    P0（1i）=P0（1i-1）+Pxy（1i）; 
end 
for i=2:m 
    P0（i1）=P0（i-11）+Pxy（i1）; 
end 
for i=2:m 
    for j=2:n
        P0（ij）=P0（ij-1）+P0（i-1j）-P0（i-1j-1）+Pxy（ij）; %累积分布概率
    end 
end 
P1=ones（mn）-P0; 
Ui（11）=0; 
for i=2:m
    Ui（1i）=Ui（1i-1）+（1-1）*Pxy（1i）; 
%       Ui（1i）=Ui（1i-1）+（i-1）*Pxy（1i）; 
end 
for i=2:m 
    Ui（i1）=Ui（i-11）+（i-1）*Pxy（i1）; 
end 
for i=2:m 
    for j=2:n
        Ui（ij）=Ui（ij-1）+Ui（i-1j）-Ui（i-1j-1）+（i-1）*Pxy（ij）; 
    end 
end 
Uj（11）=0; 
for i=2:m
    Uj（1i）=Uj（1i-1）+（i-1）*Pxy（1i）; 
end 
for i=2:m 
    Uj（i1）=Uj（i-11）+（1-1）*Pxy（i1）; 
%       Uj（i1）=Uj（i-11）+（i-1）*Pxy（i1）; 
end 
for i=2:n 
    for j=2:n
        Uj（ij）=Uj（ij-1）+Uj（i-1j）-Uj（i-1j-1）+（j-1）*Pxy（ij）; 
    end 
end 
uti=0; 
utj=0; 
for i=1:m 
    for j=1:n 
        uti=uti+（i-1）*Pxy（ij）;  
        utj=utj+（j-1）*Pxy（ij）; 
    end 
end 
for i=1:m 
    for j=1:n 
      Ui1（ij）=Ui（mn）-Ui（ij）; 
      Uj1（ij）=Uj（mn）-Uj（ij）; 
    end 
end 

c1=2;
c2=2;
d=1;
N=20;%N为粒子数目  
Wmax = 0.9;             %最大惯性因子
Wmin = 0.1;             %最小惯性因子
cmin =0;
cmax =4;
Itermax=100;                 %迭代次数
gbest1=[0 0];
gbest2=0;
pbest=zeros（N3）;
X=zeros（N2）;
V=zeros（N2）;
gprebest=0;
ts=zeros（imgmax+2imgmax+2）;      %禁忌表
% X=unidrnd（double（imgmax）-20N2）+1;
x1=unidrnd（double（imgmax-10）N1）+1;
x2=unidrnd（double（gradmax）