張 弛，吳 鑫，張建奇，李俊儒

（西安電子科技大學(xué)物理與光電工程學(xué)院，陜西西安710071）

1 引 言

霧是一種發(fā)生在大氣邊界底層的自然現(xiàn)象，它對(duì)能見(jiàn)度有重大作用，從而對(duì)人類(lèi)的正常生活產(chǎn)生一定的影響。另一方面，霧會(huì)對(duì)輻射產(chǎn)生散射和吸收作用，造成輻射衰減，嚴(yán)重影響光電類(lèi)探測(cè)設(shè)備的性能，大大降低了目標(biāo)探測(cè)的精準(zhǔn)度，容易引起虛警和漏警。因此，研究輻射在霧中的傳輸特性具有重要的意義，可為霧天氣下的光電探測(cè)與識(shí)別提供理論依據(jù)。在以往的研究中［1-3］，通常采用Beer-Lambert定律計(jì)算輻射在霧中的衰減，但該方法僅適用于能見(jiàn)度較高的天氣狀況。對(duì)于能見(jiàn)度較低的情況，霧滴粒子對(duì)輻射的多次散射作用明顯，Beer-Lambert定律的計(jì)算結(jié)果與實(shí)際情況有較大的誤差。為了實(shí)現(xiàn)輻射在霧中的透過(guò)率計(jì)算以及傳輸特性分析，必須考慮多次散射效應(yīng)的影響，準(zhǔn)確地模擬輻射在霧中的傳輸過(guò)程。

目前，常用的介質(zhì)中輻射傳輸?shù)挠?jì)算方法有：輻射傳輸方程求解法［4-5］和蒙特卡羅法［4，6-14］。輻射傳輸方程的求解過(guò)程十分復(fù)雜，需要進(jìn)行大量的近似簡(jiǎn)化才能計(jì)算出結(jié)果，因此適用范圍有限［6］。蒙特卡羅法作為一種基于統(tǒng)計(jì)學(xué)的方法，不需要求解輻射傳輸方程，而是將輻射在介質(zhì)中的傳輸看作大量光子與介質(zhì)粒子進(jìn)行碰撞的過(guò)程，光子在每次碰撞的過(guò)程中會(huì)被吸收或被散射，通過(guò)統(tǒng)計(jì)大量光子被吸收和散射的過(guò)程，可以得到較為穩(wěn)定的統(tǒng)計(jì)結(jié)果。與輻射傳輸方程求解法相比，蒙特卡羅法不僅能解決輻射在平面分層介質(zhì)中的傳輸問(wèn)題，還能解決球面介質(zhì)中的輻射傳輸問(wèn)題。由于其自身的靈活性、易實(shí)施性和高效性，蒙特卡羅法已被廣泛地應(yīng)用于氣溶膠中的輻射傳輸計(jì)算中。

Szymanski［4］等通過(guò)蒙特卡羅法和兩種基于輻射傳輸方程的近似算法分析了多次散射效應(yīng)對(duì)氣溶膠透過(guò)特性的影響。Hassan［7］等和Grabner［8］等分別采用蒙特卡羅法模擬了光在雨和霧中的傳輸，分析了雨霧的多次散射效應(yīng)對(duì)光通信的影響。Gierstad［9］等通過(guò)蒙特卡羅法研究了輻射在海洋氣溶膠中的傳輸。Van Der Laan［10］等利用蒙特卡羅法研究了偏振光在霧中的傳輸過(guò)程，指出與線偏振光相比，圓偏振光在霧中有更大的傳輸優(yōu)勢(shì)。國(guó)內(nèi)也有許多學(xué)者采用蒙特卡羅法研究了輻射在霧中的傳輸特性。王蓉蓉［11］等利用蒙特卡羅法分析了霧天氣中太赫茲波和紅外波的傳輸特性，指出相同傳輸距離下，太赫茲信號(hào)比紅外信號(hào)具有更好的透射能力。王婷［12］等人采用蒙特卡羅法研究了紅外輻射在霧中的傳輸，得到了透過(guò)率與傳輸距離以及能見(jiàn)度的關(guān)系。張澤宇［13］等通過(guò)蒙特卡羅法對(duì)3.8μm和1.55μm輻射在霧中的傳輸進(jìn)行了模擬，得到了透過(guò)率與傳輸距離的關(guān)系，研究結(jié)果表明，3.8μm輻射比1.55μm輻射更適合在霧中進(jìn)行水平鏈路的自由空間光通信。張碩［14］等利用蒙特卡羅法對(duì)水霧的多次散射進(jìn)行了模擬，分析了消光系數(shù)，水霧濃度以及不對(duì)稱(chēng)因子等參數(shù)對(duì)水霧消光性能的影響，指出提高水霧的濃度并減小粒子半徑可有效提高水霧的衰減作用。

以上研究［4-14］均基于Mie散射理論得到粒子群的平均光學(xué)參數(shù)，并作為蒙特卡羅法的仿真參量，因此光子被散射或被吸收的概率，以及介質(zhì)粒子的散射相函數(shù)在每次碰撞過(guò)程中是相同的，這與輻射在介質(zhì)中的實(shí)際傳輸過(guò)程不符。實(shí)際上，由于霧滴粒子的尺度大小并不完全相同，研究輻射在霧中的傳輸問(wèn)題時(shí)，采用動(dòng)態(tài)粒子參數(shù)選擇的研究方法，可以得到較為可靠的結(jié)果。另一方面，蒙特卡羅法需要發(fā)射大量的光子來(lái)保證計(jì)算精度，光子數(shù)的增加勢(shì)必會(huì)導(dǎo)致計(jì)算時(shí)間的增加。研究人員采用Henyey-Greenstein相函數(shù)來(lái)選擇光子的散射方向［6，9-10，13-14］，這會(huì)增加蒙特卡羅法的計(jì)算時(shí)間，降低計(jì)算效率。因此，建立一種既符合輻射傳輸物理過(guò)程，又能在不影響計(jì)算精度的情況下提高計(jì)算效率的蒙特卡羅法具有重要意義。

本文基于傳輸近似光學(xué)參數(shù)，提出一種基于動(dòng)態(tài)粒子參數(shù)選擇的蒙特卡羅法。與傳統(tǒng)的蒙特卡羅方法不同的是，該方法采用傳輸近似光學(xué)參數(shù)計(jì)算光子的碰撞步長(zhǎng)、光子被散射或吸收的概率，以及介質(zhì)粒子的散射相函數(shù)；根據(jù)介質(zhì)粒子的尺度分布計(jì)算介質(zhì)粒子與光子的碰撞概率，再根據(jù)碰撞概率對(duì)傳輸近似光學(xué)參數(shù)進(jìn)行動(dòng)態(tài)選擇。選取波長(zhǎng)分別為0.85μm和4μm的輻射在4種模式霧（重度平流霧、中度平流霧、重度輻射霧和中度輻射霧）中的傳輸進(jìn)行仿真，分析了透過(guò)率與傳輸距離之間的關(guān)系，比較了不同波長(zhǎng)下透過(guò)率的差異，為霧天氣下的目標(biāo)探測(cè)與識(shí)別等研究提供理論依據(jù)。比較了兩種蒙特卡羅法的計(jì)算時(shí)間，為相關(guān)工程應(yīng)用提供參考。

2 多粒子霧層及其光學(xué)參數(shù)

根據(jù)自然界中霧的形成機(jī)理和地域，通常將它分為兩大類(lèi)：平流霧和輻射霧。輻射在霧中的傳輸計(jì)算涉及霧滴尺度分布的建模以及光學(xué)參數(shù)的計(jì)算。實(shí)際上，作為一種大氣水凝物，霧主要由不同半徑的水滴粒子組成［6］。因此，本文將霧滴近似為球形水滴粒子。由于復(fù)雜的形成和消散過(guò)程，不同時(shí)間和不同地點(diǎn)的霧滴尺度分布差異很大。本文采用被廣泛應(yīng)用的Gammar模型作為多粒子霧層的霧滴尺度分布模型，同時(shí)研究輻射在4種霧（包括重度平流霧、中度平流霧、重度輻射霧和中度輻射霧）中的傳輸。Gammar霧滴尺度分布模型的表達(dá)式為［1］：

式中：n(r)是單位體積單位半徑內(nèi)霧滴粒子的數(shù)量，α，A和b是描述霧滴尺度分布的參數(shù)。4種霧的霧滴尺度分布參數(shù)如表1所示［15］。

表1 四種霧的霧滴尺度分布參數(shù)Tab.1 Parameters for particle distributions of four fogs

圖1 為4種霧的霧滴尺度分布。從圖中可以看出，輻射霧的霧滴數(shù)高于平流霧的霧滴數(shù)；與輻射霧相比，平流霧的霧滴數(shù)峰值向半徑較大的粒子偏移。輻射霧的霧滴數(shù)峰值主要集中在5μm附近，而平流霧的霧滴數(shù)峰值主要分布在10μm左右。

霧層作為一種半透明介質(zhì)，霧滴粒子對(duì)輻射的多次散射作用不可忽視。對(duì)于散射效應(yīng)較強(qiáng)的介質(zhì)，即使光學(xué)厚度足夠小，輻射傳輸計(jì)算仍然非常耗時(shí)。因此，為了縮短計(jì)算時(shí)間，提高計(jì)算效率，本文在輻射傳輸計(jì)算中采用傳輸近似法計(jì)算霧層的光學(xué)參數(shù)，包括傳輸近似吸收率αtr，傳輸近似衰減系數(shù)Atr和傳輸近似散射系數(shù)μtrs。各參數(shù)的計(jì)算公式如下［16-17］：

式中：μa和μs分別為粒子的吸收系數(shù)和散射系數(shù)；g為不對(duì)稱(chēng)因子，其取值在[-1，1]。g＞0時(shí)表示前向散射，g＜0時(shí)表示后向散射，g=0時(shí)表示各向同性散射。

輻射在傳輸過(guò)程中會(huì)與不同半徑大小的霧滴粒子進(jìn)行碰撞，在每次碰撞過(guò)程中，輻射被吸收或被散射的概率以及霧滴粒子的散射相函數(shù)會(huì)隨著與其碰撞的霧滴粒子的改變而變化。為了更加真實(shí)地模擬輻射在霧中的傳輸過(guò)程，本文將霧滴粒子的半徑劃分為多個(gè)區(qū)間，分別計(jì)算不同半徑區(qū)間內(nèi)霧滴粒子的光學(xué)參數(shù)。由于劃分的半徑區(qū)間越多，輻射傳輸計(jì)算所需要的時(shí)間越多。出于計(jì)算速度的考慮，本文以8個(gè)半徑區(qū)間為例，將平流霧的霧滴粒子半徑進(jìn)行如下劃分：0.01～1μm，1～2μm，2～5μm，5～10μm，10～15μm，15～18μm，18～30μm和30～60μm；對(duì)輻射霧的霧滴粒子半徑區(qū)間的劃為：0.01～1μm，1～2μm，2～3μm，3～4μm，4～5μm，5～8μm，8～10μm和10～60μm。第i個(gè)半徑區(qū)間霧滴粒子的αtr，i，Atr，i和μtr，is可以表示為：

圖1 四種霧的霧滴尺度分布Fig.1 Particle distributions of four fogs

其中：

式中：μa，μs和g的定義與式（2）中一致，上標(biāo)i表示第i個(gè)半徑區(qū)間；rmin，i和rmax，i分別為第i個(gè)半徑區(qū)間霧滴粒子的最小半徑和最大半徑。單個(gè)霧滴粒子的散射效率因子qsca，吸收效率因子qabs和不對(duì)稱(chēng)因子g均由Mie散射理論［18］計(jì)算得出。圖2為λ=0.85μm和λ=4μm時(shí)第i個(gè)半徑區(qū)間霧滴粒子的Mie吸收率［18］，αi=μia/(μia+μis)。從圖中可以看出，4種霧中的霧滴粒子對(duì)0.85μm輻射波的吸收極弱，對(duì)輻射的衰減基本是散射效應(yīng)的作用；與0.85μm輻射波相比，霧滴粒子對(duì)4μm輻射的吸收較大，吸收的作用增強(qiáng)，重度平流霧中半徑在30～60μm的霧滴粒子對(duì)輻射的吸收最大，但吸收率仍然小于0.3，這說(shuō)明散射的作用仍然比較明顯?？偟膩?lái)說(shuō)，無(wú)論是平流霧還是輻射霧，對(duì)于0.85μm輻射和4μm輻射，半徑為0.01～60μm的霧滴粒子的散射作用比較明顯。

圖2 λ=0.85μm和λ=4μm時(shí)第i個(gè)半徑區(qū)間霧滴粒子的Mie吸收率αiFig.2 Mie absorption probability of fog particles in ith radius region forλ=0.85μm andλ=4μm

式中m表示霧滴粒子的半徑區(qū)間數(shù)。表2為入射波長(zhǎng)分別為λ=0.85μm和λ=4μm時(shí)4種模式霧的總傳輸近似衰減系數(shù)從表2中可以看出，λ=4μm時(shí)，4種霧的衰減系數(shù)高于λ=0.85μm時(shí)的衰減系數(shù)；與同等程度的輻射霧相比，平流霧的衰減系數(shù)較大；同類(lèi)型霧中，重度霧的衰減系數(shù)高于中度霧的衰減系數(shù)。為了比較傳輸近似衰減系數(shù)與傳統(tǒng)的Mie衰減系數(shù)的差異，表3列出了波長(zhǎng)分別為λ=0.85μm和λ=4μm時(shí)4種霧的Mie總衰減系數(shù)比較表2和表3可知，對(duì)于λ=0.85μm和λ=4μm的入射波，4種霧的Mie總衰減系數(shù)均大于總傳輸近似衰減系數(shù)

表2 λ=0.85μm和λ=4μm時(shí)4種霧的總傳輸近似衰減系數(shù)Tab.2 Total transport attenuation coefficient of four fogs forλ=0.85μm andλ=4μm （km-1）

表3λ=0.85μm和λ=4μm時(shí)4種霧的Mie總衰減系數(shù)Tab.3 Total Mie attenuation coefficient of four fogs forλ=0.85μm andλ=4μm （km-1）

表3λ=0.85μm和λ=4μm時(shí)4種霧的Mie總衰減系數(shù)Tab.3 Total Mie attenuation coefficient of four fogs forλ=0.85μm andλ=4μm （km-1）

Wavelength/μm 0.85 4 Advection fog Heavy 29.02 31.11 Moderate 18.69 20.28 Radiation fog Heavy 16.98 21.18 Moderate 8.93 11.93

3 散射相函數(shù)模型

在輻射傳輸計(jì)算過(guò)程中，輻射的單次散射方向由介質(zhì)中粒子的散射相函數(shù)決定。散射相函數(shù)的精確度和復(fù)雜度直接影響輻射傳輸?shù)挠?jì)算結(jié)果和計(jì)算效率［19］。散射相函數(shù)p(θ，φ)用于描述入射能量經(jīng)粒子散射后的能量角分布，其中θ為散射角，φ為方位角。p(θ，φ)在所有方向上的積分滿(mǎn)足［19］：

常見(jiàn)的散射相函數(shù)計(jì)算方法有Legendre多項(xiàng)式展開(kāi)法，Henyey-Greenstein（H-G）相函數(shù)近似法［20］，以及Mie散射理論計(jì)算方法［18，20-21］。其中，Legendre多項(xiàng)式展開(kāi)法和H-G相函數(shù)近似法均是對(duì)Mie散射相函數(shù)的近似。Mie散射理論是均勻介質(zhì)球形粒子在單色平面波照射下遠(yuǎn)場(chǎng)散射的精確解，通過(guò)Mie散射計(jì)算單個(gè)球形粒子的散射相函數(shù)，可以保證計(jì)算的精度。單個(gè)粒子的Mie相函數(shù)的表達(dá)式為［18，22］：

式中：S1(θ)和S2(θ)分別表示垂直極化和平行極化方向的散射振幅函數(shù)，其具體表達(dá)式參照文獻(xiàn)［18］和［22］；χ=2πrn0λ，是粒子的尺度參數(shù)，其中r是粒子的半徑，λ是入射光波長(zhǎng)，n0是介質(zhì)的折射率。由于球形粒子的對(duì)稱(chēng)性，p(θ，φ，r)不會(huì)隨著方位角的變化而變化。

對(duì)于服從一定分布的多粒子霧層，需要結(jié)合粒子群的譜分布計(jì)算整個(gè)霧滴粒子群體系的Mie相函數(shù)［19］：

由于Mie散射理論在計(jì)算散射相函數(shù)時(shí)涉及復(fù)變量Bessel函數(shù)的計(jì)算，存在計(jì)算時(shí)間較長(zhǎng)、計(jì)算效率低下的問(wèn)題，因此在輻射傳輸計(jì)算中，通常采用相函數(shù)近似法來(lái)計(jì)算散射相函數(shù)。H-G近似相函數(shù)由于具有簡(jiǎn)單的表達(dá)式，并且能夠求解出散射角隨機(jī)抽樣的概率密度函數(shù)，因而被普遍地應(yīng)用于蒙特卡羅法中。傳輸近似相函數(shù)是一種常用于描述強(qiáng)前向散射特性的近似相函數(shù)，由描述各向同性散射分量piso和描述前向散射峰分量pf組成，其表達(dá)式為［16-17］：

式中：cosθ為散射角余弦，δ是狄拉克函數(shù)。圖3給出了不同波長(zhǎng)下r=8μm的霧滴粒子的Mie相函數(shù)、H-G近似相函數(shù)和傳輸近似相函數(shù)的分布結(jié)果。

圖3 不同波長(zhǎng)下r=8μm的霧滴粒子的Mie相函數(shù)，HG近似相函數(shù)與傳輸近似相函數(shù)的分布Fig.3 Distribution of Mie，H-G approximation and transport approximation phase functions of fog particles with r=8μm for different wavelengths

由圖3可知：H-G近似相函數(shù)與Mie相函數(shù)的分布曲線的變化趨勢(shì)大致相同，但在某些角度下的數(shù)值差異較大，表現(xiàn)在兩者前向散射峰值均較強(qiáng)，后向散射峰值均較弱，但H-G近似相函數(shù)的分布曲線相對(duì)平滑，出現(xiàn)差異的原因是Mie相函數(shù)是由麥克斯韋方程推導(dǎo)出的解析解，而H-G近似相函數(shù)將Mie相函數(shù)近似表達(dá)為與不對(duì)稱(chēng)因子有關(guān)的多項(xiàng)式。傳輸近似相函數(shù)與Mie相函數(shù)以及H-G近似相函數(shù)在數(shù)值上存在較大的差異，表現(xiàn)在-1≤cosθ＜1的范圍內(nèi)，傳輸近似相函數(shù)的數(shù)值為常數(shù)，分布為一條直線，但在cosθ=1的前向散射方向有明顯的散射峰，例如當(dāng)入射波長(zhǎng)為0.85μm時(shí)，在-1≤cosθ＜1內(nèi)，ptr恒為0.164；cosθ=1時(shí)，ptr=1.836，在θ=0°處的散射能量是在其余各方向上的11.2倍。從式（11）可以看出，傳輸近似相函數(shù)在-1≤cosθ＜1內(nèi)將散射看作各向同性散射，但在cosθ=1時(shí)存在一個(gè)散射峰。因此，即使傳輸近似相函數(shù)與Mie相函數(shù)存在較大差異，但仍可以用來(lái)合理地描述霧滴粒子較強(qiáng)的前向散射特性。文獻(xiàn)［23］指出粒子的單次散射分布特性對(duì)輻射在多次散射效應(yīng)較強(qiáng)的介質(zhì)中的傳輸特性影響較小。因此，考慮到霧是一種多次散射效應(yīng)較強(qiáng)的半透明介質(zhì)，即使傳輸近似相函數(shù)與Mie相函數(shù)差異較大，但由于其表達(dá)式簡(jiǎn)單，仍然可以將它應(yīng)用在蒙特卡羅法中。由于采用傳統(tǒng)蒙特卡羅法解決強(qiáng)散射介質(zhì)中的輻射傳輸問(wèn)題存在計(jì)算效率較低、耗時(shí)嚴(yán)重的缺點(diǎn)，本文為解決計(jì)算耗時(shí)問(wèn)題，將不同半徑區(qū)間內(nèi)的霧滴粒子的傳輸近似相函數(shù)應(yīng)用于蒙特卡羅法中來(lái)選擇粒子單次散射后的方向。第i個(gè)半徑區(qū)間的霧滴粒子的傳輸近似相函數(shù)ptr，i可以表示為：

4 動(dòng)態(tài)粒子參數(shù)選擇模型

基于動(dòng)態(tài)粒子參數(shù)選擇的蒙特卡羅法計(jì)算輻射在霧中的透過(guò)率的主要過(guò)程如下：以一定的方向發(fā)射光子，根據(jù)霧的衰減系數(shù)確定光子在兩個(gè)相鄰碰撞點(diǎn)的距離；通過(guò)霧中粒子的尺度分布和幾何尺寸選擇與光子進(jìn)行碰撞的粒子；根據(jù)粒子的光學(xué)參數(shù)確定光子被吸收的概率和光子被散射后的方向（若光子沒(méi)有被吸收）；根據(jù)光子的出射方向和數(shù)量計(jì)算透過(guò)率。該方法與傳統(tǒng)計(jì)算方法的區(qū)別在于：在每一次碰撞過(guò)程中，用于計(jì)算光子被吸收的概率和光子散射后的方向的光學(xué)參數(shù)并不完全相同，需在選擇與光子進(jìn)行碰撞的粒子后，再來(lái)確定用于計(jì)算的光學(xué)參數(shù)。

圖4 輻射在霧層中的傳輸示意圖Fig.4 Schematic of radiation transfer in fog layer

圖4 為輻射在霧層中的傳輸示意圖。霧滴的尺度分布在霧層各處均勻一致；輻射的入射天頂角為θin，入射方位角為φin，d為傳輸距離；被前向探測(cè)平面（z=d）收集到的光子即為透過(guò)的光子。基于動(dòng)態(tài)粒子參數(shù)選擇的蒙特卡羅法的具體計(jì)算步驟如下：

（1）光子的發(fā)射參數(shù)設(shè)定

每個(gè)光子的初始位置為O(0，0，0)，入射天頂角和入射方位角分別為θin和φin。

（2）光子碰撞步長(zhǎng)確定

光子在兩個(gè)相鄰的碰撞點(diǎn)之間的隨機(jī)步長(zhǎng)為［6］：

式中：ξ1是0～1之間均勻分布的隨機(jī)數(shù)而在傳統(tǒng)的蒙特卡羅法中［6-14］從式（13）可以看出，光子的碰撞步長(zhǎng)和霧的總衰減系數(shù)有關(guān)，總衰減系數(shù)越小，光子的碰撞步長(zhǎng)越大。

（3）光子與霧滴粒子的碰撞概率計(jì)算

光子與霧滴粒子的碰撞概率與霧滴粒子的尺度分布和幾何半徑有關(guān)。霧滴粒子的濃度越高，半徑越大，與光子的碰撞概率越大。因此光子與第i個(gè)半徑區(qū)間的霧滴粒子碰撞的概率可表示為：

從第1個(gè)半徑區(qū)間到第i個(gè)半徑區(qū)間的霧滴粒子累積碰撞概率為［24］：

當(dāng)光子移動(dòng)一個(gè)碰撞步長(zhǎng)時(shí)，通過(guò)在0～1均勻分布的隨機(jī)數(shù)ξ2來(lái)確定與光子進(jìn)行碰撞的霧滴粒子。若0＜ξ2≤PC，1，則光子與第1個(gè)半徑區(qū)間的霧滴粒子進(jìn)行碰撞；若PC，i-1＜ξ2≤PC，i，則光子與第i個(gè)半徑區(qū)間的霧滴粒子進(jìn)行碰撞。

（4）光子被吸收的概率分析

當(dāng)確定光子與第i個(gè)半徑區(qū)間的霧滴粒子進(jìn)行碰撞后，通過(guò)在0～1之間均勻分布的隨機(jī)數(shù)ξ3來(lái)確定光子是否被吸收。若ξ3≤αtr，i，則光子被第i個(gè)半徑區(qū)間的霧滴粒子吸收，停止跟蹤該光子；否則，光子會(huì)以新的散射方向繼續(xù)在霧層中傳輸。與傳統(tǒng)蒙特卡羅法不同的是，在基于動(dòng)態(tài)粒子參數(shù)選擇的蒙特卡羅法中，αtr，i會(huì)根據(jù)步驟（3）中選擇的第i個(gè)半徑區(qū)間的霧滴粒子進(jìn)行更新，每一次碰撞事件中的αtr，i不完全相同。

（5）散射方向確定

若光子沒(méi)有被霧滴粒子吸收，則會(huì)被散射到新的方向。將傳輸近似散射相函數(shù)轉(zhuǎn)化為累積概率密度分布函數(shù)CPDF(μ)（μ=cosθ），然后采用0～1之間均勻分布的隨機(jī)數(shù)ξ4對(duì)CPDF（μ）抽樣來(lái)獲得散射角θ。CPDF(μ)由下式得出［23］：

式中：F(g，μ)是不對(duì)稱(chēng)因子g和散射角余弦μ的函數(shù)μ，C是傳輸近似相函數(shù)在-1≤μ≤1內(nèi)積分求得的數(shù)值。CPDF(μ)是0～1之間的分布，在步驟（3）中選擇第i個(gè)半徑區(qū)間的霧滴粒子與光子進(jìn)行碰撞后，g和C已知，給定一個(gè)0～1之間均勻分布的隨機(jī)數(shù)ξ4，可以通過(guò)式（16）得到對(duì)應(yīng)的μ，進(jìn)而獲得散射角θ。方位角φ均勻分布在0～2π之間，其抽樣值為φ=2πξ5，ξ5是0～1之間均勻分布的隨機(jī)數(shù)。

（6）透過(guò)率計(jì)算

從霧層中逃逸并被前向探測(cè)器接收的光子即為透射的光子。通過(guò)統(tǒng)計(jì)透射的光子數(shù)量，可得到透過(guò)率為［25］：

式中：θout和φout分別是出射天頂角和出射方位角，Nin和Ntra分別是入射光子數(shù)和透射光子數(shù)。

5 計(jì)算結(jié)果及分析

設(shè)定入射波長(zhǎng)分別為λ=0.85μm和λ=4μm，光子數(shù)為108個(gè)，所有光子均為垂直入射（入射天頂角θin=0°，入射方位角φin=0°）。采用兩種蒙特卡羅法仿真不同入射輻射在4種霧中的傳輸：（1）傳統(tǒng)的蒙特卡羅方法［6，9-10，13-14］采用Mie光學(xué)參數(shù)以及H-G散射相函數(shù)，以下簡(jiǎn)稱(chēng)為“HG MC”；（2）基于動(dòng)態(tài)粒子參數(shù)選擇的蒙特卡羅法，以下簡(jiǎn)稱(chēng)為“DPS MC”。通過(guò)計(jì)算不同傳輸距離下的透過(guò)率和計(jì)算時(shí)間來(lái)分析輻射在霧中的傳輸特性以及蒙特卡羅法的計(jì)算效率。

圖5 為重度平流霧和中度平流霧透過(guò)率隨傳輸距離的變化曲線。由圖5可知：隨著傳輸距離的增加，兩種輻射波的透過(guò)率均逐漸減小，在相同的傳輸距離下，0.85μm輻射的透過(guò)率相對(duì)于4μm輻射的透過(guò)率更大，表現(xiàn)出良好的透射特性。這與圖2中的結(jié)果一致，由圖2可知，相對(duì)于0.85μm輻射，重度平流霧和中度平流霧對(duì)4μm輻射的吸收較大，更多的光子被霧滴粒子吸收，使得透過(guò)的光子減少；相同傳輸距離下，兩種輻射波在重度平流霧中的透過(guò)率更低，這與表2和表3中的結(jié)果一致，與中度平流霧相比，重度平流霧對(duì)兩種輻射的衰減更大。

兩種方法對(duì)透過(guò)率的計(jì)算結(jié)果具有相同的變化趨勢(shì)，說(shuō)明本文提出的DPS MC是可靠的。但由于DPS MC和H-G MC中采用的單次散射相函數(shù)不同，以及DPS MC考慮了不同半徑區(qū)間的霧滴粒子散射特性的差異性，兩種方法在數(shù)值上不完全一致。兩種方法對(duì)0.85μm輻射在平流霧中的透過(guò)率計(jì)算結(jié)果差異較明顯，對(duì)4μm輻射的計(jì)算結(jié)果則較為接近。當(dāng)傳輸距離為150 m時(shí)，兩種方法對(duì)0.85μm和4μm輻射在重度平流霧中的透過(guò)率的絕對(duì)誤差分別為0.048 3和0.001 5。兩種方法之間的絕對(duì)誤差在波長(zhǎng)上的差異也是合理的，從表2～表3以及式（13）可知，在重度平流霧中4μm輻射的衰減大于0.85μm輻射的衰減，使得光子的碰撞步長(zhǎng)減小，導(dǎo)致透射的光子經(jīng)歷的散射次數(shù)更多，而散射次數(shù)的增多，減少了散射相函數(shù)不同造成的透過(guò)率差異，兩種蒙特卡羅法的計(jì)算結(jié)果最為接近。

圖5 重度平流霧和中度平流霧的透過(guò)率隨傳輸距離的變化曲線Fig.5 Transmittances in heavy advection fog and moderate advection fog as function of propagation distance

圖6 為重度輻射霧和中度輻射霧透過(guò)率隨傳輸距離的變化曲線。由圖6可知，隨著傳輸距離的增加，兩種輻射波的透過(guò)率均逐漸減小。由圖2可知，相對(duì)于4μm輻射，重度輻射霧和中度輻射霧對(duì)0.85μm輻射的弱吸收作用，使得0.85μm輻射比4μm輻射在霧中具有更強(qiáng)的透射特性。相同傳輸距離下，兩種輻射波在重度輻射霧中的透過(guò)率更低，產(chǎn)生這種現(xiàn)象的原因與衰減系數(shù)的大小有關(guān)。相比于中度輻射霧，兩種輻射波在重度輻射霧中的衰減系數(shù)更大，使得更多的光子被吸收或者被散射到其他方向，更少的光子透過(guò)霧層被探測(cè)器收集。此外，對(duì)比圖5和圖6可知，在相同傳輸距離下，與同等程度的平流霧相比，兩種輻射波在輻射霧中的透過(guò)率更高，穿透能力更強(qiáng)。

圖6 重度輻射霧和中度輻射霧透射率隨傳輸距離的變化曲線Fig.6 Transmittance of heavy radiation fog and moderate radiation fog as functions of propagation distance

圖7 兩種蒙特卡羅法的計(jì)算時(shí)間隨傳輸距離的變化曲線Fig.7 Computational time of two Monte Carlo methods as function of propagation distance

圖7 為兩種蒙特卡羅法的計(jì)算時(shí)間隨傳輸距離的變化曲線。由圖7可知：兩種蒙特卡羅方法的計(jì)算時(shí)間均與傳輸距離成正比，隨著傳輸距離的增加，計(jì)算時(shí)間逐漸增大；相比于0.85μm輻射，兩種蒙特卡羅法對(duì)4μm輻射的計(jì)算時(shí)間較長(zhǎng)。在相同傳輸距離下，不同波長(zhǎng)計(jì)算時(shí)間的差異主要由衰減系數(shù)決定，衰減系數(shù)越大，透過(guò)的光子經(jīng)歷的散射次數(shù)越多，耗時(shí)越長(zhǎng)，這也是導(dǎo)致兩種方法對(duì)平流霧的計(jì)算時(shí)間更長(zhǎng)的原因；與H-G MC相比，DPS MC在較長(zhǎng)的傳輸距離下的計(jì)算速度有明顯的提升，隨著傳輸距離的增加，DPS MC在計(jì)算速度上的優(yōu)勢(shì)越來(lái)越明顯。入射波長(zhǎng)為4μm時(shí)，對(duì)d=100 m和d=400 m的重度平流霧H-G MC的計(jì)算時(shí)間分別為650 s和1 650 s，DPS MC的計(jì)算時(shí)間分別為610 s和920 s，分別縮短了6%和44%；相同傳輸距離下，DPS MC對(duì)重度平流霧在計(jì)算速度上的優(yōu)勢(shì)最為明顯，當(dāng)入射波長(zhǎng)為4μm，傳輸距離為400 m時(shí)，DPS MC較H-G MC對(duì)重度平流霧、中度平流霧、重度輻射霧和中度輻射霧的計(jì)算時(shí)間分別縮短了44%，33%，37%和23%。

6 結(jié) 論

本文建立了一種多粒子霧層模型，該模型中霧滴粒子為半徑不同的球形水滴粒子。霧滴粒子半徑被分為多個(gè)區(qū)間，計(jì)算了不同半徑區(qū)間內(nèi)霧滴粒子的傳輸近似光學(xué)參數(shù)。選擇0.85μm和4μm兩種輻射波以及重度平流霧、中度平流霧、重度輻射霧和中度輻射霧4種霧，采用基于動(dòng)態(tài)粒子參數(shù)選擇的蒙特卡羅法（DPS MC），模擬了兩種輻射波在4種霧中的傳輸過(guò)程，分析了它們?cè)陟F中的衰減特性，并將DPS MC的計(jì)算結(jié)果與傳統(tǒng)的蒙特卡羅法（H-G MC）的計(jì)算結(jié)果進(jìn)行了比較。結(jié)果表明：相比于4μm輻射，0.85μm輻射在霧中的衰減更小，表現(xiàn)出了良好的傳輸性；相比于同等程度的輻射霧，兩種輻射在平流霧中的衰減更加明顯；兩種蒙特卡羅法的計(jì)算結(jié)果隨傳輸距離變化的趨勢(shì)一致，但在數(shù)值上不完全相同，多次散射效應(yīng)越強(qiáng)，兩種方法在計(jì)算結(jié)果上的差異越??；與H-G MC相比，DPS MC可以有效地縮短計(jì)算時(shí)間，其計(jì)算優(yōu)勢(shì)在遠(yuǎn)距離傳輸以及多次散射效應(yīng)較強(qiáng)的介質(zhì)中傳輸?shù)那闆r下更加明顯。

DPS MC更加真實(shí)地反映了輻射在霧中的傳輸過(guò)程，并且能有效地縮短計(jì)算時(shí)間，為霧中的輻射傳輸計(jì)算提供了新的思路。在后續(xù)的研究工作中，將開(kāi)展不同入射輻射波在不同霧環(huán)境下的傳輸特性試驗(yàn)，進(jìn)一步探討霧滴粒子半徑區(qū)間的劃分與DPS MC的計(jì)算精度和計(jì)算效率之間的關(guān)系，以及該方法在相關(guān)工程應(yīng)用上的適用性。