基于X波段航海雷達的海面風場反演技術研究

王海霞，徐 進，王慶名，陳 麗，龐璽斌

(1.大連海事大學 航海學院，遼寧 大連 116026；2.大連海事大學 校友事務與合作處，遼寧 大連 116026)

0 引言

海洋風場是海洋表層運動的主要推力，它管控能源和動力的轉移，能夠將空氣帶到海洋，在氣候變化中占有非常重要的作用。海面風場常被應用于各類海洋科學研究，如海洋環流、熱帶風暴、大氣與海洋的耦合效應等[1]。因此，監測海面風場變化有著極其重要的意義。

目前，國外在航海雷達海上風場的反演方面已經存在一些研究成果。Ziemer等人[2]首次應用二維傅里葉變換，分析2張帶有海雜波的雷達圖像，計算的海浪譜與傳統浮標獲取的海浪譜值較為接近。Alpers等[3]和Vesecky等[4]對雷達海浪成像作了詳細描述，但他們均未考慮流的影響。Lee等人[5]基于電磁波海面散射機理，指出海面的粗糙度與風速之間，呈現正相關特性。Hatten等人[6]發現航海雷達散射截面風向及電磁波照射方向存在一定關系，運用這種關系可以判斷出風向，但此方法存在一定的局限性，往往不能滿足所需的全方位航海雷達散射截面。Dankert等人[7]提出了一種基于光流運動估計的技術，從X波段航海雷達圖像中獲得風速和風向，該方法相對簡單，不需對雷達系統的相位進行校正，也不存在角度模糊問題，所獲得的海面風場具有較高的時空分辨率。Liu等人[8]針對未被遮擋的低散射強度雷達圖像，提出了一種散射強度值選取法，并用對偶曲線擬合法減小了風向反演誤差。

國內應用雷達監測海面風場研究的起步較晚，楊勁松[9]用合成孔徑雷達(SAR)，采用二維傅里葉變換波數譜的方法，反演出了風向，并用COMD4模式函數算出了風速。金麗潔[10]利用浪高反演算法[11]，采用多重信號分類算法完成定向，實現了浪高場的反演。張璇等人[12]分析目前交叉極化對風場反演的優勢，總結了多極化SAR圖像反演海面風場的最新研究成果，并對多極化SAR圖像反演海面風場的業務應用發展進行了初步展望。

本文基于X波段航海雷達反演海面風場的基本原理，以采集到的雷達原始圖像為數據基礎，基于梯度算法求得主風向，并與大連氣象臺的實測值進行對比，以期為海面風場的實時監測提供有效的技術手段。

1 海面風場反演方法

1.1 數據來源及設備參數

X波段船載雷達信號采集到的原始灰度圖像(設備參數詳見表1)監測范圍為0.75海里，圖像大小為1 024×1 024(相同地點連續獲取的3幅圖像如圖1所示)。監測數據采集位置：大連港附近海域(121.830°E，38.913°N)。船載平臺：大連海事大學“育鯤”輪號，采集時間為2015年8月14日20時4分。數據中含有大量的含雜波信息，用于反演海面海浪與風向等信息的研究。

表1 育鯤輪航海雷達設備參數

參數項參數值工作極化方式水平工作頻段X波段雷達天線長度/ft8天線旋轉速度/RPM28脈沖工作頻率/Hz3 000/1 800 /785雷達控制臺23Inch液晶工控機工作量程/海里1～12

圖1 雷達圖像原圖

1.2 實驗流程

海面風場反演過程包括雷達信號的實時采集和海面風場信息提取2部分。航海雷達風場反演流程如圖2所示。

圖2 航海雷達風場反演流程

1.3 圖像預處理

1.3.1 圖像求和

圖像求和即將圖像對應像素點的值進行求和，從而增強海浪信息所在位置的灰度強度。n幅圖像中，每個像素點f(x，y)求和公式為

(1)

1.3.2 中值濾波法

中值濾波基于排序統計理論，將每一像素的灰度值設置為該點某鄰域窗口內的所有像素點灰度值的中值，公式為

g(x，y)=median({f(x-i，y-j)，(i，j)∈W)}，

(2)

式中，g(x，y)為輸出像素值；f(x-i，y-i)表示輸入像素值；W是含有基數點的模板窗口，一般為3×3窗口或5×5窗口。

1.3.3 梯度法反演風向

風向反演主要是基于海面風引起的紋理可在雷達回波圖像中成像，Dankertd等人[13]研究發現雷達圖像中存在與風向平行的紋理信息，分析雷達圖像中的海浪紋理可以獲取風向[14]。梯度是用于描述圖像灰度變化情況的二維列向量，能反映變化的劇烈程度[15]。通過將圖像進行分割并統計梯度方向與大小得到該區域的梯度特征分布的統計結果[16]，梯度的角度為圖像像素變化的方向，梯度的模為梯度方向的變化率。若求出海浪紋理變化最大的梯度方向，即可得到與之垂直的主風向。算法實施步驟如下。

(1)圖像重采樣

對雷達原始圖像采用含有一定濾波作用的算子進行計算，平滑圖像的同時，壓縮圖像像素數量。可用算子公式如下：

R=B2S2B4，

(3)

(4)

(5)

式中，S2表示2×2線性插值平均，先用四階算子平滑，然后縮減圖像，最后用二階算子再次平滑圖像。

(2)局部梯度的計算

用最優化的Sobel算子進行局部梯度計算：

(6)

Dy=DxT，

(7)

式中，T表示轉置矩陣。此時圖像的梯度G為

G=(Dx+iDy)*(A)，

(8)

式中，A表示圖像；*表示卷積，G′與G″的計算方式如下所示：

G′=R*(G2)，

(9)

G″=R*(|G2|)。

(10)

由文獻[17]可知，一致性參數C為：

(11)

(3)確定風向

應用一致性參數和梯度的模，做梯度方向角度的加權直方圖。因為風向與梯度變化最大的方向垂直，在直方圖中最大尖峰值加上90°，即為主風向。

2 實驗結果與分析

2.1 雷達原始圖像求和處理結果

在獲取的雷達圖像中包含許多信息，但在獲取過程中由于各種原因會存在不同形式的干擾。因此在反演前應當先進行圖像處理。首先將連續獲取的雷達原始圖像進行圖像求和，增強海浪信息所在位置的灰度強度，如圖3所示。

圖3 雷達原始圖像求和

2.2 中值濾波處理結果

采用5×5窗口，對雷達原始圖像求和處理后的圖像進行中值濾波處理，降低同頻干擾，如圖4所示。在進行了預處理后，為消除圖像存在的尾跡效應，截取了沒有尾跡效應的部分截取了沒有尾跡效應的部分，即(275：402，445：572)中的圖像窗口進行風向反演。

圖4 中值濾波處理

2.3 局部梯度法反演風向

采用Matlab編寫程序對圖像進行梯度計算，并做出相應的梯度加權直方圖，如圖5所示。由圖5可以看出，在梯度角度為88°時出現的頻率最大。因為主風向與梯度角度直方圖的峰值方向垂直，應當在獲得的梯度角度加上90°，即為所求的主風向178°。

圖5 梯度加權直方圖

3 驗證

根據大連氣象臺提供的2015年8月14日平均風向為183°，與實測結果相比風向偏差5°，反演的風向與預報結果是比較吻合的。因本文所用的雷達數據極化方式為水平極化，可以很好地消除180°模糊問題，因此，反演結果不存在180°模糊。

4 結束語

以大連海事大學教學實習船“育鯤”輪采集的航海雷達原始數據為基礎，基于梯度算法進行風場反演。在實驗過程中，存在多種因素可能會影響計算的結果。首先，在雷達成像過程中，獲取到的風場信息包含了雷達重力波和毛細波的影響。此外，在處理同頻干擾和尾跡效應問題時，雖然采用了中值濾波，但干擾并未完全消除。本文提出的方法還未考慮其他因素，如沒有考慮到船舶運動對于雷達回波的影響；遇到遮擋物時，在遮擋物后方向形成盲區，不能準確顯示出該區域的真實海面狀況。這些對于風場反演的結果都將存在影響，使結果產生誤差。而且，計算結果為風場的主風向，并不是該區域內的所有風向，各區域點之間也會存在差異。上述難點都需要在未來的工作中，繼續研究。