基于多尺度融合和相關性分析的全方向M型心動圖優(yōu)化研究

黃立勤 李 偉 林 強

（福州大學物理與信息工程學院，福州 350108）

引言

LEJ-2型全方向M型心動圖系統(tǒng)是福州大學生物醫(yī)學工程研究所研制開發(fā)的產(chǎn)品［1］，該系統(tǒng)是從超聲心動圖序列圖像的研究入手，得到了從一般視頻輸出序列圖像中提取動態(tài)信息的一種方法，即全方向灰度（位置）-時間波形圖法。該方法應用于心臟超聲圖像中，稱為全方向M型心動圖。它能在心臟室壁的任意結構、任意位置的某一條方向線上重建灰度（位置）—時間波形圖即組織運動曲線，通過該運動曲線可以體現(xiàn)心臟結構某個位置的運動情況。對組織運動曲線圖進行一階微分，可以得到組織瞬時速度，二階微分則可以得到瞬時組織加速度。在此基礎上作進一步研究，挖掘出速度場、加速度場、功率譜等直接和心臟結構機能和組織機體彈性素質(zhì)相關的重要信息，并讓這些信息可以以心電為時間基準互相同步地與時相相比較顯示出來。這些數(shù)據(jù)對于心血管疾病和血液動力學的研究有著非常重要意義［2］。

由于全方向M型心動圖的信息間接來源于二維超聲序列圖像，圖像質(zhì)量必然受超聲圖像質(zhì)量的影響。而超聲圖像本質(zhì)上具有模糊性和不均勻性特點，由此產(chǎn)生的圖像質(zhì)量嚴重阻礙了全方向M型心動圖運動曲線的準確提取，如圖1所示。如何去除噪聲、其他組織干擾等造成的虛假邊緣，提取出準確的組織運動曲線信息是全方向M型心動圖動態(tài)信息分析的基礎和前提，它的準確與否將嚴重影響到動態(tài)信息檢測的精度。

圖1 全方向M型心動圖Fig.1 Onmi-directional M-mode Cardiography

幾種廣泛使用的經(jīng)典的邊緣檢測方法有Sobel、Roberts和Laplacian等，這些算法的核心思想是假設邊緣點對應于原始圖像灰度級梯度的局部極值點。但是，全方向M型心動圖含有較大噪聲，這些算法對噪聲非常敏感，常常會把噪聲當作邊緣點檢測出來，而真正的邊緣由于噪聲的干擾可能被漏檢。目前全方向M型心動圖系統(tǒng)所設計使用檢測算法是針對全方向M型心動圖的構造特點，層進式地在有效搜索范圍內(nèi)運行線狀搜索模板［3］，實現(xiàn)目標運動曲線的檢測。但該方法對于目標運動曲線的檢測中的虛假的邊緣點會誤檢，仍需要較多的人工干預。

上述檢測算子對運動曲線檢測都在固定的尺度空間進行。固定尺度的邊緣檢測算子對定位精度和噪聲抑制等性能指標難以兼顧。小尺度算子雖然對噪聲敏感，但是有利于邊緣定位；大尺度算子對噪聲的抑制能力強，但對于邊緣定位精度較差，有時會丟失某些局部的細節(jié)［4］。Marr通過對神經(jīng)生理學和心理物理學的分析，指出人的視覺前期處理中有多個分辨率的邊緣算子在對圖像作卷積，如果各邊緣檢測算子輸出進行組合就可以提高定位精度，從而減少噪聲干擾［5］。由于小波算子具有良好的時頻局域化特性，特別適合于多尺度分析，因此考慮利用小波來構造多尺度邊緣檢測算子，通過對不同尺度下的多尺度運動曲線進行融合，從而實現(xiàn)運動曲線的檢測。

全方向M型心動圖是心臟二維超聲序列圖像中取樣線上各像素點因心臟結構在血液動力與心肌彈力的共同作用下引起運動產(chǎn)生的灰度（位置）隨時間變化掃描展開而重建出的波型圖，由于心臟是人體器官，它的運動前后具有時間、位置、能量相關性的特點。因此我們還可以分析這些相關性信息來結合前面經(jīng)過多尺度融合得到的數(shù)據(jù)，進一步限制運動曲線的搜索條件和搜索范圍以獲得更好的目標運動曲線。

1 算法描述

根據(jù)以上思路，處理方法可以分為3個步驟。步驟1：首先通過構建小波函數(shù)，在不同的固定尺度空間下進行運動曲線檢測。

步驟2：對不同尺度下的運動曲線進行融合。

步驟3：結合心臟運動的相關性信息生成最后需要的組織運動曲線。

1.1 固定尺度空間的運動曲線檢測

對任意的二維圖像函數(shù)f（x，y）∈L2（R2）進行小波變換，它具有兩個方向分量

寫成矢量形式為

式中，fk（x，y）是f（x，y）被θk（x，y）平滑后的圖像。檢測時把連續(xù)小波變換改成二進制形式，即取k=2j（j∈Z）。式（4）和式（5）分別代表了圖像f（x，y）沿需要x和y方向上的偏導，因此二維小波變換矢量代表的就是梯度［7］。

小波變換在尺度k=2j上梯度矢量的模以及梯度矢量與水平方向的夾角可以分別計算

平滑后的圖像f*θ2j（x，y）的突變點對應的是局部極大值的點，所以沿梯度方向檢測矢量模的極大值點可以獲得全方向M型心動圖的運動曲線（邊緣）。圖像在（x，y）點處的邊緣強度可以用式（7）來表示，圖像在（x，y）點處的方向矢量可以用式（8）來表示。具體步驟如下：

步驟1：要確定運動曲線（邊緣）必須刪除在局部范圍里面梯度不是最大的點，只保留局部范圍梯度最大的點［8］。首先對點（x，y）在梯度圖根據(jù)式（8）的方向歸并到圖2所示的4個方向，然后與該方向上的相鄰點的梯度比較，如果不是則表明不是邊緣點。經(jīng)過處理得到初步邊緣點圖fb1（x，y）。

圖2 方向圖Fig.2 Direction figure

步驟2：對初步邊緣點圖fb1（x，y）進行進一步分析。首先設定兩個上下限閾值Ta和Tb。這兩個閾值的設置可以有不同的算法，在本系統(tǒng)實驗測試中，首先確定Tb的值，Tb的取值為初步邊緣點圖fb1（x，y）中梯度直方圖中頻率最大的梯度，在實驗中發(fā)現(xiàn)Tb≈2Ta時，實驗效果較好。若邊緣點梯度大于Tb認為該點一定是邊緣點，可以得到確定邊緣點圖fb2（x，y）；若邊緣點梯度小于Tb同時大于Ta，則認為該點可能是邊緣點，可以得到可能邊緣點圖fb3（x，y）。fb2（x，y）的點使用較高的閾值計算出來，有比較多的間斷，因此通過在fb2（x，y）曲線邊緣處的點在fb3（x，y）的8個相鄰像素尋找可以鄰接的點，如果存在鄰接點，則把該點添加到fb2（x，y），一直不斷進行循環(huán)計算，最后得到fb2（x，y）是在一個固定尺度下的全方向M型心動圖的運動曲線圖。

在小尺度下檢測獲得的運動曲線比較準確，但是有較多的噪聲干繞；而在大尺度空間下檢測運動曲線的噪聲得到抑制，但許多邊緣信息被平滑掉了，而且得到的是位移后的結果。因此可以在不同的尺度空間下進行不同的閾值選擇。對于較大尺度空間，閾值可以選高一些用來減少噪聲產(chǎn)生的影響；在較小尺度空間，閾值選小一些盡可能地得到完整邊緣信息。在較大尺度空間下獲得的邊緣不完整的范圍較大。

1.2 不同尺度下的運動曲線融合

在不同尺度空間下檢測得到一系列不同尺度下的運動曲線。因為在相鄰尺度下的運行曲線最接近，因此運動曲線融合是對相鄰的尺度下的運動曲線進行，并且按照尺度大小順序逐步進行融合。在選擇尺度時是先設定最大尺度和最小尺度，再選擇尺度的間隔。雖然尺度越小檢測的運動曲線越接近真實，但噪聲較大，如果噪聲淹沒了運動曲線，則沒有意義。同樣尺度越大，越會抑制噪聲，但如果邊緣失真嚴重也沒有意義［9］。由于不同尺度下的檢測算子對全方向M型心動圖的運動曲線同一位置的響應并不相同，因此融合不能進行簡單的疊加。本方案利用了不同尺度下的運動曲線在位置、強度和方向上的聯(lián)系分析，依次采用傳遞、繼承、生長來處理運動曲線的融合。

設k1，k2是兩個相鄰尺度，k1＞k2，定義Fk1，k2（i，j）是尺度k1的3像素×3像素鄰域中的像素是尺度k2上局部模極大值點（i，j）的關聯(lián)域。定義集合Hk1是尺度k1上的局部模極大值點。則相關性可以表示為

式（9）可以用來表示在k2尺度上的點（i，j）與k1尺度的相關性。其中φk1（i，j）和φk2（i，j）代表尺度k1，k2上極大值點（i，j）的梯度方向，α表示它們方向差所設定的閾值。

1）運動曲線的傳遞

由式（9）我們可以進一步得到由k1尺度傳遞到k2尺度運動曲線。

尺度k2上的極大值點（i，j）如果與k1相關，則認為點（i，j）由k1傳承下來。Tk1，k2代表的是在k1和k2不同尺度空間反映全方向M性心動圖中同一的運動曲線。

2）運動曲線的繼承

由式（9）還可以得到由k1尺度繼承到k2尺度運動曲線。

Pk1，k2保留了k1尺度空間上的不與k2尺度存在相關性的點，代表的是k1到k2的繼承關系。

定義Dk1，k2和φk1，k2為通過上述兩個步驟對k1和k2尺度進行融合后的邊緣增強圖和梯度方向圖。

3）運動曲線的生長

前面介紹的運動曲線的傳遞僅在3像素×3像素鄰域中小窗口內(nèi)進行，這使得運動曲線信息無法傳遞到較遠的地方。如果增大了窗口尺寸，有可能有相關性的局部模極大值并不對應于同一運動曲線［10］。另外，噪聲在大尺度空間存在較少，即使大尺度鄰域內(nèi)存在噪聲，根據(jù)算法順序，優(yōu)先采用運動曲線的傳遞，其次才是運動曲線的繼承。因此算法是優(yōu)先采用傳遞的方法選擇小尺度空間上正確曲線上的點，不會采用繼承的方法得到噪聲點。得到運動曲線的傳遞和繼承后，還需要進行運動曲線的生長。本方案采用小窗口迭代的方法實現(xiàn)運動曲線的生長，具體方法如下：

不同尺度下的運動曲線融合的整個過程是在從大尺度下運動曲線附近尋找小尺度下運動曲線的相關性，若存在相關性，則用小尺度下運動曲線代替大尺度下運動曲線即完成傳遞；否則大尺度下運動曲線保留下來即繼承；然后在小尺度下完成運動曲線的迭代擴展。因為擴展是在3像素×3像素的鄰域內(nèi)逐步進行，而且對方向有限制，所以不容易擴展到噪聲里面。整個流程如圖3所示，圖3中尺度k0＞k1＞k2＞…＞kn-1＞kn。

圖3 運動曲線融合流程Fig.3 Flow chart of motion curves fusion

為去除各種干擾對全方向M型運動曲線融合結果的影響，本方案按照8鄰域細化算法對融合得到的運動曲線進行細化處理［11］，對邊緣輪廓進行8鄰域跟蹤搜索，判斷邊緣輪廓是否長度大于閾值。如果邊緣輪廓長度小于閾值，則認為是孤立點噪聲造成的假邊緣，在圖像中將其剔除掉。由此可以看出，該方法可以很好地消除超聲圖像的斑點噪聲影響，且預處理算法容易丟失精度信息，因此本實驗對要分析的全方向M型心動圖不進行降噪預處理。

1.3 全方向M型心動圖基于相關性分析的運動曲線檢測

全方向M型心動圖系統(tǒng)經(jīng)過了前面所闡述的多個尺度下進行融合檢測后，可能包含需要的和不需要的運動曲線，并且仍然可能出現(xiàn)斷點情況。根據(jù)全方向M型心動圖的具體特點，對其進行進一步的相關性分析，才能最終獲得符合要求的運動曲線圖。

1）全方向M型心動圖信息相關性

全方向M型心動圖是心臟二維超聲序列圖像中取樣線上各像素點由于心臟結構在血液動力與心肌彈力的共同作用下，引起運動產(chǎn)生的灰度（位置）隨時間變化掃描展開而重建出的波型圖如圖4所示。全方向M型心動圖相鄰兩列像素代表二維超聲序列圖像相鄰兩幀取樣線（投射）部位心臟結構的灰度像素點集合。由于心臟結構運動具有速度有限性，一般認為最大峰值速度不超過15cm/s，而相鄰兩幀二維超聲圖像的時間間隔一般在20～40ms，因此，反映在全方向M型心動圖的列像素集合之間的灰度信息就具有非常豐富的相關性。具體體現(xiàn)為：（1）時間（水平坐標）相關性 各列像素集合代表不同時刻取樣線上的灰度信息，且從左至右對應時間上的連續(xù)性；（2）位置（垂直坐標）相關性 心臟結構運動速度的有限性與幀時間間隔決定了取樣線上所代表的心臟結構的位置是密切相關的，反映在運動曲線的檢測上即相鄰列運動曲線點不可能出現(xiàn)幅度跳變；（3）灰度相關性又稱為能量相關性 待檢測的運動曲線在一定程度上近似代表心臟結構內(nèi)同一質(zhì)點的運動狀態(tài)，在未發(fā)生回波失落的情況下，反映在運動曲線的檢測上即相鄰列運動曲線點不可能出現(xiàn)灰度跳變［12］。

圖4 二維超聲序列圖像與全方向M型心動圖的聯(lián)系Fig.4 The relationship between two-dimensional ultrasound image sequences and onmi-directional M-mode cardiography

全方向M型心動圖的構造特點及其信息相關性的存在，為檢測運動曲線的搜索條件及其搜索范圍的設計提供了依據(jù)。

2）基于相關性的運動曲線檢測算法

心臟在幅度上的運動具有有限性，即在單位時間內(nèi)（按抽樣間隔算）超聲圖像在相鄰幀之間心臟各結構的運動幅度是有限的，而且是和前一幀該結構的運動位置密切相關的，這一重要信息反映在邊緣的提取上即為相鄰時間點的邊界不可能有大幅度跳變，表示如下：

（1）相鄰時間點的邊界位置是相關的。

（2）邊界搜索的區(qū)域是有限的。

因此，在邊緣檢測的算法中也應用了這一特點：

（1）找到當前時間點的邊界后，下一時間點對應的邊界的搜索范圍在該點的附近（向上和向下若干個點）進行。

（2）下一時間點邊界搜索的范圍是動態(tài)改變的，即首先在最小的區(qū)域（例如向上3個像素點和向下3個像素點）進行搜索，如果搜索不到滿足條件的邊界點，再逐級擴大搜索區(qū)域，進行下一級搜索。每一級搜索區(qū)域垂直方向上的步進值為3。

（3）最大搜索范圍是有限的，這個限度的設置也很重要。如果設置太大可能會找到一些干擾點，而且也是系統(tǒng)資源的浪費；設置太小，可能找不到滿足條件的邊界點。全方向M型心動圖是125幀/5s的序列圖像進行分析得到的，幀間隔為40ms，估計心臟結構的最大速度為30cm/s，一般心臟結構的速度都小于這個值。根據(jù)公式：最大運動范圍=最大速度×幀間隔時間，可得心臟各結構在這個時間間隔內(nèi)的最大運動范圍為1.2cm左右，一般心臟結構在這個時間間隔（40ms）內(nèi)的運動范圍均小于該距離。因此，設定最大搜索范圍為上下30個像素點，這個距離相當于1～1.5cm。實踐驗證，這個最大搜索邊界設置得比較合理。

（4）如果在這個區(qū)域內(nèi)仍然找不到滿足條件的點，即認為該處邊界是斷點或缺口（連續(xù)多個斷點）。斷點，是指一兩個沒有找到邊界點的時間點，直接用其前后兩個時間點的邊界位置的平均值替代，即采用直線插值法。缺口處理，可在缺口位置附近依據(jù)位置相關性，運用前、后時間段的邊界位置信息判斷缺口邊界所在的大致位置。

（5）根據(jù)全方向M型心動圖系統(tǒng)檢測上邊緣或下邊緣的需求，從上到下搜索或從下到上搜索獲得上邊緣或下邊緣。

（把此后的內(nèi)容按給的標題進行“填空”）

2 系統(tǒng)實驗

本實驗是在福州大學生物醫(yī)學工程研究所開發(fā)的LEJ-2型全方向M型心動圖系統(tǒng)上進行改進實驗。LEJ-2型全方向M型心動圖系統(tǒng)是用Visual C++6.0編寫，運行環(huán)境是Windows操作系統(tǒng)，CPU Intel Pentium Dual E2180 2.0G，RAM 2G。實驗分兩個部分。

2.1 多尺度融合的全方向M型運動曲線檢測

對一幅全方向M型心動圖圖像，經(jīng)過尺度分別為k=1、2、3、4的4級小波變換后，運用本研究提出的多尺度融合算法得到邊緣與傳統(tǒng)算子邊緣檢測的結果進行比較。

2.2 結合多尺度融合和相關性分析的全方向M型心動圖檢測

通過兩個實驗例子將結合多尺度融合和相關性分析的全方向M型心動圖檢測實驗結果，與目前已在使用的產(chǎn)品LEJ-2型全方向M型心動圖系統(tǒng)中使用的基于線狀模板搜索算法得到的檢測結果進行比較。

3 結果

圖5中（a）是一幅全方向M型心動圖圖像，對圖像進行小波變換，（c）～（f）是經(jīng)過4級小波變換后獲得的高于給定門限局部模極大值點位置圖，即各尺度下的邊緣圖，尺度分別為k=1、2、3、4。從邊緣圖可以看出，隨著尺度的增大，噪聲逐漸減少，邊緣逐漸平滑。在尺度1時，噪聲的影響非常大，邊界比較破碎；尺度k=4為最大尺度，提取的邊緣體現(xiàn)了原圖中的主要邊緣，基本不受噪聲的影響，但是邊緣失真比較嚴重，且提取的邊緣不完整，對于某些連續(xù)的邊緣，只檢測出其中的一段，但是在小尺度空間可以較完整地檢測出來。因此需要利用最大尺度空間提供的位置信息，融合各尺度的信息，合成精確的邊緣.運用本研究提出的多尺度融合算法，結果如圖5（b）所示。通過逐層融合。將小輪廓長度剔除原來斷裂的邊緣連接起來，而且邊緣位置越來越貼近實際邊緣位置，但是仍然可能有斷點存在。

圖6中（a）為Roberts算子（b）為Sobel算子（c）為Prewitt算子（d）為Kirsch算子檢測出的邊界。通過與圖5（b）比較可以看出：對于噪聲多、邊緣模糊的全方向M型心動圖圖像，用傳統(tǒng)算子檢測出的邊緣較模糊，去噪效果差，在定位精度、精確檢測等方面都不如本研究中所采用的方法。與其他邊緣檢測算子的比較可以獲得類似結果。

圖5 各尺度下邊緣檢測結果。（a）原圖；（b）融合之后圖；（c）k=1；（d）k=2；（e）k=3；（f）k=4Fig.5 The results of edge detection on different scale.（a）Original image；（b）Image after fusion；（c）k=1；（d）k=2；（e）k=3；（f）k=4

圖6 傳統(tǒng)算子邊緣檢測結果。（a）Roberts；（b）Sobel；（c）Prewitt；（d）KirschFig.6 The results of the traditional edge detection operator.（a）Roberts；（b）Sobel；（c）Prewitt；（d）Kirsch

圖7顯示了通過實驗1采用所提出算法與LEJ-2型全方向M型心動圖系統(tǒng)中使用的基于線狀模板搜索算法的比較結果。對比圖7中（b）和（c）可以看出，結合多尺度融合和相關性分析的全方向M型心動圖檢測在質(zhì)量較差的情況下，可以檢測到的全方向M型心動圖清晰的運動邊界，同時濾掉了第1、2波峰周圍虛假的邊緣，而傳統(tǒng)的LEJ-2全方向M型自動檢測則需要人工修正。

圖7 實驗1檢測結果與原系統(tǒng)檢測結果比較。（a）全方向M型心動圖1；（b）基于線狀模板搜索LEJ-2全方向M型系統(tǒng)自動分析結果；（c）結合多尺度融合和相關性分析的全方向M型心動圖檢測自動分析結果Fig.7 Experimental results 1 compared with the original system.（a）Onmi-directional M-mode Cardiography 1；（b）LEJ-2 Omnidirectional M-mode Echocardiography System analysis result using linear template；（c）Result using multi-scale integration and correlation analysis

圖8顯示了通過實驗2采用所提出算法與LEJ-2型全方向M型心動圖系統(tǒng)中使用的基于線狀模板搜索算法的比較結果。對比圖8中（b）和（c）可以看出，圖8（b）波峰上升的邊緣檢測到的是錯誤的邊緣，而圖8（c）檢測到正確的邊緣。

4 討論與結論

多尺度融合的全方向M型運動曲線檢測的實驗表明由于小波變換有多尺度的特點，可以利用多尺度特性，通過細節(jié)和粗節(jié)進行逼近，強于其他經(jīng)典算法；在邊緣和噪聲的取舍中，由于二者均為高頻信號，很難用頻帶劃分。使用小波變換的方法，使得可在大尺度下抑制噪聲，小尺度下得到邊緣的真實位置；而傳統(tǒng)的和經(jīng)典的邊緣檢測算法則在此問題上不能提供有效的解決辦法。不論選用怎樣的小波函數(shù)，都可以利用上述算法進行多尺度邊緣融合。因此該方法可以有效抑制噪聲的干擾，同時保證融合邊界定位的準確性。

圖8 實驗2檢測結果與原系統(tǒng)檢測結果比較。（a）全方向M型心動圖2；（b）基于線狀模板搜索LEJ-2全方向M型系統(tǒng)自動分析結果；（c）結合多尺度融合和相關性分析的全方向M型心動圖檢測自動分析結果Fig.8 Experimental results 2 compared with the original system.（a）Onmi-directional M-mode Cardiography 2；（b）LEJ-2 Omnidirectional M-mode Echocardiography System analysis result using linear template；（c）Result using multi-scale integration and correlation analysis

結合多尺度融合和相關性分析的全方向M型心動圖檢測的實驗表明：本算法對于全方向M型心動圖的虛假邊緣和錯誤邊緣的檢測效果較好。另外，檢測的運動曲線比較平滑，能較好地濾除噪聲干擾，而傳統(tǒng)的LEJ-2全方向M型自動檢測對于噪聲比較敏感；經(jīng)實驗測試，計算出運動曲線的時間小于1s，可以符合系統(tǒng)應用實際要求。

由于全方向M型心動圖的信息間接來源于二維超聲序列圖像，全方向M型圖像質(zhì)量必然受超聲圖像質(zhì)量的影響。而超聲圖像本質(zhì)上具有模糊性和不均勻性特點嚴重阻礙了邊緣的準確提取。本研究設計了一種基于多尺度融合和相關性分析的全方向M型心動圖檢測算法。通過該算法能自動去除全方向M型心動圖的虛假邊緣、無關噪聲等干擾，大大減輕了LEJ-2全方向M型系統(tǒng)的人工干預程度。本方案試驗分析是在福州大學生物醫(yī)學工程研究所開發(fā)的LEJ-2型全方向M型心動圖系統(tǒng)基礎上進行改進實驗，并提供詳細實驗對比數(shù)據(jù)和結果分析。該設計實現(xiàn)了對國家發(fā)明專利“全方向M型心動圖方法及其系統(tǒng)98 125713.5”系統(tǒng)的優(yōu)化設計。全方向M型心動圖系統(tǒng)具有完全的自主知識產(chǎn)權，改進后較大程度地減少了誤差，提高了精確度，這是加快我國醫(yī)療衛(wèi)生行業(yè)信息化程的非常價值的嘗試。

［1］林強，石江宏.超聲心動圖的一種動態(tài)信息—全方向M型心動圖［J］.儀器儀表學報，2005，26（4）：437-440.

［2］郭薇，陳斌，葉振盛.全方向M型超聲心動圖左室徑向心肌運動速度梯度評價局部心肌功能的應用價值［J］.中華超聲影像學雜志，2010，19（7）：565-568.

［3］Lin Qiang，Wu Wenji，Huang Liqin，et al.An omnidirectional M-mode echocardiography system and its clinical application［J］.Computerized Medical Imaging and Graphics.2006，30：333-338.

［4］阮松，陳松燦，王敏.采用多尺度多級組合分類器快速定位乳腺X片中的感興趣區(qū)域［J］.中國生物醫(yī)學工程學報，2008，28（5）：675-685.

［5］Marr.視覺計算理論［M］.姚國正，劉磊，汪云九，譯.北京：科學出版社，1988：256-260.

［6］張玉葉，周曉東，王春歆.小波估計圖像棱邊分布的邊緣保持規(guī)整化復原［J］.計算機輔助設計與圖形學學報，2009，21（1）：125-135.

［7］Yitzhaky Y，Peli E.A method for objective edge detection evaluation and detector parameter selection［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence.2003，25（8）：1027-1033.

［8］Canny J.A computational approach to edge detection［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1986，8（6）：679-698.

［9］Konishi S，Yuile AL，Coughlan JM，et al.Statistical edge detection：Learning and evaluating edge cues［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2003，25（1）：54-57.

［10］汪榮貴，朱靜，楊萬挺，等.基于照度分割的局部多尺度Retinex算法［J］.電子學報，2010，38（5）：1181-1186.

［11］Bao Paul，Zhang Lei.Noise reduction for magnetic resonance image via adaptive multiscale products thresholding［J］.IEEE Transactions on Medical Imaging，2003，22（9）：1089-1099.

［12］Chalana V，Linker DT，Haynor DR，et al.A multiple active contour model for cardiac boundary detection on echocardiographic sequences［J］.IEEE Transactions on Medical Imaging，1996，15（3）：290-298.