耿興寧,劉 政,李武周,許 宏,蔡 軍,陳科亦
(1.電磁空間安全全國重點實驗室,天津 300308;2.光電對抗測試評估技術重點實驗室,河南 洛陽 471003;3.中國人民解放軍93046部隊,遼寧 沈陽 110000)
大氣湍流對在其中傳播的光束的直接作用是基于空氣密度的不均勻分布而產生的折射率起伏變化。折射率的空間和時間起伏使得波前在傳播過程中不能保持一致的相位關系,出現相位變化,對光線產生了扭曲、發散、延遲等效應,進而影響光波信號的傳播,造成成像質量、功率密度以及穩定性下降等問題。理論分析建模是研究大氣湍流的一種重要手段[1]。
在對光波在大氣湍流中傳輸的仿真中,模擬大氣湍流最簡便的方法是將真空傳播與介質折射率起伏視為兩個獨立過程分別處理。即在光波傳播的路徑上劃分一定厚度的平行平板,入射光波處于平板的前表面,制造符合大氣湍流特性的相位屏,置于平板的后表面,將相位起伏疊加在真空傳播的相位之上,如圖1所示。

圖1 相位屏的作用
圖1中(a)為初始波前經過一段湍流Δz后,波前發生了畸變,(b)為初始波前在真空中傳播同樣的Δz后,再經過相位屏疊加上相位變化后的波前。適當的調整相位屏的結構可使(a)、(b)的出射波前相同。
產生相位屏的方法主要有基于傅里葉變換的譜反演法[2-8]、澤尼克多項式法[9-11]、分形方法等[12-14],澤尼克多項式方法是以單位圓內連續正交的多項式為基函數進行模擬獲得波的相位分布,但其在湍流的高空間頻率分量上的模擬有一定的局限性,并且對外尺度效應不能很好地表現出來;分形法是在相位結構函數的基礎上逐級插值得到屏上各個點的相位值,但隨著插值迭代的次數增加,準確性會降低[15]。
傅里葉變換譜反演法的基本思想是用一個復高斯矩陣對大氣湍流的功率譜進行濾波,然后通過逆傅里葉變換得到大氣導致的畸變相位,運算速度快,應用相對廣泛。但是這種方法的顯著缺點是低頻和高頻空間分量的欠采樣。其中高頻分量的誤差影響較小,可以忽略,但低頻分量對應的功率譜幾乎都受到影響,因此需對低頻分量進行補償[16]。廣泛使用的補償方法為次諧波法[17]。
以上對于大氣湍流的研究大多集中于靜態湍流,對于動態湍流以及激光束在動態湍流中傳輸的研究還相對較少。在此基礎上,本文采用傅里葉變換的方法得到大氣湍流模型,并通過疊加低頻次諧波來改善相位屏的低頻特性,之后再模擬激光在動態大氣湍流中的傳輸,研究不同湍流強度下光斑的畸變以及接收到的激光功率波動。
大氣湍流產生的折射率空間隨機分布可用一個隨機場表示。隨機場空間統計特性一般用結構函數表示。即:
Dn(r,r′)=〈[n(r)-n(r′)]2〉
(1)
式中,r、r′為空間兩點的矢量坐標;n為空間的折射率分布函數。
公式(1)表示湍流的結構函數是空間中各點的折射率方差。也就是說,結構函數表示了空間中各點的相關性。
在大多數情況下的湍流可作為Kolmogorov湍流看待。該理論認為存在大氣湍流的內尺度l0和外尺度L0,在這個區域之間的湍流是各向同性的,即:
Dn(r,r′)=Dn(|r-r′|)
(2)
公式(2)的意思是在Kolmogorov尺度內,大氣中任意兩個距離相同的點之間的結構函數值是相同的。則結構函數Dn(r,r′)可以簡化為一個固定的原點與矢量r之間的關系Dn(r),如圖2所示。

圖2 光傳播截面湍流尺度
圖2中Z軸表示光波傳播方向,X-Y平面為垂直于傳播方向的一個截面,在此截面中,以半徑為|r|的同心圓上各點與0點之間的結構函數都可以用Dn(r)表示,根據Kolmogorov理論,在該區間內的結構函數滿足2/3冪律,它的折射率結構函數可寫為:
(3)

+2.7×10-16exp(-h/1500)+Aexp(-h/100)
式中,h單位為m;A通常設為1.7×10-14m-2/5;v取為21 m/s,稱為HV-21模型。
湍流折射率起伏的三維譜密度為:
(4)


本文采用傅里葉變換譜反演法生成大氣相位屏。即利用湍流空間譜模型產生相空間隨機場并進行傅里葉變換獲得二維相位的空間分布。
首先按照上節的介紹,確定相位的頻譜函數Φn(K)。據此構造出一個二維復隨機場:
(5)

(6)

將(5)式作傅里葉變換可獲得一個二維隨機相位場:
(7)
(7)式即為在二維平面上的相位屏。
在構造中要注意一些參數的范圍,如相位屏所代表的平板厚度Δz范圍應有:
(8)
只有相位屏的間距滿足公式(8)式,才能保證建立的相位屏各自獨立,并且使得在其間的光傳播滿足幾何光學近似,這有這樣才能在公式(6)中用平板厚度代表對厚度的積分。
將一個要放置相位屏的空間平面分割為網格,設每一格寬度Δx,每邊共有N格,從采樣定理可知,最低和最高采樣頻率必須滿足:
(9)
為簡便起見,可用2π/(NΔx)和π/Δx分別作為Kmin和Kmax,而波數間隔ΔK設為Kmin。這樣可將公式(7)用離散化的形式寫出:
exp[-2πi(mp+nq)/N]
(10)
其中,
(11)
通過傅里葉變換譜反演法生成相位屏相較于其他方法更加便捷,但其空間頻率的最小和最大值分別是fmin=1/L,fmax=1/2Δx,因此相位屏沒有[-Δfx/2,Δfx/2]、[-Δfy/2,Δfy/2]這個頻段內的分量相應的功率譜[20],使得對相位屏模擬的準確性降低。因此我們需要對相位屏在低頻功率譜上進行補償。
本文采用次諧波法對大氣湍流相位屏進行低頻校正,它的主要方法是在低頻區域使用較小的采樣間隔,以實現更為精確的功率譜近似,其表達式為:
(12)
式中,p次諧波級數;f(m′,n′)這諧波函數。
(13)
將式(13)與式(10)合并,就構成了總的相位屏相位分布,從而實現了低頻的相位補償。
對于計算得到的相位屏的可信度,可用大氣湍流的相位結構函數來驗證,它是對二階相位統計特性的描述。
對于Kolmogorov譜,相位結構函數為:
(14)
式中,r0為大氣相干長度,相干長度的定義為:
(15)
式中,k=2π/λ;R為路徑長度;θz為天頂角。
現實中的大氣湍流是在持續變化的,折射率在時間上和空間上的分布都是隨機的,因此模擬大氣湍流時需要在靜態分布的基礎上加入湍流隨時間的變化。目前對于動態相位屏的構造主要有兩種方法,一種方法是樣條插值法,它的基本原理是根據傅里葉變換的平移特性,在構造相位屏時直接將時間變化的特征加進去;另一種方法是湍流凍結法,它是基于Taylor凍結假設,認為在較短的時間內,湍流在不改變折射率分布的情況下運動,整個相位屏在觀察孔徑內剛性移動,而大氣的內部運動可以被忽略,因此這種方法需要構造很大的相位屏,再隨時間平移截取一系列相位屏。
本文采用湍流凍結法實現動態湍流仿真。所考慮的波前相位可表示為S(r,t),其中t表示時間。考慮到風速對相位屏的漂移作用,引入風速v,在Taylor假設下,時間間隔τ內:
S(r,t+τ)=S(r-vτ,t)
(16)
則結構函數可以表示為:
Dn(vτ)=〈[S(r,t)-S(r-vτ,t)]2〉
(17)
表征湍流時間特性的一個參數是Greenwood頻率,它是與風速、大氣相干長度有關的物理量,簡化后的計算式為[21-22]:
(18)
按照上述的計算方法,設大氣湍流內尺度2 cm,外尺度50 m,相位屏尺寸0.4×0.4 m,分為100×100的網格,按照式(12)計算,在大氣路徑以45°仰角傳播1 km,激光波長1.064 μm時按照HV模型構造了強湍流相位屏,這是一個二維復數場,此時r0=0.03 m。疊加和未疊加次級諧波補償的相位屏分別如圖3(a),(b)所示。

圖3 譜反演法生成的相位屏
圖3中圖(b)為疊加了3級次諧波補償的相位屏,相較于未疊加次諧波補償的相位屏圖3(a),圖像更為平滑,表明其低頻成分能較好地顯示出來。
圖4為低頻補償前后的相位結構函數與理論值的對比,可以更好地看出,未疊加次諧波時的相位結構函數與理論值的誤差隨著距離的增加而逐漸增大,說明其存在著低頻成分的缺失,而補償后的相位結構函數則與理論值的符合較好,低頻分量得到改善。次諧波的級數更高會使得結構函數符合更好,但計算量也會相應增大,因此疊加的次諧波級數不宜過高。

圖4 r0=0.03 m的相位屏結構函數
動態相位屏在10~20 ms的時間尺度內滿足Taylor凍結假設,并且風速需要滿足v=0.314r0/τ0,其中τ0為大氣相干時間,只有在τ0內的大氣湍流滿足Kolmogorov統計特性[23-24]。由此,設置每秒能接收到500張圖像,此時的風速約為0.47 m/s。
按照S型方式旋轉選取相位屏,子相位屏尺寸為0.08×0.08m,記錄下連續變化的子屏圖像即可得到動態相位屏。其中的連續8幀相位屏如圖5所示。

圖5 連續變化動態相位屏
由圖5可以看出,連續變化的相位屏之間可以較為平緩地變化,具有一定的相關性。
設置光束為基模高斯光束,首先根據激光傳輸模型計算其通過相位屏前的光強分布,設1064 nm激光發射口徑100 mm,功率100 W,發散角0.05 mrad,在沒有大氣湍流干擾的情況下,在3 km處的接收平面上產生的光斑如圖6所示。

圖6 無湍流擾動時的光強分布
光斑分別經過動態弱湍流、中湍流和強湍流后的連續8幀能量分布分別如圖7~圖9所示。

圖7 弱湍流下的動態光斑功率分布變化

圖8 中湍流下的動態光斑功率分布變化

圖9 強湍流下的動態光斑功率分布變化
可以看出,隨著大氣湍流強度的增加,激光傳輸經過后更加分散,其畸變和抖動也更加嚴重。弱湍流時光斑基本沒有變化,中湍流時光斑出現一定程度的畸變,而在強湍流時的光斑畸變最大,已經沒有了本來的形狀。
為了更好地描述不同強度的湍流對激光束傳播的影響,本文計算了在連續400幀中的相同激光經過不同強度的湍流后落在接收面內的功率波動,如圖10所示。可以看出,弱湍流中的激光功率起伏較小,基本穩定在98 W左右;中湍流的激光功率有一定的起伏,大部分在94~97 W范圍內波動;強湍流中的激光起伏較大,在20~80 W的范圍內均有分布,功率起伏的方差分別為0.016、0.594、332.517。同時在整體上,接收到的激光功率隨著湍流強度的增加而減小,說明高強度的湍流對激光的衰減更大,造成的功率發散效應更嚴重。

圖10 接收面內的激光功率
本文基于傅里葉變換的譜反演法建立了靜態大氣湍流模型并利用次諧波方法進行低頻校正,在此基礎上通過湍流凍結法建立了動態大氣湍流模型,模擬實際情況下的大氣湍流變化。模擬了激光在不同強度的動態大氣湍流中的傳輸,仿真了激光通過湍流后的光斑分布情況,并計算了接收面內的激光功率。計算結果表明:功率為100 W,發散角0.05 mrad的1064 nm激光束,在相同的作用時間內經過不同湍流強度的大氣湍流相位屏,湍流越強,激光受到的影響就越明顯,光束質量就越差,并且接收到的功率降低,起伏增加;弱湍流時光斑基本無變化,接收到的功率穩定在98 W左右;中湍流的光斑有一定程度的畸變,功率在94~97 W范圍內波動;強湍流的光斑畸變最大,功率在20~80 W的范圍內波動。綜上,本文對研究動態大氣湍流建模仿真,以及激光在真實情況下的大氣湍流中的傳輸提供了一定的參考作用。