一種高頻雷達(dá)回波譜反演海浪的新算法?

常廣弘, 黎 明??, 楚曉亮, , 紀(jì)永剛, 于長軍

(1.中國海洋大學(xué), 山東 青島 266100；2.國家海洋局第一海洋研究所, 山東 青島 266061；3.哈爾濱工業(yè)大學(xué)(威海), 山東 威海 264209)

?

一種高頻雷達(dá)回波譜反演海浪的新算法?

常廣弘1, 黎 明1??, 楚曉亮1, 2, 紀(jì)永剛2, 于長軍3

(1.中國海洋大學(xué), 山東 青島 266100；2.國家海洋局第一海洋研究所, 山東 青島 266061；3.哈爾濱工業(yè)大學(xué)(威海), 山東 威海 264209)

在高頻地波雷達(dá)海浪譜反演問題中, 廣泛采用的Barrick后向散射公式屬于第一類非線性Fredholm積分方程。此類積分方程的解在本質(zhì)上是不適定的，加之高頻雷達(dá)二階回波信號信噪比較低，使得反演海浪譜存在解不穩(wěn)定的問題。本文提出一種穩(wěn)定且低復(fù)雜度的反演算法，此算法首先根據(jù)高頻雷達(dá)一階回波譜測量海浪方向，并將其引入積分方程求解過程，減少求解變量的個數(shù)，降低反演算法的復(fù)雜度。為解決反演結(jié)果不穩(wěn)定的問題，使用Tikhonov正則化方法并利用廣義交叉驗(yàn)證法(GCV)確定其正則化系數(shù)。通過在不同測試條件下對反演算法的仿真測試分析，表明此方法具有運(yùn)算量小、穩(wěn)定性好的特點(diǎn)。

高頻地波雷達(dá)； 海浪譜； 反演； 正則化

高頻地波雷達(dá)發(fā)射垂直極化電磁波(3～30MHz)繞海面?zhèn)鞑サ倪^程中與海浪相互作用產(chǎn)生后向散射，雷達(dá)回波譜中包含豐富的海表面狀況信息。Hasselmann[1]提出了二階水動力和電磁作用產(chǎn)生二階峰的概念，認(rèn)為一階峰周圍的二階峰正比于海浪頻譜。在此基礎(chǔ)上，Barrick[2-3]根據(jù)電磁波的微擾展開原理，并假設(shè)海表面為良導(dǎo)體，推導(dǎo)了后向散射公式，此公式屬于第一類非線性Fredholm積分方程。

基于Barrick后向散射公式，Barrick[4]、Lipa[5]、Wyatt[6]、Gurgel[7]和Howell[8]等人都提出了通過雷達(dá)回波譜提取海浪信息的算法。其中，Barrick法對二階峰的積分進(jìn)行處理，可測量有效波高；Lipa將回波譜分段處理，分別得到方向分布系數(shù)和長浪區(qū)域的海浪譜；Gurgel采用計(jì)算回歸系數(shù)的方法反演海浪有效波高；Howell將海浪譜展開成有限項(xiàng)傅里葉級數(shù)的形式，該方法存在截?cái)嗾`差。吳雄斌等[9]修正了Barrick反演模型，提高了模型抑制噪聲和抗干擾的能力，并且使模型應(yīng)用于寬波束雷達(dá)，通過與浮標(biāo)測量結(jié)果對比，證明了有效波高測量的有效性。

Lipa和Howell的海浪譜反演算法屬于對第一類非線性Fredholm積分方程的求解，此類方程的解在本質(zhì)上是不適定的[10](不適定問題的概念為：解存在、解唯一和解穩(wěn)定[11])，而且高頻雷達(dá)二階回波信號信噪比較低，加劇了反演結(jié)果的不穩(wěn)定性，因此，采用適當(dāng)?shù)恼齽t化方法是非常必要的。李倫等[12]對積分方程離散化處理后，采用正則化方法，得到方程的正則逼近解，通過仿真實(shí)驗(yàn)證明了算法的有效性，但此種方法需要對大型矩陣進(jìn)行操作，運(yùn)算量較大。

本文提出的反演算法使用雷達(dá)一階回波來測量海浪方向，并將其引入積分方程求解過程，減少了反演中求解變量的個數(shù)，為保證反演結(jié)果的穩(wěn)定性，使用Tikhonov正則化方法并利用廣義交叉驗(yàn)證法(GCV)確定其正則化系數(shù)。本文在詳細(xì)介紹反演算法的基礎(chǔ)上，通過不同條件下的仿真測試，分析了算法的特點(diǎn)。

1 反演算法

本節(jié)首先介紹Barrick后向散射公式，在此基礎(chǔ)上，使用兩種策略來降低反演算法的復(fù)雜度：第一，引入海浪方向信息，減小反演過程中未知變量個數(shù)；第二，對離散化處理后的海浪譜進(jìn)行近似處理，降低反演過程中的矩陣維數(shù)。

1.1 后向散射公式

無海流情況下的一階后向散射公式為[1]：

(1)

其中:m′=±1表示多普勒頻移的正負(fù)；k0為雷達(dá)電磁波波矢，以發(fā)射電磁波方向?yàn)檎?k0為其波數(shù)；S(·)為海浪方向譜；δ(·)為狄拉克函數(shù)；波數(shù)為雷達(dá)電磁波波數(shù)2倍的海浪產(chǎn)生最強(qiáng)的后向散射，并根據(jù)深水行進(jìn)波色散關(guān)系ωB與k0之間滿足：

(2)

其中:g為重力加速度。

無海流情況下的二階后向散射公式為[3]:

(3)

其中：m,m′=±1表示產(chǎn)生后向散射的4種海浪組合；k和k′為產(chǎn)生后向散射的兩列海浪波矢；k和k′分別表示兩列海浪的波數(shù);Γ為耦合系數(shù)[5]。

1.2 反演算法

對式(3)進(jìn)行歸一化處理并根據(jù)狄拉克函數(shù)的性質(zhì)進(jìn)行化簡：

(4)

其中:η=ω/ωB為歸一化多普勒角頻率;K=k/(2k0)、K′=k′/(2k0)為歸一化的海浪波數(shù)，為消除雷達(dá)增益及路徑損失等對回波的影響，將二階峰除以一階峰能量，得到歸一化的二階譜：

(5)

為了降低算法的復(fù)雜度，減少反演參數(shù)，將海浪方向譜表示為海浪譜與方向分布因子乘積的形式：

S(k)=S(k)·G(θw)

(6)

海浪方向可以通過正負(fù)一階峰的幅值差[14]來求得:

|θw|=2tan-1(10R/(10s))

(7)

其中:R=10log(B+/B-)，B+與B-分別為正負(fù)一階峰的幅值。

根據(jù)狄拉克函數(shù)及海浪的特性，將(6)帶入(5)，并引入S(k′)進(jìn)行線性化處理，在反演過程中，為了避免過多的引入誤差，只對能量最強(qiáng)的二階峰進(jìn)行操作，在此以m=m′=1為例進(jìn)行介紹：

(8)

σ=Kn×S

(9)

i=1,2…n,

其中：σ∈Rn×1為已知的歸一化二階回波譜頻點(diǎn)；Kn∈Rn×(n·m)為系數(shù)矩陣，S∈R(n·m)×1為海浪譜。由于方程數(shù)量要遠(yuǎn)小于未知數(shù)個數(shù)，且離散后的海浪譜中含有大量的重復(fù)及近似元素。因此，對離散后的S矩陣進(jìn)行排序，隔m點(diǎn)取值作為S的近似值，將S矩陣進(jìn)行降維，并將核矩陣Kn中對應(yīng)于相同S的元素求和得到C矩陣。于是式(9)可表示為：

σ=C×S

(10)

其中：C∈Rn×n為系數(shù)矩陣；S∈Rn×1為海浪譜，通過求解此方程組便可以得到海浪譜S。

2 正則化方法及算法流程

實(shí)際處理中，式(10)中系數(shù)矩陣C存在不滿秩的情況，此時C的逆不存在，而基于廣義逆方法求得的海浪譜S使殘差‖CS-σ‖2取最小值，即S=C+×σ，但此解的穩(wěn)定性得不到保證，因此采用適當(dāng)?shù)恼齽t化方法是必要的。正則化方法可以分為兩類：Tikhonov正則化方法和Landweber迭代方法。當(dāng)矩陣維數(shù)較小時，Tikhonov方法運(yùn)算量要小于Landweber迭代算法。由于已經(jīng)對矩陣進(jìn)行了降維處理，選擇Tikhonov正則化方法，來降低整個算法的運(yùn)算量。此方法中涉及正則化系數(shù)的選擇問題，為了得到最優(yōu)的正則化系數(shù)，引入GCV正則化系數(shù)選擇方法。

在此，根據(jù)Tikhonov正則化方法原理及所求海浪譜S為有界且光滑的特性，在解中增加罰項(xiàng)，若解的最初估計(jì)為S*，則增加的罰項(xiàng)為：

Ω(S)=‖S-S*‖2

(11)

增加罰項(xiàng)后必須考慮‖CS-σ‖2與Ω(S)之間的關(guān)系，Tikhonov通過引入正則化系數(shù)[15-16]來權(quán)衡兩者的大小，即

(12)

其中：λ為正則化系數(shù)，控制殘差范數(shù)及罰項(xiàng)Ω(S)使得解范數(shù)最小化，將求解不適定問題轉(zhuǎn)化為最優(yōu)化參數(shù)選擇的問題。很明顯，當(dāng)λ較大時，可使海浪譜比較平滑，但會丟失海浪譜的細(xì)節(jié)信息；當(dāng)λ較小時，所求得的海浪譜使得殘差較小，但海浪譜受雷達(dá)回波譜中噪聲的影響較大。因此，正則化系數(shù)的選擇是一個非常關(guān)鍵的問題。

計(jì)算正則化系數(shù)一般有L曲線法和GCV法，在處理離散問題時，L曲線法要進(jìn)行曲線擬合和插值處理，導(dǎo)致計(jì)算量增大。在此采用GCV法[17]，此方法的計(jì)算公式可以表示為：

(13)

其中：Cλ為正則化矩陣；tr表示矩陣的跡。GCV法通過改變λ的值，使得GCV(λ)取得最小值，來求取最優(yōu)化的正則化系數(shù)，即得到最終適定的正則化解：

Sλ=Cλσ

(14)

經(jīng)過離散化處理與正則化處理后，最終得到了后向散射方程的正則逼近解Sλ，算法流程見圖1。

3 仿真驗(yàn)證及分析

設(shè)計(jì)4組仿真實(shí)驗(yàn)來驗(yàn)證本文算法的正確性。實(shí)驗(yàn)的基本思路是首先利用式(1)和(3)，將雷達(dá)頻率、風(fēng)向、風(fēng)速和信噪比作為輸入變量，計(jì)算雷達(dá)回波譜[18]，然后通過改變不同輸入來檢驗(yàn)反演算法的效果。仿真中采用Pierson-Moskowitz 海浪譜與Longuet-Higgins方向分布函數(shù)之積來表示海浪方向譜，公式如下：

(15)

其中:α=8.1×10-3;β=0.74;kc=g/U2，U為海平面上方19.5m處的風(fēng)速;s為方向分布系數(shù);θw為風(fēng)向與參考方向之間的夾角。仿真中使用的海浪譜與海浪頻譜之間存在如下關(guān)系：

(16)

通過海浪方向譜可以計(jì)算海浪的有效波高[7]：

(17)

其中:f為海浪頻率;θ為海浪方向。

實(shí)驗(yàn)一為信噪比實(shí)驗(yàn)，雷達(dá)頻率為25MHz，風(fēng)向與雷達(dá)法線方向的夾角θw=135°，風(fēng)速U為14m/s，為了驗(yàn)證算法的穩(wěn)定性，在仿真產(chǎn)生的回波譜中增加加性高斯白噪聲，公式如下：

(18)

其中:σt為仿真產(chǎn)生回波譜;P(σt)為σt能量強(qiáng)度;SNR為信噪比;er為服從標(biāo)準(zhǔn)正態(tài)分布的隨機(jī)數(shù)序列。

圖2中黑線為仿真產(chǎn)生的雷達(dá)回波譜，紅色為回波譜中增加噪聲的結(jié)果，當(dāng)信噪比為30dB時，雷達(dá)回波譜中左側(cè)的二階峰幾乎被噪聲淹沒，在此情況下只能選擇能量較強(qiáng)的二階峰來進(jìn)行計(jì)算。圖3為對應(yīng)圖2中雷達(dá)回波譜計(jì)算出的反演結(jié)果，其中黑線為仿真產(chǎn)生的P-M海浪譜，紅線為反演結(jié)果，隨著信噪比的降低，反演結(jié)果的擾動逐漸增大，使得海浪譜中各頻率分量與真值相比產(chǎn)生的偏差逐漸增大，但整個海浪譜能量的分布及主波頻率仍較清晰。

圖1 海浪反演算法流程圖Fig.1 The algorithm flowchart of ocean wave inverse

實(shí)驗(yàn)二為風(fēng)向?qū)嶒?yàn)，雷達(dá)頻率為25MHz，風(fēng)速U為14m/s，風(fēng)向與雷達(dá)法線方向的夾角θw依次為90°和135°，為了防止噪聲將較弱的二階峰淹沒，回波譜仿真的信噪比設(shè)置為30dB，同時，為了驗(yàn)證算法對低信噪比的反演精度，反演中使用的雷達(dá)回波信噪比設(shè)置為15dB。

圖4中a圖為仿真產(chǎn)生的雷達(dá)回波數(shù)據(jù)，其中黑線對應(yīng)的風(fēng)向θw為90°的情況，紅線對應(yīng)風(fēng)向θw為135°的情況。b和c圖采用本文算法計(jì)算得到的海浪譜數(shù)據(jù)，由于算法中首先測量出海浪方向，因此，風(fēng)向的變化對海浪譜的反演產(chǎn)生的影響較小。

實(shí)驗(yàn)三為海況實(shí)驗(yàn)，雷達(dá)頻率為25MHz，風(fēng)向與雷達(dá)法線方向的夾角θw為135°，信噪比為雷達(dá)法線方向的夾角θw為135°，信噪比為15dB，風(fēng)速U為依次為8、14和20m/s。

圖5中，當(dāng)海況較低時，雷達(dá)回波譜中二階峰的幅值較低，受到噪聲的影響較大，使得反演出的結(jié)果有較明顯的抖動，但海浪譜中各頻率能量的分布趨勢仍較為明顯；隨著海況的增加，噪聲對反演結(jié)果的影響會相對減小，由于本文算法基于Barrick后向散射公式，此公式又基于電磁波的微擾展開原理，因此，反演波高的范圍存在限制k0h<4[19]，其中k0為電磁波波數(shù)，h為均方根波高。

圖2 不同信噪比下的仿真回波譜Fig.2 Simulated radar spectra with different SNR

圖3 不同信噪比下的反演海浪譜

圖4 不同風(fēng)向下仿真及反演結(jié)果

圖5 不同海況下反演海浪譜

實(shí)驗(yàn)四為雷達(dá)頻率實(shí)驗(yàn)，風(fēng)向與雷達(dá)法線方向的夾角θw為135°，信噪比為15dB，風(fēng)速U為8m/s，雷達(dá)頻率分別為13、18和25MHz。

圖6中，黑線為海浪譜真值，紅線為反演結(jié)果，為觀察不同頻率雷達(dá)所能反演海浪譜的頻率范圍，真值海浪譜只顯示反演頻帶內(nèi)的幅值。在反演過程中要進(jìn)行線性化處理，因此二階峰的選擇區(qū)域是有限制的，一般選擇的區(qū)域?yàn)?.1～1.4fB[8]，同時對海況的高低也是有限制的,較低的雷達(dá)頻率將不能測量低海況[20]。當(dāng)雷達(dá)頻率降低時，反演海浪譜的頻率上限會逐漸降低，使反演海浪譜的頻率范圍存在截?cái)喱F(xiàn)象。

表1為4組實(shí)驗(yàn)測量的有效波高誤差結(jié)果，使用本算法反演海浪譜信噪比在高于15dB時是比較可靠的，例如，雷達(dá)頻率為25MHz時，有效波高的相對誤差為8.03%。在實(shí)際系統(tǒng)中可以剔除信噪比低于15dB的數(shù)據(jù)，來增加反演海浪的精度。

圖6 不同雷達(dá)頻率下反演海浪譜

實(shí)驗(yàn)一Experiment1實(shí)驗(yàn)二Experiment2實(shí)驗(yàn)三Experiment3實(shí)驗(yàn)四Experiment450dB30dB15dB90°135°8m·s-114m·s-120m·s-113MHz18MHz25MHz有效波高①/m4.184.184.184.184.181.364.188.541.361.361.36絕對誤差②/m0.030.030.020.030.040.110.020.000.370.220.11相對誤差③/%0.640.600.470.730.988.350.400.0427.1116.228.03均方根誤差④/m0.100.110.220.290.240.450.210.561.430.860.43

Note:①Signifigent wave high; ②Absolute error; ③Relative error; ④Root-Mean-Square error

在低海況及低頻段的反演中，海浪譜的反演結(jié)果出現(xiàn)較大的擾動現(xiàn)象，主要是由兩方面原因造成的：第一，較低的海況激起的二階峰幅值太小，反演的結(jié)果受噪聲的影響大，表現(xiàn)為海浪譜的抖動；第二，在海浪譜的高頻段區(qū)域，由于雷達(dá)截面方程中等頻線未能完全覆蓋，導(dǎo)致反演時信息量小，產(chǎn)生抖動現(xiàn)象[12]。對于第一個原因，可以通過增加雷達(dá)頻率來解決；對于第二個原因，可以使用模型或f-5衰減因子來近似擬合海浪譜的抖動及截?cái)嗖糠帧?/p>

4 實(shí)驗(yàn)驗(yàn)證初步結(jié)果

為進(jìn)一步驗(yàn)證算法的有效性，將反演算法應(yīng)用于實(shí)測數(shù)據(jù)，雷達(dá)頻率為13.021MHz。圖7為實(shí)測某波束的距離-多普勒回波譜，其中橫坐標(biāo)為頻率，縱坐標(biāo)為距離單元，正負(fù)一階峰均向右產(chǎn)生偏移，此偏移由海流產(chǎn)生，正一階峰右側(cè)存在明顯的二階峰，雖然受海流的影響產(chǎn)生了嚴(yán)重的頻偏，但可以直接通過平移整個頻譜來移除海流的影響[21]。

雷達(dá)回波數(shù)據(jù)為多普勒頻率、距離和波束的三維數(shù)據(jù)，對于某個波束和距離的雷達(dá)回波譜運(yùn)用反演算法，便可以得到此波束和距離的海浪譜。在此，提取圖7中的第9距離單元的回波頻譜(見圖8)。本算法自動選擇一階峰與二階峰之間的低谷作為二階峰范圍的起始點(diǎn)，1.4fB作為終止點(diǎn)。

圖9為反演得到的海浪譜結(jié)果，有效波高為2.25m，與當(dāng)時國家海洋局預(yù)報(bào)的有效波高2m吻合，頻率范圍為0.08～0.24Hz，此頻率范圍已包含了大部分海浪的能量。

圖7 實(shí)測距離-多普勒回波譜

圖8 雷達(dá)多普勒回波譜

圖9 反演海浪譜

圖10為2013年1月9日下午采集的高頻雷達(dá)數(shù)據(jù)采用不同算法得到的有效波高測量結(jié)果，時長3.5h，其中黑線為未經(jīng)過正則化處理結(jié)果，紅線為采用正則化處理結(jié)果。很明顯黑線有3處較大的突變，是雷達(dá)回波中噪聲引起的求逆結(jié)果不穩(wěn)定造成的，而紅線更加平滑穩(wěn)定，符合海浪有效波高的變化規(guī)律。

圖10 M-P逆與正則化處理實(shí)測有效波高結(jié)果比較

5 結(jié)語

本文提出的一種低復(fù)雜度的高頻雷達(dá)海浪反演算法，通過仿真驗(yàn)證，在信噪比高于15dB時，能較準(zhǔn)確的反演海浪譜，并且本算法復(fù)雜度較低，反演的結(jié)果穩(wěn)定，如果只用來反演有效波高，則所需的信噪比還可降低。方法中海浪譜的方向分布系數(shù)s采用了經(jīng)驗(yàn)值，在后續(xù)工作中可通過回波譜數(shù)據(jù)計(jì)算得到s的具體數(shù)值，來提高反演的精度。

[1] Hasselmann K. Determination of ocean wave spectra from Doppler radio return from the sea surface [J]. Nature, 1971, 229:16-17.

[2] Barrick D E. First-order theory and analysis of MF/HF/VHF scatter from the sea [J]. IEEE Trans Antennas Propag， 1972, AP-20: 2-10.

[3] Barrick D E. Remote Sensing of Sea State by Radar [C].// Derr V E. Chapter 12 of Remote Sensing of the Troposphere, NOAA/Environmental Research Laboratories Boulder: 1972: 1-6.

[4] Barrick D E, Weber B L. On the nonlinear theory for gravity waves on the ocean’s surface. Part II: Interpretation and applications [J]. J Phys Oceanogr, 1977, 7: 11-21.

[5] Lipa B J, Barrick D E. Extraction of sea state from HF radar sea echo: Mathematical theory and modeling [J]. Radio Sci, 1986, 21: 81-100.

[6] Wyatt L R. HF radar measurements of the ocean wave directional spectrum [J]. IEEE J Oceanic Eng, 1991, 16: 163-169.

[7] Essen H H, Gurgel K W, Schlick T. Measurement of ocean wave height and direction by means of HF radar: An empirical approach [J]. Dt Hydrogr Z, 1999, 51: 369-383.

[8] Howell R, Walsh J. Measurement of ocean wave spectra using narrow beam HF radar [J]. IEEE J Oceanic Eng, 1993, 18: 296-305.

[9] WU X B, LI L, SHAO Y X, et al. Experimental Determination of Significant Waveheight by OSMAR071: Comparison with Results from Buoy [J]. Wuhan University Journal of Natural Sciences, 2009, 14( 6): 499-504.

[10] Delves L M, Mohamed J L. Computational Method for Integral Equations [M]. Cambridge: Cambridge University Press, 1985: 376.

[11] Hadamard J. Lectures on Cauchy’s Problem in Linear Partial Differential Equations [M]. New Haven: Yale University Press, 1923.

[12] Li L, Wu X B, Long C, et al. Regularization inversion method for extracting ocean wave spectra from HFSWR sea echo [J]. Chinese J Geophys (in Chinese), 2013, 56(1): 219-229.

[13] Moskowitz M. Estimates of the power spectrums of fully developed seas for wind speeds of 20 to 40knots [J]. J Geophys Res, 1964, 69: 5161-5179.

[14] Huang W, Gill E. Extraction of Sea Surface Wind Direction from Bistatic High-Frequency Radar Doppler Spectra [C]. Oceans’11, Hawail： IEEE Conference and Hawaii Exhibition: 2011.

[15] Tikhonov A N， Arsenin V Y. Solution of Ill-Posed problem [M]. Washington D. C: Winston & Sons, 1977.

[16] Tikhonov A N, Goncharsky A V. Ill-Posed Problems in the Natural Sciences [M]. Moscow: MIR Publishers, 1987.

[17] Golub G, Von Matt U. Generalized cross-validation for large-scale problems [J]. Journal of Computational and Graphical Statistics, 1997, 6(1): 1-34.

[18] Lipa B J, Barrick D E. Analysis Methods for Narrow-Beam High Frequency Radar Sea Echo [R]. NOAA Technical Report ERL 420-WPL, 1982, 56.

[19] Wyatt L R. An evaluation of wave parameters measured using a single HF radar system [J]. Can J Remote Sens, 2002, 28( 2): 205-218.

[20] Wyatt L R, Green J J. Measuring high and low waves with HF radar [C]. Bramen： OCEANS 2009- EUROPE, 2009: 1- 5.

[21] Vincent C E. The interaction of wind-generated sea waves with tidal currents [J]. J Phys Oceanogr, 1979, 9: 748-755.

責(zé)任編輯 陳呈超

Ocean Wave Inversion Algorithm From HF Radar Echo Spectrum

CHANG Guang-Hong1, LI Ming1, CHU Xiao-Liang1,2, JI Yong-Gang2, YU Chang-Jun3

(1. Ocean University of China, Qingdao 266100, China; 2.The First Institute of Oceanography, SOA, Qingdao 266061, China;3.Harbin Institute of Technology at Weihai, Weihai 264209, China)

In the classification of ocean wave spectrum inversions from HF radar echo, widely used Barrick′s backscatter formula belongs to the nonlinear first kind Fredholm integral equation. However, this kind of integral equation is Ill-posed, moreover, the SNR of HF radar echo is relatively low, the inversion ocean wave spectrum is unstable. In this paper, a stable and low complexity inverse algorithm is presented, which uses the first-order peak to measure the wave direction, and then introduces the direction to the progress of solving the integral function, which deduced the number of variables. In order to solve the instability problem of inversion results, we use Thikhonov Regularization methods, and we use Generalized Cross Validation method to determine the regularization coefficients. Simulation results show that this method is stable and small calculating amount.

High-Frequency Surface Wave Radar (HFSWR); ocean wave spectrum; inversion; regularization

國家自然科學(xué)基金重點(diǎn)項(xiàng)目(61032011、61132005)資助

2013-10-13；

2014-05-10

常廣弘(1990-)，男，碩士生。E-mail:changguanghong@126.com

?? 通訊作者： E-mail：limingneu@ouc.edu.cn

TN958

A

1672-5174(2015)02-127-07

10.16441/j.cnki.hdxb.20130262