一、簡介

分水嶺演算法是一種影象區域分割法，分割的過程中將圖片轉化為灰度圖，然後我會將灰度值看作是海拔，然後向較低點注水，這種基於地形學的解釋，我們著重考慮三種點：

極小值點，該點對應一個盆地的最低點，當我們在盆地裡滴一滴水的時候，由於重力作用，水最終會匯聚到該點。注意：可能存在一個最小值面，該平面內的都是極小值點。
盆地的其它位置點，該位置滴的水滴會匯聚到區域性最小點。
盆地的邊緣點，是該盆地和其它盆地交接點，在該點滴一滴水，會等概率的流向任何一個盆地。

明白上述三種點之後，我們開始往盆地的極小值點注水，然後隨著注水的深入，每一個極小值點慢慢的向外擴充套件，然後知道兩個盆地的水匯合，匯合處就是我們需要的分水嶺。

從下圖可以直觀理解一下，首先這三塊區域都含有極小值點

然後逐漸填充就能獲得分水嶺（即分界線）

得到分界線就能完成影象分割：

二、原始碼

clear;
close all;
%------------------
%程式中定義影象變數說明
%Image->原圖變數;
%Image_BW->二值化圖象;
%Image_BW_medfilt->中值濾波後的二值化影象;
%Optimized_Image_BW-〉通過“初次二值化影象”與“中值濾波後的二值化影象”進行“或”運算優化影象效果；
%Reverse_Image_BW-〉優化後二值化圖象取反；
%Filled_Image_BW-〉已填充背景色的二進位制影象；
%Open_Image_BW-〉開運算後的影象。
%------------------
%--------------------------------------
%-------圖片前期處理-------------------
%--------------------------------------
%第一步：讀取原圖,並顯示
a= imread('b1.bmp');
subplot(331)
imshow(a);
title('原圖');
Image=rgb2gray(a);
%第二步：進行二值化
Theshold = graythresh(Image);%取得圖象的全域性域值
Image_BW = im2bw(Image,Theshold);%二值化圖象
subplot(332)
imshow(Image_BW);
title('初次二值化影象');

%第三步二值化影象進行
Image_BW_medfilt= medfilt2(Image_BW,[13 13]);
subplot(333)
imshow(Image_BW_medfilt);
title('中值濾波後的二值化影象');

%第四步：通過“初次二值化影象”與“中值濾波後的二值化影象”進行“或”運算優化影象效果
Optimized_Image_BW = Image_BW_medfilt|Image_BW;
subplot(334)
imshow(Optimized_Image_BW);
title('進行“或”運算優化影象效果');
%第五步：優化後二值化圖象取反，保證：‘1’-〉‘白色’，‘0’-〉‘黑色’
%方便下面的操作
Reverse_Image_BW = ~Optimized_Image_BW;
subplot(335)
imshow(Reverse_Image_BW);
title('優化後二值化圖象取反');

%第六步：填充二進位制影象的背景色，去掉細胞內的黑色空隙

Filled_Image_BW = bwfill(Reverse_Image_BW,'holes');
subplot(336) 
imshow(Filled_Image_BW);
title('已填充背景色的二進位制影象');
%第七步：對影象進行開運算，去掉細胞與細胞之間相粘連的部分
SE = strel('disk',10);
BW = imopen(Filled_Image_BW,SE);
subplot(337)
imshow(BW);
title('開運算後的影象');

subplot(338),imhist(Image);

%-----------------------------------------------

%-------------開始計算細胞數--------------------

%-----------------------------------------------

[Label Number]=bwlabel(BW,8)%初步取得細胞個數
Array = bwlabel(BW,8);%取得貼標籤處理後的影象
Sum = [];
%依次統計貼標籤後陣列
for i=1:Number
[r,c] = find(Array==i);%獲取相同標籤號的位置，將位置資訊存入[r,c]
rc = [r c];
Num = length(rc);%取得vc陣列的元素的個數
Sum([i])=Num;%將元素個數存入Sum陣列
end
Sum
N = 0;
%假如Sum陣列中的元素大於了1500，表示有兩個細胞相連，畫素點較多，即分為兩個細胞數
ss=0;
for i=1:length(Sum)
    ss=ss+Sum([i]); 

end
dd=ss/length(Sum);
dd
for i=1:length(Sum)
if(Sum([i]))>dd
N = N+1;
end
if(Sum([i])) >(dd*2)
N = N+2;
end
function [cen,copy]=Kmeans(I,k)       % K均值聚類函式
mm=max(max(I));                   % mm為I的最大灰度值
copy=I;   
cen=(1:k)*mm/k;                    % 求出首次中心值
cen1=zeros(1,k);
d=ones(1,k)./255;
[m,n]=size(I);
J=zeros(m,n);
while(1)
    for i=1:m
        for j=1:n
            c=abs(I(i,j)-cen);          % c為I中元素與各中心點的差值陣列
            cc=find(min(c)==c);       % cc為距中心值最近的元素下標
            J(i,j)=cc;
        end
    end
    % 求新的中心點
    for r=1:k
        h=0;
        [J1,J2]=find(J==r);
        p=length(J1);
        for s=1:p
            h=h+I(J1(s),J2(s));
            copy(J1(s),J2(s))=cen(r);    % 按中心值，把原影象灰度分為K類
        end
        cen1(r)=h/p;
    end
    if (abs(cen1-cen)<=d)
        break;
    else cen=cen1;
    end
end


I=imread('b.bmp');       % 讀入處理影象
R=I(:,:,1);                 % 分離R分量
G=I(:,:,2);
B=I(:,:,3);
clc;                      % 清屏
close all;                 % 關閉影象視窗
figure,                   % 開闢影象顯示視窗
subplot(2,3,1),imshow(R),title('R分量');
subplot(2,3,2),imshow(G),title('G分量');
subplot(2,3,3),imshow(B),title('B分量');
K=rgb2hsv(I);            % RGB影象轉換成HSV影象
H=K(:,:,1);
S=K(:,:,2);
V=K(:,:,3);
subplot(2,3,4),imshow(H),title('H分量');
subplot(2,3,5),imshow(S),title('S分量');
subplot(2,3,6),imshow(V),title('V分量');

%****************************** 白細胞計數 ******************************
% 分割出只有白細胞的二值影象 
[cenwhite,copywhite]=Kmeans(H,3);      % K均值聚類函式呼叫       
cenwhite                     % 檢視聚類後的中心值
figure,
subplot(1,3,1),imshow(copywhite,[ ]),title('H分量K均值聚類結果');
A=copywhite;
[m,n]=size(A);                % 把聚類後圖像二值化
for i=1:m
    for j=1:n
        if  A(i,j)==cenwhite(2)
           A(i,j)=1;
        else A(i,j)=0;
        end
    end
end
subplot(1,3,2),imshow(A,[ ]),title('白細胞二值化結果');
C1=imdilate(A,ones(5));        % 對二值化影象腐蝕膨脹處理達到去噪效果
C2=imerode(C1,ones(11));
C=imdilate(C2,ones(6));
subplot(1,3,3),imshow(C,[ ]),title('白細胞去噪結果');

% 下面兩種方法中任選一種作為白細胞計數輸出：
% 方法一：用4連通分割影象計算白細胞個數 
L4=bwlabel(C,4);          % 4連通運算
figure,subplot(1,2,1),imshow(L4,[ ]),title('法一：白細胞4連通結果');
Numwhite=max(max(L4));   % 計算4連通後圖像中的最大值，即為白細胞個數
Numwhite                % Numwhite為4連通計數結果，顯示結果

% 方法二：用分水嶺演算法分割影象計算白細胞個數 %
E=bwdist(~C);
E=-E;
E(~C)=-Inf;
L0=watershed(E);
rgb=label2rgb(L0,'jet');
subplot(1,2,2),imshow(rgb,[ ]),title('法二：白細胞分水嶺結果');
NumWhite=max(max(L0))-1;               
NumWhite             % NumWhite為分水嶺演算法計數結果，顯示結果

%****************************** 紅細胞計數 ******************************
% 分割出只有紅細胞的二值影象，並對其做必需處理
I1(:,:,1)=double(I(:,:,1)).*C;        % 把二值影象還原為彩色影象
I1(:,:,2)=double(I(:,:,2)).*C;
I1(:,:,3)=double(I(:,:,3)).*C;
figure,subplot(1,2,1),imshow(uint8(I1)),title('白細胞');
J=double(I)-double(I1);           % 原圖減去白細胞圖得紅細胞彩圖 
subplot(1,2,2),imshow(uint8(J)),title('紅細胞');
J=rgb2hsv(J);                   % 紅細胞圖轉換成HSV圖
HH=J(:,:,1);
SS=J(:,:,2);                     % 提取S分量
VV=J(:,:,3);
figure,
subplot(1,3,1),imshow(HH),title('紅細胞圖H分量');
subplot(1,3,2),imshow(SS),title('紅細胞圖S分量');
subplot(1,3,3),imshow(VV),title('紅細胞圖V分量');
[cenred,copyred]=Kmeans(SS,2);  % 對S分量圖均值聚類
cenred                        % 檢視聚類後的中心值
figure,subplot(2,2,1),imshow(copyred,[ ]),title('S分量K均值聚類結果');
H=copyred;
[m,n]=size(H);                 % 把聚類後圖像二值化
for i=1:m
    for j=1:n
        if  H(i,j)==cenred(2)
           H(i,j)=1;
        else H(i,j)=0;
        end
    end
end
subplot(2,2,2),imshow(H,[ ]),title('紅細胞二值化結果');
G=imfill(H,'holes');         % 紅細胞圖填充孔洞
subplot(2,2,3),imshow(G,[ ]),title('紅細胞填洞結果');
C1=imerode(G,ones(9));     % 開運算去除噪聲
G=imdilate(C1,ones(9));
subplot(2,2,4),imshow(G,[ ]),title('紅細胞去噪結果');

% 方法一：通過判斷面積大小計算紅細胞個數 
[J,num]=bwlabel(G,8);    % 8連通運算，儲存區域個數num
figure,subplot(1,2,1),imshow(J,[ ]),title('紅細胞8連通結果');
sum=0;              % 定義sum為總面積變數
numred1=0;          % 定義numred1為紅細胞個數變數
for i=1:num
    [J1,J2]=find(J==i);    % 找出J中值為i的畫素座標分別存入J1,J2中
    area(i)=length(J1);    % area(i)為第i個紅細胞的面積
    sum=sum+area(i);    % 面積累加
end
ave=sum/num;           % 求紅細胞的平均面積ave
for i=1:num             % 通過面積比較，判斷個數的增加值
    if 0.5*ave<area(i) && area(i)<=1.5*ave
        numred1=numred1+1;    % 如果細胞區域的面積在0.5～1.5倍之間，視為一個細胞
    end
    if 1.5*ave<area(i) && area(i)<=2.5*ave
        numred1=numred1+2;    % 如果細胞區域的面積在1.5～2.5倍之間，視為兩個細胞
    end
    if 2.5*ave<area(i) && area(i)<=3.5*ave
        numred1=numred1+3;    % 如果細胞區域的面積在2.5～3.5倍之間，視為三個細胞
    end
 

三、執行結果

四、備註

版本：2014a
完整程式碼或代寫加1564658423