劉 兵，李 聰，向磊磊，孫玉秋 （長江大學(xué)信息與數(shù)學(xué)學(xué)院，湖北荊州434023）

Norden Huang于1988年提出一種新的用于分析非線性和非平穩(wěn)的信號處理方法——經(jīng)驗?zāi)Ｊ椒纸?（Empirical mode decomposition，EMD）［1］，這種方法的關(guān)鍵在于：任何復(fù)雜的數(shù)據(jù)都可以分解為一系列有限且是少量的固有模式函數(shù) （Intrinc mode function，IMF），其中IMF滿足以下2個條件：①在整個數(shù)據(jù)序列中，極值點(diǎn)的數(shù)量與過零點(diǎn)的數(shù)量必須相差不超過1個；②在任一時間點(diǎn)上，信號局部極大值確定的上包絡(luò)線和局部極小值確定的下包絡(luò)線的均值為零，即滿足傳統(tǒng)平穩(wěn)高斯過程。這種分解是自適應(yīng)的，因此具有很高的效率，將這種分解方式應(yīng)用于信號分解時能獲得比較好的IMF，并且對非線性和非平穩(wěn)的信號進(jìn)行處理時，能夠取得小波分析等分析方法難以取得的效果。

EMD方法在處理一維信號時，對極值的選取一般只需與左右像素比較，在擬合時可以得到足夠的擬合點(diǎn)。但是在二維圖像中進(jìn)行選取極值的過程，一般需要比較每個像素的3鄰域，5鄰域或8鄰域，然后在擬合時對每列或每行進(jìn)行擬合，可能會出現(xiàn)某行或某列的極值點(diǎn)低于2個，這樣算法就不具有普遍性。應(yīng)用形態(tài)學(xué)的算法可以很好的解決此類問題，并且所需的存儲空間遠(yuǎn)低于原算法。在信號中，對極值點(diǎn)的擬合一般采用三次樣條擬合，這種擬合方式會出現(xiàn)擬合過沖的情況，且在端點(diǎn)處也有可能會發(fā)生端點(diǎn)飛翼?？墒沁@些變異點(diǎn)在信號處理中不是很多，所以沒有引起研究者的重視。但是在圖像處理時，這些變異點(diǎn)會隨著計算次數(shù)的增加來影響整個圖像［2］。針對這些問題，筆者提出1種EMD的改進(jìn)算法：先運(yùn)用形態(tài)學(xué)算子［3］提取出圖像的極大值和極小值點(diǎn)，對極值點(diǎn)采用Delaunay三角擬合［4-5］算法進(jìn)行曲面擬合，得到上下包絡(luò)面。

1 算法改進(jìn)

1．1 EMD算法

一維EMD算法原理如圖1所示。把二維EMD算法［6］運(yùn)用到圖像處理方面，將一幅圖像分解為若干固有模式函數(shù)和一個余量函數(shù)的集合，和一維情況類似，二維EMD定義如下：

1）把原始數(shù)據(jù)X作為待處理數(shù)據(jù)，確定該圖像的所有局部極值點(diǎn) （包括極大值點(diǎn)和極小值點(diǎn)）；

2）對所有極大值和極小值分別進(jìn)行曲面擬合得到包絡(luò)曲面emax，emin；

4）計算余量h＝X-eave；

5）判斷h是否是一個固有模式分量 （IMF），若不是，將h返回到第一步。

經(jīng)過第5）步的判定，假如條件滿足，就可以得到一個固有模式函數(shù)，將原始圖像減去所求得的IMF1，即得到一個余量數(shù)據(jù)r，將r作為原始數(shù)據(jù)再重復(fù)前面5個步驟，得到IMF2，IMF3，…，IMFn以及一個余量R。余量R滿足預(yù)先設(shè)定的停止準(zhǔn)則后即可停止，最后剩下原始數(shù)據(jù)的余項R。

將這些固有模式函數(shù)和余量相加，就可以實現(xiàn)對信號的重構(gòu)。

圖1 一維EMD算法流程

1．2 利用形態(tài)學(xué)提取極值點(diǎn)

用原始的方式進(jìn)行極值點(diǎn)的提取，對于邊界上的灰度值要單獨(dú)處理，耗費(fèi)大量時間且尋找的極值點(diǎn)在某一行 （列）只有一個或零個，對后續(xù)操作有難以克服的影響。因此，考慮用形態(tài)學(xué)的方法來求極值。

形態(tài)學(xué)有以下4種基本運(yùn)算：膨脹，腐蝕，開運(yùn)算，閉運(yùn)算，用形態(tài)學(xué)求極值主要運(yùn)用的是膨脹和腐蝕運(yùn)算。所謂在灰度級圖像上的膨脹操作就是指輸入圖像在結(jié)構(gòu)元素 （本身可以看作是一個子圖像函數(shù)）下的一種灰度擴(kuò)張。原圖像滑過結(jié)構(gòu)元素時，圖像有部分與結(jié)構(gòu)元素重合，此時對應(yīng)點(diǎn)的灰度值根據(jù)結(jié)構(gòu)元素的取值進(jìn)行取最大值處理。處理過后的原圖像主要取決于結(jié)構(gòu)元素的值和形狀，如果所有結(jié)構(gòu)元素的值為正，則輸出圖像會趨于比輸入圖像更亮并且暗的部分全部減少或者被消除掉。同樣腐蝕操作是指輸入圖像在結(jié)構(gòu)元素下的一種灰度降低。與膨脹運(yùn)算相比，僅在對應(yīng)點(diǎn)灰度值取值時，所采取的運(yùn)算方式不一樣，腐蝕操作是根據(jù)結(jié)構(gòu)元素采取最小值操作。

選取一個3階元素全為1的矩陣作為結(jié)構(gòu)元素，運(yùn)用形態(tài)學(xué)的膨脹 （腐蝕）算子，分別得到的是每一個點(diǎn)在其3鄰域內(nèi)的最大值 （最小值）為這一點(diǎn)的值。然后與原圖像進(jìn)行比較，具有相同數(shù)值的點(diǎn)就是原圖像的極大 （小）值點(diǎn)。

1．3 利用Delaunay三角插值擬合

一維的插值算法在二維圖像中直接應(yīng)用時會產(chǎn)生較大的誤差，因此，采用Delaunay三角插值算法，直接對極值點(diǎn)進(jìn)行插值擬合，這樣可以極大地減小圖像因為擬合而引起的誤差。

二維點(diǎn)集三角剖分是指將二維平面上的點(diǎn)集用不相交的直線段連接起來，使得所形成的凸包內(nèi)每一個區(qū)域都是三角形。Delaunay三角插值算法是先用Voronoi多邊形將離散的極值點(diǎn)分開，使每一個極值點(diǎn)屬于一個Voronoi多邊形，而連接3個共點(diǎn)的Voronoi多邊形內(nèi)的極值點(diǎn)則形成一個Delaunay三角形，所有這些Delaunay三角形的集合構(gòu)成Delaunay三角剖分。在每個三角形剖分內(nèi)進(jìn)行三次樣條插值，就得到了點(diǎn)集的擬合面。

Delaunay準(zhǔn)則實現(xiàn)的條件是凸包內(nèi)每一個三角形外接圓中不包含點(diǎn)集中的其他任何點(diǎn)，這使得每個三角形都盡可能接近于等邊三角形，避免產(chǎn)生狹長的三角形也就是使各離散點(diǎn)對整個三角形有限元網(wǎng)格的影響僅限于局部。

但是用Delaunay插值擬合時，在邊界處會產(chǎn)生一些奇異點(diǎn)；而這與之后進(jìn)行的形態(tài)學(xué)操作有一些矛盾，因此，要對這些奇異點(diǎn)進(jìn)行如下處理：用離奇異點(diǎn)相對近的點(diǎn)的灰度值來代替奇異點(diǎn)的灰度值，這一做法符合圖像的基本特性且處理之后的效果顯著。

2 試驗結(jié)果

采用Matlab［7］編程，筆者提出的改進(jìn)算法運(yùn)行的時間為1．797s，直接用傳統(tǒng)的經(jīng)驗?zāi)Ｊ椒纸膺\(yùn)行的時間是6．797s，改進(jìn)算法運(yùn)行比傳統(tǒng)方法快了5s；速度提高了73．56%，降低了時間復(fù)雜度，提高了運(yùn)算效率。圖像分解效果如圖2所示。

3 結(jié) 語

對Norden Huang提出的EMD算法進(jìn)行了改進(jìn)，改進(jìn)算法大大地降低了時間復(fù)雜度和空間復(fù)雜度。試驗結(jié)果表明，該算法能大大的減小對計算機(jī)存儲容量的要求，能有效的進(jìn)行分解處理，在圖像分解中有良好的效果，進(jìn)一步研究可以應(yīng)用在圖像去噪、圖像合成、圖像復(fù)原等方面。

圖2 試驗結(jié)果

［1］Norden E H．A new method for nonlinear and nonstationary time series analysis：empirical mode decomposition and hillbert spectral analysis［J］．Proc SPIE，2000，4056：197-205．

［2］徐曉剛，徐冠雷，王曉通，等 ．經(jīng)驗?zāi)Ｊ椒纸?（EMD）及其應(yīng)用 ［J］．電子學(xué)報，2009，03：151-153．

［3］岡薩雷斯 ．?dāng)?shù)字圖象處理 ［M］．北京：電子工業(yè)出版社，2007：420-454．

［4］張合勇，任德明，趙衛(wèi)疆，等 ．圖像處理中二維經(jīng)驗?zāi)Ｊ椒纸獾母倪M(jìn)算法 ［J］．光學(xué)學(xué)報，2009，29（5）：1248-1253．

［5］Waston D F．Computing the n-Dimensional Delaunay Tessellation with Application to Voronoi Polytops［J］．The Computer Jour-nal，1981，24 （2）：16-20．

［6］Nunes J C，Guyot S，Delechelle E．Texture analysis based on local analysis of the bidimensional empirical mode decomposition ［J］．Machine Vision and Applications，2005，16 （3）：177-188．

［7］劉衛(wèi)國．Matlab程序設(shè)計教程 ［M］．北京：中國水利水電出版社，2005．