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

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

0 引言

海洋風場是海洋表層運動的主要推力，它管控能源和動力的轉(zhuǎn)移，能夠?qū)⒖諝鈳У胶Ｑ?，在氣候變化中占有非常重要的作用。海面風場常被應(yīng)用于各類海洋科學研究，如海洋環(huán)流、熱帶風暴、大氣與海洋的耦合效應(yīng)等[1]。因此，監(jiān)測海面風場變化有著極其重要的意義。

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

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

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

1 海面風場反演方法

1.1 數(shù)據(jù)來源及設(shè)備參數(shù)

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

表1 育鯤輪航海雷達設(shè)備參數(shù)

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

圖1 雷達圖像原圖

1.2 實驗流程

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

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

1.3 圖像預(yù)處理

1.3.1 圖像求和

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

(1)

1.3.2 中值濾波法

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

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

(2)

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

1.3.3 梯度法反演風向

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

(1)圖像重采樣

對雷達原始圖像采用含有一定濾波作用的算子進行計算，平滑圖像的同時，壓縮圖像像素數(shù)量?？捎盟阕庸饺缦拢?/p>

R=B2S2B4，

(3)

(4)

(5)

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

(2)局部梯度的計算

用最優(yōu)化的Sobel算子進行局部梯度計算：

(6)

Dy=DxT，

(7)

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

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

(8)

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

G′=R*(G2)，

(9)

G″=R*(|G2|)。

(10)

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

(11)

(3)確定風向

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

2 實驗結(jié)果與分析

2.1 雷達原始圖像求和處理結(jié)果

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

圖3 雷達原始圖像求和

2.2 中值濾波處理結(jié)果

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

圖4 中值濾波處理

2.3 局部梯度法反演風向

采用Matlab編寫程序?qū)D像進行梯度計算，并做出相應(yīng)的梯度加權(quán)直方圖，如圖5所示。由圖5可以看出，在梯度角度為88°時出現(xiàn)的頻率最大。因為主風向與梯度角度直方圖的峰值方向垂直，應(yīng)當在獲得的梯度角度加上90°，即為所求的主風向178°。

圖5 梯度加權(quán)直方圖

3 驗證

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

4 結(jié)束語

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