李焜 方世良 安良

(東南大學水聲信號處理教育部重點實驗室,南京 210096)

1 引言

確定水下聲源的位置是水聲領域中的一個關鍵問題.考慮到水聲環(huán)境的復雜性,為了能夠更準確地對水下目標實施定位,多數(shù)文獻從聲傳播的角度出發(fā),使用匹配場定位技術來確定聲源的位置[1-4].匹配場定位技術是聲場傳播規(guī)律與水聲信號處理相結合的技術,它充分考慮了水聲傳播規(guī)律的特點,利用聲場傳播模型計算預測的拷貝聲場,并與實際測量聲場進行相關,從而確定水下目標的位置.傳統(tǒng)的匹配場定位技術,一般多采用陣列的處理方式,具有大的孔徑,以獲得良好的陣增益和分辨性能.但是采用多陣元的大陣列,一方面增加了系統(tǒng)的開銷,給基陣的設計帶來不便;另一方面,在實際海水中布放時會受到諸如陣傾斜以及陣元失效等問題,增加了對水下目標定位的難度.此外,在某些應用方面,由于受安裝平臺尺寸的限制,也使得多陣元的布放無法實現(xiàn).因此,可否采用較少的水聽器尤其是單個水聽器來對目標實施定位,激發(fā)了研究人員的興趣,不斷激勵著相關研究人員為此進行探索.

使用單水聽器進行定位的一個難點在于可利用的信息量太少,主要是空間信息的缺乏.多數(shù)文獻借助寬帶信號的多頻點特性,采用“頻點換孔徑”的思想,對寬帶目標信號實施定位.文獻[5]將范數(shù)的不同表達式引入匹配場定位函數(shù)中,研究了利用單個水聽器來對目標進行定位;文獻[6]采用寬帶相干匹配場處理的方式使用單個水聽器來對目標進行定位;文獻[7]針對寬帶信號提出了基于時延匹配的單水聽器定位方法;文獻[8]針對超低頻信號采用模態(tài)濾波的方法來實現(xiàn)單水聽器的聲源定位;文獻[9]利用單水聽器采用直方圖濾波的方法對淺海中的移動目標進行定位;文獻[10]討論了在合作方式下借助波導不變量的性質(zhì)來對目標實施定位;文獻[11]利用單水聽器使用迭代優(yōu)化算法反演聲源的位置并進行地聲參數(shù)估計;文獻[12]則研究了利用單水聽器對海洋生物實施定位的方法.

對于淺海中傳播的信號而言,受海洋媒質(zhì)的影響,會產(chǎn)生頻散效應.頻散現(xiàn)象是信號的傳播速度與信號的頻率之間存在一定的關系,使得不同的頻率分量會以不同的速度傳播,同時受海水吸收和衰減的影響,造成接收波形的失真.就頻散本身而言,它一方面對發(fā)射的信號進行了更為復雜的變化,引入了具有非線性時頻形狀的多分量結構,增加了接收信號的復雜性,使得分量之間不容易進行辨識和分離;但另一方面,頻散本身蘊含了關于海洋環(huán)境和信號的相關信息,通過分析頻散波導中所接收到的信號,可有助于我們提取目標信號的特征并獲得目標的位置.

對海洋波導中頻散現(xiàn)象的研究,多集中于時頻分析的方法,通過時頻域2維結構來進行聯(lián)合表征[13-15].根據(jù)簡正波理論[16],接收信號是一系列傳播模式疊加所組成的,不同的模式以不同的群速度傳播,因而將以不同的時間間隔到達接收機,其在時頻面上的分布蘊含了傳播所包含的位置信息,聲源的位置不同就會表現(xiàn)出不同的頻散結構.利用時頻表征來分析頻散現(xiàn)象,最關鍵的是所采用的時頻分析方法可準確地反映所分析的信號,并通過相應的處理方法能在時頻面上孤立或是分離出每個傳播的模式.Bonnel等[17-19]提出了基于單水聽器的時頻翹曲算法(warping).此算法通過對接收信號進行酉變換操作,將每個模式轉變成近似Dirac函數(shù),從而補償波導的頻散效應,使得每個模式在短時傅里葉變換(STFT)后能夠在時頻域上實現(xiàn)分離,更好地獲取每個模式的特征,從而進行聲源位置估計[20]或地聲參數(shù)反演[21].

傳統(tǒng)的時頻分析方法如基于線性時頻表示的STFT受不確定性原理的影響,其時頻分辨率較低;而基于二次型時頻表示的Wigner-Ville分布、Choi-Williams分布等存在嚴重的交叉項干擾問題.另外,傳統(tǒng)的時頻分析方法都是基于固定核函數(shù)的設計方法,并不能很好地反映水聲脈沖信號短時瞬態(tài)的非平穩(wěn)特性.針對這一不足,Baraniuk等[22]提出了基于信號的自適應徑向高斯核函數(shù)(ARGK)的時頻分析方法,其核函數(shù)的形狀根據(jù)所分析的信號自適應地變化,提高了對于非平穩(wěn)信號的分辨能力.Li等[23-25]通過應用與ARGK時頻分布相類似的自適應最優(yōu)核函數(shù)的時頻分析方法來獲取簡正波的頻散特征,對海底聲速和密度等參數(shù)的反演進行了研究.

本文針對淺海傳播的低頻寬帶水聲脈沖信號,借鑒Li等反演時采用自適應時頻分布提取簡正波群延遲的思路,對利用單水聽器進行定位的問題進行了研究.分析了經(jīng)典Pekeris波導模型中的頻散現(xiàn)象,采用自適應徑向高斯核函數(shù)的時頻分析方法來表征接收信號的頻散特征;通過采用自適應徑向高斯核函數(shù)的時頻分布來提取頻散關系曲線中傳播模式的到達時間差,利用模式的到達時間差估計聲源的距離;通過二值掩模濾波,估計模式的能量,并采用多模式能量聯(lián)合匹配的方式,確定聲源的深度.

2 海洋中的模式傳播及頻散

在與距離無關的環(huán)境中,位于深度為zs的一個脈沖聲源,經(jīng)過海洋波導傳播后,在距離為r,深度為z處所接收到的聲場可以表示為[16]

其中,S(f)為聲源信號的頻譜,Ψm(z)為第m階與深度有關的模式函數(shù),krm(f)為第m階模式的水平波數(shù),ρ(zs)為聲源深度處的海水密度,M為總的傳播模式數(shù).

(1)式可進一步寫為[21]

從(2)式可以看出,接收點處的聲場是各階簡正波模式疊加所組成的.對于每一個模式,定義如下的相速度和群速度分別為

相速度和群速度刻畫了頻率與波數(shù)之間的關系,其中相速度是以某一特定相位傳播的速度,它代表了等相位面的傳播速度,而群速度則表示信號不同頻率分量傳播的速度,它反映了信號水平傳播的速度,也是能量傳播的速度,它是頻散關系中最為重要的量.從以上的分析可以看出,頻散現(xiàn)象是信號的傳播速度與信號的頻率之間存在一定的關系.對于在淺海中傳播的低頻脈沖信號而言,其頻散效應尤為明顯.為說明波導中的頻散效應,下面針對典型淺海環(huán)境模型Pekeris波導進行分析.

3 典型海洋環(huán)境的頻散曲線

Pekeris波導是一個具有兩層的分層結構的海洋波導,它與實際的海洋環(huán)境更為接近,如圖1所示,相關參數(shù)設置為海底聲速c2=1800 m/s,海底密度ρ2=1800 kg/m3,海水聲速c1=1500 m/s,海水密度ρ1=1000 kg/m3.Pekeris波導中,海底處于液態(tài)海底,海底聲速大于海水中的聲速,同時海底的密度也大于海水中的密度.

圖1 Pekeris波導模型

Pekeris波導中的頻散關系,可以通過如下的特征方程來表示[16]:

其中,kzm1為第m階模式在海水中的垂直波數(shù);kzm2為海底的垂直波數(shù).

Pekeris波導中各模式的最低頻率稱為截止頻率,可表示為

模式的截止頻率定義了波導中各階模式所能傳播的下限頻率.當頻率低于截止頻率時,不能激發(fā)有效的模式,因而對聲傳播沒有貢獻.

通過(3)式和(4)式可以獲得相應的相速度和群速度.圖2給出了Pekeris波導中D=120 m的前4個模式的相速度和群速度隨頻率變化的曲線.

圖2 Pekeris波導中相速度和群速度隨頻率的變化曲線

從圖中可以看出,Pekeris波導中相速度和群速度在高頻段均趨近于海水中的聲速c1;而在截止頻率處,這兩種速度均趨近于海底聲速c2.相速度隨頻率單調(diào)下降,而群速度則在某一頻率上存在極小值,此極小值稱為艾里相.

4 頻散效應的時頻表征

從以上模式的頻散曲線中可以看出,同一頻率處,不同的傳播模式具有不同的傳播速度,描述的是模式之間的頻散關系,稱為模態(tài)間頻散;而同一模式在不同的頻率處具有不同的傳播速度,則反映了單一模式的頻散現(xiàn)象,稱為模態(tài)內(nèi)頻散[21].對于寬帶脈沖而言,既存在模態(tài)間頻散也包含模態(tài)內(nèi)頻散,每個模式以不同的群速度傳播,因而以不同的時間到達接收水聽器.各模式的到達時間定義為[18]

由此可以看出,各模式的到達時間通過聲源與接收機之間的距離以及各模式的群速度相聯(lián)系.圖3給出了Pekeris波導中傳播距離為15km的前4個傳播模式到達時間隨頻率的變化關系.

圖3 模式到達時間隨頻率的變化

從圖中可以看出,受激發(fā)頻率的不同,每個模式從低頻的艾里相開始一直延伸到最高頻率,各模式都具有不同的到達時間.從頻散曲線上看,高頻分量的傳播要快于低頻分量,在高頻段,頻散的作用漸漸減弱,頻散曲線逐漸變成直線,各模式間傳播的時間間隔變小;而在低頻段,頻散曲線隨著頻率的降低,各模式間傳播的時間間隔變大,頻散效應更加地明顯.

通過以上分析可知,受模式群速度的影響,各模式在不同頻率處具有不同的到達時間.因此,可通過時頻分析的方法對接收信號的頻散效應進行表征.相應的時頻分布可以表示為

其中,(Am(f)代表)接收信號中各模式所對應的幅度,則反映了每個模式在時頻面上的位置.

由于每個模式在時頻面上的分布不同,具有不同的形狀,因此,時頻分析最關鍵的是所采用的時頻分析方法可準確地反映所分析的信號,且能夠更好地在時頻域辨識和分離各自的模式.采用STFT的方法,對Pekeris波導深度為120 m,距離分別為5 km和15 km以及距離為15 km,波導深度分別為40 m和80 m的4個傳播模式進行仿真,相應的結果如圖4所示,其中的虛線表示理論的頻散曲線.

從圖中可以看出,在同一波導深度下,隨著距離的增加,模式之間的分離性變大,模式低頻段的頻散效應更加明顯;而隨著波導深度的減小,各模式的截止頻率升高,模式的頻散效應更加的明顯,頻散向更高的頻段擴散,在時頻面上的分離性更大.但從整體上看,STFT的方法受時頻分辨率的限制,在時頻表征上并不理想.

為了提高時頻分辨率,采用ARGK的時頻分析方法,相應的結果如圖5所示.可以看出,ARGK較STFT能夠獲得更高的分辨率,可明顯地反映出所分析信號的特征.

5 定位原理

5.1 距離估計

由以上的分析可知,頻散信道中接收到的信號是由多個分量所組成的,每個分量根據(jù)傳播模式的不同,在時頻面上會有不同的形狀.各模式的到達時間蘊含了聲源的距離信息,因此可通過測定模式間的到達時間差來估計聲源的距離.

圖4 STFT時頻分析的結果 (a)D=120 m,R=5 km;(b)D=120 m,R=15 km;(c)R=15 km,D=40 m;(d)R=15 km,D=80 m

圖5 ARGK時頻分析的結果 (a)D=120 m,R=5 km;(b)D=120 m,R=15 km;(c)R=15 km,D=40 m;(d)R=15 km,D=80 m

對于時間隨頻率變化的非平穩(wěn)信號,在頻率 f處,第m階模式傳播的時間可表示為

其中,te(f)為發(fā)射信號的時間隨頻率的變化,也即是發(fā)射信號的調(diào)制率.

則兩個模式的到達時間差可表示為

其中,Kmn(f)稱為群慢度之差.

模式的群速度可通過理論計算求出,而模式的到達時間差可通過ARGK時頻分析來獲取,從而可確定出聲源的距離.

5.2 深度估計

對于聲源深度的估計,本文借鑒文獻[26]利用水平陣列采用頻率波數(shù)域變換來定深的匹配模式的思想,對提取出的各模式采用模式能量匹配的方法進行深度估計.令am,real為從實際接收信號提取出的第m階模式的能量,am,replica為拷貝場信號中提取出的第m階模式的能量,則構造如下的代價函數(shù):

其中,Nm為估計所用的模式數(shù).

通過(12)式所確定的代價函數(shù),在聲源深度范圍內(nèi)進行峰值搜索,則可確定聲源的深度

5.3 二值掩模濾波

由以上的分析可以看出,使用單水聽器進行聲源的位置估計需要提取所需的傳播模式,估計出模式的到達結構以及模式的能量.由于各模式在時頻面上的出現(xiàn)時間不同且能量分布不均勻,為準確提取定位所需模式,在對接收信號進行ARGK時頻分析的基礎上,本文采用類似文獻[26]中的二值掩模濾波的方法進行模式的特征提取.

二值掩模濾波是一種時頻濾波的方法,通過二值掩模濾波,可以提取所需的模式.所謂的二值,是指其在時頻面某個位置上的取值只有0和1兩種,即

由于各階模式在時頻面上的分布并非是理想的曲線,而是具有一定寬度的帶狀區(qū)域,因此,需要對所設計的二值掩模濾波器中的頻散曲線進行加寬處理.采用圖像形態(tài)學中的膨脹方法[27],實現(xiàn)對各模式的擴充.圖6給出了圖5(d)表示的第4階模式進行膨脹前后的掩模濾波輸出.

通過ARGK時頻分析的方法,確定模式在時頻面上的位置,設計二值掩模濾波矩陣,將二值掩模濾波矩陣與接收信號的時頻矩陣相乘,即可提取出相應的模式.

6 仿真實驗

6.1 仿真環(huán)境及參數(shù)設置

本文選用經(jīng)典的Pekeris波導模型,波導環(huán)境參數(shù)如圖1所示,其中海水深度D=120 m.發(fā)射信號為線性調(diào)頻信號,頻帶范圍為40-120 Hz,脈沖寬度為0.05 s.聲源位于水下30 m,距離接收機為15 km.單一接收水聽器放置于海底.此波導中頻率為120 Hz所激發(fā)的9個模式形狀函數(shù)隨深度的變化如圖7所示.從圖中可以看出,各階模式隨深度的變化類似正弦曲線的形狀,第m階模式有m個過零點.

圖8分別給出了脈沖響應各階模式的理論到達時間隨頻率變化的頻散曲線以及發(fā)射信號經(jīng)過此波導后頻散曲線的變化.

圖6 二值掩模濾波器膨脹前后的輸出 (a)未經(jīng)膨脹的二值掩模濾波器;(b)經(jīng)膨脹之后的二值掩模濾波器

圖7 模式形狀函數(shù)隨深度的變化

圖8 各階模式到達時間隨頻率變化的頻散曲線 (a)脈沖響應的理論頻散曲線;(b)接收信號的理論頻散曲線

從圖8可以看出,傳播的不同模式隨頻率的變化具有不同的到達時間,高頻分量的傳播要快于低頻分量,各模式的到達時間隨著模式數(shù)的增加而變慢.受激發(fā)頻率的不同,每個模式從低頻的艾里相開始一直延伸至最高頻率,各模式的頻散曲線結構各不相同.1階模式的頻散效應持續(xù)時間最短,頻散曲線主要集中在低頻段;隨著模式數(shù)的增加,各模式頻散曲線的持續(xù)時間逐漸增大.在高頻段,頻散的作用漸漸減弱,各模式間傳播的時間間隔變小,頻散曲線逐漸變成直線;隨著頻率的降低,各模式間傳播的時間間隔變大,頻散效應更加地明顯.對于接收信號的頻散曲線而言,其頻散結構與圖8(a)的情形相類似,但受發(fā)射信號調(diào)制率的影響,頻散曲線將以發(fā)射信號的調(diào)制率作為斜率而發(fā)生傾斜.

6.2 時頻表征結果

發(fā)射信號通過此波導后,單一接收水聽器上所接收到的信號如圖9(a)所示,傳播的各階模式的出現(xiàn)范圍已在圖中標出,相應的時頻結構如圖9(b)所示,其中的虛線代表理論的頻散曲線.

從時域波形圖中可以看出,接收信號是由傳播的各階模式疊加所組成的,信號波形在時間上被展寬,產(chǎn)生波形失真,且受波導頻散效應的影響,各模式傳播所出現(xiàn)的時間不同.從時頻圖中可以看出,不同的模式具有不同的時頻形狀,每個模式在時頻面上的分布具有不同的能量,能量的變化反映了模式形狀函數(shù)隨深度的變化.其中,模式6和模式7由于靠近模式函數(shù)波峰(谷)的位置,因而其能量最強,相應的時頻表征具有很強的亮度;而模式1的持續(xù)時間最短,所引起的頻散效應最弱,到達接收機呈短時脈沖信號,因而其能量較小,在時頻面上未被顯示;而模式4由于處在模式函數(shù)的過零點位置,因此不被激發(fā),從而在時頻面上不出現(xiàn)第4階模式.另外,由于峰值能量在模式8和模式9的激發(fā)頻率之下,因此模式8和模式9的能量較小,相應的時頻表征的亮度較弱.

6.3 模式提取結果

在ARGK時頻面上,對每個模式使用經(jīng)膨脹之后的二值掩模濾波方法進行提取,相應的結果如圖10所示.

圖9 單水聽器上的接收信號和相應的時頻顯示 (a)接收信號的時域波形及各階模式的顯示范圍;(b)接收信號的時頻結構

圖10 7個模式提取后的綜合顯示

6.4 定位結果分析

在各模式的出現(xiàn)頻段內(nèi),提取出各頻點對應的到達時間差獲得聲源的距離.由于各模式所對應的截止頻率不同,因此本文在各模式最大相同的頻帶范圍來衡量測距結果.表1給出了在相同頻帶內(nèi)選用不同模式組合情況下的距離估計結果.

表1 單水聽器的距離估計結果

從表中的結果可以看出,不同模式組合下的距離估計結果稍有不同.其中,模式6和模式7因在時頻面上的能量較高,易于辨識和分離,因此,估計的結果較為準確,其距離估計相對誤差小于2%.而模式3因其所具有的能量較低,時頻表征的效果較差,使得估計得到的相對到達時間的誤差較大,從而與其他模式組合之后對測距的結果會有一定的影響.

圖11 各階模式的能量分布

采用模式能量匹配的方式進行深度估計,各模式的能量分布如圖11所示.選用能量較高的模式5、模式6以及模式7三個模式,以步距為5 m,在深度搜索范圍為5-105 m的范圍內(nèi)進行多模式能量匹配確定聲源的深度,相應的結果如圖12所示.從圖中可以看出,對所提取的模式采用模式能量匹配的方式進行搜索,所確定出的峰值較為明顯,但同時,不同模式的組合方式在定深效果上會有所差別.其中,模式5和模式7的模式能量匹配效果相對于模式6和模式7存在更多的偽峰,同時參與模式匹配的個數(shù)越多,匹配所能獲得的判別信息增加,使得主峰值更加地尖銳,同時具有低的偽峰,從而深度估計的性能有所提升.

圖12 深度估計的結果

7 結論

利用單水聽器來對水下聲源實施定位一直是國內(nèi)外研究的熱點和難點.本文針對淺海環(huán)境中低頻寬帶水聲脈沖信號,通過理論分析和仿真計算對基于頻散特征的單水聽器定位問題進行了研究.分析了經(jīng)典Pekeris波導環(huán)境下的頻散現(xiàn)象,采用自適應徑向高斯核函數(shù)的時頻分析方法來表征接收信號的頻散特征,提取頻散關系曲線中傳播模式的到達時間差估計聲源的距離.通過二值掩模濾波的時頻方法,估計模式的能量,采用多模式聯(lián)合匹配的方式,確定聲源的深度.

通過對基于Pekeris波導模型的淺海環(huán)境進行仿真驗證,結果表明:自適應徑向高斯核函數(shù)的時頻分析方法能夠很好地反映信號本身的頻散特征,具有較高的時頻分辨率,克服了傳統(tǒng)短時傅里葉變換時頻表征的限制,使得模式在時頻域更加容易辨識和分離,模式頻散曲線的提取更加地準確;從測距結果來看,不同模式組合下的距離估計結果不同,采用在時頻面上具有較高能量的模式,可得到較為準確的距離估計;采用二值掩模濾波的方法能夠獲得可靠的模式的特征參數(shù);在定深方面,參與聯(lián)合匹配的模式個數(shù)越多,深度估計的性能會進一步的提升.

本文對經(jīng)典的Pekeris波導模型中利用模式頻散特征的定位問題進行了研究,證明了此方法具有一定的可行性.但由于海洋環(huán)境存在的時變性以及不確定性,如何對更為復雜的海洋環(huán)境利用頻散特征來進行定位的問題,有待進一步的研究.

[1]Porter M B,Tolstoy A 1994 J.Acoust.Soc.Am.2 161

[2]Xu W,Xiao Z,Yu L 2011 IEEE J.Ocean.Eng.36 273

[3]Wu K M,Ling Q,Wu L X 2011 IEEE International Conference on Signal Processing,Communications and Computing(ICSPCC),Xi’an,September14-16,2011 p1

[4]Wang Q,Jiang Q 2010 EURASIP J.Advances.Signal Processing.483524 1

[5]Frazer L N,Pecholcs P I 1990 J.Acoust.Soc.Am.88 995

[6]Lee Y P 1998 IEEE Oceans’98 Conference Proceedings Nice,September 28-October 1 1998 p1074

[7]Jesus S M,Porter M B,Y St′ephan,D′emoulin X,Rodriguez O C,Ferreira Coelho E M M 2000 IEEE J.Ocean.Eng.25 337

[8]Touz′e G L,Torras J,Nicolas B,Mars J 2008 IEEE Oceans,Quebec City,September 15-18,2008 p1

[9]Jemmott C W,Culver R L,Bose N K 2008 IEEE Conference on Signal,Systems and Computers,Paci fic Grove,October 26-29,2008 p283

[10]Tao H L,Hickman G,Krolik J L,Kemp M 2007 IEEE Oceans Aberdeen,June 18-21,2007 p1

[11]Gac L J C,Asch Mark,St′ephan Y,Demoulin X 2003 IEEE J.Ocean.Eng.28 479

[12]Tiemann C O,Thode A M,Straley J,O’Connell Victoria,Folkert K 2006J.Acoust.Soc.Am.120 2355

[13]Chen C S,Miller J H,Boudreaux G F,Potty G R,Lazauski C J 2003 IEEE Oceans 2003.proceedings,San Diego,September 22-26,2003 p2903

[14]Touz′e G L,Nicolas B,Lacoume J L,Mars J,Fattaccioli D 2005 IEEE Europe Oceans Brest,June 20-23,2005 p725

[15]Ioana C,Jarrot A,Gervaise C,St′ephan Y,Quinquis A 2010 IEEE Trans.Signal Processing.58 4093

[16]Jensen F B,Kuperman W A,Porter M B,Schmidt H 1994 Computational Ocean Acoustics(New York:American Institute of Physics)p337

[17]Bonnel J,Nicolas B,Mars J I,Fattaccioli D 2009 IEEE Oceans’2009 Biloxi,October 26-29,2009 p1

[18]Bonnel J,Gervaise C,Roux P,Nicolas B,Mars J I 2011 J.Acoust.Soc.Am.130 61

[19]Bonnel J,Gervaise C,Nicolas B,Mars J I 2012 J.Acoust.Soc.Am.131 119

[20]Lopatka M,Touz′e G L,Nicolas B,Cristol X,Mars J I,Fattaccioli D 2010 EURASIP J.Advances.Signal Processing.304103 1

[21]Bonnel J,Nicolas B,Mars J I,Walker S C 2010 J.Acoust.Soc.Am.128 719

[22]Baraniuk R G,Jones D L 1993 Signal Processing 32 263

[23]Li Z L,Zhang R H 2007 Chin.Phys.Lett.24 471

[24]Zhang D M,Li Z L,Zhang R H 2005 Acta Acoustica 30 415(in Chinese)[張德明,李整林,張仁和2005聲學學報30 415]

[25]Zhang X L,Li ZL,Huang XD 2009 Acta Acoustica 3454(in Chinese)[張學磊,李整林,黃曉砥2009聲學學報30 54]

[26]Nicolas B,Mars J I,Lacoume J L 2006 EURASIP J.Applied Signal Processing 65901 1

[27]Rein van den B,Richard van B 1992 Computer Vision,Graphics,And Image Processing:Graphical Models And Image Processing 54 252