t_start = clock

pnorm=2; 
merge_thresh = 0.8;         % the threshold for cluster merging.
c=6;  c0=4;                % the range of the number of clusters.
cluster_thresh = 0.2;     % the threshold for clustering [0.01, 0.0001].
g = 0;
m = 2;
r = 1;                      % the number of linearly independent vectors
DB_Old = realmax;

ImgFile = '..\test images\baboon.jpg';

[Img,map] = imread(ImgFile);
if (isind(Img))
    Img = ind2rgb(Img, map);
end
Img = im2double(Img);
Img = imresize(Img, [128 128]);
LAB = rgb2LAB(Img);
[rows,cols,dim] = size(LAB);

leng = rows*cols;
feature = reshape(LAB, leng, dim);
feature = feature(:,2:3); dim = dim-1;

U = U_initialize(c,leng);
V = V_initialize(U,feature,m,dim);

for k = 1:100
    [U,V] = ufc_pr(feature,U,V,c,r,m,g,pnorm,cluster_thresh);
    [FR,DB] = FuzzyDBIndex(feature,U,V,c,m);
    FRmax = max(max(FR))
    if FRmax < merge_thresh
        break;
    end
    if DB > DB_Old | c <= c0
        break;
    end
    [i j] = find(FR == FRmax);
    c = c - 1
    U(i,:) = U(i,:) + U(j,:);
    V(i,:) = (V(i,:)+V(j,:)) / 2;
    U(j(1),:) = []; V(j(1),:) = [];
end 

a = eye(c) - tril(ones(c,c),-1);
Umax = (U == repmat(max(U),c,1));
[Labels,j] = find(a*Umax == 1);

Labels = reshape(Labels, rows, cols);
figure, imshow(mat2gray(Labels));


t_process = etime(clock,t_start)