基于SAR圖像速度聚束調制的海浪反演研究

許蕎暉，張彥敏*，王運華,2

( 1. 中國海洋大學 信息科學與工程學院，山東 青島 266100；2. 青島海洋科學與技術試點國家實驗室 區域海洋動力學與數值模擬功能實驗室，山東 青島 266237)

1 引言

海浪是海洋動力學的重要組成部分，海浪信息對于海洋工程、近海結構設計和航海以及理解和預測惡劣的海洋天氣都是至關重要的。由于合成孔徑雷達（SAR）成像過程不受天氣和光照的影響，基于SAR圖像的海浪反演技術已成為全天時、全天候、大面積海浪觀測的主要途徑之一。合成孔徑雷達的海浪成像機制主要包括3種調制作用：傾斜調制、流體力學調制以及速度聚束調制[1-3]。其中，傾斜調制是由于大尺度海浪斜率改變了雷達局地入射角度，從而引起了回波強度的變化；流體力學調制則是由于大尺度海浪軌道速度所引起了Bragg共振波振幅的變化，進而導致雷達回波強度的變化；速度聚束調制則是SAR所特有的，是由于海浪的雷達視向速度導致海面散射面元沿SAR圖像方位向發生位置偏移，進而引起回波強度沿著雷達方位向產生輻聚輻散效應[4-8]。以上3種調制是SAR圖像海浪譜反演算法的基礎，1991年，Hasselmann和Hasselmann[9]基于傾斜調制、流體力學調制和速度聚束調制，推導給出了SAR圖像譜和海浪譜之間的非線性映射關系，并基于初猜譜與SAR圖像譜構建價值函數，通過不斷迭代計算使價值函數最小，從而獲取最優海浪方向譜。1994年，Brüning等[10]為提高效率，對迭代求逆過程進行了改進。1996年，Hasselmann等[11]進一步改進了成像過程中非線性映射關系。2000年，Mastenbroek和De Valk[12]提出了基于Hasselmann所得非線性映射關系的半參數化海浪譜反演方法。2005年，Schulz-Stellenfleth等[13]則在MPI方法基礎上，應用交叉譜提出了PARSA （Partition Rescaling and Shift Algorithm）海浪譜反演算法，該算法可有效解決海浪譜180°模糊問題。2010年，Zhang等[14]在He等[15-16]提出的海浪斜率譜極化調制算法基礎上，利用全極化的Radarsat-2數據分別反演了海浪沿SAR圖像方位向和距離向的海浪斜率譜，進而求取了海浪參數，該方法有效消除了流體力學調制的影響，但是，當入射角較小時，水平極化圖像和垂直極化圖像之間的差異較小，從而導致新構建圖像中的海浪紋理變弱，從而不利于海浪反演。

近年來，基于SAR圖像回波散射系數及截斷波長等參數的海浪有效波高、主波波向等的經驗化反演算法也引起了廣泛關注。2007年，Schulz-Stellenfleth 等[17]提出了CWAVE經驗化算法，由于這種方法主要是針對ERS SAR數據開發的，所以也稱作CWAVEERS算法，該算法可直接反演得到海浪有效波高和平均波周期。在CWAVE基礎上，后來不同學者又相繼提出了針對其他衛星SAR數據的經驗算法，例如適用于ENVISAT ASAR數據的CWAVE_ENV算法[18]，適用于Sentinel-1 SAR數據的CWAVE_S1A算法[19]和適用于Terra SAR-X數據的XWAVE算法[20]等。2015年，Romeiser等[21]則針對颶風極端海況條件下，建立了一種海浪有效波高與NRCS之間的經驗函數來反演有效波高。2016年，Shao等[22]和Grieco等[23]則分別提出了關于有效波高和截止波長[24-25]、雷達入射角、海浪傳播方向之間的半經驗化關系。需要說明的是：基于海浪參數與SAR參數之間經驗關系的反演算法，并不能得到海浪譜。

本文首先對比分析了傾斜調制、流體力學調制和速度聚束調制對SAR海浪圖像的影響，結果表明，當海浪傳播方向在偏離SAR距離向的過程中，3種調制因素中速度聚束調制對海浪SAR圖像的影響始終最為顯著。基于此，文中建立了基于速度聚束調制海浪方位向斜率反演算法，進而基于經驗關系提取海浪有效波高。通過將Radarsat-2數據海浪反演結果與浮標數據相比較，本文反演算法較MPI算法和同極化調制函數算法反演結果的誤差更小。

2 影響SAR圖像的調制因素的比較

2.1 調制因素影響的理論分析

在線性調制理論的框架下，海面高度 ζ(r,t)和后向散射截面 σ(r,t)隨大尺度海浪的變化均可表示為線性疊加的形式[9]，即

對于垂直極化（VV），傾斜調制函數為

式中，θ為雷達入射角；kl是雷達視向上的大尺度海浪波數分量。流體力學調制傳遞函數可根據短波與長波相互作用的雙尺度模型推導得出[28]

式中，μ為阻尼系數；Yr+iYi為 負反饋因子；kx是海浪波數距離向分量。

對于實孔徑雷達，定義歸一化圖像強度為

于是結合式（2）和式（8）可得IR(r)圖像的傅里葉變換系數為

進而可得，IR(r) 圖 像譜與海浪譜Fk的關系為[9]

然而，對于海浪SAR圖像而言，海浪運動所引起的速度聚束效應非常顯著，在線性條件下，考慮到速度聚束效應的影響后海浪歸一化SAR圖像的傅里葉變換系數為[9]

式中，β=R/V，其中R為斜距；v為SAR平臺飛行速度；ky是大尺度海浪波數沿雷達方位向的分量。

由于SAR圖像同時存在傾斜調制、流體力學調制、速度聚束調制3種調制因素，基于式（11），可將海浪歸一化SAR圖像的傅里葉變換系數分解為為僅考慮實孔徑雷達調制的傅里葉變換系數為僅考慮速度聚束調制的傅里葉變換系數，二者分別為

圖1 仿真海面SAR圖像（a），僅考慮速度聚束調制的海面SAR圖像（b）和僅考慮傾斜調制和流體力學調制的海面SAR圖像（c）Fig. 1 Simulated SAR image of sea surface (a), simulated SAR image of sea surface with the velocity bunching modulation (b) and simulated SAR image of sea surface with the tilt modulation and the hydrodynamic modulation (c)

圖2 a至圖2c則給出了圖1中3幅海面SAR圖像所對應的圖像譜。通過圖像譜可以看出，僅考慮速度聚束效應影響的圖像譜（圖2b）與考慮3種調制的海面仿真SAR圖像譜（圖2a）的譜值在同一數量級，而僅考慮傾斜調制和流體力學調制影響的SAR圖像譜（圖2c）的譜值則小得多。另外，從譜的形態上，僅考慮速度聚束效應影響的圖像譜與考慮3種調制的海面仿真SAR圖像譜更為相近。

2.2 實測海面SAR圖像調制因素影響分析

為了進一步驗證圖1和圖2所示的結果，我們采用真實的SAR衛星數據進行分析。本文后續使用的實驗數據為Radarsat-2 SAR數據，Radarsat-2 SAR具有高分辨率多種極化通道的成像能力，本文中使用的Radarsat-2 SAR數據為超精細全極化模式數據，方位向分辨率為4.73 m，距離向分辨率約為4.9 m，表1為本文所用到的所有數據目錄，圖3給出了一幅Radarsat-2 數據2 SAR圖像的總體位置示例。

圖3 Radarsat-2數據2的SAR圖像Fig. 3 SAR image of Radarsat-2 data 2

圖4 和圖5中基于Radarsat-2真實測量SAR數據進了不同調制因素的影響分析，這里選用SAR圖像來自表1中的數據2。圖4a為截取的512×512個像素大小的HH極化歸一化SAR圖像，圖4b為僅考慮速度聚束調制時的歸一化SAR圖像，圖4c則是僅考慮傾斜調制和流體力學調制影響的歸一化SAR圖像。由圖可見，對于真實SAR海面圖像而言，所得結果與圖1仿真結果所得結論相一致，即：僅考慮速度聚束調制影響的歸一化SAR圖像與考慮3種調制的歸一化SAR圖像的強度和海浪紋理非常接近，而只考慮傾斜調制和流體力學調制影響的SAR海面圖像的強度值明顯偏低且海浪紋理并不明顯。圖4中3幅海面SAR圖像對應的圖像譜分別如圖5a至圖5c所示，圖5中的結果與圖2仿真SAR圖像的結果一致，僅含有速度聚束調制影響的SAR圖像譜（圖5b）與原SAR圖像譜（圖5a）的譜形狀和譜密度都非常接近，而僅包含傾斜調制和流體力學調制影響的圖像譜（圖5c）的譜值就小得多，且峰值位置與原SAR圖像譜相比有一定的偏離。

圖2 仿真海面SAR圖像譜（a），僅考慮速度聚束調制的仿真海面SAR圖像譜（b）和僅考慮傾斜調制和流體力學調制的仿真海面SAR圖像譜（c）Fig. 2 The spectrum of the simulated SAR image of sea surface (a), the spectrum of the simulated SAR image of sea surface with the velocity bunching modulation (b) and the spectrum of the simulated SAR image of sea surface with the tilt modulation and the hydrodynamic modulation (c)

圖4 Radarsat-2數據2中選取的512×512像素SAR圖像（a），僅考慮速度聚束調制的SAR圖像（b）和僅考慮傾斜調制和流體力學調制的SAR圖像（c）Fig. 4 SAR image of 512×512 size selected in Radarsat-2 data 2 (a), SAR image that only the velocity bunching modulation is considered(b) and SAR image that only the tilt modulation and the hydrodynamic modulation are considered (c)

圖5 SAR圖像譜（a），僅考慮速度聚束調制的SAR圖像譜（b）和僅考慮傾斜調制和流體力學調制的SAR圖像譜（c）Fig. 5 The spectrum of the SAR image (a), the spectrum of the SAR image that only the velocity bunching modulation is considered (b)and the spectrum of the SAR image that only the tilt modulation and the hydrodynamic modulation are considered (c)

表1 Radarsat-2全極化SAR數據信息Table 1 Radarsat-2 full polarized SAR data information

2.3 海浪不同傳播方向時的調制因素比較

為了進一步分析不同海浪傳播方向時各調制函數對SAR圖像的影響，我們仿真不同涌浪傳播方向角（與距離向的夾角分別為0°、15°、30°、45°、60°、75°、90°），其他參數同圖1情況下的海面SAR圖像，進而將考慮3種調制效應的歸一化SAR圖像譜、僅考慮速度聚束調制效應的歸一化SAR圖像譜與僅考慮傾斜調制和流體力學調制時的歸一化SAR圖像譜等3種譜的積分能量進行了比較。圖6a為給出考慮不同調制因素時圖像譜積分能量隨海浪傳播方向的變化，圖6b則為僅考慮速度聚束調制時的SAR圖像譜積分和僅考慮傾斜調制和流體力學調制影響時的SAR圖像的譜積分能量分別占考慮3種調制效應的圖像譜能量的比重。由圖可見：在海浪傳播方向角遠離距離向的變化過程中，僅考慮速度聚束調制的SAR圖像譜積分能量和考慮3種調制效應的圖像譜積分能量呈現上升趨勢，而且通過圖6b可以看出，僅考慮速度聚束調制影響的SAR圖像譜積分能量占考慮3種調制效應的SAR圖像譜積分的比重不斷增大，而且總保持在70%以上。實孔徑調制（傾斜調制和流體力學調制）對SAR圖像的影響隨著海浪傳播方向與距離向的夾角的增加呈現出減小趨勢，而且傾斜調制與流體力學調制影響下的圖像譜積分能量較考慮3種調制效應的SAR圖像譜積分能量來說所占比重較小，基本都在20%以下。以上結果進一步表明速度聚束調制對SAR圖像紋理的影響遠大于傾斜調制和流體力學調制的影響。

圖6 仿真SAR圖像譜積分能量隨涌浪傳播方向角的變化（a），在不同傳播方向角下，考慮不同調制的譜積分能量占總能量的比重（b）Fig. 6 The spectral integration value variation of the simulated SAR images with different wave propagation direction angle (a), the proportion of the spectral integration variation under the different modulation conditions with different wave propagation direction angles (b)

通過上述結果我們發現在SAR成像過程中，速度聚束效應對SAR圖像中海浪紋理的影響遠大于傾斜調制和流體力學調制的影響。通過式（5）、式（7）、式（12）可知，海浪沿方位向的波數只要稍微有所增大，此時，由于β的值較大，從而導致影響非常顯著，相比較而言，傾斜調制和流體力學調制的影響則非常小。

3 相干斑噪聲影響

在海浪反演過程中相干斑噪聲的影響也是不可忽視的。Arsenault 和 April[29]證明了相干斑噪聲是乘性獨立同分布的。因此在空間域時，一般采用的相干斑噪聲模型可表示為[30-32]

式中，I是觀測到的SAR圖像中像素點的強度；x是相應的真實反射強度；乘性噪聲n滿足Gamma分布，其均值為1，并且對于均勻區域它的標準差σn為[28]

式中，var( )表示求隨機變量的方差；E( )表示求隨機變量的均值；對于單視圖像來說，方差σn為1。據此對模擬的二維海面回波信號（圖7a）添加乘性噪聲如圖7b所示。

圖7 未添加（a）和添加（b）乘性噪聲的二維海面回波信號Fig. 7 Two-dimensional echo signal of ocean wave without (a) and with multiplicative noise (b)

為分析乘性噪聲對圖像譜的影響，將未添加乘性噪聲和添加乘性噪聲的二維信號的圖像譜進行比較，圖8中給出沿著方位向波數平均的平均功率譜密度，可以看出添加乘性噪聲后的所得功率譜密度要高于未添加乘性噪聲的平均功率譜密度。對于基于SAR圖像的海浪反演算法而言，尤其值得關注的是：在圖像譜的低波數區域，乘性噪聲引起了功率譜的增大，此時基于調制函數直接進行海浪反演時，根據Hasselmann海浪譜與圖像譜的映射關系可知，乘性噪聲在低波數區必將引起顯著反演誤差[9]。

圖8 沿方位向的平均圖像譜密度Fig. 8 Mean density of image spectrum along azimuth direction

4 基于速度聚束調制的SAR圖像海浪反演

目前，在Hasselmann海浪譜反演算法中，由于Bragg理論并不能準確反映海面回波隨大尺度海浪斜率的變化，因此，基于Bragg理論所得海浪傾斜調制函數并不準確；此外，流體力學調制至今也尚不成熟。然而，通過前面的分析可見，相對于傾斜調制和流體力學調制而言，速度聚束調制對SAR圖像紋理特征的影響起主導作用，這使得單純基于速度聚束調制的海浪譜反演算法所得海浪參數可能更為準確。另外，通過上一節討論可知，相干斑會給SAR圖像譜帶來影響，尤其是低波數范圍內的噪聲會導致海浪譜的反演產生顯著誤差，為減小該因素的影響，我們基于速度聚束調制反演海浪方位向斜率譜。

由式（1）可得，海面沿方位向斜率可表示為

式中，ky代表方位向波數。根據SAR圖像海浪調制理論，僅考慮速度聚束調制影響時的歸一化SAR圖像可表示為

式中，

根據斜率譜可求主波波長為

式中，kp為主波波數，主波周期為

當海浪并非嚴格地沿距離向傳播時，海浪斜率均方根

式中，σas為方位向斜率均方根；?表示海浪傳播方向與距離向的夾角；有效波高HS可以通過主波波長和斜率均方根求出[14]

圖9 中分別給出了基于數據2和數據3 SAR圖像選取的一塊子圖像反演所得海浪方位向斜率譜。圖9a中為數據2 SAR圖像基于速度聚束調制反演所得斜率譜，計算得到的有效波高為3.39 m，主波波長為256.74 m，海浪傳播方向與方位向夾角為63.23 °。圖9b中，反演所得到的有效波高為2.18 m，主波波長216.85 m，傳播方向角為126.69 °。一般情況下，對于開闊海域中的海浪而言，在統計意義上相近海域的海浪是一樣的，因此，在整幅SAR圖像中選取一塊區域反演所得海浪譜與SAR圖像中其他區域反演所得海浪譜具有相似性。為了更好地說明這一情況， 作為示例，我們基于數據2和數據3兩幅SAR圖像的整幅圖像對海浪譜進行了反演，所得結果為圖10a和圖10b，將圖10反演結果與圖9中選取SAR子圖像反演結果進行對比可見：整幅圖像與局域圖像反演所得海浪譜的形狀、有效波高、主波波長和傳播方向等趨于一致。

圖9 Radarsat-2數據2（a）和數據3（b）中 SAR子圖像反演所得方位向斜率譜Fig. 9 The azimuth slope spectrum retrieved from the sub-image of SAR in Radarsat-2 data 2 (a) and data 3 (b)

圖10 Radarsat-2數據2（a）和數據3（b）整幅SAR圖像反演所得方位向斜率譜Fig. 10 The azimuth slope spectrum retrieved from the SAR in Radarsat-2 data 2 (a) and data 3 (b)

為驗證本文方法反演結果的有效性，我們將表1中10景Radarsat-2 SAR數據反演所得結果分別與浮標數據和ECMWF再分析數據進行了對比，結果如表2所示，作為比較，表2中還給出了MPI法和同極化調制函數法反演所得結果，為了減少對外部數據的依賴，在MPI法中嘗試將通過SAR圖像譜和海浪譜的一階線性映射關系得到的海浪譜作為初猜譜進行輸入。3種不同反演方法所得有效波高與浮標觀測結果對比的散點圖如圖11所示，與ECMWF再分析有效波高數據的對比散點圖見圖12。在與浮標實測數據對比的圖11中，紅色散點和紅色實線分別代表本文方法反演所得有效波高及其擬合線，反演所得波高的均方根誤差為0.79 m，擬合線斜率為1.13；綠色散點和綠色實線分別為MPI方法反演的有效波高結果及其擬合線，其均方根誤差為1.47 m，擬合線的斜率為0.44；藍色散點和藍色實線代表同極化調制函數法反演的有效波高結果及其擬合線，該方法所得結果的均方根誤差為1.46 m，擬合直線斜率為0.62。通過與ECMWF再分析數據對比的圖10可以看出，其結果的整體趨勢與圖9中一致，并且本文方法反演所得波高的均方根誤差0.92 m也為3種方法中最小。通過圖11和圖12不同方法反演結果的對比分析，可見：與MPI方法和同極化調制函數法反演的有效波高相比，利用本文算法反演斜率譜進而計算得到的有效波高的誤差更小，與浮標實測數據對比的擬合直線斜率更接近1。但是需要說明的是：當海浪沿距離向傳播時，由式（24）可得，本文方法不再適用，而同極化調制函數法則仍適用。

圖11 不同方法反演有效波高與浮標觀測結果對比散點圖Fig. 11 The scatter plots of retrieved significant wave height by different methods with buoy data

圖12 不同方法反演有效波高與ECMWF再分析數據對比散點圖Fig. 12 The scatter plots of retrieved significant wave height by different methods with ECMWF reanalysis data

表2 Radarsat-2數據不同反演方法的有效波高結果Table 2 The retrieved significant wave height of different inversion methods for Radarsat-2 data

5 結論

本文首先分析了海面SAR圖像中傾斜調制、流體力學調制和速度聚束調制3種調制函數的影響，將考慮不同調制函數的SAR圖像分別提取對比，并在不同海浪傳播方向情況下統計分析比較，結合理論分析，發現速度聚束調制在3種調制中對SAR圖像的影響尤為顯著，遠大于傾斜調制和流體力學調制的影響。通過對相干斑噪聲的影響分析發現相干斑噪聲對海浪SAR圖像反演有一定影響，特別是在低波數區域相干斑噪聲對反演結果的影響顯著，基于此，我們利用只考慮速度聚束調制的SAR圖像直接反演方位向斜率譜的方法以降低相干斑噪聲的影響。通過比較，可得：本文算法反演所得有效波高結果較MPI方法和同極化調制函數法所得結果與浮標數據和ECMWF再分析數據吻合更好。與Zhang等[14]提出的全極化斜率譜反演算法不同的是，本文算法是基于單極化SAR數據，從中提取出僅考慮速度聚束調制影響的SAR圖像反演海浪方位向斜率譜。然而，需要說明的是：當海浪沿距離向傳播時，本文反演海浪有效波高的方法不再適用。