基于蒙特卡羅的光散射傳輸三維動態仿真*

張陽陽，高 雋，范之國，吳良海，高麗娟

(合肥工業大學 計算機與信息學院，安徽 合肥 230009)

天空中存在許多不同偏振度、不同偏振化方向的散射光，便形成了特定的包含大量偏振信息的偏振分布模式[1]，這種偏振分布模式蘊含豐富的方位信息和大氣參數信息。這種偏振特性而主要是由于光在傳輸過程中與介質中粒子散射碰撞而產生的，研究光在散射介質中的傳輸和散射過程對偏振導航、探測都有重要的意義[2-3]。

光在散射介質中的傳輸是一個非常復雜的過程。典型的傳輸理論為輻射傳輸理論，該理論一般難以得到解析解。對散射介質中的光傳輸問題也可以用統計方法進行描述。目前解決光線在大氣中的輻射傳輸問題比較有效的方法是蒙特卡羅方法，其思路清晰，能方便地應用到實際計算中。WILSON和ADAM首次將蒙特卡羅思想引入激光同組織相互作用方面[4]，隨后，蒙特卡羅模擬方法就被廣泛地應用到光在組織中的傳輸模擬。SEBASTIAN B應用斯托克斯-穆勒(Stokes-Mueller)形式模擬了偏振光在半無限介質中的傳輸，將蒙特卡羅模擬方法的應用擴展到了偏振領域[5]。JESSICA總結并給出了半無限介質中的偏振光傳輸的三種蒙特卡羅模擬方法[6]，蒙特卡羅方法得到的結果依賴于對大量光子的統計，由于模擬得到的結果與實驗結果比較符合，因此蒙特卡羅方法被越來越多的科研工作者所采用。

本文采用蒙特卡羅統計方法模擬光束在散射介質中的傳輸和散射過程。同時，利用開源圖形函數庫OpenGL工具在VC++環境下設計實現了三維動態仿真，提供了三維空間中光散射傳輸的可視化系統，并對光散射傳輸后的偏振特性進行了分析。

1 偏振蒙特卡羅模型

英國物理學家Stokes于1852年提出利用一種四維向量(即Stokes矢量)，來描述光的偏振狀態，由于其具有可加性，因此非常適用于對光線的分析工作。用Stokes矢量來描述光的偏振狀態：S=[I Q U V]T,其中，I表示光術的總強度，±Q反映線偏振光水平或垂直的程度；±U反映線偏振光與Q成 ±45°的程度；±V反映光束中圓偏振光的左旋或右旋的程度。其偏振度可以表示為：

入射光Sin在介質層中傳輸時，輸出光為Sout，可以用一個4×4的穆勒矩陣來描述其中的變化過程，有：

其中,πn和 τn是散射角度 α 的勒讓德多項式，系數 an和bn是復折射率m=n+n′i(球形粒子相對于外界環境的折射系數)以及尺寸參數的函數。

偏振蒙特卡羅模型中粒子運動示意圖如圖1所示。其描述的是光子在半無限單一介質層中的散射傳輸過程，并用一定數量的光子依次循環來模擬現實光束散射的效果。

散射介質層可以看作是一個無限大但厚度有限的面板，在本文中厚度取4/us(us為散射系數)，介質層中為同一介質。如圖1所示，沿著Z軸正方向從Z=0(即X-Y面)到Z=4/us面為介質層，可以認為在發射點 O有一個光子發射器沿著Z軸正方向發射光子，在Z=4/us面上放置一個正方形接收器(圖中灰色區域所示)，本文中接收器的邊長取14/us。為了較真實地反應一束光穿過介質層之后偏振態的變化情況，在本文中，如果沒有特殊說明，取1e6個光子連續發射。接收器取100×100個方格點采樣接收。

蒙特卡羅方法的算法流程圖如圖2所示，其中灰色步驟只有在偏振蒙特卡羅模型中需要。

圖2 蒙特卡羅算法流程圖

下面對蒙特卡羅方法的關鍵步驟做簡要解釋（可參考JESSICA總結的偏振蒙特卡羅子午面模型[6]）：

(1)初始參考平面為X-Z平面;

(2)啟動時刻光的E矢量是在X-Z平面內定義的;

(3)根據已知的散射物質和拒絕方法的相位函數，隨機確定角度α(散射角度)和角度β(入射平面和散射平面的夾角)。

(4)根據散射參數確定穆勒矩陣計算出射光;

(5)在光子散射出介質層后，若打在接收器上，進行累計至所有光子散射完。

另外，當光子其權重低于閾值時，光子運動終止。

2 光散射傳輸的三維動態仿真系統設計實現

按照上述的偏振蒙特卡羅模型，利用OpenGL圖形庫和VC++工具設計實現光在介質中散射傳輸的三維動態仿真系統，圖3所示為仿真系統的功能框圖。

仿真系統采用MFC框架的單文檔結構，使用雙線程來處理散射傳輸計算與可視化交互，界面線程負責光子散射傳輸過程的三維動態顯示，結合OpenGL圖形庫，將設備描述符DC和渲染描述符RC聯系起來，以完成顯示功能；后臺工作線程負責光子散射傳輸過程的路經計算，根據所述的偏振蒙特卡羅模型實現光子能量、散射方向、偏振態等信息的計算，兩個線程通過消息隊列相互通信，將后臺計算的光子路徑和接收器偏振信息顯示在三維視圖中。圖4所示為仿真系統演示圖。

利用此三維動態仿真系統可以直觀地在三維空間中查看光子在介質中的散射傳輸過程和接收器的累積過程，且操作方便，可擴展性強。

3 系統仿真結果分析

利用所述仿真系統來模擬光束散射傳輸的過程，查看分析接收器上接收累積的光子偏振特性，即前向散射的偏振特性。參數為：入射光是自然光[1,0,0,0]，粒子半徑 1.0 μm，光子波長 0.6 μm，粒子復折射系數(1.59+0i)。圖5所示分別為接收器上 Stokes矢量（I,Q,U,V）的各個參數、偏振度DOP和Q分量的退偏振度DeP-Q。

從圖5中可以發現，當入射光為自然光[1,0,0,0]時，經過介質散射傳輸后，強度I主要集中在中心區域且從中心向外擴散，這是符合一般常理的；對于Q分量和U分量，在相應區域均有一個線偏振較大的區域，其他區域基本趨于零；圓偏振強度V分量的分布比較散亂，無規律，且強度較低；對于偏振度OPP和Q分量的退偏振度DeP-Q分布基本上是以中心點對稱的，且退偏振度較大的區域退偏現象較明顯，故相應區域偏振度較低。

U分量和Q分量的分布結果呈一定的規律性，可以將其利用于偏振應用的相關方面。對于V分量，含量較少，且獲取有一定的難度，故一般在偏振信息獲取時可以將其忽略。結果中多數圖像都是對稱的，這些對稱性包含了基于蒙特卡羅的光散射傳輸的一些偏振特性，值得后期更詳細地研究。

文中采用斯托克斯-穆勒(Stokes-Mueller)形式描述偏振光散射傳輸過程，介紹了偏振蒙特卡羅模型，并以

此為基礎，利用VC和OpenGL工具設計實現了三維

動態仿真系統，可以直觀地在三維空間中查看光子在介質中的散射傳輸過程和接收器的累積過程，為光散射傳輸過程中的偏振特性分析提供了可視化工具。同時，對系統仿真的存儲結果進行了分析，數據結果表現了偏振信息的對稱性，具有一定的后續研究價值。

[1]COULSON K.Polarization and intensity of light in the atmosphere[M].Hampton VA:A.Deepak Pub，1988.

[2]黃旭鋒，步揚，王向朝，等.基于米氏散射理論的太陽光散射偏振特性[J].中國激光，2010，37(12):3002-3006.

[3]吳良海，高雋，范之國，等.大氣粒子散射特性及其對空間偏振分布的影響[J].光學學報，2011，31(7):27-33.

[4]WILSON B C，ADAM G.A Monte Carlo model for the absorption and flux distributions of light in tissue[J].Med Phys，1983，10(6):824-830.

[5]SEBASTIAN B，HIELSCHER A H.Monte Carlo simulations of the diffuse backscattering mueller matrix for highly scattering media[J].Appl Opt，2000，39(10):1580-1588.

[6]JESSICA R R，SCOTT P，STEVE J.Three Monte Carlo programs of polarized light transport into scattering media:part I[J].Opt.Express，2005，13(12):4420-4438.