區域生長算法原理及MATLAB實現


1. 基於區域生長算法的圖像分割原理

數字圖像分割算法一般是基於灰度值的兩個基本特性之一:不連續性和相似性。前一種性質的應用途徑是基於圖像灰度的不連續變化分割圖像,比如圖像的邊緣。第二種性質的主要應用途徑是依據實現指定的准則將圖像分割為相似的區域。區域生長算法就是基於圖像的第二種性質,即圖像灰度值的相似性。

1.1 基本公式

令R表示整幅圖像區域,那么分割可以看成將區域R划分為n個子區域R1,,R2,......Rn的過程,並需要滿足以下條件:
a: U(Ri) = R;
b: Ri是一個連通區域,i=1,2,3,......n;
c: Ri ∩ Rj = 空集,對於任何的i,j;都有i≠j;
d: P(Ri) = Ture, 對i=1,2,......n;
e: R(Pi U Rj) = False, i≠j;
正如“區域生長”的名字所暗示的:區域生長是根據一種事先定義的准則將像素或者子區域聚合成更大區域的過程,並且要充分保證分割后的區域滿足a~e的條件。

1.2 區域生長算法設計思路

區域生長算法的設計主要由以下三點:生長種子點的確定,區域生長的條件,區域生長停止的條件。種子點的個數根據具體的問題可以選擇一個或者多個,並且根據具體的問題不同可以采用完全自動確定或者人機交互確定。
區域生長的條件實際上就是根據像素灰度間的連續性而定義的一些相似性准則,而區域生長停止的條件定義了一個終止規則,基本上,在沒有像素滿足加入某個區域的條件的時候,區域生長就會停止。在算法里面,定義一個變量,最大像素灰度值距離reg_maxdist.當待加入像素點的灰度值和已經分割好的區域所有像素點的平均灰度值的差的絕對值小於或等於reg_maxdist時,該像素點加入到已經分割到的區域。相反,則區域生長算法停止。
在種子店1的4鄰域連通像素中,即2、3、4、5點,像素點5的灰度值與種子點的灰度值最接近,所以像素點5被加入到分割區域中,並且像素點5會作為新的種子點執行后面的過程。在第二次循環過程中,由於待分析圖像中,即2、3、4、6、7、8,像素7的灰度值和已經分割的區域(由1和5組成)的灰度均值10.5最接近,所以像素點7被加入到分割區域中。最后一幅圖,示意了區域生長的方向(由淺入深)。
從上面的分析中,我們可以看出,在區域生長過程中,需要知道待分析像素點的編號(通過像素點的x和y坐標值來表示),同時還要知道這些待分析點的像素的灰度值。

1.3 區域生長算法偽代碼

begin
初始化變量 pixdist = 0;
交互式選擇一個種子點,並且初始化區域的灰度均值變量reg_mean為該種子點的灰度值。
while (pixdist < reg_maxdist)
將當前種子點的四鄰域像素點加入到鏈表neg_list中;
分別計算neg_list中所有元素的灰度值和reg_mean差的絕對值,並且得到最小值的元素i(x,y).令pixdist = abs(neg_list(i,3)-reg_mean);
更新 reg_mean = (reg_mean*reg_size + neg_list(i,3))/(reg_size + 1);(注:reg_size表示分割好區域內像素點的數量)
將舊的種子點標記為已經分割好的區域像素點;
將i(x,y)當做新的種子點,並將新的種子點i(x,y)從鏈表neg_list中移除
end
end

1.4 基於MATLAB的區域生長算法設計

 1 % Segment based on area, Region Growing;  2 clear all; close all; clc  3 [fileName,pathName] = uigetfile('*.*','Please select an image');%文件筐,選擇文件  4 if(fileName)  5     fileName = strcat(pathName,fileName);  6     fileName = lower(fileName);%一致的小寫字母形式  7 else 
 8     J = 0;%記錄區域生長所分割得到的區域  9     msgbox('Please select an image');  10     return; %退出程序  11 end  12  
 13 I = imread(fileName);  14 if( ~( size(I,3)-3 ))  15     I = rgb2gray(I);%轉化為單通道灰度圖  16 end  17 I = im2double(I); %圖像灰度值歸一化到[0,1]之間  18 Ireshape = imresize(I,[600,800]);  19 I = Ireshape(51:475,200:699);  20 gausFilter = fspecial('gaussian',[5 5],0.5);  21 I = imfilter(I,gausFilter,'replicate');  22  
 23 %種子點的交互式選擇  24 if( exist('x','var') == 0 && exist('y','var') == 0)  25     subplot(2,2,1),imshow(I,[]);  26  hold on;  27     [y,x] = getpts;%鼠標取點 回車確定  28     x = round(x(1));%選擇種子點  29     y = round(y(1));  30 end  31  
 32 if( nargin == 0)  33     reg_maxdist = 0.1;  34     %nargin是matlab代碼編寫中常用的一個技巧,主要用於計算當前主函數的輸入參數個  35     %數,一般可以根據nargin的返回值來確定主函數輸入參數的缺省值。在實現中,如果  36     %用戶輸入的參數個數為零,那么默認為0.2
 37 end  38 J = zeros(size(I)); % 主函數的返回值,記錄區域生長所得到的區域  39 Isizes = size(I);  40 reg_mean = I(x,y);%表示分割好的區域內的平均值,初始化為種子點的灰度值  41 reg_size = 1;%分割的到的區域,初始化只有種子點一個  42 neg_free = 10000; %動態分配內存的時候每次申請的連續空間大小  43 neg_list = zeros(neg_free,3);  44 %定義鄰域列表,並且預先分配用於儲存待分析的像素點的坐標值和灰度值的空間,加速  45 %如果圖像比較大,需要結合neg_free來實現matlab內存的動態分配  46 neg_pos = 0;%用於記錄neg_list中的待分析的像素點的個數  47 pixdist = 0;  48 %記錄最新像素點增加到分割區域后的距離測度  49 %下一次待分析的四個鄰域像素點和當前種子點的距離  50 %如果當前坐標為(x,y)那么通過neigb我們可以得到其四個鄰域像素的位置  51 neigb = [ -1 0;  52           1  0;  53           0 -1;  54           0  1];  55  %開始進行區域生長,當所有待分析的鄰域像素點和已經分割好的區域像素點的灰度值距離  56  %大於reg_maxdis,區域生長結束  57  
 58  while (pixdist < 0.06 && reg_size < numel(I))  59      %增加新的鄰域像素到neg_list中  60      for j=1:4
 61          xn = x + neigb(j,1);  62          yn = y + neigb(j,2);  63          %檢查鄰域像素是否超過了圖像的邊界  64          ins = (xn>=1)&&(yn>=1)&&(xn<=Isizes(1))&&(yn<=Isizes(1));  65          %如果鄰域像素在圖像內部,並且尚未分割好;那么將它添加到鄰域列表中  66          if( ins && J(xn,yn)==0)  67              neg_pos = neg_pos+1;  68              neg_list(neg_pos,:) =[ xn, yn, I(xn,yn)];%存儲對應點的灰度值  69              J(xn,yn) = 1;%標注該鄰域像素點已經被訪問過 並不意味着,他在分割區域內  70  end  71  end  72     %如果分配的內存空問不夠,申請新的內存空間  73     if (neg_pos+10>neg_free)  74         neg_free = neg_free + 100000;  75         neg_list((neg_pos +1):neg_free,:) = 0;  76  end  77     %從所有待分析的像素點中選擇一個像素點,該點的灰度值和已經分割好區域灰度均值的  78     %差的絕對值時所待分析像素中最小的  79     dist = abs(neg_list(1:neg_pos,3)-reg_mean);  80     [pixdist,index] = min(dist);  81     %計算區域的新的均值  82     reg_mean = (reg_mean * reg_size +neg_list(index,3))/(reg_size + 1);  83     reg_size = reg_size + 1;  84     %將舊的種子點標記為已經分割好的區域像素點  85     J(x,y)=2;%標志該像素點已經是分割好的像素點  86     x = neg_list(index,1);  87     y = neg_list(index,2);  88 %     pause(0.0005);%動態繪制  89 %     if(J(x,y)==2)  90 %     plot(x,y,'r.');  91 % end  92     %將新的種子點從待分析的鄰域像素列表中移除  93     neg_list(index,:) = neg_list(neg_pos,:);  94     neg_pos = neg_pos -1;  95  end  96  
 97  J = (J==2);%我們之前將分割好的像素點標記為2  98  hold off;  99  subplot(2,2,2),imshow(J); 100  J = bwmorph(J,'dilate');%補充空洞 101  subplot(2,2,3),imshow(J); 102  subplot(2,2,4),imshow(I+J);

2. 算法實驗分析

2.1 CT圖像實驗-分割肝臟實質


2.2 US圖像實驗-分割肝臟血管

          分割不完整
                                                                                                                                        過分割
      對有明顯邊界的分割效果好
 
 

3.算法應用總結

1.對CT圖像、MR圖像以及所有的具有邊界效應,或者是該區域與外界區域有明顯差距的圖像,分割效果很好。
2.對感興趣區域與外接區域存在邊緣連通現象的圖像,分割效果很差。如超聲圖像,肝臟對超聲的反應就是“均勻性”散點回聲。這造成“基於鄰域像素相似性”分割很難應用。閾值設置的小,造成分割不完整;閾值設置得太大,造成過分割現象。


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM