張榮國(guó)，劉小君，黨偉超，劉焜

(1.太原科技大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，山西太原 030024;2.合肥工業(yè)大學(xué)機(jī)械與汽車工程學(xué)院，安徽合肥 230009)

多目標(biāo)輪廓Mumford－Shah水平集提取

張榮國(guó)1，2，劉小君2，黨偉超1，劉焜2

(1.太原科技大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，山西太原 030024;2.合肥工業(yè)大學(xué)機(jī)械與汽車工程學(xué)院，安徽合肥 230009)

目標(biāo)輪廓的快速檢測(cè)進(jìn)而提取其幾何形狀，在圖形圖像處理中有著重要的作用.提出了一種多目標(biāo)輪廓的水平集提取方法，對(duì)基于Mumford－Shah模型的C－V方法從兩方面進(jìn)行了改進(jìn):增加梯度矢量場(chǎng)和曲線法方向的融合作為邊界吸引場(chǎng)，生成可以驅(qū)動(dòng)主動(dòng)輪廓向邊緣進(jìn)化的雙向幾何變形流，保留原圖像分布信息作為區(qū)域進(jìn)化能，解決未考慮局部幾何信息造成的區(qū)域能量捕捉信息不全，或邊緣梯度場(chǎng)和演化曲線法線方向正交時(shí)無(wú)法實(shí)現(xiàn)拓?fù)浣Y(jié)構(gòu)變化的缺陷;對(duì)水平集函數(shù)進(jìn)行修正，使得它在收斂過(guò)程中能自動(dòng)進(jìn)行調(diào)整，確保其滿足符號(hào)距離函數(shù)的要求，擴(kuò)大初始化前迭代搜索區(qū)域，減少初始化次數(shù)，提高收斂效率;最后給出所提方法的數(shù)字化求解方案.實(shí)驗(yàn)表明該方法可行且具有較好的魯棒性.

Mumford－Shah模型;水平集方法;多目標(biāo)輪廓;能量方程

主動(dòng)輪廓模型使用能量最小化方法對(duì)數(shù)字圖像中的目標(biāo)對(duì)象進(jìn)行處理，為圖像分析和計(jì)算機(jī)視覺(jué)中各種問(wèn)題的解決提供了一種高效便捷的方案［1］.該模型既可以用于圖像分割和理解，也可以用于動(dòng)態(tài)圖像數(shù)據(jù)分析或3－D圖像數(shù)據(jù)的重建，還適用于3－D圖形建模以及任意形狀變形輪廓的提取，因而在計(jì)算機(jī)視覺(jué)和圖形圖像處理中得到廣泛的應(yīng)用.目前，對(duì)主動(dòng)輪廓模型的研究主要從2個(gè)方面進(jìn)行［2－3］:參數(shù)主動(dòng)輪廓模型和幾何主動(dòng)輪廓模型.參數(shù)主動(dòng)輪廓模型適用于圖像中單目標(biāo)對(duì)象輪廓的提?。?－5］，它計(jì)算簡(jiǎn)單、速度快，輪廓曲線動(dòng)態(tài)演化過(guò)程中拓?fù)浣Y(jié)構(gòu)不能發(fā)生變化，無(wú)法直接對(duì)多目標(biāo)輪廓邊緣進(jìn)行提?。?］.

幾何主動(dòng)輪廓模型，利用曲線進(jìn)化理論和水平集方法，將低維平面上的曲線演化問(wèn)題轉(zhuǎn)化為高維空間的曲面演化問(wèn)題，由對(duì)參數(shù)方程的直接求解轉(zhuǎn)化為水平集函數(shù)的一個(gè)零水平集隱含方式求解，通過(guò)控制曲面演化的偏微分方程進(jìn)行數(shù)值計(jì)算，最終得到所希望輪廓［7］.幾何主動(dòng)輪廓在曲線演化過(guò)程中，能自適應(yīng)所搜索對(duì)象的幾何拓?fù)浣Y(jié)構(gòu)變化，對(duì)多目標(biāo)圖形對(duì)象的輪廓提取有著較好的效果.傳統(tǒng)的幾何主動(dòng)輪廓模型，采用圖像中像素分布的均值曲率作為演化曲線運(yùn)動(dòng)的幾何變形流，通過(guò)水平集的數(shù)值計(jì)算，在對(duì)比度好的圖像中，能獲得較為滿意的目標(biāo)輪廓提取效果.隨后Caselles等人對(duì)傳統(tǒng)方法進(jìn)行了改進(jìn)，提出了測(cè)地線主動(dòng)輪廓模型［8］，他們通過(guò)使用較大尺度的濾波器，對(duì)原圖像進(jìn)行濾波，并使用較大的膨脹力，使得演化曲線在沿測(cè)地線的運(yùn)動(dòng)過(guò)程中，能越過(guò)這些局部干擾，對(duì)具有部分不連續(xù)的邊緣檢測(cè)具有一定效果.其后，Kuhne和Paragio等人對(duì)測(cè)地線模型進(jìn)行了一些改進(jìn)和融合［9－11］，擴(kuò)大了模型的適用范圍，改善了模型的檢測(cè)效果.這些模型都是基于梯度信息的圖像邊緣檢測(cè)，對(duì)有裂口、縫隙等的圖像具有較好的分割結(jié)果.若圖像中含有大量噪聲，或者所提取目標(biāo)對(duì)象邊緣模糊，這些方法就很難奏效，甚至達(dá)不到所要提取的目標(biāo)要求.為此，Chan和 Vese等人提出了基于 Mumford－Shah模型的圖像分割方法［12］，利用勻質(zhì)區(qū)域的圖像分布信息，通過(guò)求解能量函數(shù)的極小值來(lái)得到感興趣目標(biāo)的檢測(cè)，該方法將圖像去噪、圖像目標(biāo)檢測(cè)與目標(biāo)對(duì)象圖形形狀重建這些問(wèn)題用一個(gè)廣義能量泛函的最小值問(wèn)題描述出來(lái)，曲線的演化和目標(biāo)的搜索不依賴于圖像的邊緣信息，適用于非常弱的邊緣圖像和有噪圖像的邊緣形狀提取.Gao和龔等人通過(guò)對(duì)圖像中的目標(biāo)對(duì)象像素分布情況進(jìn)行分析，反復(fù)多次地利用單一的水平集函數(shù)，對(duì)感興趣目標(biāo)對(duì)象進(jìn)行分級(jí)別、分層次的圖像劃分［13－14］.為了加快水平集的收斂速度，李等人采用源點(diǎn)掃描方法來(lái)快速計(jì)算符號(hào)距離函數(shù)［15］，而Li等人則提出了無(wú)需重新初始化的變分方法［16］，隨后楊等人利用小波變換對(duì)圖像進(jìn)行多分辨率分解，在對(duì)各子區(qū)域劃分的基礎(chǔ)上進(jìn)行水平集演化，實(shí)現(xiàn)多目標(biāo)輪廓的提?。?7］.Lie和Tai等人把圖像處理為不連續(xù)的二值分布狀態(tài)，用分段常數(shù)的 Mumford－Shah 模型對(duì)圖像進(jìn)行分割［18－19］，葉等人則是通過(guò)最小生成樹(shù)的方式，先對(duì)區(qū)域進(jìn)行劃分，然后再對(duì)圖像進(jìn)行分割［20］.這些改進(jìn)主要針對(duì)圖像全局特征的不足，沒(méi)有考慮原能量函數(shù)中不含有局部特征信息，何等人把測(cè)地線主動(dòng)輪廓模型中的邊緣檢測(cè)函數(shù)引入到水平集中完成對(duì)圖像的分割［21］，但梯度信息只存在于邊緣附近，無(wú)法利用圖像的局部信息對(duì)演化過(guò)程進(jìn)行調(diào)控.為克服這些缺點(diǎn)，本文提出了一種新的能量模型，它在Mumford－Shah模型的基礎(chǔ)上，增加了梯度矢量場(chǎng)和曲線法方向融合的雙向幾何變形流，既考慮了目標(biāo)對(duì)象局部幾何信息，也使用了目標(biāo)所圍區(qū)域內(nèi)外灰度圖像的整體分布信息，充分利用數(shù)字圖像所提供的完整信息;同時(shí)，通過(guò)調(diào)節(jié)項(xiàng)在迭代過(guò)程中的自適應(yīng)變化，減少了求解過(guò)程中為滿足符號(hào)距離函數(shù)要求對(duì)水平集函數(shù)重新初始化的次數(shù)，提高了收斂效率.

1 Mumford－Shah 模型的 Level Set方法

1.1 Level Set方法

給定一個(gè)區(qū)域Ω，設(shè)它為R2上的一個(gè)有界開(kāi)區(qū)間，C為所考慮的封閉演化曲線，它由Ω內(nèi)光滑曲線的有限集構(gòu)成，C將Ω分成i個(gè)子區(qū)域.若二維圖像定義為u0:Ω→k，那么區(qū)域Ω就是平面圖像上像素點(diǎn)信息所構(gòu)成的二維長(zhǎng)方形網(wǎng)格，C可看作為所檢測(cè)目標(biāo)在其上演化的運(yùn)動(dòng)曲線.Level Set方法把隨時(shí)間運(yùn)動(dòng)的曲線C看作是某個(gè)函數(shù)φ(x，y，t)的零等值面，要得到某個(gè)時(shí)刻t處的運(yùn)動(dòng)曲線，只需要通過(guò)函數(shù)φ求出該時(shí)刻的零等值面位置即可，也就是:

式中:signdist表示(x，y)到C(x，y，0)的符號(hào)距離，在運(yùn)動(dòng)曲線C封閉區(qū)域內(nèi)取正值，在區(qū)域外取負(fù)值，在曲線上取零值.在任意時(shí)刻t，要使水平集函數(shù)φ的零等值面就是主動(dòng)輪廓演化曲線C，φ的控制方程要滿足以下的偏微分方程:

式中:函數(shù) φ(x，y，t)∈Ω 是水平集函數(shù)，F(xiàn)表示零等值面曲線上點(diǎn)沿法線方向的演化速度，它通常和運(yùn)動(dòng)曲線的位置、輪廓曲線的幾何形狀、當(dāng)前運(yùn)動(dòng)時(shí)刻、以及其他一些外部物理特性相關(guān);緊跟F的|▽?duì)諀項(xiàng)表示水平集函數(shù)的梯度范數(shù).

水平集函數(shù)φ的單位外法向以及演化曲線在法向上的速度滿足:

1.2 Mumford－Shah 模型

設(shè)圖像u0(x，y)為有界開(kāi)區(qū)間Ω上含有多個(gè)目標(biāo)的含噪圖像，圖像中目標(biāo)輪廓C可以用Ω上不連續(xù)閉集，尋找圖像中感興趣的目標(biāo)對(duì)象，就必須對(duì)區(qū)域Ω用子區(qū)域Ωi進(jìn)行有效地分割，使Ωi成為Ω的分片光滑區(qū)，并最終成為圖像u0中的目標(biāo)圖像u，也就是要尋找原始圖像u0的分片光滑區(qū)域圖像u，目標(biāo)輪廓的邊緣為C.求解這樣的問(wèn)題，可以采用以下的Mumford－Shah能量最小化模型的式子來(lái)表示:

如果在每一個(gè)分片區(qū)域Ωi中所得到的逼近圖像u的分布是一個(gè)常量ci，該常量ci一般取區(qū)域Ωi內(nèi)圖像分布的均值，那么上面的式子就可以進(jìn)行簡(jiǎn)化.特別地，當(dāng)所研究的圖像u0被封閉輪廓曲線C分割成2個(gè)勻質(zhì)區(qū)域Ω1和Ω2，Ω1為感興趣圖形目標(biāo)對(duì)象，位于封閉輪廓內(nèi)部，其像素分布的均值為c1;Ω2為圖形對(duì)象所處的背景區(qū)域，位于封閉輪廓外部，其像素分布的均值為c2;輪廓C所圍曲線的長(zhǎng)度和面積用Heaviside函數(shù)和Dirac函數(shù)來(lái)表示，能量公式可以重新寫為如下:

要使得式(1)得到最小值，它需要滿足Euler－Lagrange方程:

Mumford－Shah方法是一種變差圖像分割方法，它將圖像按像素分布的灰度相似性來(lái)劃分，劃分φ的粗細(xì)通過(guò)演化曲線中長(zhǎng)度項(xiàng)權(quán)值來(lái)進(jìn)行控制.當(dāng)所檢測(cè)對(duì)象目標(biāo)灰度與背景灰度呈階躍特性，灰度差異明顯，此時(shí)按灰度相似性進(jìn)行目標(biāo)區(qū)分，效果較好.若灰度差異不大，且各自內(nèi)部灰度分布不均勻，用此方法就不能很好地對(duì)目標(biāo)輪廓加以檢測(cè)，為此，提出以下的改進(jìn)方法.

2 梯度場(chǎng)的擴(kuò)散及幾何流的生成

2.1 梯度場(chǎng)的擴(kuò)散

要使主動(dòng)輪廓曲線的運(yùn)動(dòng)隨著水平集函數(shù)φ的變化而演化，就需要加一個(gè)和φ相關(guān)的擴(kuò)散速度場(chǎng)，由于原始方程中F只在零水平集附近有意義，因此需要將和零水平集相關(guān)的幾何流擴(kuò)散到整個(gè)函數(shù)曲面.傳統(tǒng)的水平集幾何主動(dòng)輪廓是使得演化曲線沿法線方向以曲率K為演化速度向目標(biāo)輪廓收斂，但演化曲線的法線方向不一定就是輪廓曲線的最速下降方向.梯度矢量場(chǎng)可把邊緣信息在整個(gè)零水平集平面內(nèi)進(jìn)行擴(kuò)散［17］，以獲取輪廓曲線進(jìn)行演化時(shí)所需要的幾何矢量場(chǎng)，該矢量場(chǎng)中的矢量在封閉輪廓內(nèi)部時(shí)由里向外指向邊緣，當(dāng)在封閉輪廓外部時(shí)由外向內(nèi)指向邊緣，即無(wú)論當(dāng)前進(jìn)化點(diǎn)處于什么位置，它始終指向目標(biāo)輪廓邊緣，并且該矢量場(chǎng)是輪廓曲線的最速下降方向，演化曲線沿此方向進(jìn)化，可以到達(dá)所期望的輪廓邊緣.為了充分利用這一信息，把擴(kuò)散后的梯度矢量場(chǎng)和法線方向通過(guò)一個(gè)內(nèi)積來(lái)進(jìn)行融合，得到式(2).

式中:(f1，f2)(x，y)為圖像輪廓邊緣的梯度擴(kuò)散場(chǎng)，N(x，y)為其法線方向，可以通過(guò)極小化以下的能量函數(shù)獲得:

從式(2)看到，當(dāng)梯度矢量場(chǎng)和外法線方向一致時(shí)，對(duì)曲線進(jìn)行膨脹，當(dāng)梯度矢量場(chǎng)和外法線方向相反時(shí)，使演化曲線收縮，當(dāng)梯度矢量場(chǎng)和外法線方向正交時(shí)，對(duì)演化曲線的運(yùn)動(dòng)不產(chǎn)生影響，此時(shí)可以通過(guò)區(qū)域力使得演化曲線向目標(biāo)邊緣輪廓靠近.該擴(kuò)散場(chǎng)可以從主動(dòng)輪廓的任一邊驅(qū)使輪廓向目標(biāo)邊緣靠近，當(dāng)主動(dòng)輪廓曲線接近目標(biāo)輪廓邊緣時(shí)，該梯度矢量場(chǎng)不再對(duì)主動(dòng)輪廓施加任何影響，從而使得主動(dòng)輪廓停留在邊緣輪廓上.可以將此影響稱為邊緣吸引項(xiàng)，其演化方程可以定義如式(3).

2.2 調(diào)節(jié)項(xiàng)的建立

為了運(yùn)算方便，在主動(dòng)輪廓曲線的演化過(guò)程中，保持水平集函數(shù)在任意時(shí)刻t滿足符號(hào)距離函數(shù)是非常重要的，對(duì)φ初始值的設(shè)定是滿足這一要求的，但在微分方程數(shù)值求解過(guò)程中，經(jīng)過(guò)幾個(gè)時(shí)間步長(zhǎng)的迭代以后，φ(x，y，t)將不再滿足所定義的符號(hào)距離式子，為了保證φ盡量保持這一良好特性，需要對(duì)φ重新初始化，也就是要重新構(gòu)造φ函數(shù)，使得它既要滿足符號(hào)距離函數(shù)，又要與φ0具有相同的零等值面.通過(guò)求解初值問(wèn)題的穩(wěn)定解，就可以得到滿足要求的函數(shù)，但這一過(guò)程計(jì)算量大，計(jì)算相當(dāng)耗時(shí)，因此希望能尋找一種辦法能減少初始化函數(shù)的調(diào)用次數(shù).由于在0水平集處|▽?duì)諀=1，可以類似于文獻(xiàn)［17］，在能量函數(shù)中引入一調(diào)節(jié)項(xiàng):

來(lái)衡量水平集函數(shù)偏離符號(hào)距離函數(shù)的程度.當(dāng)|▽?duì)諀接近于1時(shí)，該式的值很小，不會(huì)對(duì)能量函數(shù)的極值產(chǎn)生較大影響;當(dāng)|▽?duì)諀偏離距離函數(shù)時(shí)，該項(xiàng)的值會(huì)變得越來(lái)越大，對(duì)能量函數(shù)極小化過(guò)程會(huì)產(chǎn)生較大影響，通過(guò)對(duì)主動(dòng)輪廓曲線的演化，使得它接近于0水平集，從而使該項(xiàng)的值變小.由此可以看到，通過(guò)加入此調(diào)節(jié)項(xiàng)，可以增加微分方程PDF數(shù)值求解過(guò)程中的迭代次數(shù)，保持水平集函數(shù)和符號(hào)距離函數(shù)盡量接近，以使演化曲線盡量靠近目標(biāo)輪廓，這樣可大大減少重新初始化的次數(shù)，提高能量方程的求解效率.要使得上式得到最小值，它仍要滿足Euler－Lagrange 方程:

從式(5)可以看到，當(dāng)|▽?duì)諀大于1時(shí)，擴(kuò)散系數(shù)大于0，產(chǎn)生正向擴(kuò)散，降低|▽?duì)諀的取值可以達(dá)到對(duì)圖像進(jìn)行平滑的作用;當(dāng)|▽?duì)諀小于1時(shí)，擴(kuò)散系數(shù)小于0，產(chǎn)生逆向擴(kuò)散，降低|▽?duì)諀的取值可以對(duì)圖像進(jìn)行銳化，因此，選用該項(xiàng)作為調(diào)節(jié)項(xiàng)，能保證其滿足Hamilton－Jacobi偏微分方程對(duì)函數(shù)進(jìn)行演化的要求.

至此，曲線演化的幾何流就可以看作為圖像分布的均值情況、圖像邊緣梯度場(chǎng)的擴(kuò)散情況和調(diào)節(jié)項(xiàng)共同作用下的變形流，可以定義和能量函數(shù)對(duì)應(yīng)的新的演化函數(shù)如下:

其中的3個(gè)數(shù)據(jù)項(xiàng)分別是區(qū)域項(xiàng)、擴(kuò)散項(xiàng)和調(diào)節(jié)項(xiàng)對(duì)應(yīng)的演化方程.

3 數(shù)值計(jì)算和實(shí)現(xiàn)

為了便于計(jì)算，需要對(duì)前面所用到的公式進(jìn)行離散化.首先來(lái)看重新初始化方程的離散化求解.前面已經(jīng)提到，對(duì)輪廓曲線進(jìn)行演化迭代后，需要對(duì)φ函數(shù)重新進(jìn)行初始化，以使其滿足符號(hào)距離函數(shù)的要求.初始化過(guò)程可以通過(guò)求解初值問(wèn)題的式(4)來(lái)進(jìn)行，從式子中可以看到，符號(hào)距離函數(shù)的初始值是通過(guò)t時(shí)刻對(duì)演化曲線的求解結(jié)果，其離散化過(guò)程如下:

式中:

式中:a1、a2、a3、a4的值可以通過(guò)式(6)進(jìn)行計(jì)算.

通常的Dirac函數(shù)δ(x)定義的范圍過(guò)于狹窄，使得用它對(duì)圖像輪廓邊緣進(jìn)行檢測(cè)時(shí)受到一些限制.為了使得在圖像定義域范圍內(nèi)，Dirac函數(shù)能保證在每一個(gè)演化曲線點(diǎn)值都是趨于零的正值，需要對(duì)該函數(shù)進(jìn)行正則化，參考C－V方法［12］，采用了以下的正則化式子，對(duì)Heaviside函數(shù)和Dirac函數(shù)進(jìn)行正則化.

此正則式保證了用前述的能量方程檢測(cè)出帶內(nèi)部空洞的圖像邊緣，但若內(nèi)部空洞相對(duì)于檢測(cè)區(qū)域尺寸太小，則此函數(shù)對(duì)原主動(dòng)輪廓的邊緣檢測(cè)產(chǎn)生抑制作用，導(dǎo)致不能穩(wěn)定地檢測(cè)出內(nèi)部區(qū)域邊緣，需要對(duì)此函數(shù)進(jìn)行適當(dāng)?shù)恼{(diào)整.

區(qū)域相似項(xiàng)也需要滿足Euler－Lagrange方程，如果設(shè)離散間隔步長(zhǎng)為h=△x=△y，時(shí)間步長(zhǎng)為△t，離散點(diǎn)(xi，yi)的值為(ih，jh)，使用半隱含的有限差分方法，由此可以得到如下的離散化演化方程:

式中:系數(shù)b1、b2、b3、b4可以用以下式子計(jì)算:

系數(shù)m和C可以用以式(7)計(jì)算.

分片常數(shù)c1、c2可以用式(8)計(jì)算.

梯度擴(kuò)散場(chǎng)的離散化方法，在先前所做的工作中已有介紹，詳細(xì)過(guò)程可查閱參考文獻(xiàn)［22］，此處不再贅述.

4 實(shí)驗(yàn)結(jié)果和分析

為了進(jìn)一步驗(yàn)證本文所敘述方法，在P4 2.99 GHz、內(nèi)存為512 MB的個(gè)人計(jì)算機(jī)上進(jìn)行了實(shí)驗(yàn)，所用的開(kāi)發(fā)工具為Matlab 7.0.首先對(duì)掃描電鏡下所獲得的磨粒圖片進(jìn)行幾何輪廓提取，實(shí)驗(yàn)結(jié)果如圖1所示.從圖1(a)的原始圖片中可以看到，圖中有一大兩小3個(gè)磨粒組成，且輪廓邊緣呈現(xiàn)不規(guī)則的復(fù)雜形狀，磨粒的影像特征和背景對(duì)比度并不十分強(qiáng).取時(shí)間步長(zhǎng)△t=0.1，離散間隔步長(zhǎng)h=1，控制長(zhǎng)度項(xiàng)的系數(shù)μ實(shí)際上是一個(gè)尺度參數(shù)，在進(jìn)化過(guò)程中檢測(cè)物體的大小，μ值大，可以檢測(cè)到較大物體，μ值小，則可以檢測(cè)任意小的物體，此處取μ=0.001 5×255×255.面積控制項(xiàng)v用來(lái)控制演化曲線所圍區(qū)域面積對(duì)能量方程的影響，此處取v=1.λ1用來(lái)控制演化曲線所圍區(qū)域內(nèi)部圖像分布信息對(duì)能量控制方程的影響，λ2用來(lái)控制演化曲線所圍區(qū)域外部圖像分布信息對(duì)能量控制方程的影響，此處取λ1=λ2=1，所得結(jié)果如圖1所示.圖1(b)為中間演化過(guò)程，圖1(c)為演化結(jié)果圖，圖1(d)為所提取的邊緣輪廓圖.

圖2是一個(gè)人工合成圖像的輪廓邊緣提取，它由4個(gè)不同形狀的目標(biāo)對(duì)象組成，圖2(b)和圖2(d)分別為C－V方法［12］和本文方法所提取的邊緣輪廓圖，從結(jié)果圖中可以看到，由于圖2中四邊形的圖像像素統(tǒng)計(jì)分布和背景圖的區(qū)別不大，且多個(gè)目標(biāo)之間的梯度矢量場(chǎng)和進(jìn)化曲線的法線方向有正交的情況出現(xiàn)，因此用C－V方法未能檢測(cè)出其邊緣輪廓，而本文方法則成功地將4個(gè)目標(biāo)邊緣合理地提取出來(lái).

圖1 磨粒圖輪廓提取過(guò)程Fig.1 Wearing particulate boundaries extraction

圖2 人工合成圖輪廓提取過(guò)程Fig.2 Artificial synthetic image edges detection

圖3是對(duì)醫(yī)學(xué)CT圖片中人體胸腔的輪廓邊緣提取，圖3(a)和圖3(b)為用Li的方法［16］的中間演化過(guò)程和輪廓提取結(jié)果，從中可以看到這種方法只能提取其外部輪廓，而內(nèi)部輪廓?jiǎng)t無(wú)法提取;圖3(c)和圖3(d)為本文方法的中間演化過(guò)程和輪廓提取結(jié)果，除了對(duì)2個(gè)內(nèi)部目標(biāo)的外邊緣輪廓成功提取外，還對(duì)胸腔的內(nèi)部邊緣輪廓也得以成功提取，顯示了本文方法不僅和先驗(yàn)的初始輪廓相關(guān)性很小，而且在演化過(guò)程中能自適應(yīng)地進(jìn)行幾何結(jié)構(gòu)的拓?fù)渥兓瓤梢赃M(jìn)行邊緣輪廓的自動(dòng)分裂，也可以有效地進(jìn)行多邊緣輪廓的自動(dòng)合并.從圖中可以看到，盡管3個(gè)圖像目標(biāo)圖形的邊緣都是有凹有凸的不規(guī)則形狀，但用本文方法卻能很好地提取它們的形狀輪廓.

圖3 CT胸腔圖輪廓提取過(guò)程Fig.3 Thorax CT image contours extraction

用 C－V 方法［12］、Li方法［16］和本文所提的改進(jìn)方法對(duì)圖1、圖2和圖3內(nèi)容作了對(duì)比實(shí)驗(yàn)，詳細(xì)參數(shù)及運(yùn)行結(jié)果如表1所示.表1中最小誤差Dmin、最大誤差Dmax、均方誤差Eerr指的是在多目標(biāo)輪廓提取過(guò)程中，演化曲線相鄰2次迭代輪廓線之間的最小誤差、最大誤差和均方誤差.從表中可以看出，在得到相近收斂效果的情況下，本文方法較其他2種方法迭代次數(shù)更少、收斂速度更快.實(shí)驗(yàn)表明，圖像的區(qū)域分布統(tǒng)計(jì)信息和圖像邊緣的幾何特征信息，對(duì)多目標(biāo)輪廓的提取有著不同程度的影響，合理選擇相關(guān)的控制項(xiàng)，能加快收斂速度;同時(shí)，有效地減少對(duì)水平集的重新初始化次數(shù)，仍是加快輪廓演化過(guò)程的一個(gè)有效途徑.

表1 不同方法性能比較Table 1 Performance comparison of various methods

5 結(jié)束語(yǔ)

基于Mumford－Shah的水平集方法在輪廓邊緣提取中有著廣泛的應(yīng)用，本文提出的快速提取方法，不僅利用了原圖像在整個(gè)區(qū)域的像素分布信息，還將圖像的邊緣梯度場(chǎng)擴(kuò)散后結(jié)合曲線法線方向信息來(lái)作為曲線演化的幾何變形流，克服了片面強(qiáng)調(diào)區(qū)域分布信息或邊緣梯度信息的缺陷.同時(shí)，由于引入對(duì)水平集函數(shù)的修正，使得其盡量滿足符號(hào)距離函數(shù)要求，曲線在演化過(guò)程中基本保持在0水平集附近，減少了對(duì)水平集函數(shù)重新初始化的次數(shù)，從而提高了收斂效率.實(shí)驗(yàn)結(jié)果也表明，本文所述方法能有效地對(duì)多目標(biāo)輪廓邊緣的幾何形狀進(jìn)行提取，但由于要對(duì)原始圖像的梯度矢量場(chǎng)進(jìn)行擴(kuò)散，需要耗費(fèi)一定的時(shí)間，特別是灰度分布信息豐富的圖像，預(yù)處理需要花費(fèi)較多的時(shí)間，這是后期將要研究解決的一個(gè)問(wèn)題.

［1］KASS M，WITKIN A，TERZOPOULS D.Snake:active contour models［J］.International Journal of Computer Vision，1987，1(4):321－331.

［2］XU C，PRINCE J L.Snake，shapes，and gradient vector flow［J］.IEEE Trans on Image Processing，1998，7:359－369.

［3］NING Jifang，WU Chengke，LIU Shigang，et al.NGVF:an improved external force field for active contour model［J］.Pattern Recognition Letters，2007，28(1):58－63.

［4］SAKALLI M，LAM K M，YAN H.A faster converging snake algorithm to locate object boundaries［J］.IEEE Transactions on Image Processing，2006，15(5):1182－1191.

［5］SUM K W，CHEUNG P Y S.Boundary vector for parametric active contours［J］.Pattern Recognition，2007，40(6):1635－1645.

［6］XIE X H，MIRMEHDI M.MAC:magnetostatic active contour model［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2008，30(4):632－646.

［7］OSHER S，SETHIAN J A.Fronts propagating with curvature－dependent speed:algorithms based on Hamilton－Jacobi formulation［J］.Journal of Computational Physics，1988，79(1):12－49.

［8］CASSELLES V，KIMMEL R，SAPIRO G.Geodesic active contours［J］.International Journal of Computer Vision，1997，22:61－79.

［9］OSHER S，PARAGIOS N.Geometric level set methods in imaging，vision，and graphics［M］.Berlin:Springer Verlag，2003:43－57.

［10］PARAGIOS N，DERICH R.Geodesic active regions for supervised texture segmentation［C］//IEEE International Conference on Computer Vision.Kerkyra，Greece，1999，2:926－932.

［11］PARAGIOS N，MELLINA G O，RAMESH V.Gradient vector flow fast geometric active contours［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2004，26(3):402－407.

［12］CHAN T F，VESE L A.Active contours without edges［J］.IEEE Transactions on Image Processing，2001，10(2):266－277.

［13］龔永義，羅笑南，黃輝，等.基于單水平集的多目標(biāo)輪廓提取［J］.計(jì)算機(jī)學(xué)報(bào)，2007，30(1):120－128.

GONG Yongyi，LUO Xiaonan，HUANG Hui，et al.Multiobjects extracted based on single level set［J］.Chinese Journal of Computers，2007，30(1):120－128.

［14］GAO S，TIEN D.Image segmentation and selective smoothing by using Mumford－Shah mode［J］.IEEE Transact－ions on Image Processing，2005，14(10):1537－1549.

［15］李俊，楊新，施鵬飛.基于Mumford－Shah模型的快速水平集圖像分割方法［J］.計(jì)算機(jī)學(xué)報(bào)，2002，25(11):1175－1183.

LI Jun，YANG Xin，SHI Pengfei.A fast level set approach to image segmentation based on Mumford－Shah model［J］.Chinese Journal of Computers，2002，25(11):1175－1183.

［16］LI Chunming，XU Chenyang，GUI Changfeng，et al.Level set evolution without re－initialization:a new variational formulation［C］//IEEE International Conference on Computer Vision and Pattern Recognition.San Diego，USA，2005，1:430－436.

［17］楊莉，楊新.基于區(qū)域劃分的曲線演化多目標(biāo)分割［J］.計(jì)算機(jī)學(xué)報(bào)，2004，27(3):420－425.

YANG Li，YANG Xin.Multi－object segmentation based on curve evolving and region division［J］.Chinese Journal of Computers，2004，27(3):420－425.

［18］LIE J，LYSAKER M，TAI X C.A binary level set model and some applications for Mumford－Shah image segmentation［J］.IEEE Transactions on Image Processing，2006，15(5):1171－1181.

［19］TAI Xuecheng，LI Hongwei.A piecewise constant level set methods for elliptic inverse problems［J］.Applied Numerical Mthematics，2007，57(5/6/7):686－696.

［20］葉偉，王遠(yuǎn)軍.基于Mumford－Shah理論的最小生成樹(shù)圖像分割方法［J］.計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)，2009，21(8):1127－1133.

YE Wei，WANG Yuanjun.MST image segmentation based on Mumford－shah theory method based on boundary and region information［J］.Journal of Computer—Aided Design＆ Computer Graphics，2009，21(8):1127－1133.

［21］何寧，張朋.基于邊緣和區(qū)域信息相結(jié)合的變分水平集圖像分割方法［J］.電子學(xué)報(bào)，2009，37(10):2215－2219.

HE Ning，ZHANG Peng.Variational level set image segmentation method based on boundary and region information［J］.Acta Electronic Sinica，2009，37(10):2215－2219.

［22］張榮國(guó)，劉小君，王蓉，劉焜.自適應(yīng)梯度矢量流輪廓提取方法研究［J］.模式識(shí)別與人工智能，2008，21(6):799－805.

ZHANG Rongguo，LIU Xiaojun，WANG Rong，LIU Kun.Adaptive gradient vector flow algorithm for boundary extraction［J］.Pattern Recognition and Artificial Intelligence，2008，21(6):799－805.

張榮國(guó)，男，1964年生，教授，博士，主要研究方向?yàn)閳D形圖像處理、CAD/CG和計(jì)算機(jī)支持的協(xié)同設(shè)計(jì)等.

劉小君，女，1965年生，教授，博士，主要研究方向?yàn)閿?shù)字化設(shè)計(jì)和圖像處理.

黨偉超，男，1974年生，副教授，主要研究方向?yàn)閳D像處理與信息系統(tǒng).

Mumford-Shah level set method for multi-objective contour extraction

ZHANG Rongguo1，2，LIU Xiaojun2，DANG Weichao1，LIU Kun2
(1.School of Computer Science and Technology，Taiyuan University of Science＆ Technology，Taiyuan 030024，China;2.School of Mechanical and Automotive Engineering，Hefei University of Technology，Hefei 230009，China)

Fast detection of objective contours and extraction of its geometric shape have important roles in graphics and image processing.Based on the Mumford－Shah model，a novel level set method for multi－objective contour extraction was presented.First，the gradient vector field was combined with normal direction of the curves as boundary abstracted fields，so as to generate a bi－directional geometric deformable flow field which can drive active contours evolving towards the boundary from inside or outside edges.Furthermore，the distributed information of the image would be left as area evolution energy.This method can solve problems that arise when area energy information is lost because local geometric information isn＇t considered，or when topological structure should not be changed because the gradient vector field is orthogonal with normal direction.Then the level set function was modified so that it could change adaptively in curve convergence.Other reasons for this modification were to make sure that the level set changes could maintain signal distance function，the search area could be covered sufficiently before reinitialization，and the iterative number could be decreased.The convergence efficiency was also raised.Finally，a numerical solving scheme was given.Experimental results illustrate that the method proposed in this paper is feasible and robust.

Mumford－Shah model;level set method;multi－objective contours;energy equation

TP391

A

1673－4785(2011)04－0360－07

10.3969/j.issn.1673－4785.2011.04.014

2010－07－15.

國(guó)家自然科學(xué)基金資助項(xiàng)目(51075113).

張榮國(guó).E－mail:rg_zh@163.com.