腦外科CT圖像的綜合邊緣提取算法

寧 祎，高紅偉

（河南工業(yè)大學(xué) 機(jī)器人研究所，河南 鄭州 450001）

眾所周知，腦外科手術(shù)是最復(fù)雜而危險的外科手術(shù)，只要一個微小的錯誤動作，就足以造成中樞神經(jīng)系統(tǒng)的永久損害。因此若需進(jìn)行腦外科手術(shù)，事先必須縝密的規(guī)劃，并做好充分的準(zhǔn)備工作，包括分毫不差地確認(rèn)腫瘤位置、決定最佳的切除方式以及了解腦細(xì)胞與其他重要大腦部位的相對位置等。 隨 著計算機(jī) 技術(shù) 、CT （Computed Tomography）、MRI（Magnetic Resonance Imaging）等醫(yī)學(xué)影像技術(shù)的不斷發(fā)展，虛擬現(xiàn)實(shí)技術(shù)也越來越多地應(yīng)用到現(xiàn)代醫(yī)療領(lǐng)域[1]。利用計算機(jī)圖像處理和數(shù)據(jù)可視化技術(shù)，根據(jù)醫(yī)學(xué)影像設(shè)備提供的二維斷層圖像，進(jìn)行人體器官的三維重建已是現(xiàn)代醫(yī)學(xué)重要發(fā)展方向之一。圖像最基本的特征是邊緣，邊緣是區(qū)域?qū)傩园l(fā)生突變的地方，也是圖像中不確定性最大的地方和圖像信息最集中的地方，包含著圖像的豐富信息[2－3]。圖像的邊緣對大腦的三維重建具有重要意義，是三維重建的重要基礎(chǔ)，圖像邊緣提取的正確性和可靠性將直接影響到三維重建后的效果。CT二維圖像的邊緣提取作為器官三維重建的第一步，如何快速、準(zhǔn)確地提取圖像邊緣信息一直受到國內(nèi)外學(xué)者的關(guān)注。1959年，Julez[4]最早提出了邊緣檢測；1965年，Roberts[5]最早開始系統(tǒng)地研究邊緣檢測；1986年，Canny提出基于最優(yōu)化算法的邊緣檢測算子。經(jīng)過50多年的研究，已經(jīng)出現(xiàn)了許多不同的邊緣檢測方法，各自有其特點(diǎn)和局限性，但針對腦外科CT圖像的研究報道并不多。因此，本研究針對腦外科CT圖像紋理復(fù)雜、細(xì)節(jié)豐富等特點(diǎn)，設(shè)計了一種采用圖像增強(qiáng)、Canny算子與數(shù)學(xué)形態(tài)學(xué)相結(jié)合的綜合邊緣提取算法，達(dá)到了理想的邊緣提取效果。

1 本文綜合算法的提出與圖像增強(qiáng)

1.1 本文綜合算法的提出

經(jīng)典的邊緣檢測算法是對原始圖像中像素的某小鄰域來構(gòu)造邊緣檢測算子，常用的邊緣檢測算子有Prewitt、Log、Roberts和Canny算子等。近些年，隨著數(shù)學(xué)形態(tài)學(xué)理論不斷完善與發(fā)展，數(shù)學(xué)形態(tài)學(xué)在圖像邊緣檢測中也得到了廣泛的應(yīng)用[6]。Prewitt算子的特點(diǎn)是抑制噪聲效果好，但對邊緣的定位不準(zhǔn)確；Log算子的特點(diǎn)是邊緣定位較準(zhǔn)確，但對噪聲較敏感，對圖像中的某些邊緣容易產(chǎn)生雙重響應(yīng)；Roberts算子的特點(diǎn)是計算時采用較少的像素，邊緣定位準(zhǔn)確，但是對噪聲高度敏感。Canny算子是在圖像處理中人們公認(rèn)的功能較為強(qiáng)大的檢測算子之一，但有時對于灰度分布不均勻的圖像，容易將非邊緣點(diǎn)檢測為邊緣。數(shù)學(xué)形態(tài)學(xué)是利用一個結(jié)構(gòu)元素去探測一個圖像，看是否能夠?qū)⑦@個結(jié)構(gòu)元素很好地填放在圖像的內(nèi)部，同時驗(yàn)證填放結(jié)構(gòu)元素的方法是否有效。數(shù)學(xué)形態(tài)學(xué)的腐蝕算法會導(dǎo)致圖像區(qū)域的分裂，而膨脹算法容易產(chǎn)生零星的像素點(diǎn)。

從目前的研究情況來看，單純采用某種邊緣檢測算法很難實(shí)現(xiàn)腦外科CT圖像準(zhǔn)確的邊緣提取。因而本文針對腦外科圖像紋理復(fù)雜、細(xì)節(jié)豐富等特點(diǎn)，考慮到醫(yī)學(xué)圖像數(shù)據(jù)多層切片間的相關(guān)性，提出了一種綜合的圖像邊緣檢測算法，該算法集成了用平均值濾波器對圖像增強(qiáng)，Canny邊緣檢測算子與數(shù)學(xué)形態(tài)學(xué)相結(jié)合的算法。

該算法的實(shí)現(xiàn)過程是在Matlab 7.5軟件上對圖像進(jìn)行處理分析，運(yùn)用了Matlab圖像處理工具箱 （Image Processing Toolbox）中的大量函數(shù)。本文邊緣提取算法的主要步驟為：

1）對原始圖像運(yùn)用平均值濾波器來對圖像增強(qiáng)；

2）對增強(qiáng)后的圖像用Canny算子進(jìn)行邊緣檢測；

3）利用數(shù)學(xué)形態(tài)學(xué)中的膨脹算法進(jìn)行邊緣處理。

1.2 圖像的增強(qiáng)

圖像增強(qiáng)是數(shù)字圖像處理過程中經(jīng)常采用的一種方法，這種根據(jù)圖像的特點(diǎn)或存在的問題采取的簡單改善方法或者加強(qiáng)特征的措施稱為圖像增強(qiáng)（Image Enhancement）。圖像增強(qiáng)的目標(biāo)有2個：第一是改善圖像的視覺效果，提高圖像的清晰度；第二是使圖像變得更利于計算機(jī)處理。所有的圖像處理就是對這些組成圖像的象素灰度值進(jìn)行增或減的計算處理。從增強(qiáng)處理的作用域出發(fā)，圖像增強(qiáng)的方法分為兩大類：空間域方法和頻域方法。由于CT圖像存在大量的噪聲，本文采用了平滑線性濾波器增強(qiáng)的方法，它是空間域的一種。平滑線性濾波器的輸出是包含在濾波掩模鄰域內(nèi)像素的簡單平均值，因此也稱為均值濾波器[7]。

平均值濾波器定義如下，對于給定的圖像 f（i，j）中的每個像素點(diǎn)（m，n），取其鄰域S。設(shè)S含有M個像素，取其平均值作為處理后所得圖像像素點(diǎn)（m，n）處的灰度值，即用一像素鄰域內(nèi)各像素灰度平均值來代替該像素原來的灰度值。鄰域S的形狀和大小根據(jù)圖像特點(diǎn)確定，一般取的形狀是正方形、矩形及十字形等，鄰域S的形狀和大小可以在全圖處理過程中保持不變，也可根據(jù)圖像的局部統(tǒng)計特性而變化，點(diǎn)（m，n）一般位于鄰域 S的中心。

文中根據(jù)腦外科圖像的特點(diǎn)和大量的實(shí)驗(yàn)比較，得出比較理想的效果是采用鄰域S為正方形5×5鄰域，像素點(diǎn)（m，n）位于鄰域S的中心，共25個像素點(diǎn)，通過公式（1）計算該點(diǎn)的灰度平均值：

式中，fˉ（m，n）為像素灰度平均值；m，n分別為像素點(diǎn)的橫、縱坐標(biāo)值；i，j分別為鄰域S的橫、縱坐標(biāo)值。

原始圖像如圖1所示。在Matlab命令窗口中，首先應(yīng)用imread（filename）函數(shù)將原始圖像導(dǎo)入到Matlab中，然后應(yīng)用filter2（fspecial（‘a(chǎn)verage’，5），I）/255 語句對原始圖像進(jìn)行增強(qiáng)，增強(qiáng)后圖像如圖2所示。

圖1 原始圖像Fig.1 Original image

圖2 圖像增強(qiáng)后Fig.2 Image enhancement

從圖2可以看出，經(jīng)過平滑線性濾波器后的圖像，去除了CT圖像中因掃描存在的噪聲、平滑了圖像信號，為后面能夠更好地進(jìn)行邊緣檢測打下了基礎(chǔ)。

2 Canny算子的邊緣檢測

邊緣檢測是對圖像的邊緣進(jìn)行處理獲得需要的邊緣信息，常用的邊緣檢測算子有Prewitt、Log、Roberts和Canny算子等。每種邊緣檢測算子都有各自的優(yōu)缺點(diǎn)，Prewitt算子抑制噪聲效果好，但對邊緣的定位不準(zhǔn)確，Log和Roberts算子都是邊緣定位較準(zhǔn)確，但對噪聲較敏感，因此本文選用了具有高信噪比、高定位精度和單邊緣響應(yīng)等優(yōu)良性能的Canny算子。

Canny算子采用二維高斯函數(shù)的任意方向上的一階方向?qū)?shù)為噪聲濾波器，通過與圖像卷積進(jìn)行濾波，然后對濾波后的圖像尋找圖像梯度的局部最大值，以此來確定圖像邊緣[8]。由于噪聲的影響，一個閾值的判斷往往是不夠的，邊緣信號的響應(yīng)只有近一半是大于這個閾值的，由此造成了邊緣斷裂。如果降低這個閾值，又會發(fā)現(xiàn)一些錯誤的“邊緣”。為了解決這個問題，Canny提出了雙閾值方法[9]。Canny算子進(jìn)行邊緣檢測的具體步驟如下：

1）用二維高斯濾波器對圖像進(jìn)行濾波，圖像通過二維高斯濾波函數(shù)。

進(jìn)行平滑，抑制圖像噪聲，其中σ為高斯濾波器參數(shù)，它控制著平滑程度。

2）用導(dǎo)數(shù)算子找到圖像灰度沿著兩個方向的偏導(dǎo)數(shù)（Gx，Gy），并計算出梯度的大小和方向。

3）對梯度幅值進(jìn)行非極大值抑制。僅僅得到全局的梯度并不足以確定邊緣，因此為確定邊緣，必須保留局部梯度最大的點(diǎn)，抑制非極大值。若某個像素的灰度值與其梯度方向上前后兩個像素的灰度值相比不是最大的，那么該像素值置為0，即不是邊緣，這個過程稱為“非極大值抑制”。

4）用雙閾值法檢測和連接邊緣。對梯度取兩次閾值Tl和T2，T1=T2，0＜α＜1，可根據(jù)噪聲的大小適當(dāng)調(diào)整。 用 T1、T2兩個閾值對非極大值抑制圖像進(jìn)行雙閾值化后，可得兩個檢測結(jié)果，分別記為圖像1和圖像2。圖像l閾值較低，則保留了較多信息；圖像2閾值較高，所以噪聲較少，但會造成邊緣信息的損失。于是以圖像2為基礎(chǔ)，以圖像l為補(bǔ)充，連接圖像的邊緣。

應(yīng)用 Matlab 圖像處理工具箱中 edge（I，’ ’，thresh）函數(shù)，幾種經(jīng)典的邊緣檢測算子對增強(qiáng)后的圖像進(jìn)行邊緣檢測，效果如圖3所示。

通過對圖3的對比分析，可以清楚地看到：圖 3（a）Prewitt算子對噪聲有抑制作用，檢測的圖像沒有較多孤立的邊緣像素點(diǎn)，但對邊緣的定位不是很準(zhǔn)確，圖像的邊界損失較多；圖3（b）Log算子檢測的圖像邊緣信息損失較少，但對圖像中的某些邊緣產(chǎn)生雙重響應(yīng)；圖3（c）Roberts算子對邊緣定位比較準(zhǔn)確，但存在較多孤立的像素點(diǎn)，抑制噪聲能力差；圖3（d）用Canny算子檢測后的效果最好，由于Canny算子具有非極大值抑制和雙閾值法連接邊緣等優(yōu)點(diǎn)，比前3種提取的邊緣清晰而且連續(xù)，虛假邊緣少，但是提取出的邊緣仍然存在不光滑和斷裂的現(xiàn)象。為了實(shí)現(xiàn)更好的邊緣檢測效果，下一步就要用本文提出的數(shù)學(xué)形態(tài)學(xué)算法來具體解決該問題。

3 數(shù)學(xué)形態(tài)學(xué)邊緣檢測

3.1 數(shù)學(xué)形態(tài)學(xué)的基本原理

數(shù)學(xué)形態(tài)學(xué)是一門綜合了多學(xué)科知識的交叉科學(xué)，它建立在嚴(yán)格的數(shù)學(xué)理論基礎(chǔ)上。其基本思想是用具有一定形態(tài)的結(jié)構(gòu)元素 SE（structuring element）來度量和提取圖像中的對應(yīng)形狀，以達(dá)到對圖像的分析和識別的目的[10]。

膨脹和腐蝕是形態(tài)學(xué)的兩個基本操作，膨脹過程是在圖像中對圖象的邊界添加像素點(diǎn)，向外擴(kuò)張，而腐蝕是它的逆過程。膨脹和腐蝕這兩種變換都對灰度值變化明顯的圖像邊緣較為敏感，膨脹可以填充圖像中的小孔及圖像邊緣上較小的凹陷部分，腐蝕可以消除圖像中細(xì)小的成分。下面介紹2個基本運(yùn)算的定義：

膨脹運(yùn)算定義為：

設(shè)A是原始圖像，B是結(jié)構(gòu)元素，式中A、B為Z2中的集合，φ為空集，B?為B的映射，A被B膨脹，記為A⊕B，⊕為膨脹算子。該式表明膨脹過程為B首先做關(guān)于原點(diǎn)的映射，然后平移x。A被B膨脹是被所有x平移后與A至少有一個非零的公共元素。

腐蝕運(yùn)算定義為：

圖3 邊緣檢測Fig.3 Edge detection

設(shè)A是原始圖像，B是結(jié)構(gòu)元素，式中A、B為Z2中的集合，A被B腐蝕，記為AΘB。該式表明A被B腐蝕的結(jié)果是所有使B被x平移后包含于A的點(diǎn)x的集合。

3.2 數(shù)學(xué)形態(tài)學(xué)方法在圖像處理中的運(yùn)用

應(yīng)用 Matlab圖像處理工具箱中 imdilate（I，SE）函數(shù)，對經(jīng)過Canny算子邊緣檢測后的圖像進(jìn)行數(shù)學(xué)形態(tài)法中的膨脹運(yùn)算。利用數(shù)學(xué)形態(tài)學(xué)進(jìn)行圖像處理時，選擇簡單、表現(xiàn)力強(qiáng)的結(jié)構(gòu)元素是關(guān)鍵，是形態(tài)變換中最重要的參數(shù)；其次，還要綜合考慮目標(biāo)體的清晰度和噪聲的大小來選取結(jié)構(gòu)元素的大小。一般目標(biāo)體輪廓不清晰時，選擇較小的結(jié)構(gòu)元素；噪聲顆粒較大時，選擇較大的結(jié)構(gòu)元素。本文根據(jù)腦外科圖像的特點(diǎn)和大量的實(shí)驗(yàn)比較，采用平面結(jié)構(gòu)元素為5×5的正方形，處理后的效果如圖4所示。

圖4 膨脹后效果圖Fig.4 Dilated image

從圖4中可以看到，經(jīng)過膨脹后得到的邊緣清晰而且光滑連續(xù)，基本上消除了圖3（d）中所存在的邊緣模糊和斷裂現(xiàn)象，取得了較好的效果。

4 結(jié) 論

從以上的實(shí)驗(yàn)結(jié)果和分析可以看出：針對腦外科CT圖像的紋理復(fù)雜、細(xì)節(jié)豐富等特點(diǎn)，本文提出的用平均值濾波器對圖像增強(qiáng)、Canny邊緣檢測算子與數(shù)學(xué)形態(tài)學(xué)相結(jié)合的綜合邊緣提取算法，不僅具有簡單快捷的處理能力，而且具有很好的邊緣提取、較強(qiáng)的抗噪性能，提取出的邊緣清晰連續(xù)，能夠比較準(zhǔn)確的提取腦外科CT圖像的邊緣，有效地保留了圖像中的重要信息，為腦外科CT圖像三維重建和虛擬手術(shù)的研究等打下了堅實(shí)的基礎(chǔ)。

[1]M_engen S B， Marjunath B S， KenneyC.Image segmentation using multi-region stability and edge strength[C]//IEEE International Conference on Image Processing.Barcelona：Spain，2003：429-432.

[2]黃金國，戴志明，周培源.一種改進(jìn)的基于梯度的圖像邊緣檢測算法[J].武漢科技學(xué)院學(xué)報，2010，23（3）：33-35.

HUANG Jin-guo， DAI Zhi-ming， ZHOU Pei-yuan.An image edge detection algorithm based on gradient[J].Journal of Wuhan University of Science and Engineering，2010，23（3）：33-35.

[3]蘇恒陽，袁先珍.一種改進(jìn)的Canny的圖像邊緣檢測算法[J].計算機(jī)仿真，2010，27（10）：242-245.

SU Heng-yang， YUAN Xian-zhen.An improved image edge detectionalgorithm ofCannyoperator[J].Computersimulation，2010，27（10）：242-245.

[4]Julez B.A method of coding TV signals based on edge detection[J].BellSystem Tech，1959，38（4）:1001-1020.

[5]Roberts L D.Machine perception of three-dimension solids in optica and electro-optimal information processing[C]//Massachusetts Institute of Technology Press,1966:157-197.

[6]Huang C P，Wang R Z.An intergrated edge detection method using mathematical morphology[J].Pattern Recgnition and Image Analysis，2006，16（3）：406-412.

[7]孫兆林.MATLAB 6.x圖像處理[M].北京：清華大學(xué)出版社，2002.

[8]張斌,賀賽先.基于Canny算子的邊緣提取改善方法[J].紅外技術(shù)，2006，28（3）：165-169.

ZHANG Bin,HE Sai-xian.An improved method of edge detection based on Canny operator[J].Journal of Infrared technology，2006，28（3）：165-169.

[9]王家文，曹宇.MATLAB 6.5圖形圖像處理[M].北京：國防工業(yè)出版社，2004.

[10]張小琴,楊闊,郝葆青.基于Canny算子與數(shù)學(xué)形態(tài)學(xué)的細(xì)胞邊緣提取[J].生物醫(yī)學(xué)研究，2009，28（4）：275-279.

ZHANG Xiao-qin,YANG Kuo,HAO Bao-qing.Cell edge detection based on Canny operatorand Mathematical morphology[J].Journal of Biomedical Research，2009，28（4）：275-279.