王 珂 洪 峻 明 峰 丁 巖

(1.中國科學(xué)院空間信息處理與應(yīng)用系統(tǒng)技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京 100190； 2.中國科學(xué)院電子學(xué)研究所，北京 100190；3.中國科學(xué)院研究生院，北京 100049)

1． 引 言

雷達(dá)可遙感探測淺海水下地形，可實(shí)現(xiàn)對淺海地質(zhì)環(huán)境的快速、動態(tài)、大范圍監(jiān)測，對此國外已開展了大量的理論和實(shí)驗(yàn)研究[1-3]，但其主要是基于布拉格(Bragg)共振散射模型的。Hennings I[4-6]對基于高斯分布的海面鏡面散射模型的雷達(dá)探測淺海水下地形機(jī)制，進(jìn)行了理論研究，并利用船載雷達(dá)進(jìn)行了實(shí)驗(yàn)驗(yàn)證。小入射角時雷達(dá)海洋遙感，需要掌握海面鏡面散射建模及其應(yīng)用方法。對于海面鏡面散射建模，國外學(xué)者已有了一些研究[7-9]，但有關(guān)Gram-Charlier分布在鏡面散射建模中的性能及其在海洋遙感中的應(yīng)用效果，卻少有論述。首先建立了基于Gram-Charlier分布的海面鏡面散射模型，深入研究其性能表現(xiàn)，為由海面鏡面散射產(chǎn)生的雷達(dá)后向散射截面，提供了定量計算方法。其次，驗(yàn)證了基于Gram-Charlier分布的海面鏡面散射模型在中低風(fēng)速時性能可靠后，將其應(yīng)用于雷達(dá)探測淺海水下地形。為基于海面鏡面散射模型的雷達(dá)海洋遙感，提供了理論驗(yàn)證和應(yīng)用方法。

2. 理論與模型

2.1 海面鏡面散射理論

海面電磁波的散射機(jī)制，按入射角范圍，主要分為[7,10]：鏡面散射和Bragg共振散射。在海面風(fēng)場作用下，波浪可發(fā)生相對于水平面大至25°的傾斜[11]，則雷達(dá)入射角小于20°時，主要是由鏡面散射產(chǎn)生了雷達(dá)后向散射截面[7]。

Valenzuela G R[7]提出當(dāng)海面斜率為各向同性高斯分布時，可計算海面鏡面散射產(chǎn)生的雷達(dá)后向散射截面為

(1)

式中：θ為雷達(dá)入射角；R(0)為雷達(dá)波垂直入射時(θ=0°)的菲涅耳(Fresnel)反射系數(shù)；s2為總的海面均方斜率。

經(jīng)典理論認(rèn)為：鏡面散射正比于海面斜率的概率密度函數(shù)[7]。目前海面斜率分布有三種：高斯分布、Gram-Charlier分布、Liu Y分布[12]。Cox C和Munk W用航空拍攝海面太陽耀斑的方法測量了海面斜率的實(shí)際分布[11]，提出海面斜率可服從Gram-Charlier分布。Liu Y分布在高風(fēng)速時比Gram-Charlier分布性能更好。但Liu Y分布沒有給出其中峰度系數(shù)n和偏度系數(shù)λsk隨連續(xù)風(fēng)速的取值方法[12]，嚴(yán)重限制了其實(shí)際應(yīng)用。為了計算在不同雷達(dá)和海況條件連續(xù)輸入下，由海面鏡面散射產(chǎn)生的σ0，建立了基于Gram-Charlier分布的海面鏡面散射模型。

2.2 基于Gram-Charlier分布的海面鏡面散射模型

Barrick D E[13]提出的隨機(jī)粗糙表面鏡面散射模型，可用于計算由海面鏡面散射產(chǎn)生的雷達(dá)后向散射截面[7],即

σ0=πsec4θ|R(0)|2p(ζx,ζy)

(2)

式中：θ為雷達(dá)入射角；R(0)為雷達(dá)波垂直入射時(θ=0°)的Fresnel反射系數(shù)；p(ξx,ξy)為海面x方向斜率ζx和y方向斜率ζy的聯(lián)合概率密度函數(shù)。

假設(shè)在海面上存在確定坐標(biāo)系(x,y)，y軸的正向指向逆風(fēng)方向；x軸的正向指向?yàn)閥軸正向順時針旋轉(zhuǎn)90°，即側(cè)風(fēng)方向。根據(jù)入射雷達(dá)波在海面傾斜波面上發(fā)生鏡面散射時的幾何光學(xué)關(guān)系[11,13]，對于產(chǎn)生鏡面散射的海面斜率，則有

ζx=tanθsinφ

(3)

ζy=tanθcosφ

(4)

式中：θ為雷達(dá)入射角；φ為雷達(dá)波與逆風(fēng)方向之間的夾角，φ=0°時雷達(dá)波指向逆風(fēng)方向。

對于海面斜率的Gram-Charlier分布[11]，即有

c04(η4-6η2+3)]}

(5)

其中，逆風(fēng)(up-wind)和側(cè)風(fēng)(cross-wind)方向上的海面均方斜率規(guī)范化為

(6)

(7)

Cox C和Munk W由實(shí)測數(shù)據(jù)，給出了逆風(fēng)和側(cè)風(fēng)方向上的海面均方斜率為

(8)

(9)

式中,u12.5為海面上12.5 m高度的風(fēng)速。

總的海面均方斜率為

(10)

Wu J[14]對Cox C和Munk W的實(shí)測數(shù)據(jù)，重新進(jìn)行了統(tǒng)計，給出了更為準(zhǔn)確的表達(dá)

s2= (0.90+1.20·Inu10)×10-2

(u10<7 m/s)

(11)

s2= (-8.40+6.00·Inu10)×10-2

(u10>7 m/s)

(12)

且有如下近似關(guān)系

(13)

偏度系數(shù)為

c21=m-nu12.5

(14)

c03=p-qu12.5

(15)

而c40、c22、c04為峰度系數(shù)。式(2)～(15)中的參數(shù)，取清潔水體(clean surface)時[11]，如表1所示。

表1 Gram-Charlier分布中的參數(shù)取值

對于電磁波的鏡面散射可有[13]

(16)

(17)

RVH(0)=RHV(0)→0

(18)

式中：εr為海水的相對介電常數(shù)，與雷達(dá)波段有關(guān)[15]；μr為海水的相對磁導(dǎo)率，可取μr=1，代入式(16)、(17)，則有RVV(0)=RHH(0)； 下標(biāo)VV、HH等表示極化。所以可取

|R(0)|2=|RVV(0)|2=|RHH(0)|2

(19)

將式(16)～(19)代入式(2)可見

(20)

2.3 風(fēng)速的統(tǒng)一

風(fēng)速在實(shí)際使用時，通常是以海面上10 m高度的u10為基準(zhǔn)，所以需將u12.5換算為u10。在大氣中性穩(wěn)定的情況下，摩擦風(fēng)速u*與風(fēng)速的關(guān)系為

(21)

式中：uz為海面上z(m)高度處的風(fēng)速；κ為馮·卡門(von Karman)常數(shù)；z0為粗糙長度，可由下式計算[3]

(22)

c10=(0.8+0.065u10)×10-3

(23)

綜合式(21)～(23)，在模型中風(fēng)速將統(tǒng)一由u10表示，以減少選取風(fēng)速時帶來的誤差。

2.4 雷達(dá)探測模型

根據(jù)Hennings I[6]的理論探討，在u10≤8 m/s時應(yīng)用Gram-Charlier分布的雷達(dá)探測淺海水下地形應(yīng)具有可行性。典型淺海水下地形如圖1所示，為鋸齒狀沙坡，設(shè)地形梯度方向?yàn)閤perp方向。

圖1 典型淺海水下地形示意圖

能量譜F(k)與作用量譜N(k)和波高譜Ψ(k)的關(guān)系為[16]

(24)

式中：g為重力加速度；k為波數(shù)，固有頻率為

ω′=(gk)1/2

(25)

可取波高譜的形式為[17]

Ψ(k)=apk-4

(26)

式中，ap是菲利普斯(Phillips)常數(shù)。Stolte S[18]經(jīng)測量給出了ap與風(fēng)速(u10≤8 m/s時)的關(guān)系為

log10ap= -2.90+3.06·10-1u10-

(27)

水下地形對能量譜(假設(shè)潮流流經(jīng)水下地形之前處于譜平衡狀態(tài)F0)的一階擾動調(diào)制為[1，4-5]

(28)

(29a)

(29b)

(29c)

(29d)

(29e)

對于真實(shí)的粗糙海面，有效入射角可重寫為

θ=θplane+θrough

(30)

式中,θplane為平坦海面時的入射角，粗糙海面引起的入射角可表示為

tan2θrough=s2

(31)

被水下地形擾動的入射角，可由下式得到

(32a)

(32b)

能量譜與海面均方斜率的關(guān)系為[17]

(33)

所以有

(34)

(35)

式中： 波數(shù)的積分區(qū)間為k0=2π/λ0≤k≤kc=2π/λradar,λ0為雷達(dá)空間分辨率，λradar為雷達(dá)波長；α是逆風(fēng)方向與xperp方向上潮流分量之間的夾角。

由式(3)、(4)，被水下地形擾動的海面斜率為

δζx=tanδθ·sinφ

(36)

δζy=tanδθ·cosφ

(37)

由式(6)、(7)，被水下地形擾動后，規(guī)范化的海面均方斜率可重寫為

(38)

(39)

綜合式(24)～(39)代入式(2)，可得水下地形擾動引起的雷達(dá)后向散射截面的變化為

δσ=σ-σ0

π·|R(0)|2·sec4(θ)·p(ζx,ζy)

(40)

則水下地形擾動對雷達(dá)后向散射截面的調(diào)制為

(41)

3. 實(shí)驗(yàn)與分析

3.1 海面鏡面散射模型的實(shí)驗(yàn)分析

為了驗(yàn)證基于Gram-Charlier分布的海面鏡面散射模型的正確性，并研究其性能，給出了其與基于高斯分布時(即式(1))的實(shí)驗(yàn)比較，如圖2～5所示。模型計算的σ0，與Plant W J[8]的Ku波段實(shí)測數(shù)據(jù)和CMOD4模型[19]數(shù)據(jù)(可表示C波段實(shí)測的σ0[12])相比較。

圖2是Ku波段，u10=8 m/s，θ=10°時基于Gram-Charlier分布和高斯分布，其中采用Wu J方法估計海面均方斜率，分別計算的σ0，與Plant W J實(shí)測數(shù)據(jù)相比較?？梢?，基于Gram-Charlier分布計算的σ0，總體上能合理地反映鏡面散射隨風(fēng)向的變化，更接近實(shí)測數(shù)據(jù)；相比基于高斯分布時具有一定優(yōu)越性。

圖2 Ku波段，u10=8 m/s，Wu J方法

圖3是C波段，u10=8 m/s，θ=16°時基于Gram-Charlier分布和高斯分布，采用Wu J方法估計海面均方斜率，分別計算的σ0，與CMOD4模型數(shù)據(jù)相比較?？色@得與圖2中同樣的結(jié)論。

圖3 C波段，u10=8 m/s，Wu J方法

圖4是C波段，u10=2 m/s，θ=16°時基于Gram-Charlier分布和高斯分布，采用Cox C和Munk W方法估計海面均方斜率，分別計算的σ0，與CMOD4模型數(shù)據(jù)相比較。可見，兩種鏡面散射模型計算的σ0，均有較大誤差；特別是基于Gram-Charlier分布時，在側(cè)風(fēng)方向(φ=90°)計算的σ0大于逆風(fēng)方向(φ=0°)。還將Ku波段，u10=2 m/s，θ=10°時基于Gram-Charlier分布，采用Cox C和Munk W方法估計海面均方斜率，計算的σ0，與Plant W J實(shí)測數(shù)據(jù)進(jìn)行了比較，出現(xiàn)了同樣的問題。這明顯不符合Plant W J由實(shí)測數(shù)據(jù)，獲得的σ0隨風(fēng)向的“類余弦”變化規(guī)律。

圖5是C波段，u10=2 m/s，θ=16°時，基于Gram-Charlier分布和高斯分布，采用Wu J方法估計海面均方斜率，分別計算的σ0，與CMOD4模型數(shù)據(jù)相比較?？梢姡捎肳u J方法估計海面均方斜率，改進(jìn)了兩種鏡面散射模型的性能表現(xiàn)，計算的σ0更接近實(shí)測數(shù)據(jù)。還將Ku波段，u10=2 m/s，θ=10°時基于Gram-Charlier分布，采用Wu J方法估計海面均方斜率，計算的σ0，與Plant W J實(shí)測數(shù)據(jù)進(jìn)行了比較，出現(xiàn)了同樣的改進(jìn)效果。

圖4 C波段，u10=2m/s，Cox C、Munk W方法

圖5 C波段，u10=2 m/s，Wu J方法

研究基于Gram-Charlier分布的海面鏡面散射模型的性能時，發(fā)現(xiàn)導(dǎo)致其偏差的主要原因是Cox C和Munk W估計的海面均方斜率、峰度系數(shù)、偏度系數(shù)有一定的誤差范圍，可用“±標(biāo)準(zhǔn)差”來表示[11]。通過實(shí)驗(yàn)分析，認(rèn)為采用Wu J方法估計海面均方斜率， 可改進(jìn)Gram-Charlier分布在中低風(fēng)速時的性能，以減少σ0的計算誤差。但發(fā)現(xiàn)改進(jìn)后模

型的準(zhǔn)確性仍與風(fēng)速條件有關(guān)， 過高風(fēng)速時(u10>14 m/s)會引起基于Gram-Charlier分布的建模偏差增大，這與Liu Y的分析相符合[12]。而通過與CMOD4模型數(shù)據(jù)的比較，發(fā)現(xiàn)過低風(fēng)速時(u10<2 m/s)也有較大的建模偏差。所以認(rèn)為模型需要首先進(jìn)行風(fēng)速條件下的驗(yàn)證，以避免應(yīng)用時出現(xiàn)錯誤。通過與Plant W J實(shí)測數(shù)據(jù)和CMOD4模型數(shù)據(jù)的比較，采用Wu J方法估計海面均方斜率，將風(fēng)速統(tǒng)一為u10，在中低風(fēng)速(u10=(2～8)m/s)范圍內(nèi)，基于Gram-Charlier分布的海面鏡面散射模型性能可靠。

3.2 雷達(dá)探測的實(shí)驗(yàn)分析

圖6 淺海水下地形對雷達(dá)后向散射截面的調(diào)制

4. 結(jié) 論

深入研究了基于Gram-Charlier分布的海面鏡面散射建模方法，分析了其極化特性，將風(fēng)速統(tǒng)一以海面上10 m高度為基準(zhǔn)，以減少選取風(fēng)速時帶來的誤差。通過與實(shí)測數(shù)據(jù)的比較，證明了基于Gram-Charlier分布的海面鏡面散射模型的正確性，及相比基于高斯分布時具有優(yōu)越性，特別是明確地給出了σ0隨風(fēng)向的變化。提出了采用Wu J方法估計海面均方斜率，可改進(jìn)Gram-Charlier分布在中低風(fēng)速時的性能，以減少σ0的計算誤差。指出了過高(u10>14 m/s)或過低(u10<2 m/s)的風(fēng)速均會引起建模偏差的增大，并分析了其原因。提出了基于Gram-Charlier分布的海面鏡面散射模型在應(yīng)用時，需要首先進(jìn)行風(fēng)速條件下的驗(yàn)證，以避免出現(xiàn)錯誤，這是Hennings I[6]未曾考慮到的問題。進(jìn)而將基于Gram-Charlier分布的海面鏡面散射模型應(yīng)用于雷達(dá)探測淺海水下地形，通過與實(shí)驗(yàn)結(jié)果的比較，證明了其可行性和優(yōu)越性。對Gram-Charlier分布在海面鏡面散射建模中的性能，及其在雷達(dá)探測淺海水下地形中應(yīng)用方法的研究，一方面為海面鏡面散射提供了定量計算方法；另一方面為基于海面鏡面散射模型的遙感探測提供了理論驗(yàn)證和應(yīng)用方法，獲得的結(jié)論對實(shí)際應(yīng)用有指導(dǎo)作用。由于Hennings I理論是一階模型[4-6]，其準(zhǔn)確性有限，適用條件也有待更多的實(shí)驗(yàn)驗(yàn)證，這將是以后的研究方向。

[1] ALPERS W, HENNINGS I. A theory of the imaging mechanism of under water bottom topography by real and synthetic aperture radar[J]. Journal of Geophysical Research, 1984, 89(6): 10529-10546.

[2] HENNINGS I. A historical overview of radar imaging mechanism of sea bottom topography[J]. International Journal of Remote Sensing, 1998, 19(7): 1447-1454.

[3] VOGELZANG J. Mapping submarine sand waves with multiband imaging radar 1.Model development and sensitivity analysis[J]. Journal of Geophysical Research, 1997, 102(C1): 1163-1181.

[4] HENNINGS I, HERBERS D. Radar imaging mechanism of marine sand waves at very low grazing angle illumination caused by unique hydrodynamic interactions[J]. Journal of Geophysical Research, 2006, 111(C10): C10008. doi:10.1029/2005JC003302.

[5] HENNINGS I, HERBERS D. Radar imaging mechanism of marine sand waves at very low grazing angle illumination[C]//EARSeL eProceedings, 2005, 4: 226-240.

[6] HENNINGS I, HERBERS D. The expected potential of TerraSAR-X high resolution spotlight mode data for shallow sea bottom topography imaging:a preview[C]//EARSeL eProceedings, 2007, 6: 67-81.

[7] VALENZUELA G R. Theories of the interaction of electromagnetic and ocean waves-a review[J]. Boundary-Layer Meteorology, 1978, 13: 61-85.

[8] PLANT W J. A stochastic, multiscale model of microwave backscatter from the ocean[J]. Journal of Geophysical Research, 2002, 107 (C9): 3120. doi:10.1029/2001JC000909.

[9] FREILICH M H, VANHOFF B A. The relationship between winds,surface roughness,and radar backscatter at low incidence angles from TRMM precipitation radar measurements[J]. Journal of Atmospheric and Oceanic Technology, 2003, 20(4): 549-562.

[10] 邵連軍, 張訓(xùn)械, 孫 強(qiáng), 等. GPS信號海面散射特性和初步實(shí)驗(yàn)結(jié)果[J]. 電波科學(xué)學(xué)報, 2008, 23(4): 699-703.

SHAO Lianjun, ZHANG Xunxie, SUN Qiang, et al. Property of GPS signal scattered from sea surface and preliminary experimental results[J]. Chinese Journal of Radio Science, 2008, 23(4): 699-703. (in Chinese)

[11] COX C, MUNK W. Measurement of the roughness of the sea surface from photographs of the sun's glitter[J]. Journal of the Optical Society of America, 1954, 44(11): 838-850.

[12] LIU Y, YAN X H, LIU W T, et al. The probability density function of ocean surface slopes and its effects on radar backscatter[J]. Journal of Physical Oceanography, 1997, 27(5): 782-797.

[13] BARRICK D E. Rough surface scattering based on the specular point theory[J]. IEEE Transactions on Antennas and Propagation, 1968, AP16(4): 449-454.

[14] WU J. Mean square slopes of the wind-disturbed water surface,their magnitude, directionality,and composition[J]. Radio Science, 1990, 25(1): 37-48.

[15] ULABY F T, MOORE R K, FUNG A K. Microwave Remote Sensing[M]. MA: Artech House, 1986: Vol.Ⅲ.

[16] HOLLIDAY D, ST-CYR G, WOODS N E. A radar ocean imaging model for small to moderate incidence angles[J]. International Journal of Remote Sensing, 1986, 7(12): 1809 -1834.

[17] PHILLIPS O M. The dynamics of the upper ocean(second edition)[M]. Cambridge University Press, 1977.

[18] STOLTE S. Dynamics of short waves and wave breaking[R]. (Federal Armed Forces Underwater Acoustics and Marine Geophysics Research Institute,Kiel, Germany)Report1990-4, 1990.

[19] STOFFELEN A C M, ANDERSON D L T. Scatterometer data interpretation: derivation of the transfer function CMOD4[J]. Journal of Geophysical Research, 1997, 102(C3): 5767-5780.