warning off;

clear all;
close all;
clc;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Image Preprocessing

% Image acquisition

ImageName=uigetfile('*.jpg;.tif;*.bmp','my pic');
Image=imread(ImageName);
InImage=Image;
figure,imshow(Image);
title('Input Image');
Imo=Image;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[R C D]=size(Image);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Image Restoration using Hybrid Median Filter (Proposed System)

HMF = HybridMedianFilter(Image); %Image restoration for both color and gray image (in case of color
image
% no need to convert from rgb2gray.. (Approximation of color image is not
% required ... We can give color image directly)

figure,imshow(HMF);
title('Hybrid Median Filtered Image');

Diff=imsubtract(Image,HMF);
figure,imshow(Diff*10);
title('Noisy Pixels using Hybrid Median Filter');

[psn1,mse] = Psnr(Image,HMF);
disp('PSNR for Hybrid Median Filter')
disp(psn1);
disp('MSE for Hybrid Median Filter');
disp(mse);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%
% LOG COLOR FILTER

[re ce de]=size(Image);
if de<2
Img1x(:,:,1)=HMF;
Img1x(:,:,2)=HMF;
Img1x(:,:,3)=HMF;
 Image=Img1x;
end

Img1=double(Image)/255;
w = 5;
sigma = [3 0.1];

B = LogColorFilter(Img1,w,sigma(1),sigma(2));
B=ceil(B.*255);
figure,imshow(uint8(B));
title('Log Color Filter (Gabor Method)');
[psn2,mse] = Psnr(Image,B);
disp('PSNR for Log Color Filter')
disp(psn2);
disp('MSE for Log Color Filter');
disp(mse);

Imgdiff=imsubtract(Image,uint8(B));

figure,imshow(Imgdiff*10,[])
title('Noisy Pixels using Log Color Filter');

ImageDenoising=uint8(B);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%

% Negative Abundance-Oriented Hyperspectral Unmixing Algorithm

Image=ImageDenoising;
Clusters=2;
[C1,Rh] = imhist(Image(:,:,1));
figure,
subplot(1,3,1);
bar(C1);
title('Histogram of Red Channel');
[C2,Gh] = imhist(Image(:,:,2));
subplot(1,3,2);
bar(C2);
title('Histogram of Green Channel');
[C3,Bh] = imhist(Image(:,:,3));
subplot(1,3,3);
bar(C3);
title('Histogram of Blue Channel');
S1=Image(:,:,1);
figure,subplot(1,3,1);
imshow(S1);
title('Red Channel');
S2=Image(:,:,2);
subplot(1,3,2);
imshow(S2);
title('Green Channel');
S3=Image(:,:,3);
subplot(1,3,3);
imshow(S3);
title('Blue Channel');
Lmax=256;
Nt=size(Image,1)*size(Image,2);

% CLUSTERING
for i=1:Lmax
NormalizedRed(i)=C1(i)/Nt;
NormalizedGreen(i)=C2(i)/Nt;
NormalizedBlue(i)=C3(i)/Nt;
end
Centroid = 150;
ClustParameter = Clusters-1;
CentroidParameter = 150;
Phase1 = 0.8;
Phase2 = 0.8;
Weight = 1.2;
Minima=-5;
Maxima=5;
Redvector=zeros(Centroid,ClustParameter);
Greenvector=zeros(Centroid,ClustParameter);
Bluevector=zeros(Centroid,ClustParameter);
Xred = Lmax*ones(1,ClustParameter);
Xgreen = ones(1,ClustParameter);
XxGreen = Lmax*ones(1,ClustParameter);
XyGreen = ones(1,ClustParameter);
XxBlue = Lmax*ones(1,ClustParameter);
X_MINB = ones(1,ClustParameter);
Sub1 = zeros(1,ClustParameter);
Value = -10000;
Gaussian = ones(Centroid,1);
xBestR=zeros(Centroid,ClustParameter);
fitBestR=zeros(Centroid,1);
AvRed = zeros(Centroid,1);
Sub2 = zeros(1,ClustParameter);
gbestvalueG = -10000;
gauxG = ones(Centroid,1);
xBestG=zeros(Centroid,ClustParameter);
fitBestG=zeros(Centroid,1);
AvGreen = zeros(Centroid,1);
Sub3 = zeros(1,ClustParameter);
gbestvalueB = -10000;
gauxB = ones(Centroid,1);
xBestB=zeros(Centroid,ClustParameter);
fitBestB=zeros(Centroid,1);
AvBlue = zeros(Centroid,1);
RedCluster = zeros(Centroid,ClustParameter);

%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i = 1: Centroid
for j = 1: ClustParameter
RedCluster(i,j) = fix(rand(1,1) * ( Xred(j)-Xgreen(j) ) + Xgreen(j));
end
end
GreenCluster = zeros(Centroid,ClustParameter);
for i = 1: Centroid
for j = 1: ClustParameter
GreenCluster(i,j) = fix(rand(1,1) * ( XxGreen(j)-XyGreen(j) ) + XyGreen(j));
end
end
BlueCluster = zeros(Centroid,ClustParameter);
for i = 1: Centroid
for j = 1: ClustParameter
BlueCluster(i,j) = fix(rand(1,1) * ( XxBlue(j)-X_MINB(j) ) + X_MINB(j));
end
end
for si=1:length(RedCluster)
RedCluster(si,:)=sort(RedCluster(si,:));
end
for si=1:length(GreenCluster)
GreenCluster(si,:)=sort(GreenCluster(si,:));
end
for si=1:length(BlueCluster)
BlueCluster(si,:)=sort(BlueCluster(si,:));
end
nger=1;

for j=1:Centroid

AvRed(j)=sum(NormalizedRed(1:RedCluster(j,1)))*(sum((1:RedCluster(j,1)).*NormalizedRed(1:RedClu
ster(j,1))/sum(NormalizedRed(1:RedCluster(j,1)))) - sum((1:Lmax).*NormalizedRed(1:Lmax)) )^2;
for jlevel=2:Clusters-1
AvRed(j)=AvRed(j)+sum(NormalizedRed(RedCluster(j,jlevel-
1)+1:RedCluster(j,jlevel)))*(sum((RedCluster(j,jlevel-
1)+1:RedCluster(j,jlevel)).*NormalizedRed(RedCluster(j,jlevel-1)+1:RedCluster(j,jlevel))/
sum(NormalizedRed(RedCluster(j,jlevel-1)+1:RedCluster(j,jlevel))))-
sum((1:Lmax).*NormalizedRed(1:Lmax)))^2;
end
AvRed(j)=AvRed(j)+sum(NormalizedRed(RedCluster(j,Clusters-
1)+1:Lmax))*(sum((RedCluster(j,Clusters-1)+1:Lmax).*NormalizedRed(RedCluster(j,Clusters-
1)+1:Lmax)/sum(NormalizedRed(RedCluster(j,Clusters-1)+1:Lmax)))-
sum((1:Lmax).*NormalizedRed(1:Lmax)))^2;
fitBestR(j)=AvRed(j);

AvGreen(j)=sum(NormalizedGreen(1:GreenCluster(j,1)))*(sum((1:GreenCluster(j,1)).*NormalizedGreen
(1:GreenCluster(j,1))/sum(NormalizedGreen(1:GreenCluster(j,1)))) -
sum((1:Lmax).*NormalizedGreen(1:Lmax)) )^2;
 for jlevel=2:Clusters-1
AvGreen(j)=AvGreen(j)+sum(NormalizedGreen(GreenCluster(j,jlevel-
1)+1:GreenCluster(j,jlevel)))*(sum((GreenCluster(j,jlevel-
1)+1:GreenCluster(j,jlevel)).*NormalizedGreen(GreenCluster(j,jlevel-1)+1:GreenCluster(j,jlevel))/
sum(NormalizedGreen(GreenCluster(j,jlevel-1)+1:GreenCluster(j,jlevel))))-
sum((1:Lmax).*NormalizedGreen(1:Lmax)))^2;
end
AvGreen(j)=AvGreen(j)+sum(NormalizedGreen(GreenCluster(j,Clusters-
1)+1:Lmax))*(sum((GreenCluster(j,Clusters-1)+1:Lmax).*NormalizedGreen(GreenCluster(j,Clusters-
1)+1:Lmax)/sum(NormalizedGreen(GreenCluster(j,Clusters-1)+1:Lmax)))-
sum((1:Lmax).*NormalizedGreen(1:Lmax)))^2;
fitBestG(j)=AvGreen(j);

AvBlue(j)=sum(NormalizedBlue(1:BlueCluster(j,1)))*(sum((1:BlueCluster(j,1)).*NormalizedBlue(1:Blu
eCluster(j,1))/sum(NormalizedBlue(1:BlueCluster(j,1)))) -
sum((1:Lmax).*NormalizedBlue(1:Lmax)) )^2;
for jlevel=2:Clusters-1
AvBlue(j)=AvBlue(j)+sum(NormalizedBlue(BlueCluster(j,jlevel-
1)+1:BlueCluster(j,jlevel)))*(sum((BlueCluster(j,jlevel-
1)+1:BlueCluster(j,jlevel)).*NormalizedBlue(BlueCluster(j,jlevel-1)+1:BlueCluster(j,jlevel))/
sum(NormalizedBlue(BlueCluster(j,jlevel-1)+1:BlueCluster(j,jlevel))))-
sum((1:Lmax).*NormalizedBlue(1:Lmax)))^2;
end
AvBlue(j)=AvBlue(j)+sum(NormalizedBlue(BlueCluster(j,Clusters-
1)+1:Lmax))*(sum((BlueCluster(j,Clusters-1)+1:Lmax).*NormalizedBlue(BlueCluster(j,Clusters-
1)+1:Lmax)/sum(NormalizedBlue(BlueCluster(j,Clusters-1)+1:Lmax)))-
sum((1:Lmax).*NormalizedBlue(1:Lmax)))^2;
fitBestB(j)=AvBlue(j);
[aR,bR]=max(AvRed);
Sub1=RedCluster(bR,:);
Value = AvRed(bR);
[aG,bG]=max(AvGreen);
Sub2=GreenCluster(bG,:);
gbestvalueG = AvGreen(bG);
[aB,bB]=max(AvBlue);
Sub3=RedCluster(bB,:);
gbestvalueB = AvBlue(bB);
xBestR = RedCluster;
xBestG = GreenCluster;
xBestB = BlueCluster;

while(nger<=CentroidParameter)
i=1;
randnum1 = rand ([Centroid, ClustParameter]);
randnum2 = rand ([Centroid, ClustParameter]);
Redvector = fix(Weight.*Redvector + randnum1.*(Phase1.*(xBestR-RedCluster)) +
randnum2.*(Phase2.*(Gaussian*Sub1-RedCluster)));
Redvector = ( (Redvector <= Minima).*Minima ) + ( (Redvector > Minima).*Redvector );
Redvector = ( (Redvector >= Maxima).*Maxima ) + ( (Redvector < Maxima).*Redvector );
RedCluster = RedCluster+Redvector;
 Greenvector = fix(Weight.*Greenvector + randnum1.*(Phase1.*(xBestG-GreenCluster)) +
randnum2.*(Phase2.*(gauxG*Sub2-GreenCluster)));
Greenvector = ( (Greenvector <= Minima).*Minima ) + ( (Greenvector > Minima).*Greenvector );
Greenvector = ( (Greenvector >= Maxima).*Maxima ) + ( (Greenvector < Maxima).*Greenvector );
GreenCluster = GreenCluster+Greenvector;
Bluevector = fix(Weight.*Bluevector + randnum1.*(Phase1.*(xBestB-BlueCluster)) +
randnum2.*(Phase2.*(gauxB*Sub3-BlueCluster)));
Bluevector = ( (Bluevector <= Minima).*Minima ) + ( (Bluevector > Minima).*Bluevector );
Bluevector = ( (Bluevector >= Maxima).*Maxima ) + ( (Bluevector < Maxima).*Bluevector );
BlueCluster = BlueCluster+Bluevector;
RedCluster = ( (RedCluster <= Xgreen(1)).*Xgreen(1) ) + ( (RedCluster >
Xgreen(1)).*RedCluster );
RedCluster = ( (RedCluster >= Xred(1)).*Xred(1) ) + ( (RedCluster < Xred(1)).*RedCluster );
GreenCluster = ( (GreenCluster <= XyGreen(1)).*XyGreen(1) ) + ( (GreenCluster >
XyGreen(1)).*GreenCluster );
GreenCluster = ( (GreenCluster >= XxGreen(1)).*XxGreen(1) ) + ( (GreenCluster <
XxGreen(1)).*GreenCluster );
BlueCluster = ( (BlueCluster <= X_MINB(1)).*X_MINB(1) ) + ( (BlueCluster >
X_MINB(1)).*BlueCluster );
BlueCluster = ( (BlueCluster >= XxBlue(1)).*XxBlue(1) ) + ( (BlueCluster <
XxBlue(1)).*BlueCluster );
for j = 1:Centroid
for k = 1:ClustParameter
if ((k>1)&&(k<ClustParameter))
if RedCluster(j,k) < RedCluster(j,k-1)
RedCluster(j,k) = RedCluster(j,k-1);
elseif RedCluster(j,k) > RedCluster(j,k+1)
RedCluster(j,k) = RedCluster(j,k+1);
end
end
if (k==ClustParameter)&&(k~=1)
if RedCluster(j,k) < RedCluster(j,k-1)
RedCluster(j,k) = RedCluster(j,k-1);
elseif RedCluster(j,k) > Xred(k)
RedCluster(j,k) = Xred(k);
end
end
if (k==1)&&(k==ClustParameter)
if RedCluster(j,k) < Xgreen(k)
RedCluster(j,k) = Xgreen(k);
elseif RedCluster(j,k) > Xred(k)
RedCluster(j,k) = Xred(k);
end
end
if (k==1)&&(k~=ClustParameter)
if RedCluster(j,k) < Xgreen(k)
RedCluster(j,k) = Xgreen(k);
elseif RedCluster(j,k) > RedCluster(j,k+1)
RedCluster(j,k) = RedCluster(j,k+1);
end
end
if ((k>1)&&(k<ClustParameter))
if RedCluster(j,k) < RedCluster(j,k-1)
RedCluster(j,k) = RedCluster(j,k-1);
elseif RedCluster(j,k) > RedCluster(j,k+1)
RedCluster(j,k) = RedCluster(j,k+1);
end
end
if (k==ClustParameter)&&(k~=1)
if RedCluster(j,k) < RedCluster(j,k-1)
RedCluster(j,k) = RedCluster(j,k-1);
elseif RedCluster(j,k) > Xred(k)
RedCluster(j,k) = Xred(k);
end
end
if (k==1)&&(k==ClustParameter)
if RedCluster(j,k) < Xgreen(k)
RedCluster(j,k) = Xgreen(k);
elseif RedCluster(j,k) > Xred(k)
RedCluster(j,k) = Xred(k);
end
end
if (k==1)&&(k~=ClustParameter)
if GreenCluster(j,k) < XyGreen(k)
GreenCluster(j,k) = XyGreen(k);
elseif GreenCluster(j,k) > GreenCluster(j,k+1)
GreenCluster(j,k) = GreenCluster(j,k+1);
end
end
if ((k>1)&&(k<ClustParameter))
if GreenCluster(j,k) < GreenCluster(j,k-1)
GreenCluster(j,k) = GreenCluster(j,k-1);
elseif GreenCluster(j,k) > GreenCluster(j,k+1)
GreenCluster(j,k) = GreenCluster(j,k+1);
end
end
if (k==ClustParameter)&&(k~=1)
if GreenCluster(j,k) < GreenCluster(j,k-1)
GreenCluster(j,k) = GreenCluster(j,k-1);
elseif GreenCluster(j,k) > XxGreen(k)
GreenCluster(j,k) = XxGreen(k);
end
end
if (k==1)&&(k==ClustParameter)
if GreenCluster(j,k) < XyGreen(k)
GreenCluster(j,k) = XyGreen(k);
elseif GreenCluster(j,k) > XxGreen(k)
GreenCluster(j,k) = XxGreen(k);
end
end
if (k==1)&&(k~=ClustParameter)
if BlueCluster(j,k) < X_MINB(k)
 BlueCluster(j,k) = X_MINB(k);
elseif BlueCluster(j,k) > BlueCluster(j,k+1)
BlueCluster(j,k) = BlueCluster(j,k+1);
end
end
if ((k>1)&&(k<ClustParameter))
if BlueCluster(j,k) < BlueCluster(j,k-1)
BlueCluster(j,k) = BlueCluster(j,k-1);
elseif BlueCluster(j,k) > BlueCluster(j,k+1)
BlueCluster(j,k) = BlueCluster(j,k+1);
end
end
if (k==ClustParameter)&&(k~=1)
if BlueCluster(j,k) < BlueCluster(j,k-1)
BlueCluster(j,k) = BlueCluster(j,k-1);
elseif BlueCluster(j,k) > XxBlue(k)
BlueCluster(j,k) = XxBlue(k);
end
end
if (k==1)&&(k==ClustParameter)
if BlueCluster(j,k) < X_MINB(k)
BlueCluster(j,k) = X_MINB(k);
elseif BlueCluster(j,k) > XxBlue(k)
BlueCluster(j,k) = XxBlue(k);
end
end
end
end
while(i<=Centroid)
if(i==Centroid)
for j=1:Centroid

AvRed(j)=sum(NormalizedRed(1:RedCluster(j,1)))*(sum((1:RedCluster(j,1)).*NormalizedRed(1:RedClu
ster(j,1))/sum(NormalizedRed(1:RedCluster(j,1)))) - sum((1:Lmax).*NormalizedRed(1:Lmax)) )^2;
for jlevel=2:Clusters-1
AvRed(j)=AvRed(j)+sum(NormalizedRed(RedCluster(j,jlevel-
1)+1:RedCluster(j,jlevel)))*(sum((RedCluster(j,jlevel-
1)+1:RedCluster(j,jlevel)).*NormalizedRed(RedCluster(j,jlevel-1)+1:RedCluster(j,jlevel))/
sum(NormalizedRed(RedCluster(j,jlevel-1)+1:RedCluster(j,jlevel))))-
sum((1:Lmax).*NormalizedRed(1:Lmax)))^2;
end
AvRed(j)=AvRed(j)+sum(NormalizedRed(RedCluster(j,Clusters-
1)+1:Lmax))*(sum((RedCluster(j,Clusters-1)+1:Lmax).*NormalizedRed(RedCluster(j,Clusters-
1)+1:Lmax)/sum(NormalizedRed(RedCluster(j,Clusters-1)+1:Lmax)))-
sum((1:Lmax).*NormalizedRed(1:Lmax)))^2;

AvGreen(j)=sum(NormalizedGreen(1:GreenCluster(j,1)))*(sum((1:GreenCluster(j,1)).*NormalizedGreen
(1:GreenCluster(j,1))/sum(NormalizedGreen(1:GreenCluster(j,1)))) -
sum((1:Lmax).*NormalizedGreen(1:Lmax)) )^2;
for jlevel=2:Clusters-1
 AvGreen(j)=AvGreen(j)+sum(NormalizedGreen(GreenCluster(j,jlevel-
1)+1:GreenCluster(j,jlevel)))*(sum((GreenCluster(j,jlevel-
1)+1:GreenCluster(j,jlevel)).*NormalizedGreen(GreenCluster(j,jlevel-1)+1:GreenCluster(j,jlevel))/
sum(NormalizedGreen(GreenCluster(j,jlevel-1)+1:GreenCluster(j,jlevel))))-
sum((1:Lmax).*NormalizedGreen(1:Lmax)))^2;
end
AvGreen(j)=AvGreen(j)+sum(NormalizedGreen(GreenCluster(j,Clusters-
1)+1:Lmax))*(sum((GreenCluster(j,Clusters-1)+1:Lmax).*NormalizedGreen(GreenCluster(j,Clusters-
1)+1:Lmax)/sum(NormalizedGreen(GreenCluster(j,Clusters-1)+1:Lmax)))-
sum((1:Lmax).*NormalizedGreen(1:Lmax)))^2;

AvBlue(j)=sum(NormalizedBlue(1:BlueCluster(j,1)))*(sum((1:BlueCluster(j,1)).*NormalizedBlue(1:Blu
eCluster(j,1))/sum(NormalizedBlue(1:BlueCluster(j,1)))) -
sum((1:Lmax).*NormalizedBlue(1:Lmax)) )^2;
for jlevel=2:Clusters-1
AvBlue(j)=AvBlue(j)+sum(NormalizedBlue(BlueCluster(j,jlevel-
1)+1:BlueCluster(j,jlevel)))*(sum((BlueCluster(j,jlevel-
1)+1:BlueCluster(j,jlevel)).*NormalizedBlue(BlueCluster(j,jlevel-1)+1:BlueCluster(j,jlevel))/
sum(NormalizedBlue(BlueCluster(j,jlevel-1)+1:BlueCluster(j,jlevel))))-
sum((1:Lmax).*NormalizedBlue(1:Lmax)))^2;
end
AvBlue(j)=AvBlue(j)+sum(NormalizedBlue(BlueCluster(j,Clusters-
1)+1:Lmax))*(sum((BlueCluster(j,Clusters-1)+1:Lmax).*NormalizedBlue(BlueCluster(j,Clusters-
1)+1:Lmax)/sum(NormalizedBlue(BlueCluster(j,Clusters-1)+1:Lmax)))-
sum((1:Lmax).*NormalizedBlue(1:Lmax)))^2;
if AvRed(j) > fitBestR(j)
fitBestR(j) = AvRed(j);
xBestR(j,:) = RedCluster(j,:);
end
if AvGreen(j) > fitBestG(j)
fitBestG(j) = AvGreen(j);
xBestG(j,:) = GreenCluster(j,:);
end
if AvBlue(j) > fitBestB(j)
fitBestB(j) = AvBlue(j);
xBestB(j,:) = BlueCluster(j,:);
end
end
[aR,bR] = max (AvRed);
[aG,bG] = max (AvGreen);
[aB,bB] = max (AvBlue);
if (AvRed(bR) > Value)
Sub1=RedCluster(bR,:)-1;
Value = AvRed(bR);
end
if (AvGreen(bG) > gbestvalueG)
Sub2=GreenCluster(bG,:)-1;
gbestvalueG = AvGreen(bG);
end
if (AvBlue(bB) > gbestvalueB)
Sub3=BlueCluster(bB,:)-1;
 gbestvalueB = AvBlue(bB);
end
nger=nger+1;
end
i=i+1;
end
end
end
Endmember=ceil((Sub1+Sub2+Sub3)/3);
Sub1=sort(Sub1);
Sub2=sort(Sub2);
Sub3=sort(Sub3);
disp('End member');
disp(Endmember);
EI=ImageReconstruction(Image,Sub1,Sub2,Sub3);
figure,imshow(EI);
title('Oriented Satellite Image based on Deep Learning');
figure,imagesc(rgb2gray(EI));
colormap('default');
title('Clustered Image using Fuzzy C Means Algorithm');
[psn,mse] = Parameters(EI,Image);
disp('Parameters of Fuzzy C Means Algorithm');
disp('PSNR');
disp(psn);
disp('Mean Square Error');
disp(mse);
disp(' ');
I=ImageDenoising;
[n_countR,x_valueR] = imhist(I(:,:,1));
[n_countG,x_valueG] = imhist(I(:,:,2));
[n_countB,x_valueB] = imhist(I(:,:,3));
Nt=size(I,1)*size(I,2);
Lmax=256;
level=4;
for i=1:Lmax
probR(i)=n_countR(i)/Nt;
probG(i)=n_countG(i)/Nt;
probB(i)=n_countB(i)/Nt;
end
N = 100;
MIN_POP = 10;
MAX_POP = 50;
POP_INICIAL = N;
N_SWARMS = 4;
N_SWARMSR = N_SWARMS;
N_SWARMSG = N_SWARMS;
N_SWARMSB = N_SWARMS;
MIN_SWARMS = 2;
MAX_SWARMS = 6;
STAGNANCY = 10;
N_PAR = level-1;
N_GER = 150;
PHI1 = 0.8;
PHI2 = 0.8;
W = 1.2;
vmin=-1.5;
vmax=1.5;
vR=zeros(N*N_SWARMS,N_PAR);
vG=zeros(N*N_SWARMS,N_PAR);
vB=zeros(N*N_SWARMS,N_PAR);
X_MAXR = Lmax*ones(1,N_PAR);
X_MINR = ones(1,N_PAR);
X_MAXG = Lmax*ones(1,N_PAR);
X_MING = ones(1,N_PAR);
X_MAXB = Lmax*ones(1,N_PAR);
X_MINB = ones(1,N_PAR);
Red = zeros(N_SWARMS,N_PAR);
gbestvalueR = -1000*ones(N_SWARMS,1);
gauxR = ones(N*N_SWARMS,1);
xBestR=zeros(N*N_SWARMS,N_PAR);
fitBestR=zeros(N*N_SWARMS,1);
fitR = zeros(N*N_SWARMS,1);
Green = zeros(N_SWARMS,N_PAR);
gbestvalueG = -1000*ones(N_SWARMS,1);
gauxG = ones(N*N_SWARMS,1);
xBestG=zeros(N*N_SWARMS,N_PAR);
fitBestG=zeros(N*N_SWARMS,1);
fitG = zeros(N*N_SWARMS,1);
Blue = zeros(N_SWARMS,N_PAR);
gbestvalueB = -1000*ones(N_SWARMS,1);
gauxB = ones(N*N_SWARMS,1);
xBestB=zeros(N*N_SWARMS,N_PAR);
fitBestB=zeros(N*N_SWARMS,1);
fitB = zeros(N*N_SWARMS,1);
xR = zeros(N*N_SWARMS,N_PAR);
for i = 1: N*N_SWARMS
for j = 1: N_PAR
xR(i,j) = fix(rand(1,1) * ( X_MAXR(j)-X_MINR(j) ) + X_MINR(j));
end
end
xG = zeros(N*N_SWARMS,N_PAR);
for i = 1: N*N_SWARMS
for j = 1: N_PAR
xG(i,j) = fix(rand(1,1) * ( X_MAXG(j)-X_MING(j) ) + X_MING(j));
end
end
xB = zeros(N*N_SWARMS,N_PAR);
for i = 1: N*N_SWARMS
for j = 1: N_PAR
xB(i,j) = fix(rand(1,1) * ( X_MAXB(j)-X_MINB(j) ) + X_MINB(j));
end
end
 for si=1:length(xR)
xR(si,:)=sort(xR(si,:));
end
for si=1:length(xG)
xG(si,:)=sort(xG(si,:));
end
for si=1:length(xB)
xB(si,:)=sort(xB(si,:));
end

nger=1;
for j=1:N*N_SWARMS
fitR(j)=sum(probR(1:xR(j,1)))*(sum((1:xR(j,1)).*probR(1:xR(j,1))/sum(probR(1:xR(j,1)))) -
sum((1:Lmax).*probR(1:Lmax)) )^2;
for jlevel=2:level-1
fitR(j)=fitR(j)+sum(probR(xR(j,jlevel-1)+1:xR(j,jlevel)))*(sum((xR(j,jlevel-
1)+1:xR(j,jlevel)).*probR(xR(j,jlevel-1)+1:xR(j,jlevel))/sum(probR(xR(j,jlevel-1)+1:xR(j,jlevel))))-
sum((1:Lmax).*probR(1:Lmax)))^2;
end
fitR(j)=fitR(j)+sum(probR(xR(j,level-1)+1:Lmax))*(sum((xR(j,level-
1)+1:Lmax).*probR(xR(j,level-1)+1:Lmax)/sum(probR(xR(j,level-1)+1:Lmax)))-
sum((1:Lmax).*probR(1:Lmax)))^2;
fitBestR(j)=fitR(j);
fitG(j)=sum(probG(1:xG(j,1)))*(sum((1:xG(j,1)).*probG(1:xG(j,1))/sum(probG(1:xG(j,1)))) -
sum((1:Lmax).*probG(1:Lmax)) )^2;
for jlevel=2:level-1
fitG(j)=fitG(j)+sum(probG(xG(j,jlevel-1)+1:xG(j,jlevel)))*(sum((xG(j,jlevel-
1)+1:xG(j,jlevel)).*probG(xG(j,jlevel-1)+1:xG(j,jlevel))/sum(probG(xG(j,jlevel-1)+1:xG(j,jlevel))))-
sum((1:Lmax).*probG(1:Lmax)))^2;
end
fitG(j)=fitG(j)+sum(probG(xG(j,level-1)+1:Lmax))*(sum((xG(j,level-
1)+1:Lmax).*probG(xG(j,level-1)+1:Lmax)/sum(probG(xG(j,level-1)+1:Lmax)))-
sum((1:Lmax).*probG(1:Lmax)))^2;
fitBestG(j)=fitG(j);
fitB(j)=sum(probB(1:xB(j,1)))*(sum((1:xB(j,1)).*probB(1:xB(j,1))/sum(probB(1:xB(j,1)))) -
sum((1:Lmax).*probB(1:Lmax)) )^2;
for jlevel=2:level-1
fitB(j)=fitB(j)+sum(probB(xB(j,jlevel-1)+1:xB(j,jlevel)))*(sum((xB(j,jlevel-
1)+1:xB(j,jlevel)).*probB(xB(j,jlevel-1)+1:xB(j,jlevel))/sum(probB(xB(j,jlevel-1)+1:xB(j,jlevel))))-
sum((1:Lmax).*probB(1:Lmax)))^2;
end
fitB(j)=fitB(j)+sum(probB(xB(j,level-1)+1:Lmax))*(sum((xB(j,level-
1)+1:Lmax).*probB(xB(j,level-1)+1:Lmax)/sum(probB(xB(j,level-1)+1:Lmax)))-
sum((1:Lmax).*probB(1:Lmax)))^2;
fitBestB(j)=fitB(j);
end
for i=1:N_SWARMS

k=i*N;
j=k-N+1;
[aR,bR]=max(fitR(j:k,:));
 Red(i,:)=xR(j-1+bR,:);
gbestvalueR(i,1)= fitR(j-1+bR);

xBestR = xR;

if i==N_SWARMS
[LS3,LSC3]=max(fitR);
end

[aG,bG]=max(fitG(j:k,:));

Green(i,:)=xG(j-1+bG,:);
gbestvalueG(i,1)= fitG(j-1+bG);

xBestG = xG;

if i==N_SWARMS
[LS2,LSC2]=max(fitG);
end

[aB,bB]=max(fitB(j:k,:));

Blue(i,:)=xB(j-1+bB,:);
gbestvalueB(i,1)= fitB(j-1+bB);

xBestB = xB;

if i==N_SWARMS
[LS1,LSC1]=max(fitB);
end
end
xBestR = xR;
xBestG = xG;
xBestB = xB;
% end

NR(1:N_SWARMSR,1)=N;
NG(1:N_SWARMSG,1)=N;
NB(1:N_SWARMSB,1)=N;

stagnancy_counterR=zeros(N_SWARMSR,1);
stagnancy_counterG=zeros(N_SWARMSG,1);
stagnancy_counterB=zeros(N_SWARMSB,1);

already_deletedR=zeros(N_SWARMSR,1);
already_deletedG=zeros(N_SWARMSG,1);
already_deletedB=zeros(N_SWARMSB,1);

nger=1;
while(nger<=N_GER)
 i=1;
while ((i<=max([N_SWARMSR,N_SWARMSG, N_SWARMSB])))

if (i<=N_SWARMSR)
kR=(sum(NR(1:i)));
jR=kR-NR(i)+1;
end

if (i<=N_SWARMSG)
kG=(sum(NG(1:i)));
jG=kG-NG(i)+1;
end

if (i<=N_SWARMSB)
kB=(sum(NB(1:i)));
jB=kB-NB(i)+1;
end
if (i<=N_SWARMSR)
v_swarmR = vR(jR:kR,:);
gBest_swarmR = Red(i,:);
gbestvalue_swarmR = gbestvalueR(i,1);
x_swarmR = xR(jR:kR,:);
xBest_swarmR = xBestR(jR:kR,:);
fit_swarmR = fitR(jR:kR,:);
fitBest_swarmR = fitBestR(jR:kR,:);
gaux_swarmR = ones(NR(i),1);
end

if (i<=N_SWARMSG)
v_swarmG = vG(jG:kG,:);
gBest_swarmG = Green(i,:);
gbestvalue_swarmG = gbestvalueG(i,1);
x_swarmG = xG(jG:kG,:);
xBest_swarmG = xBestG(jG:kG,:);
fit_swarmG = fitG(jG:kG,:);
fitBest_swarmG = fitBestG(jG:kG,:);
gaux_swarmG = ones(NG(i),1);
end

if (i<=N_SWARMSB)
v_swarmB = vB(jB:kB,:);
gBest_swarmB = Blue(i,:);
gbestvalue_swarmB = gbestvalueB(i,1);
x_swarmB = xB(jB:kB,:);
xBest_swarmB = xBestB(jB:kB,:);
fit_swarmB = fitB(jB:kB,:);
fitBest_swarmB = fitBestB(jB:kB,:);
gaux_swarmB = ones(NB(i),1);
end
% end
 if (i<=N_SWARMSR)
randnum1R = rand ([NR(i), N_PAR]);
randnum2R = rand ([NR(i), N_PAR]);
end

if (i<=N_SWARMSG)
randnum1G = rand ([NG(i), N_PAR]);
randnum2G = rand ([NG(i), N_PAR]);
end

if (i<=N_SWARMSB)
randnum1B = rand ([NB(i), N_PAR]);
randnum2B = rand ([NB(i), N_PAR]);
end
if (i<=N_SWARMSR)
v_swarmR = fix(W.*v_swarmR + randnum1R.*(PHI1.*(xBest_swarmR-x_swarmR)) +
randnum2R.*(PHI2.*(gaux_swarmR*gBest_swarmR-x_swarmR)));
v_swarmR = ( (v_swarmR <= vmin).*vmin ) + ( (v_swarmR > vmin).*v_swarmR );
v_swarmR = ( (v_swarmR >= vmax).*vmax ) + ( (v_swarmR < vmax).*v_swarmR );
x_swarmR = round(x_swarmR+v_swarmR);
end

if (i<=N_SWARMSG)
v_swarmG = fix(W.*v_swarmG + randnum1G.*(PHI1.*(xBest_swarmG-x_swarmG)) +
randnum2G.*(PHI2.*(gaux_swarmG*gBest_swarmG-x_swarmG)));
v_swarmG = ( (v_swarmG <= vmin).*vmin ) + ( (v_swarmG > vmin).*v_swarmG );
v_swarmG = ( (v_swarmG >= vmax).*vmax ) + ( (v_swarmG < vmax).*v_swarmG );
x_swarmG = round(x_swarmG+v_swarmG);
end

if (i<=N_SWARMSB)
v_swarmB = fix(W.*v_swarmB + randnum1B.*(PHI1.*(xBest_swarmB-x_swarmB)) +
randnum2B.*(PHI2.*(gaux_swarmB*gBest_swarmB-x_swarmB)));
v_swarmB = ( (v_swarmB <= vmin).*vmin ) + ( (v_swarmB > vmin).*v_swarmB );
v_swarmB = ( (v_swarmB >= vmax).*vmax ) + ( (v_swarmB < vmax).*v_swarmB );
x_swarmB = round(x_swarmB+v_swarmB);
end
if (i<=N_SWARMSR)
x_swarmR = ( (x_swarmR <= X_MINR(1)).*X_MINR(1) ) + ( (x_swarmR >
X_MINR(1)).*x_swarmR );
x_swarmR = ( (x_swarmR >= X_MAXR(1)).*X_MAXR(1) ) + ( (x_swarmR <
X_MAXR(1)).*x_swarmR );
end

if (i<=N_SWARMSG)
x_swarmG = ( (x_swarmG <= X_MING(1)).*X_MING(1) ) + ( (x_swarmG >
X_MING(1)).*x_swarmG );
x_swarmG = ( (x_swarmG >= X_MAXG(1)).*X_MAXG(1) ) + ( (x_swarmG <
X_MAXG(1)).*x_swarmG );
end
 if (i<=N_SWARMSB)
x_swarmB = ( (x_swarmB <= X_MINB(1)).*X_MINB(1) ) + ( (x_swarmB >
X_MINB(1)).*x_swarmB );
x_swarmB = ( (x_swarmB >= X_MAXB(1)).*X_MAXB(1) ) + ( (x_swarmB <
X_MAXB(1)).*x_swarmB );
end
if (i<=N_SWARMSR)
for jj = 1:NR(i)
for kk = 1:N_PAR
if (kk==1)&&(kk~=N_PAR)
if x_swarmR(jj,kk) < X_MINR(kk)
x_swarmR(jj,kk) = X_MINR(kk);
elseif x_swarmR(jj,kk) > x_swarmR(jj,kk+1)
x_swarmR(jj,kk) = x_swarmR(jj,kk+1);
end
end
if ((kk>1)&&(kk<N_PAR))
if x_swarmR(jj,kk) < x_swarmR(jj,kk-1)
x_swarmR(jj,kk) = x_swarmR(jj,kk-1);
elseif x_swarmR(jj,kk) > x_swarmR(jj,kk+1)
x_swarmR(jj,kk) = x_swarmR(jj,kk+1);
end
end
if (kk==N_PAR)&&(kk~=1)
if x_swarmR(jj,kk) < x_swarmR(jj,kk-1)
x_swarmR(jj,kk) = x_swarmR(jj,kk-1);
elseif x_swarmR(jj,kk) > X_MAXR(kk)
x_swarmR(jj,kk) = X_MAXR(kk);
end
end
if (kk==1)&&(kk==N_PAR)
if x_swarmR(jj,kk) < X_MINR(kk)
x_swarmR(jj,kk) = X_MINR(kk);
elseif x_swarmR(jj,kk) > X_MAXR(kk)
x_swarmR(jj,kk) = X_MAXR(kk);
end
end
end
end
end

if (i<=N_SWARMSG)
for jj = 1:NG(i)
for kk = 1:N_PAR
if (kk==1)&&(kk~=N_PAR)
if x_swarmG(jj,kk) < X_MING(kk)
x_swarmG(jj,kk) = X_MING(kk);
elseif x_swarmG(jj,kk) > x_swarmG(jj,kk+1)
x_swarmG(jj,kk) = x_swarmG(jj,kk+1);
end
 end
if ((kk>1)&&(kk<N_PAR))
if x_swarmG(jj,kk) < x_swarmG(jj,kk-1)
x_swarmG(jj,kk) = x_swarmG(jj,kk-1);
elseif x_swarmG(jj,kk) > x_swarmG(jj,kk+1)
x_swarmG(jj,kk) = x_swarmG(jj,kk+1);
end
end
if (kk==N_PAR)&&(kk~=1)
if x_swarmG(jj,kk) < x_swarmG(jj,kk-1)
x_swarmG(jj,kk) = x_swarmG(jj,kk-1);
elseif x_swarmG(jj,kk) > X_MAXG(kk)
x_swarmG(jj,kk) = X_MAXG(kk);
end
end
if (kk==1)&&(kk==N_PAR)
if x_swarmG(jj,kk) < X_MING(kk)
x_swarmG(jj,kk) = X_MING(kk);
elseif x_swarmG(jj,kk) > X_MAXG(kk)
x_swarmG(jj,kk) = X_MAXG(kk);
end
end
end
end

if (i<=N_SWARMSB)
for jj = 1:NB(i)
for kk = 1:N_PAR
if (kk==1)&&(kk~=N_PAR)
if x_swarmB(jj,kk) < X_MINB(kk)
x_swarmB(jj,kk) = X_MINB(kk);
elseif x_swarmB(jj,kk) > x_swarmB(jj,kk+1)
x_swarmB(jj,kk) = x_swarmB(jj,kk+1);
end
end
if ((kk>1)&&(kk<N_PAR))
if x_swarmB(jj,kk) < x_swarmB(jj,kk-1)
x_swarmB(jj,kk) = x_swarmB(jj,kk-1);
elseif x_swarmB(jj,kk) > x_swarmB(jj,kk+1)
x_swarmB(jj,kk) = x_swarmB(jj,kk+1);
end
end
if (kk==N_PAR)&&(kk~=1)
if x_swarmB(jj,kk) < x_swarmB(jj,kk-1)
x_swarmB(jj,kk) = x_swarmB(jj,kk-1);
elseif x_swarmB(jj,kk) > X_MAXB(kk)
x_swarmB(jj,kk) = X_MAXB(kk);
end
end
if (kk==1)&&(kk==N_PAR)
if x_swarmB(jj,kk) < X_MINB(kk)
 x_swarmB(jj,kk) = X_MINB(kk);
elseif x_swarmB(jj,kk) > X_MAXB(kk)
x_swarmB(jj,kk) = X_MAXB(kk);
end
end
end
end
end
end
if (i<=N_SWARMSR)
for jj=1:NR(i)

fit_swarmR(jj)=sum(probR(1:x_swarmR(jj,1)))*(sum((1:x_swarmR(jj,1)).*probR(1:x_swarmR(jj,1))/
sum(probR(1:x_swarmR(jj,1)))) - sum((1:Lmax).*probR(1:Lmax)) )^2;
for jlevel=2:level-1
fit_swarmR(jj)=fit_swarmR(jj)+sum(probR(x_swarmR(jj,jlevel-
1)+1:x_swarmR(jj,jlevel)))*(sum((x_swarmR(jj,jlevel-
1)+1:x_swarmR(jj,jlevel)).*probR(x_swarmR(jj,jlevel-1)+1:x_swarmR(jj,jlevel))/
sum(probR(x_swarmR(jj,jlevel-1)+1:x_swarmR(jj,jlevel))))- sum((1:Lmax).*probR(1:Lmax)))^2;
end
fit_swarmR(jj)=fit_swarmR(jj)+sum(probR(x_swarmR(jj,level-
1)+1:Lmax))*(sum((x_swarmR(jj,level-1)+1:Lmax).*probR(x_swarmR(jj,level-1)+1:Lmax)/
sum(probR(x_swarmR(jj,level-1)+1:Lmax)))- sum((1:Lmax).*probR(1:Lmax)))^2;
if fit_swarmR(jj) > fitBest_swarmR(jj)
fitBest_swarmR(jj) = fit_swarmR(jj);
xBest_swarmR(jj,:) = x_swarmR(jj,:);
end
end
end

if (i<=N_SWARMSG)
for jj=1:NG(i)

fit_swarmG(jj)=sum(probG(1:x_swarmG(jj,1)))*(sum((1:x_swarmG(jj,1)).*probG(1:x_swarmG(jj,1))/
sum(probG(1:x_swarmG(jj,1)))) - sum((1:Lmax).*probG(1:Lmax)) )^2;
for jlevel=2:level-1
fit_swarmG(jj)=fit_swarmG(jj)+sum(probG(x_swarmG(jj,jlevel-
1)+1:x_swarmG(jj,jlevel)))*(sum((x_swarmG(jj,jlevel-
1)+1:x_swarmG(jj,jlevel)).*probG(x_swarmG(jj,jlevel-1)+1:x_swarmG(jj,jlevel))/
sum(probG(x_swarmG(jj,jlevel-1)+1:x_swarmG(jj,jlevel))))- sum((1:Lmax).*probG(1:Lmax)))^2;
end
fit_swarmG(jj)=fit_swarmG(jj)+sum(probG(x_swarmG(jj,level-
1)+1:Lmax))*(sum((x_swarmG(jj,level-1)+1:Lmax).*probG(x_swarmG(jj,level-1)+1:Lmax)/
sum(probG(x_swarmG(jj,level-1)+1:Lmax)))- sum((1:Lmax).*probG(1:Lmax)))^2;
if fit_swarmG(jj) > fitBest_swarmG(jj)
fitBest_swarmG(jj) = fit_swarmG(jj);
xBest_swarmG(jj,:) = x_swarmG(jj,:);
end
end
end
 if (i<=N_SWARMSB)
for jj=1:NB(i)

fit_swarmB(jj)=sum(probB(1:x_swarmB(jj,1)))*(sum((1:x_swarmB(jj,1)).*probB(1:x_swarmB(jj,1))/
sum(probB(1:x_swarmB(jj,1)))) - sum((1:Lmax).*probB(1:Lmax)) )^2;
for jlevel=2:level-1
fit_swarmB(jj)=fit_swarmB(jj)+sum(probB(x_swarmB(jj,jlevel-
1)+1:x_swarmB(jj,jlevel)))*(sum((x_swarmB(jj,jlevel-
1)+1:x_swarmB(jj,jlevel)).*probB(x_swarmB(jj,jlevel-1)+1:x_swarmB(jj,jlevel))/
sum(probB(x_swarmB(jj,jlevel-1)+1:x_swarmB(jj,jlevel))))- sum((1:Lmax).*probB(1:Lmax)))^2;
end
fit_swarmB(jj)=fit_swarmB(jj)+sum(probB(x_swarmB(jj,level-
1)+1:Lmax))*(sum((x_swarmB(jj,level-1)+1:Lmax).*probB(x_swarmB(jj,level-1)+1:Lmax)/
sum(probB(x_swarmB(jj,level-1)+1:Lmax)))- sum((1:Lmax).*probB(1:Lmax)))^2;
if fit_swarmB(jj) > fitBest_swarmB(jj)
fitBest_swarmB(jj) = fit_swarmB(jj);
xBest_swarmB(jj,:) = x_swarmB(jj,:);
end
end
end
if (i<=N_SWARMSR)
[aR,bR] = max (fit_swarmR);
if (fit_swarmR(bR) > gbestvalue_swarmR)
gBest_swarmR=x_swarmR(bR,:)-1;
gbestvalue_swarmR = fit_swarmR(bR);
stagnancy_counterR(i)=0;
else
stagnancy_counterR(i)=stagnancy_counterR(i)+1;
end
end

if (i<=N_SWARMSG)
[aG,bG] = max (fit_swarmG);
if (fit_swarmG(bG) > gbestvalue_swarmG)
gBest_swarmG=x_swarmG(bG,:)-1;
gbestvalue_swarmG = fit_swarmG(bG);
stagnancy_counterG(i)=0;
else
stagnancy_counterG(i)=stagnancy_counterG(i)+1;
end
end

if (i<=N_SWARMSB)
[aB,bB] = max (fit_swarmB);
if (fit_swarmB(bB) > gbestvalue_swarmB)
gBest_swarmB=x_swarmB(bB,:)-1;
gbestvalue_swarmB = fit_swarmB(bB);
stagnancy_counterB(i)=0;
else
stagnancy_counterB(i)=stagnancy_counterB(i)+1;
end
 end

manteveR = 1;
manteveG = 1;
manteveB = 1;
if (i<=N_SWARMSR)
if( stagnancy_counterR(i)==0 && NR(i)<MAX_POP )
NR(i)=NR(i)+1;
x_swarmR(NR(i),:) = fix(rand(1,1) * ( X_MAXR(:)-X_MINR(:) ) + X_MINR(:));
v_swarmR(NR(i),:) = zeros (1,N_PAR);
xBest_swarmR(NR(i),:) = x_swarmR(NR(i),:);

fit_swarmR(NR(i))=sum(probR(1:x_swarmR(NR(i),1)))*(sum((1:x_swarmR(NR(i),1)).*probR(1:x_swar
mR(NR(i),1))/sum(probR(1:x_swarmR(NR(i),1)))) - sum((1:Lmax).*probR(1:Lmax)) )^2;
for jlevel=2:level-1
fit_swarmR(NR(i))=fit_swarmR(NR(i))+sum(probR(x_swarmR(NR(i),jlevel-
1)+1:x_swarmR(NR(i),jlevel)))*(sum((x_swarmR(NR(i),jlevel-
1)+1:x_swarmR(NR(i),jlevel)).*probR(x_swarmR(NR(i),jlevel-1)+1:x_swarmR(NR(i),jlevel))/
sum(probR(x_swarmR(NR(i),jlevel-1)+1:x_swarmR(NR(i),jlevel))))-
sum((1:Lmax).*probR(1:Lmax)))^2;
end
fit_swarmR(NR(i))=fit_swarmR(NR(i))+sum(probR(x_swarmR(NR(i),level-
1)+1:Lmax))*(sum((x_swarmR(NR(i),level-1)+1:Lmax).*probR(x_swarmR(NR(i),level-1)+1:Lmax)/
sum(probR(x_swarmR(NR(i),level-1)+1:Lmax)))- sum((1:Lmax).*probR(1:Lmax)))^2;
[aR,bR]=max(fit_swarmR);
gBest_swarmR=x_swarmR(bR,:);
gbestvalue_swarmR = fit_swarmR(bR);
xBest_swarmR = x_swarmR;
fitBest_swarmR(NR(i),:)=fit_swarmR(NR(i),:);
vR=Dimension(vR,v_swarmR(NR(i),:),kR);
vR(jR:kR+1,:) = v_swarmR;
xR=Dimension(xR,x_swarmR(NR(i),:),kR);
xR(jR:kR+1,:) = x_swarmR;
xBestR=Dimension(xBestR,xBest_swarmR(NR(i),:),kR);
xBestR(jR:kR+1,:) = xBest_swarmR;
fitR=Dimension(fitR,fit_swarmR(NR(i),:),kR);
fitR(jR:kR+1,:) = fit_swarmR;
fitBestR=Dimension(fitBestR,fitBest_swarmR(NR(i),:),kR);
fitBestR(jR:kR+1,:) = fitBest_swarmR;
Red(i,:)=gBest_swarmR;
gbestvalueR(i,:)=gbestvalue_swarmR;
kR=kR+1;
manteveR = 0;
end
end

if (i<=N_SWARMSG)
if( stagnancy_counterG(i)==0 && NG(i)<MAX_POP )
NG(i)=NG(i)+1;
x_swarmG(NG(i),:) = fix(rand(1,1) * ( X_MAXG(:)-X_MING(:) ) + X_MING(:));
v_swarmG(NG(i),:) = zeros (1,N_PAR);
 xBest_swarmG(NG(i),:) = x_swarmG(NG(i),:);

fit_swarmG(NG(i))=sum(probG(1:x_swarmG(NG(i),1)))*(sum((1:x_swarmG(NG(i),1)).*probG(1:x_swa
rmG(NG(i),1))/sum(probG(1:x_swarmG(NG(i),1)))) - sum((1:Lmax).*probG(1:Lmax)) )^2;
for jlevel=2:level-1
fit_swarmG(NG(i))=fit_swarmG(NG(i))+sum(probG(x_swarmG(NG(i),jlevel-
1)+1:x_swarmG(NG(i),jlevel)))*(sum((x_swarmG(NG(i),jlevel-
1)+1:x_swarmG(NG(i),jlevel)).*probG(x_swarmG(NG(i),jlevel-1)+1:x_swarmG(NG(i),jlevel))/
sum(probG(x_swarmG(NG(i),jlevel-1)+1:x_swarmG(NG(i),jlevel))))-
sum((1:Lmax).*probG(1:Lmax)))^2;
end
fit_swarmG(NG(i))=fit_swarmG(NG(i))+sum(probG(x_swarmG(NG(i),level-
1)+1:Lmax))*(sum((x_swarmG(NG(i),level-1)+1:Lmax).*probG(x_swarmG(NG(i),level-1)+1:Lmax)/
sum(probG(x_swarmG(NG(i),level-1)+1:Lmax)))- sum((1:Lmax).*probG(1:Lmax)))^2;
[aG,bG]=max(fit_swarmG);
gBest_swarmG=x_swarmG(bG,:);
gbestvalue_swarmG = fit_swarmG(bG);
xBest_swarmG = x_swarmG;
fitBest_swarmG(NG(i),:)=fit_swarmG(NG(i),:);
vG=Dimension(vG,v_swarmG(NG(i),:),kG);
vG(jG:kG+1,:) = v_swarmG;
xG=Dimension(xG,x_swarmG(NG(i),:),kG);
xG(jG:kG+1,:) = x_swarmG;
xBestG=Dimension(xBestG,xBest_swarmG(NG(i),:),kG);
xBestG(jG:kG+1,:) = xBest_swarmG;
fitG=Dimension(fitG,fit_swarmG(NG(i),:),kG);
fitG(jG:kG+1,:) = fit_swarmG;
fitBestG=Dimension(fitBestG,fitBest_swarmG(NG(i),:),kG);
fitBestG(jG:kG+1,:) = fitBest_swarmG;
Green(i,:)=gBest_swarmG;
gbestvalueG(i,:)=gbestvalue_swarmG;
kG=kG+1;
manteveG = 0;
end
end

if (i<=N_SWARMSB)
if( stagnancy_counterB(i)==0 && NB(i)<MAX_POP )
NB(i)=NB(i)+1;
x_swarmB(NB(i),:) = fix(rand(1,1) * ( X_MAXB(:)-X_MINB(:) ) + X_MINB(:));
v_swarmB(NB(i),:) = zeros (1,N_PAR);
xBest_swarmB(NB(i),:) = x_swarmB(NB(i),:);

fit_swarmB(NB(i))=sum(probB(1:x_swarmB(NB(i),1)))*(sum((1:x_swarmB(NB(i),1)).*probB(1:x_swar
mB(NB(i),1))/sum(probB(1:x_swarmB(NB(i),1)))) - sum((1:Lmax).*probB(1:Lmax)) )^2;
for jlevel=2:level-1
fit_swarmB(NB(i))=fit_swarmB(NB(i))+sum(probB(x_swarmB(NB(i),jlevel-
1)+1:x_swarmB(NB(i),jlevel)))*(sum((x_swarmB(NB(i),jlevel-
1)+1:x_swarmB(NB(i),jlevel)).*probB(x_swarmB(NB(i),jlevel-1)+1:x_swarmB(NB(i),jlevel))/
sum(probB(x_swarmB(NB(i),jlevel-1)+1:x_swarmB(NB(i),jlevel))))-
sum((1:Lmax).*probB(1:Lmax)))^2;
 end
fit_swarmB(NB(i))=fit_swarmB(NB(i))+sum(probB(x_swarmB(NB(i),level-
1)+1:Lmax))*(sum((x_swarmB(NB(i),level-1)+1:Lmax).*probB(x_swarmB(NB(i),level-1)+1:Lmax)/
sum(probB(x_swarmB(NB(i),level-1)+1:Lmax)))- sum((1:Lmax).*probB(1:Lmax)))^2;
[aB,bB]=max(fit_swarmB);
gBest_swarmB=x_swarmB(bB,:);
gbestvalue_swarmB = fit_swarmB(bB);
xBest_swarmB = x_swarmB;
fitBest_swarmB(NB(i),:)=fit_swarmB(NB(i),:);
vB=Dimension(vB,v_swarmB(NB(i),:),kB);
vB(jB:kB+1,:) = v_swarmB;
xB=Dimension(xB,x_swarmB(NB(i),:),kB);
xB(jB:kB+1,:) = x_swarmB;
xBestB=Dimension(xBestB,xBest_swarmB(NB(i),:),kB);
xBestB(jB:kB+1,:) = xBest_swarmB;
fitB=Dimension(fitB,fit_swarmB(NB(i),:),kB);
fitB(jB:kB+1,:) = fit_swarmB;
fitBestB=Dimension(fitBestB,fitBest_swarmB(NB(i),:),kB);
fitBestB(jB:kB+1,:) = fitBest_swarmB;
Blue(i,:)=gBest_swarmB;
gbestvalueB(i,:)=gbestvalue_swarmB;
kB=kB+1;
manteveB = 0;
end
end
if (i<=N_SWARMSR)
if( stagnancy_counterR(i)==0 && already_deletedR(i)==0 && N_SWARMSR <
MAX_SWARMS )
p=rand/N_SWARMSR;
cria_swarm=rand;
if (cria_swarm<=p)
for n=1:ceil( NR(i)/2 )
idx_parent = randi (NR(i),1);
x_nova(n,:) = x_swarmR(idx_parent,:);
end
swarm_aleatoria = randi(N_SWARMSR,1);
while (swarm_aleatoria == i)
swarm_aleatoria = randi(N_SWARMSR,1);
end
inicio=sum( NR(1:swarm_aleatoria-1) );
for m=1:ceil( NR(swarm_aleatoria)/2 )
idx_parent = randi ( NR(swarm_aleatoria),1 );
x_nova(n+m,:) = xR(inicio+idx_parent,:);
if (n+m >= POP_INICIAL)
break;
end
end
m = size(x_nova,1);
if (m<POP_INICIAL)
x_nova(m+1:POP_INICIAL,:) = fix(rand(POP_INICIAL-m,1) * (( X_MAXR-
X_MINR ) + X_MINR)); %new particle
 m=POP_INICIAL;
end
N_SWARMSR = N_SWARMSR+1;
NR(N_SWARMSR) = m;
n=size(xR,1);
o=n;
xR=[xR;x_nova];
xBestR=[xBestR;x_nova];
vR=[vR;zeros(m,N_PAR)];
for n=o+1:o+m

fitR(n)=sum(probR(1:xR(n,1)))*(sum((1:xR(n,1)).*probR(1:xR(n,1))/sum(probR(1:xR(n,1)))) -
sum((1:Lmax).*probR(1:Lmax)) )^2;
for jlevel=2:level-1
fitR(n)=fitR(n)+sum(probR(xR(n,jlevel-1)+1:xR(n,jlevel)))*(sum((xR(n,jlevel-
1)+1:xR(n,jlevel)).*probR(xR(n,jlevel-1)+1:xR(n,jlevel))/sum(probR(xR(n,jlevel-1)+1:xR(n,jlevel))))-
sum((1:Lmax).*probR(1:Lmax)))^2;
end
fitR(n)=fitR(n)+sum(probR(xR(n,level-1)+1:Lmax))*(sum((xR(n,level-
1)+1:Lmax).*probR(xR(n,level-1)+1:Lmax)/sum(probR(xR(n,level-1)+1:Lmax)))-
sum((1:Lmax).*probR(1:Lmax)))^2;
fitBestR(n)=fitR(n);
end
[aR,bR]=max(fitR(o+1:o+m,:));
Red(N_SWARMSR,:)=xR(o+bR,:);
gbestvalueR(N_SWARMSR,1) = fitR(o+bR);
already_deletedR(N_SWARMSR)=0;
stagnancy_counterR(N_SWARMSR)=0;
clear x_nova;
end
% end
end

if (i<=N_SWARMSG)
if( stagnancy_counterG(i)==0 && already_deletedG(i)==0 && N_SWARMSG <
MAX_SWARMS )
p=rand/N_SWARMSG;
cria_swarm=rand;
if (cria_swarm<=p)
for n=1:ceil( NG(i)/2 )
idx_parent = randi (NG(i),1);
x_nova(n,:) = x_swarmG(idx_parent,:);
end
swarm_aleatoria = randi(N_SWARMSG,1);
while (swarm_aleatoria == i)
swarm_aleatoria = randi(N_SWARMSG,1);
end
inicio=sum( NG(1:swarm_aleatoria-1) );
for m=1:ceil( NG(swarm_aleatoria)/2 )
idx_parent = randi( NG(swarm_aleatoria),1 );
x_nova(n+m,:) = xG(inicio+idx_parent,:);
 if (n+m >= POP_INICIAL)
break;
end
end
m = size(x_nova,1);
if (m<POP_INICIAL)
x_nova(m+1:POP_INICIAL,:) = fix(rand(POP_INICIAL-m,1) * (( X_MAXG-
X_MING ) + X_MING));
m=POP_INICIAL;
end
N_SWARMSG = N_SWARMSG+1;
NG(N_SWARMSG) = m;
n=size(xG,1);
o=n;
xG=[xG;x_nova];
xBestG=[xBestG;x_nova];
vG=[vG;zeros(m,N_PAR)];
for n=o+1:o+m

fitG(n)=sum(probG(1:xG(n,1)))*(sum((1:xG(n,1)).*probG(1:xG(n,1))/sum(probG(1:xG(n,1)))) -
sum((1:Lmax).*probG(1:Lmax)) )^2;
for jlevel=2:level-1
fitG(n)=fitG(n)+sum(probG(xG(n,jlevel-1)+1:xG(n,jlevel)))*(sum((xG(n,jlevel-
1)+1:xG(n,jlevel)).*probG(xG(n,jlevel-1)+1:xG(n,jlevel))/sum(probG(xG(n,jlevel-1)+1:xG(n,jlevel))))-
sum((1:Lmax).*probG(1:Lmax)))^2;
end
fitG(n)=fitG(n)+sum(probG(xG(n,level-1)+1:Lmax))*(sum((xG(n,level-
1)+1:Lmax).*probG(xG(n,level-1)+1:Lmax)/sum(probG(xG(n,level-1)+1:Lmax)))-
sum((1:Lmax).*probG(1:Lmax)))^2;
fitBestG(n)=fitG(n);
end
[aG,bG]=max(fitG(o+1:o+m,:));
Green(N_SWARMSG,:)=xG(o+bG,:);
gbestvalueG(N_SWARMSG,1) = fitG(o+bG);
already_deletedG(N_SWARMSG)=0;
stagnancy_counterG(N_SWARMSG)=0;
clear x_nova;
end
end
end

if (i<=N_SWARMSB)
if( stagnancy_counterB(i)==0 && already_deletedB(i)==0 && N_SWARMSB <
MAX_SWARMS )
p=rand/N_SWARMSB;
cria_swarm=rand;
if (cria_swarm<=p)
for n=1:ceil( NB(i)/2 )
idx_parent = randi (NB(i),1);
x_nova(n,:) = x_swarmB(idx_parent,:);
end
 swarm_aleatoria = randi(N_SWARMSB,1);
while (swarm_aleatoria == i)
swarm_aleatoria = randi(N_SWARMSB,1);
end
inicio=sum( NB(1:swarm_aleatoria-1) );
for m=1:ceil( NB(swarm_aleatoria)/2 )
idx_parent = randi ( NB(swarm_aleatoria),1 );
x_nova(n+m,:) = xB(inicio+idx_parent,:);
if (n+m >= POP_INICIAL)
break;
end
end
m = size(x_nova,1);
if (m<POP_INICIAL)
x_nova(m+1:POP_INICIAL,:) = fix(rand(POP_INICIAL-m,1) * (( X_MAXB-
X_MINB ) + X_MINB));
m=POP_INICIAL;
end
N_SWARMSB = N_SWARMSB+1;
NB(N_SWARMSB) = m;
n=size(xB,1);
o=n;
xB=[xB;x_nova];
xBestB=[xBestB;x_nova];
vB=[vB;zeros(m,N_PAR)];
for n=o+1:o+m

fitB(n)=sum(probB(1:xB(n,1)))*(sum((1:xB(n,1)).*probB(1:xB(n,1))/sum(probB(1:xB(n,1)))) -
sum((1:Lmax).*probB(1:Lmax)) )^2;
for jlevel=2:level-1
fitB(n)=fitB(n)+sum(probB(xB(n,jlevel-1)+1:xB(n,jlevel)))*(sum((xB(n,jlevel-
1)+1:xB(n,jlevel)).*probB(xB(n,jlevel-1)+1:xB(n,jlevel))/sum(probB(xB(n,jlevel-1)+1:xB(n,jlevel))))-
sum((1:Lmax).*probB(1:Lmax)))^2;
end
fitB(n)=fitB(n)+sum(probB(xB(n,level-1)+1:Lmax))*(sum((xB(n,level-
1)+1:Lmax).*probB(xB(n,level-1)+1:Lmax)/sum(probB(xB(n,level-1)+1:Lmax)))-
sum((1:Lmax).*probB(1:Lmax)))^2;
fitBestB(n)=fitB(n);
end
[aB,bB]=max(fitB(o+1:o+m,:));
Blue(N_SWARMSB,:)=xB(o+bB,:);
gbestvalueB(N_SWARMSB,1) = fitB(o+bB);
already_deletedB(N_SWARMSB)=0;
stagnancy_counterB(N_SWARMSB)=0;
end
end
end
end

if (i<=N_SWARMSR)
if (stagnancy_counterR(i) == STAGNANCY)
 [a,idx_pior]=min(fit_swarmR);
NR(i)=NR(i)-1;
x_swarmR(idx_pior,:)=[];
xBest_swarmR(idx_pior,:)=[];
v_swarmR(idx_pior,:)=[];
fit_swarmR(idx_pior,:)=[];
fitBest_swarmR(idx_pior,:)=[];
vR(jR+idx_pior-1,:)=[];
vR(jR:kR-1,:) = v_swarmR;
xR(jR+idx_pior-1,:)=[];
xR(jR:kR-1,:) = x_swarmR;
xBestR(jR+idx_pior-1,:)=[];
xBestR(jR:kR-1,:) = xBest_swarmR;
fitR(jR+idx_pior-1,:)=[];
fitR(jR:kR-1,:) = fit_swarmR;
fitBestR(jR+idx_pior-1,:)=[];
fitBestR(jR:kR-1,:) = fitBest_swarmR;
Red(i,:)=gBest_swarmR;
gbestvalueR(i,:)=gbestvalue_swarmR;
kR=kR-1;
manteveR = 0;
already_deletedR(i)=already_deletedR(i)+1;
stagnancy_counterR(i) = round( STAGNANCY * (1- ( 1/(already_deletedR(i)+1) ) ) );
end

if (i<=N_SWARMSG)
if (stagnancy_counterG(i) == STAGNANCY)
[a,idx_pior]=min(fit_swarmG);
NG(i)=NG(i)-1;
x_swarmG(idx_pior,:)=[];
xBest_swarmG(idx_pior,:)=[];
v_swarmG(idx_pior,:)=[];
fit_swarmG(idx_pior,:)=[];
fitBest_swarmG(idx_pior,:)=[];
vG(jG+idx_pior-1,:)=[];
vG(jG:kG-1,:) = v_swarmG;
xG(jG+idx_pior-1,:)=[];
xG(jG:kG-1,:) = x_swarmG;
xBestG(jG+idx_pior-1,:)=[];
xBestG(jG:kG-1,:) = xBest_swarmG;
fitG(jG+idx_pior-1,:)=[];
fitG(jG:kG-1,:) = fit_swarmG;
fitBestG(jG+idx_pior-1,:)=[];
fitBestG(jG:kG-1,:) = fitBest_swarmG;
Green(i,:)=gBest_swarmG;
gbestvalueG(i,:)=gbestvalue_swarmG;
kG=kG-1;
manteveG = 0;
already_deletedG(i)=already_deletedG(i)+1;
stagnancy_counterG(i) = round( STAGNANCY * (1- ( 1/(already_deletedG(i)+1) ) ) );
end
 end

if (i<=N_SWARMSB)
if (stagnancy_counterB(i) == STAGNANCY)
[a,idx_pior]=min(fit_swarmB);
NB(i)=NB(i)-1;
x_swarmB(idx_pior,:)=[];
xBest_swarmB(idx_pior,:)=[];
v_swarmB(idx_pior,:)=[];
fit_swarmB(idx_pior,:)=[];
fitBest_swarmB(idx_pior,:)=[];
vB(jB+idx_pior-1,:)=[];
vB(jB:kB-1,:) = v_swarmB;
xB(jB+idx_pior-1,:)=[];
xB(jB:kB-1,:) = x_swarmB;
xBestB(jB+idx_pior-1,:)=[];
xBestB(jB:kB-1,:) = xBest_swarmB;
fitB(jB+idx_pior-1,:)=[];
fitB(jB:kB-1,:) = fit_swarmB;
fitBestB(jB+idx_pior-1,:)=[];
fitBestB(jB:kB-1,:) = fitBest_swarmB;
Blue(i,:)=gBest_swarmB;
gbestvalueB(i,:)=gbestvalue_swarmB;
kB=kB-1;
manteveB = 0;
already_deletedB(i)=already_deletedB(i)+1;
stagnancy_counterB(i) = round( STAGNANCY * (1- ( 1/(already_deletedB(i)+1) ) ) );
end
end
end

if (i<=N_SWARMSR)
if (NR(i)==MIN_POP)&&(N_SWARMSR>MIN_SWARMS)
vR(jR:kR,:)=[];
xR(jR:kR,:)=[];
xBestR(jR:kR,:)=[];
fitR(jR:kR,:)=[];
fitBestR(jR:kR)=[];
Red(i,:)=[];
gbestvalueR(i,:)=[];
NR(i)=[];
already_deletedR(i)=[];
stagnancy_counterR(i)=[];
N_SWARMSR = N_SWARMSR - 1;
manteveR=0;
end
end

if (i<=N_SWARMSG)
if (NG(i)==MIN_POP)&&(N_SWARMSG>MIN_SWARMS)
vG(jG:kG,:)=[];
 xG(jG:kG,:)=[];
xBestG(jG:kG,:)=[];
fitG(jG:kG,:)=[];
fitBestG(jG:kG)=[];
Green(i,:)=[];
gbestvalueG(i,:)=[];
NG(i)=[];
already_deletedG(i)=[];
stagnancy_counterG(i)=[];
N_SWARMSG = N_SWARMSG - 1;
manteveG=0;
end
end

if (i<=N_SWARMSB)&&(N_SWARMSB>MIN_SWARMS)
if (NB(i)==MIN_POP)
vB(jB:kB,:)=[];
xB(jB:kB,:)=[];
xBestB(jB:kB,:)=[];
fitB(jB:kB,:)=[];
fitBestB(jB:kB)=[];
Blue(i,:)=[];
gbestvalueB(i,:)=[];
NB(i)=[];
already_deletedB(i)=[];
stagnancy_counterB(i)=[];
N_SWARMSB = N_SWARMSB - 1;
manteveB=0;
end
end

if (i<=N_SWARMSR)
if (manteveR==1)
vR(jR:kR,:) = v_swarmR;
Red(i,:) = gBest_swarmR;
gbestvalueR(i,1) = gbestvalue_swarmR;
xR(jR:kR,:) = x_swarmR;
xBestR(jR:kR,:) = xBest_swarmR;
fitR(jR:kR,:) = fit_swarmR;
fitBestR(jR:kR,:) = fitBest_swarmR;
end
end
if (i<=N_SWARMSG)
if (manteveG==1)
vG(jG:kG,:) = v_swarmG;
Green(i,:) = gBest_swarmG;
gbestvalueG(i,1) = gbestvalue_swarmG;
xG(jG:kG,:) = x_swarmG;
xBestG(jG:kG,:) = xBest_swarmG;
fitG(jG:kG,:) = fit_swarmG;
fitBestG(jG:kG,:) = fitBest_swarmG;
 end
end
if (i<=N_SWARMSB)
if (manteveB==1)
vB(jB:kB,:) = v_swarmB;
Blue(i,:) = gBest_swarmB;
gbestvalueB(i,1) = gbestvalue_swarmB;
xB(jB:kB,:) = x_swarmB;
xBestB(jB:kB,:) = xBest_swarmB;
fitB(jB:kB,:) = fit_swarmB;
fitBestB(jB:kB,:) = fitBest_swarmB;
end
end
i=i+1;
end
antigoR = LS3;
[LS3,LSC3]=max(gbestvalueR);
antigoG = LS2;
[LS2,LSC2]=max(gbestvalueG);
antigoB = LS1;
[LS1,LSC1]=max(gbestvalueB);
nger=nger+1;
end

Red = round(Red(LSC3,:));
gbestvalueR=LS3;
Green = round(Green(LSC2,:));
gbestvalueG=LS2;
Blue = round(Blue(LSC1,:));
gbestvalueB=LS1;
Red=sort(Red);
Green=sort(Green);
Blue=sort(Blue);
ReconImage=ImageReconstruction(I,Red,Green,Blue);
X=adapthisteq(ReconImage(:,:,1));
Y=adapthisteq(ReconImage(:,:,2));
Z=adapthisteq(ReconImage(:,:,3));
Cat=cat(3,X,Y,Z);
figure,imshow(Cat)
title('Convolutional Neural Network Algorithm for Satellite Image Segmentation');
figure,imshow(rgb2gray(Cat))
colormap('default');
title('Convolutional Neural Network Algorithm for Satellite Image Segmentation');
[psn,mse] = Parameters(Cat,Image);
disp('Parameters of Convolutional Neural Network Algorithm');
disp('PSNR');
disp(psn);
disp('Mean Square Error');
disp(mse);
I=Image;
Ref=20;
nBins = Ref;
Wavelength=1000000;
rHist = imhist(I(:,:,1), nBins)/Wavelength/2;
gHist = imhist(I(:,:,2), nBins)/Wavelength/1.5;
bHist = imhist(I(:,:,3), nBins)/Wavelength;
figure,
X=linspace(0,3,Ref);
h(1) = plot(X, rHist,'b');
hold on;
h(2) = plot(X , gHist,'g');
h(3) = plot(X , bHist,'r');
grid on;
xlabel('Wavelength');
ylabel('Reflectance');
legend('Mean Signature','Actionolite','Muscovite')
title('(a) Mean Value of Convolutional Neural Network');

disp('GLCM Feature Extraction for Segmented Image');

echo on;
% Autocorrelation
% amplitude
% Cross Correlation minima
% Auto Correlation maxima
% Cluster Prominence
% Cluster Shade
% Dissimilarity
% Energy
% Entropy
% Homogeneity
% Maximum probability
% Sum of Squares
% Sum Average
% Sum Variance
% Sum Entropy
% Difference variance
% Difference entropy
% Information measures of correlation (1)
% Information measures of correlation (2)
% Maximal correlation coefficient
% Inverse difference normalized (INN)
% Inverse difference moment normalized (IDN)
echo off;
FEA = graycomatrix(im2double(Cat(:,:,1)),'Offset',[2 0;0 2]);
statsL = TexFeatureEstm(FEA,0);
disp(statsL);