基于船舶運動的微波多普勒海面回波模擬與分析

李漢東, 陳澤宗,2*, 趙 晨, 陳 曦, 陳小莉

(1.武漢大學電子信息學院， 武漢 430072； 2.武漢大學地球空間信息技術協同創新中心， 武漢 430079)

船載微波多普勒雷達從散射回波中獲取海面動態信息，實現了走航觀測海浪的新手段[1]。目前對海面散射回波的研究大多建立在靜態雷達平臺上，而船載雷達平臺下天線接收的回波信號中不僅包含海態信息，也包含了平臺運動信息，改變了海面電磁散射的計算模型，給海面回波的研究帶來了新的挑戰。

目前海面電磁散射的計算方法主要分為數值法和近似法[2]。矩量法等數值法的計算精度高，但存在計算量大的缺點；小斜率近似法[3]是一種廣泛應用的高精度近似方法，通過保留表面斜率級數的各階展開項獲取散射振幅的近似值，適用于復合尺度的海面電磁散射問題，具有計算簡便、適用性廣且精度較高等優點。

船體在航行中容易受到海浪和海流等外部環境影響而產生橫搖、縱搖等六自由度運動，導致回波多普勒譜頻移和譜寬發生改變，進而影響到海態參數的反演。近年來，中外學者提出了一系列實時或非實時的船舶運動預測算法。Fossen等[4]提出了動態船舶運動預測的條帶理論公式；Ueng等[5]提出了實時預測的船舶運動仿真有效計算模型；Sandaruwan等[6]提出了高性能船舶姿態預測算法。

現基于船載微波多普勒雷達系統，根據船載雷達探測海面時的工作特點，研究走航狀態下微波電磁波與海表面的作用機理。通過二維時變海面建模與船舶運動建模，采用小斜率近似法計算海面電磁散射幅度，仿真得到船載雷達平臺下的海面散射回波，并分析不同運動方式(勻速運動、變速運動、搖晃運動)對多普勒譜特征的影響。根據空間海浪反演算法計算得到海浪譜并與理論值對比分析，驗證該海面回波模擬方法的正確性。

1 時變海面建模

有向海浪譜S(ω,θ)描述了海面組成諧波分量在頻域和方向上的分布，可以表示為無向海浪譜S(ω)與方向分布函數G(ω,θ)的乘積。采用Pierson-Moscowitz譜和Longuet-Higgins方向函數的表達式為

(1)

(2)

S(ω,θ)=S(ω)G(ω,θ)

(3)

式中：S(ω)為無向海浪普；α為無因次常量，α=0.008 1；β=0.74；g為重力加速度；U19.5表示海面上方19.5 m處的平均風速；ω為海浪角頻率；G(ω,θ)為方向分布函數；S(ω,θ)為有向海浪譜；s為角散系數；θ為方向角；θ0為主效浪方向；Γ(*)為伽瑪函數。

根據線性波浪理論，時變海面的波面位移能夠表示為1組不同振幅、頻率和相位的余弦波的疊加[7]。采用線性疊加法計算得到波面位移h(x,y,t)為

(4)

式(4)中：(x,y)表示海面離散點坐標；t為時間；M、N分別為頻率和方向維度的離散點數；ωi和ki分別為角頻率與波數；θj為方向角；εij為均勻分布在(0, 2π)間的隨機相位。

2 船舶運動建模

2.1 船舶運動模型

當船舶在海面上航行時，建立圖1所示的船體坐標系G-XsYsZs。由于受到海洋環境的影響，船體會沿三坐標軸分別產生平行運動與旋轉運動，統稱為六自由度運動。

圖1 船體六自由度運動Fig.1 Six degrees of freedom motion of hull

當船舶沿固定航向在海面上航行時，海浪作用于船體使姿態角改變。由于橫蕩、縱蕩和垂蕩等平行運動不影響搖晃速度，在航向固定時艏搖對船體影響較小，因此重點研究橫搖和縱搖運動。根據船舶水動力理論，該類運動可看作是受遭遇海浪激勵與流體阻力影響的阻尼運動，在平衡位置上下不斷振動，因此會在天線上疊加時變的搖晃速度。根據阿基米德定律，相對于靜止海面的波面凸起會對船底接觸面施加等效于相同體積海水質量的轉動作用力F，引起船體姿態角改變(圖2)。

圖2 海浪引起船體搖晃運動Fig.2 Shipshaking caused by waves

當船舶沿固定航向行駛時，將船體垂直投影到二維海平面上，并劃分為1 m×1 m的網格，如圖3所示。將船體形狀簡化為長(L)、寬(W)、高(D)的立方體，則投影邊界框內近似包含L×W個采樣點。根據時變海面位移計算每個采樣點當前時刻的波高(hi，j)與轉動作用力，減去海水流體阻力并求和得到瞬時合力，根據流體動力學公式計算得到船體搖晃運動加速度，進而得到搖晃運動速度與瞬時姿態角。

圖3 船體搖晃計算示意圖Fig.3 Schematic diagram for calculation of ship shaking

由圖3可知，搖晃矢量rm為轉動徑向矢量在平行于船艏(船舷)方向的分量，搖晃力矩Mm為搖晃矢量與轉動作用力的矢量積，表達式為

(5)

式(5)中：ρsea為海水密度；Aunit為網格面積;n為方向矢量。搖晃阻力(R)與作用力(F)表達式分別為

R=bIωm

(6)

(7)

式中：b為阻尼系數；ωm為搖晃角速度；K為取決于船舶類型的搖晃系數；轉動慣量(I)表達式為

(8)

根據流體動力學公式計算得

am=F/I,ωm=ωm+amΔt,θm=ωmΔt

(9)

式(9)中：αm為搖晃角加速度；θm為姿態角；Δt為時間步長。

2.2 船舶運動仿真結果分析

以2018年9月南海海試實驗為例，船載微波多普勒雷達被安裝在“實驗1號”科考船的船頭[8-9]。設定船舶運動仿真參數，海面10 m風速分別為8、10 m/s，計算得到3 min內5.14 m/s航速下的搖晃運動速度，如圖4所示。

由圖4可知，海況越高時，搖晃運動幅度越劇烈，角速度為1～2(°)/s，周期為8～12 s。

圖5給出了2018年9月11日一段3 min姿態傳感器實測數據，期間海面10 m風速為8 m/s。由于傳感器安裝位置存在傾斜角，導致橫縱搖角存在固定偏移量。傳感器測得的橫縱搖角幅度變化范圍同樣在1～2(°)/s以內，周期為8～12 s，與仿真結果基本一致，驗證了該船舶運動仿真算法的有效性。

3 海面回波模擬算法

小斜率近似法綜合了基爾霍夫近似法與微擾法，在掠入射角下也具有較高的精度，適用于復合尺度的海面電磁散射問題。當微波段電磁波照射在二維海面h(r,t)=h(x,y,t)上時，海面電磁散射幾何模型如圖6所示，其中Ki為入射波矢量，Ks為散射波矢量，θi、φi、θs、φs分別為入射角、入射方位角、散射角和散射方位角。

根據圖6中的幾何關系可知：

Ki=k0-q0z，Ks=k1+q1z

(10)

式(10)中：k0和q0分別是Ki投影在海表面與z軸上的分量；k1和q1分別是Ks投影在海表面與z軸上的分量。根據一階小斜率近似法計算得到散射振幅表達式為

圖4 縱搖、橫搖運動仿真結果Fig 4 Simulation results of pitch and roll

圖5 實測姿態角數據Fig.5 Measured attitude angle data

圖6 海面電磁散射幾何模型Fig.6 Geometrical model of EM scattering from sea surface

(11)

G(r,h)exp[-j(k1-k0)r+j(q0+q1)h(r)]

(12)

式中：Pinc為入射波功率；G(r,h)表示消除截斷誤差的錐形入射波[2]。將海水視為理想導體，計算得到不同極化方式下的散射系數矩陣B(k0,k1)。其中BHH代表水平極化，BVV代表垂直極化。

(13)

(14)

假設船載雷達平臺以徑向速度(v)向斜距為Ri的目標距離元靠近，T為海面采樣時間間隔，則第i次采樣時回波信號接收時延(τi)表達式為

(15)

式(15)中：c為電磁波傳播速度。f0為入射波頻率。則船載雷達天線接收回波的散射振幅Aship為

Aship(k0,k1,t)=S(k0,k1,t)exp(-j2πf0τ)

(16)

式(16)中：根據散射振幅可推導得到多普勒譜Sdop為

(17)

采用能量均衡矩估計法計算多普勒譜的頻移fd和譜寬Bw，其中f表示頻率。

(18)

(19)

4 回波多普勒譜仿真與分析

船載雷達平臺下不同時刻天線照射的海面區域不同，需要根據實際運動情況選取每一時刻對應的時變海面區域，并計算其電磁散射幅度。

設定入射電磁波頻率為2.85 GHz，海面10 m風速為10 m/s，雷達采樣時間間隔設定為0.003 9 s，共采樣128點，計算0.5 s內時變海面的電磁散射幅度和多普勒譜。為了減小海面隨機性的影響，方便多普勒譜的觀察與分析，仿真多普勒譜取100個隨機海面樣本求取均值。

圖7(a)為在勻速運動情況下船速分別為0、2.57、5.14 m/s時的海面回波多普勒譜，圖7(b)為多普勒譜頻移和譜寬隨船速變化。由圖7可知，隨著船速的增大，多普勒譜頻移呈線性增長趨勢，譜峰值逐漸減小，譜寬逐漸變大，且船速越大，譜寬增長越快。當船速越大時，采樣時間間隔內移動的距離越大，參與計算的相鄰海面回波的相關性越弱，導致多普勒譜峰值降低與譜寬展寬。

圖7 不同船速的海面回波仿真(勻速)Fig 7 Sea echo simulation with different ship speeds (constant speed)

在實際的走航觀測中，船舶的航行并不總是理想的勻速運動，船速可能會在相干積累時間內發生變化。假設雷達平臺的徑向速度v(t)是如式(20)所示的時變量。

v(t)=v0+at

(20)

式(20)中：v0為初始速度；a為恒定加速度。圖8(a)給出了初始速度是2.57 m/s，加速度分別是0、1、2 m/s2的海面回波多普勒譜，圖8(b)為多普勒譜頻移和譜寬隨船加速度變化。

由圖8可知，隨著船加速度的增大，多普勒譜頻移呈線性增長趨勢。相比勻速運動的情況，多普勒譜譜型對船速的變化十分敏感。當船速改變時，譜峰值迅速減小，譜寬向速度變化側嚴重展寬，且加速度越大，譜峰值越小，譜寬增長越快。在存在初始船速的情況下，正加速度的譜寬變化曲線斜率更大，表明相同加速度的勻加速運動相比勻減速運動對譜寬的影響更明顯。

圖9分別為未添加搖晃速度、添加了海面10 m風速為8、10 m/s時搖晃速度的海面回波多普勒譜。

計算得到頻移和譜寬如表1所示。

由表1可知，當雷達天線添加搖晃速度時，頻移會發生改變，多普勒譜峰值減小，譜寬展寬，且海況越大時，搖晃速度越大，頻移變化量越大，譜寬展寬程度越劇烈。由于不同時刻的船體搖晃速度不同，頻移和譜寬變化量也會隨時間改變。

圖8 不同船加速度的海面回波仿真(變速)Fig.8 Sea echo simulation with different ship accelerations (variable speed)

表1 多普勒譜參數

5 基于海浪反演的模擬算法驗證

根據二維時變海面不同距離元的相位差，擴展得到探測范圍300～900 m共80 個距離元的距離多普勒譜，提取多普勒譜中心頻率并轉換為多距離元徑向速度序列。根據空間采樣的海浪反演算法，推導得到海浪譜，并反演海浪參數。圖10給出了不同均勻船速下反演海浪譜與理論靶譜的對比。

采用矩估計法對反演海浪譜作誤差分析，如表2所示。

表2 海浪譜誤差分析(勻速)

當均勻船速增大時，有效浪高反演結果誤差增大，海浪譜相關系數減小。由圖10可知，不同船速下的反演海浪譜均在理論值附近上下波動，整體趨勢和理論靶譜高度重合，誤差在可接受的范圍內，驗證了該海面回波模擬算法的正確性。

圖10 不同船速下的反演海浪譜Fig.10 Inversion wave spectrum with different ship speeds

圖11給出了初始速度為2.57 m/s、不同船加速度下反演海浪譜與理論靶譜的對比。

圖11 不同船加速度下的反演海浪譜Fig.11 Inversion wave spectrum of different ship acceleration

同樣采用矩估計法對反演海浪譜作誤差分析，如表3所示。

表3 海浪譜誤差分析(變速)

根據相干多普勒雷達原理，當船加速度增大時，船速在采樣時間內發生改變，天線接收的海面回波信號的相關性被削弱，導致多普勒譜峰值驟降與譜寬嚴重展寬，提取的徑向速度序列起伏更劇烈，反演海浪譜能量增加。由表3可知，有效浪高反演誤差明顯增大，海浪譜相關系數減小。

6 結論

基于時變海面建立了船舶運動模型，并通過與實測數據對比驗證了模型的有效性。采用一階小斜率近似法計算海面電磁散射幅度，仿真得到了船載雷達平臺下的海面散射回波，并分析了不同運動方式(勻速運動、變速運動、搖晃運動)對多普勒譜特征的影響。結果表明，船舶運動會導致多普勒譜譜型改變，船速和船加速度越大，頻移和譜寬變化越明顯。根據空間海浪反演算法計算得到海浪譜并與理論值對比分析，驗證了該海面回波模擬算法的正確性，研究結果對船載雷達平臺下海浪反演與算法改進具有重要的指導意義。