clear all, clc;

img = imread('q.jpg');
img = im2bw(img, graythresh(img));
imgDist = -bwdist(~img, 'cityblock');
figure, imshow(imgDist, [])
% imgDist(~img)=-inf; 
% figure, imshow(imgDist, [])

%% check local minimums
Dis_img = imregionalmin(imgDist);
figure, imshow(~Dis_img);  

se = strel('disk', 5);
marker = imdilate(Dis_img, se);
figure, imshow(marker);

% marker0 = zeros(size(marker));
% [label, num] = bwlabel(marker, 8);
% for i = 1 : num
%     cen = regionprops(label == i, 'centroid');   
%     cen = cat(1, cen.Centroid);   % cen(:,1): horizontal, cen(:,2): vertical
%     cen = floor(cen);
%     marker0(cen(2), cen(1)) = 1;
% end
% figure, imshow(marker0);

imgDist = imimposemin(imgDist, marker);
figure, imshow(imgDist, []);

bw = img;
dis = marker .* bw;
bw0 = img;
bw0(find(dis == 1)) = 0;
figure, imshow(bw0);

%% watershed
imgLabel = watershed(imgDist);    
bgm = imgLabel == 0;
figure, imshow(bgm,'InitialMagnification','fit')
img(bgm) = 0;
img = bwareaopen(img,5);
figure, imshow(img);

colorimg = bwlabel(img);
figure,imshow(label2rgb(colorimg,'jet','k','shuffle'));











