粗糙面上粒子層矢量輻射傳輸方程的二階迭代解法*

姚海 何姿? 丁大志 陳如山 黨訓旺 陳勇

1)(南京理工大學電子工程與光電技術(shù)學院,南京 210018)

2)(電磁散射重點實驗室,北京 100089)

為了探究典型粗糙面上隨機粒子層中能量傳輸?shù)亩啻紊⑸錂C制,提出了一種基于矢量輻射傳輸方程的建模二階計算方法.該建模方法將建模場景(粗糙面上粒子層)在高度維(Z 軸)劃分為多個傳輸散射層,基于矢量輻射傳輸理論中的一階迭代散射解,利用典型粗糙面的半經(jīng)驗半解析方法,求解出整個場景的二階迭代散射解.同時,研究粒子層能量在粒子與粒子間的多次散射機制,以及粒子與地表粗糙面間的多次散射機制.數(shù)值結(jié)果表明,該二階迭代解法相較于矢量輻射傳輸方程的一階迭代散射解,能夠更完整地探究互作用的散射特性,且可從能量傳輸角度解譯建模場景中物體間的相干作用,從而可用于植被地物環(huán)境下的多次散射機制的解析以及散射系數(shù)變化趨勢的預估.

1 引言

與地物電磁散射特性有關(guān)的研究,在遙感、探測、反隱身等領(lǐng)域具有重要的應用價值.地物對象包括各種地形地貌,如雪地、土壤、沙漠、海面等,以及大自然中的各種植被,如草地、樹木、莊稼等.地物研究的種類之多、范圍之廣,足以可見研究地物電磁散射特性的強烈需求和重要性.

對于典型地表植被散射的建模,在Ulaby 等[1]提出的密歇根森林植被散射模型(MIMICS)的基礎(chǔ)上,文獻[2-5]等基于矢量輻射傳輸(vector radiative transfer,VRT)方程發(fā)展了粒子層散射矢量輻射傳輸模型,采用隨機取向的粒子模擬樹冠和樹干[6-9],求解VRT 方程的一階解.同時,引入由Fung 等[10]提出基于積分方程方法的隨機粗糙表面散射,可用于求解出包含粗糙面散射的Mueller矩陣解,模擬出非球型粒子體散射和粗糙面的散射耦合,建立粗糙面上混雜非球型粒子散射模型[11-13].除此之外,對于非均勻散射層中VRT方程多次散射求解的問題,也有相關(guān)的討論[14,15],但是,對于包含粗糙面散射和非球型粒子體的多次散射耦合的VRT 二階Mueller 矩陣解,目前暫沒有詳細的推導和討論.

VRT 方程是一種半解析建模方法,在分析討論建模對象的過程中,結(jié)合物體參數(shù)和環(huán)境特性在實際情況中的復雜性,用解析公式加上由大量實際的實驗數(shù)據(jù)擬合的經(jīng)驗公式,來對背景環(huán)境進行建模.VRT 描述了電磁波在傳輸?shù)倪^程中由于環(huán)境背景的散射、吸收、衰減以及增益的影響,使電磁波能量發(fā)生相應的變化,可用一定的微積分方程來描述.通常,對典型植被地面環(huán)境的研究分析,可將其建模成粗糙面上混雜粒子層,使用VRT 半解析方法進行數(shù)值仿真[16-18].

對于粗糙面上粒子層的模型,本文在VRT 一階迭代解的基礎(chǔ)上,求解出VRT 二階迭代解的表達式,探究粒子層中粒子間、粒子與粗糙面間的多次散射作用,分析粒子層整體的后向散射系數(shù)隨入射角、占空比的變化規(guī)律,以及二階迭代解的數(shù)值作用.

2 VRT 二階散射的迭代解方法

2.1 VRT 方程及一階迭代解

如圖1 所示,將典型地物場景的建模在垂直方向(Z軸)分為兩層,上層為模擬離散地物模型的粒子層,在層厚d∈(-d1,0)均勻分布隨機粒子模型.下層d=-d1處為模擬隨機粗糙面的地表層.當電磁波在粒子層和地表層之間發(fā)生傳播時,電磁波的能量隨著電磁波傳播距離的增大而改變,此時電磁波的強度變化可通過VRT 方程方法來分析.一般情況下,VRT 方程表示為

圖1 粗糙面上粒子層模型示意圖,圖中為層厚為d1 的粒子層,最下層為隨機粗糙面模型,綠色部分為粒子層中橢球粒子模型,深色部分為圓柱粒子模型Fig.1.Particle layer model on rough surface.The figure shows the particle layer with layer thickness of d1,the bottom layer is the random rough surface model,the green parts is the ellipsoid particle model in the particle layer,and the dark parts is the cylinder particle model.

為解決粗糙面上粒子層的一維VRT 問題,將(1)式改寫為關(guān)于z方向的一維偏微分方程,即有坐標變化:

其中,k為波數(shù),n0為單位體積粒子數(shù),fvv為粒子前向散射矩陣.

微積分方程通解為

其中A為迭代相關(guān)系數(shù),等式右側(cè)第二項為迭代求解需要更新的部分.(5)式可求解出(3)式的通解,即VRT 上行波和下行波迭代解的一般表達式:

通過邊界條件求A和B,代入邊界條件:

VRT 邊界迭代關(guān)系為

可求得迭代系數(shù)的一般表達式為

VRT 一階解忽略電磁波多次散射的作用,即忽略迭代更新部分置零,求解迭代系數(shù)A,推導VRT 一階散射迭代解[2],本文直接給出建模環(huán)境任意位置中VRT 一階散射解的一維表達式,上行波、下行波的具體表達式分別為

2.2 VRT 方程的二階迭代解

在建立的粗糙面上方隨機粒子層模型中,對于所建立的一維VRT 方程,可設(shè)置粒子層上下邊界條件為[2]

VRT 散射迭代二階解推導過程如下.

步驟1給出在粒子層中,VRT 二階迭代解的通解,其中上行波和下行波分別為

步驟2給出的VRT 一階解,并作為迭代的散射源,即(10)式和(11)式,一階解在任意位置的具體表達式,代入二階解通解公式中,得到代入新的散射源后的二階解通解.

步驟3將(13)式和(14)式聯(lián)立VRT 邊界條件(7),更新得到二階系數(shù)A2的表達式,如(15)式.將VRT 一階解也代入A2,得到A2完整形式,如(16)式:

步驟4將A2代入VRT 二階解通解,可得粒子層任意位置VRT 二階迭代的多次散射解.令VRT 二階解上行波z=0,即為電磁波離開粒子層(最上層表面)的散射強度:

3 VRT 二階迭代解的Mueller 矩陣

如圖2,電磁波在粒子層中傳輸,通過不同的散射過程,最終從粒子層傳輸出去.VRT 方程反映了電磁波能量的變化,從物理意義的角度,將VRT 迭代解表達式分解為多次散射項,加入描述電磁波與目標之間作用的Mueller 矩陣,最終形式表示:

圖2 VRT 方程各項散射機制示意圖.圖中最上層為粗糙面模型,綠色部分表示粒子層隨機分布的粒子模型,粒子層層厚為d1,Mi 為VRT 方程不同散射機制Muller 矩陣解Fig.2.Various scattering mechanisms of VRT equation.In the figure,the top layer is the rough surface model,and the green parts represent the particle model with random distribution of particle layers.The thickness of particle layers is d1,and Mi is the solution of Muller matrix with different scattering mechanism of VRT equation.

其中,1)當S為粒子的散射矩陣時,可表示為粒子的相函數(shù)矩陣形式,且,n0為單位體積內(nèi)粒子數(shù).

2)當S為下墊粗糙面的散射矩陣時,可表示為粗糙面的相函數(shù)矩陣形式,即.

當Mueller 矩陣描述電磁波在環(huán)境中傳播的能量變化時,Mueller 矩陣解可表示為VRT 多次迭代解之和:

二階解的多次散射項,從實際散射方式劃可分為4 種,其詳細散射過程參考本節(jié)后面內(nèi)容.對于粒子層存在的其他更多次散射項,散射機制具體表達式也可得到,由于其數(shù)值對于求解的VRT 方程的效果極小,無實際意義,故不再討論和分析.

3.1 粒子-粒子的二次散射

如圖3(a)所示,該散射過程代表電磁波首先進入粒子層,然后電磁波經(jīng)過兩次粒子散射作用,最終從粒子層散射出去.其數(shù)學表達式為

3.2 粒子-粒子-粗糙面之間的三次散射機制

如圖3(b)所示,該散射過程代表電磁波首先進入粒子層,然后電磁波首先經(jīng)過了兩次粒子散射作用,又經(jīng)過一次粗糙面散射作用,最終從粒子層散射出去.其數(shù)學表達式為

3.3 粗糙面-粒子-粒子的三次散射機制

如圖3(c)所示,該散射過程代表電磁波首先進入粒子層,然后電磁波經(jīng)過了一次地表粗糙面的散射作用,又經(jīng)過了兩次粒子散射作用,最終從粒子層散射出去.其數(shù)學表達式為

圖3 粒子層多次散射過程示意圖 (a)(粒子-粒子)二次散射機制示意圖;(b)(粒子-粒子-粗糙面)三次散射機制示意圖;(c)(粗糙面-粒子-粒子)三次散射機制示意圖;(d)(粒子-粗糙面-粒子)三次散射機制示意圖Fig.3.Diagram of multiple scattering process of particle layer:(a)(Particle-particle)secondary scattering mechanism;(b)(particleparticle-rough surface)third-order scattering mechanism;(c)(rough surface-particle-particle)third-order scattering mechanism;(d)(particle-rough surface-particle)third-order scattering mechanism.

3.4 粒子-粗糙面-粒子的三次散射機制

如圖3(d)所示,該散射過程代表電磁波首先進入粒子層,然后電磁波經(jīng)過一次粒子散射作用,再經(jīng)過地表粗糙面的散射作用,又經(jīng)過一次粒子散射作用,最終從粒子層散射出去.其數(shù)學表達式為

(26)—(30)式中的散射角在圖(3)中標出并說明,角度積分(包括θ和π),本文在程序計算的過程都采用高斯積分處理和計算,根據(jù)角度積分的范圍可選擇不同的積分點個數(shù),一般對角度2π 積分可選取最多20 個高斯積分點.

4 數(shù)值仿真和分析

4.1 粒子形狀對后向散射系數(shù)的影響

算例1首先,考慮一個下墊面上粒子層模型,粒子層高度為300 cm.分別考慮當粒子層中均勻分布3 種不同形態(tài)的粒子的情況:

1)GRG 近似橢球粒子,粒子半徑為0.1 cm,半高為5 cm,粒子介電常數(shù)為(14.1,4.5),隨機取向歐拉角為β∈(0°,30°),γ∈(0°,360°)取平均,粒子在粒子層中的占空比為0.003;

2)豎直圓柱粒子,圓柱半徑為1 cm,半高為50 cm,粒子介電常數(shù)為(19.6,8.1),粒子在粒子層中占空比為0.003;

3)傾斜圓柱粒子,圓柱尺寸與2)相同,固定取向歐拉角為β,γ=10..

其中,下墊面采用IEM 半解析粗糙面建模方法,粗糙面的均方根高度σ為0.2 cm,相關(guān)長度l為5 cm,粗糙面介電常數(shù)εr為(14.1,5.0).

分別仿真粒子層在C 波段(頻率為5.3 GHz)下,粒子層后向散射系數(shù)〈σvv〉,〈σhh〉,〈σvh〉,〈σhv〉隨入射角θi的變化規(guī)律曲線,不同形態(tài)粒子(GRG近似橢球,圓柱粒子)的散射矩陣表達式,以及相矩陣形式引用參考文獻[6-9].

圖4(a)和圖4(d)是粒子為均勻取向、均勻分布的GRG 橢球粒子時,粒子層的同極化、交叉極化的后向散射系數(shù)隨入射角度變化的曲線.圖4(b)和圖4(e)是粒子為均勻分布的豎直圓柱粒子時,粒子層的同極化、交叉極化的后向散射系數(shù)隨入射角度變化的曲線.圖4(c)和圖4(f)是粒子為均勻分布的傾斜圓柱粒子時,粒子層的同極化、交叉極化的后向散射系數(shù)隨入射角度變化的曲線.

圖4 不同形態(tài)粒子的粒子層后向散射系數(shù)對比 (a)取向角均勻橢球粒子層同極化后向散射系數(shù);(b)豎直圓柱粒子層的同極化后向散射系數(shù);(c)固定傾斜角圓柱粒子層的同極化后向散射系數(shù);(d)取向角均勻橢球粒子層的交叉極化后向散射系數(shù);(e)豎直圓柱粒子層的交叉極化后向散射系數(shù);(f)固定傾斜角粒子層的交叉極化后向散射系數(shù).Fig.4.Backscattering coefficients of particle layers with different morphology particles:(a)Co-polarization backscattering coefficient of a layer of ellipsoid particles with uniform orientation angles;(b)the co-polarization backscattering coefficient of a layer of vertical cylindrical particles;(c)the co-polarization backscattering coefficient of a layer of cylindrical particles with fixed inclination angles;(d)the cross-polarization backscattering coefficient of a layer of ellipsoid particles with uniform orientation angles;(e)the cross-polarization backscattering coefficient of the vertical cylindrical particle layer;(f)the cross-polarization backscattering coefficient of the particle layer with a fixed inclination angle.

從數(shù)值結(jié)果可以看出,當粒子層中均勻分布取向角均勻的GRG 橢球粒子時,VRT 方程的二階解和一階解量級和趨勢一致,整體區(qū)別較小,在交叉極化下二階解的散射系數(shù)略大于一階解,這是由于對于分布和取向都均勻的橢球粒子層,電磁進入粒子層后,在各個方向的散射強度都基本一致,無較強散射,粒子間的二次散射以及粒子-粗糙面的三次散射的作用都較小,因此,二階解的散射系數(shù)與一階解的數(shù)值基本一致.

當粒子層分布圓柱粒子時,由于圓柱粒子隨入射角改變有一定規(guī)律變化,例如豎直圓柱粒子隨入射角增大呈波動性上升趨勢,所以當入射角較大時,圓柱粒子間產(chǎn)生強散射互作用,在各個方向的散射增強,圓柱粒子二次散射的強度和概率增大,粒子層存在較強的多次散射,影響粒子層整體的后向散射系數(shù),VRT 二階解的數(shù)值和一階解區(qū)別較大,交叉極化及其曲線趨勢變化明顯.當圓柱粒子的取向角發(fā)生改變,產(chǎn)生強散射的角度發(fā)生變化,二階解的多次散射互作用減小.

圖4(b),(c),(e)和(f)中產(chǎn)生“多峰現(xiàn)象”,單個豎直圓柱粒子的散射特性隨入射角度變化有明顯的波動起伏,會出現(xiàn)若干個波峰波谷;對于本文出現(xiàn)的現(xiàn)象,分析將圓柱看作直線目標(特別是該參數(shù)的細長圓柱),直線目標的散射強度可以看作若干個散射中心點的共同作用[1,2],不同入射角電磁波入射圓柱散射后,圓柱上散射強度會在角度改變的過程中發(fā)生相位疊加或者相位消減,因此圓柱的散射特性曲線是隨入射角度增大而“波動形增長”

VRT 仿真中圓柱粒子的散射相矩陣(衰減置零)和VRT 粒子一次散射項有一致性,考慮到對粒子不同位置的積分(Z方向)以及消光衰減的作用,VRT 圓柱粒子散射層會出現(xiàn)“多峰現(xiàn)象”.在VRT 二階解中,不僅需要考慮粒子在高度上不同位置、粒子與粒子間相對距離,還有電磁波從一次散射到二次散射中漫散射的過程(文中用對散射角做積分處理),因此VRT 粒子二次散射“波動起伏”會依舊存在,但整體趨勢較粒子一次散射較為平緩.同時,VRT 一階存在粒子與粗糙面的二次散射相干作用,二階解存在粒子與粗糙面間三次散射作用,粗糙面的單次散射作用是隨入射角逐漸減小的,一般減小得快慢與粗糙面、譜分布有關(guān),而圓柱粒子的散射系數(shù)總體是隨入射角增大的,粒子與粗糙面的互作用也是波峰產(chǎn)生的原因,這些因素都是本文算例波峰波谷產(chǎn)生的原因.

算例2仿真含有混合粒子的粒子層后向散射系數(shù)隨入射角后向散射系數(shù),混合粒子的參數(shù)為算例1)中的橢球粒子和算例2)中的圓柱粒子按照其各自占空比混合計算,占空比混合模型計算公式為

圖5 混合粒子的同極化散射系數(shù)對比Fig.5.Comparison of co-polarization scattering coefficients of mixed particles.

圖6 混合粒子的交叉極化散射系數(shù)對比Fig.6.Comparison of cross-polarization scattering coefficients of composite particles.

從結(jié)果可以看出,混合粒子層中存在兩種形態(tài)不同的粒子,取向均勻的針狀橢球和豎直圓柱,兩種粒子的占空比一致.由之前的結(jié)果可知,橢球粒子隨入射角度增大散射系數(shù)減小,VRT 二階解的作用較小,散射系數(shù)隨角度變化趨勢基本和一階解的基本一致.

入射角增大后,圓柱粒子的散射系數(shù)逐漸增大,粒子間互作用增大,粒子層多次散射增強,VRT 二階解對粒子層的散射作用明顯.當入射角逐漸增大后,圓柱粒子的散射系數(shù)大于橢球粒子,因此在混合粒子-粒子層模型中,根據(jù)(24)式可知圓柱粒子的散射系數(shù)對粒子層整體的影響明顯,VRT 二階解對混合粒子粒子層的作用增大.在混合粒子層中,在低掠角入射范圍內(nèi),VRT 二階解的多次散射作用較大,影響后向散射系數(shù)的數(shù)值結(jié)果.可進一步研究更多種混合粒子的后向散射系數(shù)的變化規(guī)律.

4.2 粒子占空比對后向散射系數(shù)的影響

算例3為了驗證VRT 方程二階迭代散射解的正確性,以及其對VRT 數(shù)值求解后向散射系數(shù)的作用,固定粒子模型為豎直圓柱,傾斜歐拉角均固定為10°,圓柱半徑為1 cm,半長為50 cm,介電常數(shù)為(19.6,8.1),粒子層厚為300 cm,固定入射角θi=50°,仿真下墊面上一層圓柱粒子的后向散射系數(shù)〈σvv〉,〈σhh〉,〈σvh〉,〈σhv〉隨粒子密度(占空比)的變化曲線,仿真結(jié)果如圖7 和圖8 所示.

圖7 不同占空比下粒子層的同極化散射系數(shù)對比Fig.7.Comparison of the co-polarization scattering coefficients of particle layers with different duty cycles.

圖8 不同占空比下粒子層的交叉極化散射系數(shù)對比Fig.8.Comparison of cross-polarization scattering coefficients of particle layers with different duty cycles.

結(jié)果表明:隨著粒子層粒子密度(粒子占空比)的不斷增減,VRT 方程的一階解和二階解的后向散射系數(shù)都在不斷增大,但是無論是同極化還是交叉極化,VRT 二階的后向散射系數(shù)的增加幅度和趨勢明顯大于VRT 一階解.二階解同極化、交叉極化的散射系數(shù)都不同于一階解散射系數(shù)變化相對平緩,由于隨著粒子層粒子數(shù)量不斷增加,粒子間二次散射作用的概率不斷增大,使得粒子層的總散射作用逐漸增強.同時,粒子和粗糙面間的多次散射也是交叉極化散射系數(shù)逐漸增大的重要原因.在高密度粒子的粒子層中,二階解的多次散射對整體環(huán)境后向散射系數(shù)起到較大增強作用..

4.3 草地沙盤實測與VRT 模型的對比

算例4測量一個半徑0.7 m,厚0.1 m 草地沙盤的散射系數(shù),用該草盤模型的物理參數(shù)建立一個VRT 模型,如圖9 和圖10 所示.該VRT 模型的參數(shù)如下:粒子層厚10 cm,GRG 近似橢球粒子尺寸為h=5 cm,r=0.1 cm,粒子介電常數(shù)為(14.98,1.466)且占空比為0.05,歐拉角取向為,β,γ∈(0,10°).粗糙面模型為IEM,其介電常數(shù)為(26.7,11.7),粗糙度σ=0.3 cm,l=1.0 cm,電磁波頻率為10 GHz.

圖9 VRT 模型示意圖Fig.9.Grass of VRT model.

圖10 仿真結(jié)果和實測數(shù)據(jù)對比圖Fig.10.Comparison diagram of simulation results and measured data.

分別計算該模型的VRT 低階解和VRT 二階解的后向散射系數(shù),并與該草盤的實測數(shù)據(jù)進行對比.

圖10 仿真結(jié)果表明,VRT 建模仿真草地能較好地擬合所測得的實測數(shù)據(jù),VRT 高階解仿真結(jié)果曲線與實測結(jié)果基本吻合且無明顯波動,呈現(xiàn)單調(diào)遞減,且GRG 近似橢球取向角范圍較大,仿真粒子的隨機性能更好模擬實測沙盤草地的特性.

5 總結(jié)

本文將下墊面上方粒子層在Z軸方向上分成若干層,推導出了含有下墊面的VRT 方程的二階散射解.分析不同形態(tài)的粒子對散射系數(shù)的影響,比較了非球型粒子的一階、二階散射迭代解的差異.數(shù)值結(jié)果表明:取向均勻的粒子形態(tài),粒子間多次散射作用較弱,從而對粒子層散射系數(shù)的影響較小.對于在特定角度下存在較強散射作用的粒子形態(tài),粒子間多次散射作用較大,從而對粒子層散射系數(shù)的影響較大.隨著粒子數(shù)量增多,粒子層散射強度因多次散射作用影響較大,VRT 二階解對數(shù)值結(jié)果有較大幫助.該方法可以更好地描述粒子間、粒子與粗糙面之間的相干作用.