陶 法, 馬舒慶

（成都信息工程學(xué)院電子工程學(xué)院,四川成都 610225）

1 引言

隨著海洋遙感技術(shù)的不斷發(fā)展,基于Seasat上的合成孔徑雷達(dá)（SAR）獲取海浪遙感圖像和無(wú)人機(jī)航拍獲取海浪圖像,對(duì)海洋和海洋氣象發(fā)揮了越來(lái)越重要的作用。準(zhǔn)確有效地獲取海浪波長(zhǎng)和波向,不僅對(duì)海浪預(yù)報(bào)、海洋工程,而且對(duì)反演海面風(fēng)場(chǎng)同樣具有重要意義[1]。文中主要針對(duì)利用無(wú)人機(jī)航拍海浪圖像獲取海浪波長(zhǎng)和波向算法進(jìn)行了研究,無(wú)人機(jī)航拍海浪圖如圖1所示。

2 海浪模型

根據(jù)Louguet-Higgins的理論,將平穩(wěn)海況下的海浪視為平穩(wěn)的具有各態(tài)歷經(jīng)性的隨機(jī)過(guò)程,波動(dòng)可看作是無(wú)限多個(gè)振幅不等、頻率不等、初相位不等的簡(jiǎn)單余弦波的疊加而成。忽略波與波之間的相互作用,可得到三維隨機(jī)海浪波高方程[2]:

其中,m,n分別為頻率分割數(shù)和方向分割數(shù),km為波數(shù),Φn為相對(duì)于x軸的傳播方向,x,y分別表示平行和垂直于主導(dǎo)波傳播方向的笛卡爾坐標(biāo),εmn為隨機(jī)相位角,Amn為不同頻率的簡(jiǎn)單正弦波的振幅,其表達(dá)式為:

其中,S（ωm,Φn）為海浪方向譜,且 S（ω,Φ）=S（ω）G（Φ）,Δωm 、ΔΦn分別為ωm 、Φn的增量,S（ω）為頻率譜函數(shù),即不同頻率間隔內(nèi)的組成波提供的能量,代表海浪能量相對(duì)于組成波頻率的分布,又稱為能量譜,G（Φ）為方向函數(shù)。在方向函數(shù)G（Φ）一定的條件下海浪的振幅 Amn與頻率譜函數(shù)S（ω）成正比關(guān)系,根據(jù)Neumann于1952年得到的一般表達(dá)式為[3]:

圖1 截取達(dá)里諾爾湖的無(wú)人機(jī)航拍海浪圖

其中,指數(shù) p取 4～6,q為2～4,A和B中包含風(fēng)要素（風(fēng)速、風(fēng)時(shí)、風(fēng)距）或波高要素（波高、周期）。

基于Neumann理論P(yáng)ierson-Moscowitz 1964年提出[3]:

其中 α=8.1×10-3,β=0.74,α g2=0.78。從理論上講,S（ω）分布在 ω=0～ ∞,實(shí)際上其顯著部分集中在一定頻帶內(nèi),其P-M譜圖如圖2所示。

3 基于重心坐標(biāo)的改進(jìn)算法

目前從海浪圖像中提取海浪波長(zhǎng)的方法主要有:二維傅里葉變換法（2D-FFT）[4]、一維傅里葉變換法（1D-FFT）[4]、基于聚類和hough變換法[5]、以及基于轉(zhuǎn)動(dòng)慣量[6]等方法。

3.1 2D-FFT法

方法對(duì)海浪圖像進(jìn)行二維傅氏變換,尋找代表海浪主頻率的最大點(diǎn)的位置,從而求出海浪傳播方向和波長(zhǎng)。

其中L,θ代表海浪的波長(zhǎng)和方向。xm,ym是功率譜圖中代表海浪主頻率最大值點(diǎn)的坐標(biāo),N代表海浪圖像的采樣數(shù),Δx代表采樣間隔。

2D-FFT法比較直觀,但主要缺點(diǎn)有:由于噪聲的存在和主頻不明顯,極易產(chǎn)生誤差,而用平滑法去除噪聲又會(huì)引起主頻位置的偏移;在計(jì)算海浪波長(zhǎng)和方向時(shí)只用到最大值點(diǎn),大部分求出來(lái)的頻率分量都沒(méi)有用到。

圖2 P-M譜圖

3.2 基于轉(zhuǎn)動(dòng)慣量法

假設(shè)海浪功率譜圖是一個(gè)有許多離散粒子組成的薄平面,那么波浪的方向就是使這些粒子圍繞某一軸的轉(zhuǎn)動(dòng)慣量最小的軸的方向,設(shè)該軸為y=kx,通過(guò)點(diǎn)（x0,y0）到直線的距離公式可得到:

那么其轉(zhuǎn)動(dòng)慣量為:

令轉(zhuǎn)動(dòng)慣量最小即dJ/dk=0,可求出:

其中k1×k2=-1。k1,k2分別代表使轉(zhuǎn)動(dòng)慣量最大和最小的值,并舍去最大的值。再沿著這條軸的方向在某一領(lǐng)域?qū)ふ易畲笾?該最大值點(diǎn)的坐標(biāo)就是海浪主波峰的位置,從而可以計(jì)算出波浪的波長(zhǎng)和波向。

方法在求取方向值時(shí)較方便,但求取最大值點(diǎn)時(shí)同樣受到噪聲點(diǎn)和主頻不明顯的影響。

3.3 對(duì)2D-FFT方法的改進(jìn)

從圖2中可以看出在中心頻率點(diǎn)處S（ω）最大,即由頻率為中心頻率的組成波提供的能量最大,可認(rèn)為該頻率為海浪的主頻率。但在實(shí)際中海浪不是單一頻率組成而且海浪的方向也不是單一方向的,由公式（1）可知每一個(gè)頻率組成波對(duì)海浪的波長(zhǎng)和波向都起一定的作用,其作用的大小由其系數(shù)值決定。所以僅以主頻率計(jì)算海浪波長(zhǎng)和波向存在誤差。

3.3.1 算法實(shí)現(xiàn)

將波浪圖像進(jìn)行二維傅立葉變換得到頻譜圖并通過(guò)頻移使其低頻部分位于圖像中心,（為了計(jì)算方便可以保留頻譜中大于某一閾值的值,其他值設(shè)為0）。由于頻譜是中心對(duì)稱的所以只需選取頻譜圖的一半,如圖3所示。

設(shè)像素點(diǎn)為（m×n）的頻譜圖中每一點(diǎn)坐標(biāo)代表圖像的頻率分量可設(shè)為 f（x,y）,其值代表圖像中該頻率分量所占的比重可設(shè)為 A（x,y）,則由所有頻率組成波提供的海浪功率譜為:

如果用某一特征頻率點(diǎn)為 f（x0,y0）來(lái)代表所有頻率點(diǎn)的話,則由該點(diǎn)計(jì)算出海浪功率譜為:

利用合力矩定理可以計(jì)算出該點(diǎn)坐標(biāo):

圖3 提取前1024個(gè)最大值點(diǎn)的圖像頻譜圖

即將圖像的功率譜看成密度不均勻的薄片平面,對(duì)其求重心,所得重心的坐標(biāo)即是所要求的特征頻率點(diǎn)。

3.3.2 仿真結(jié)果

選取閾值為32點(diǎn)到1024點(diǎn)得出到中心點(diǎn)的距離和方向值,結(jié)果如圖4所示。

從圖中可以看到,隨著取樣點(diǎn)的增多計(jì)算出的方向角漸漸變小,距中心點(diǎn)距離漸漸變大,但變化平穩(wěn)。其原因在于,隨著采樣點(diǎn)的增多圖像的高頻分量會(huì)相應(yīng)地增多,其方向角變化趨勢(shì)與海浪紋理細(xì)節(jié)的變化趨勢(shì)相一致;距中心點(diǎn)距離變大和P-M譜中隨著海浪譜的增大其主頻率點(diǎn)向低頻方向漂移相一致。

其中閾值為32點(diǎn)所求的重心坐標(biāo)為:（x0=5,y0=13）,方向角θg=71°,其空域仿真結(jié)果如圖5、6所示。

圖4 基于重心法的距離值Dg和方向值θg

圖5 基于重心點(diǎn)的空域海浪圖

圖6 截取閾值為32點(diǎn)的空域海浪圖

對(duì)2D-FFT方法的改進(jìn)算法的優(yōu)點(diǎn)為:所求的f（x0,y0）點(diǎn)是所有頻率點(diǎn)的重心擬合。不會(huì)造成大部分頻率點(diǎn)用不上的問(wèn)題;該重心坐標(biāo)是唯一的,不會(huì)產(chǎn)生主頻率點(diǎn)不明顯等問(wèn)題;能夠根據(jù)選取的閾值合理地反映實(shí)際海浪的波長(zhǎng)和波向。

4 分析與比較

圖7 基于重心法、最大值點(diǎn)的距離值 Δ D和方向值Δ θ

圖8 基于重心法和轉(zhuǎn)動(dòng)慣量法的方向值

在對(duì)圖1加入不同比重的鹽椒噪聲條件下分別基于重心法、轉(zhuǎn)動(dòng)慣量法、最大值法獲取到中心點(diǎn)的距離（Dg、Dm）、方向值（θg、θI、θm）然后作比較,結(jié)果如圖 9、10 所示 。

從圖中可以看到:與通過(guò)最大值點(diǎn)求取的相應(yīng)值作比較,結(jié)果是在一定范圍內(nèi)方向增加6°,距離增加1個(gè)像素;與基于轉(zhuǎn)動(dòng)慣量法所測(cè)得的方向值數(shù)據(jù)基本一致;基于最大值法求取的距離值Dm、方向值θm受加入鹽椒噪聲影響明顯;加入鹽椒噪聲后基于重心法和轉(zhuǎn)動(dòng)慣量法獲取的距離值Dg、方向值（θg、θI）基本保持不變。

圖9 加入鹽椒噪聲下重心法與最大值法的距中心點(diǎn)距離值

圖10 加入鹽椒噪聲下重心法、轉(zhuǎn)動(dòng)慣量法、最大值法的方向值

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

利用重心坐標(biāo)法對(duì)基于圖像功率譜法獲取最大值點(diǎn)計(jì)算海浪波長(zhǎng)、波向算法加以改進(jìn),并通過(guò)大量的試驗(yàn)與最大值點(diǎn)和轉(zhuǎn)動(dòng)慣量法求取的數(shù)據(jù)對(duì)比,驗(yàn)證了該方法計(jì)算方便、抗噪聲能力強(qiáng)、物理概念清晰,能夠更合理地表示海浪波長(zhǎng)和波向,并有效地解決了基于二維傅里葉變換帶來(lái)的主頻不明顯和大部分求出來(lái)的頻率分量都沒(méi)有用到等缺點(diǎn)。

[1] 陳艷玲,黃城,丁曉利.ERS22 SAR反演海洋風(fēng)矢量的研究[J].地球物理學(xué)報(bào),2007,50（6）:1688-1694.

[2] 俞韋修.隨機(jī)波浪及其工程應(yīng)用[M].大連:大連理工大學(xué)出版社,2003.

[3] 文圣常,余宙文.海浪理論與計(jì)算原理[M].北京:科學(xué)出版社,1984.

[4] 孫京生,劉征凱.利用SeasatSAR遙感圖象分析海浪的一種新方法[J].電子科學(xué)學(xué)刊,1988,10（6）.

[5] 劉政凱,畢亞雷,黃潤(rùn)恒.HOUGH變換在海浪遙感圖像處理中的應(yīng)用[J].環(huán)境遙感,1989,6（1）.

[6] Tang Huiming,Xu Shengrong.Detection of Direction and Wavelength of Ocean Wave by Power Spectrum of Ocean Wave Image[M].Pattern Recognition,1992:164-166.

[7] Rafael C,Gonzalez.Digital Image Processing[M].Prentice Hall,2002.