商景誠 吳濤 何興道 楊傳音
(南昌航空大學(xué)測試與光電工程學(xué)院,江西省光電檢測技術(shù)工程實驗室,無損檢測與光電傳感技術(shù)及應(yīng)用國家地方聯(lián)合工程實驗室,南昌 330063)
目前,隨著人們對瑞利-布里淵散射的深入研究[1?3],基于瑞利-布里淵散射激光雷達對氣體參數(shù)進行遙測的技術(shù)逐漸被應(yīng)用于環(huán)境監(jiān)測、國防工業(yè)和空間工程等領(lǐng)域.同時,隨著我國航空發(fā)動機研究的自主化,對航空發(fā)動機各類參數(shù)計量測試的發(fā)展提出了新的要求.溫度、壓強等高速氣流參數(shù)是評價吸氣式發(fā)動機的重要指標(biāo),傳統(tǒng)對這些參數(shù)的測試是在流場內(nèi)放置溫度、壓力等傳感器的侵入式的測試方法,該方法中的探測器不僅會因受到高溫氣體的影響縮短壽命,還會對發(fā)動機內(nèi)部流場造成干擾,影響測試精度[4],而通過對發(fā)動機尾焰的瑞利-布里淵散射的測量就可以實時獲得高溫氣體溫度和壓力的廓線分布,實現(xiàn)無擾測量.氣體自發(fā)瑞利-布里淵散射主要產(chǎn)生于氣體分子隨機熱運動引起的密度波動,這種熱密度波動可以被認為是一種熱聲波.氣體的許多物理參數(shù)(溫度、壓強、體黏滯系數(shù)等)都可以通過分析瑞利-布里淵散射光譜獲得[5?7].Lock等[8]在1992年利用瑞利-布里淵散射估算了氣體的溫度和與氣體密度數(shù)成正比的無量綱的參數(shù);2004年,Pan等[9]利用測得的瑞利-布里淵散射光譜獲得N2,O2及CO2的體黏滯系數(shù);2014年,Witschas等[10]利用在大氣中測量的瑞利-布里淵散射光譜實現(xiàn)了大氣溫度的反演.上述工作大都是基于Tenti S6模型[11]的卷積光譜來實現(xiàn)的,雖然科研工作者也提出了一些新的模型來描述瑞利-布里淵散射,但是這些模型的準(zhǔn)確性相對于S6模型還有待進一步提高,S6模型仍然是目前公認的描述瑞利-布里淵散射最準(zhǔn)確的模型.目前文獻報道的多為利用測量的高純氣體的瑞利-布里淵散射信號和S6模型的比較來實現(xiàn)氣體溫度和體黏滯系數(shù)等參數(shù)的反演[12,13],而對于氣體壓強,特別是對于在氣溶膠存在的條件下氣體壓強的反演研究卻很少有報道.在實際大氣中不可避免地存在氣溶膠,利用瑞利-布里淵散射實現(xiàn)對含有氣溶膠氣體壓強的反演研究,不僅可以檢驗Tenti S6模型的準(zhǔn)確性,而且能夠為高壓氣體壓強無擾測量方面的實際應(yīng)用提供參考.另外,在已有報道的研究中多是利用理想S6模型與儀器傳遞函數(shù)的卷積光譜對測量光譜直接擬合對比來分析瑞利-布里淵散射,而通過對測量光譜進行反卷積處理后直接與理想S6模型對比并實現(xiàn)壓強反演的研究卻少見.通過反卷積對光譜分析,不僅可以消除儀器函數(shù)對光譜線型的影響,提高光譜分辨率,而且可以提高氣體參數(shù)反演測量的準(zhǔn)確性.
本文利用卷積和反卷積的方法分別對氣體瑞利-布里淵散射光譜進行了理論和實驗研究,測量了存在氣溶膠時,N2氣體在壓強分別為2,4和6 atm(1 atm=1.01325×105Pa),溫度為297 K(室溫更容易得到高精度的穩(wěn)定性控制,避免溫度波動對測量結(jié)果產(chǎn)生影響)條件下90°散射方向的自發(fā)瑞利-布里淵散射光譜,采用卷積和反卷積方法對測量光譜進行擬合來反演氣體壓強,并將兩種方法反演得到的壓強值與氣池壓力計記錄的壓力值進行比較,驗證反演結(jié)果的準(zhǔn)確性,最后計算了求和歸一化后的均方根誤差.
布里淵散射是由氣體分子或原子運動引起的非彈性散射,散射光的頻率相對于入射光會發(fā)生變化,即產(chǎn)生布里淵頻移ΔvB,其可由下式來表達:
為了實現(xiàn)對瑞利-布里淵散射更為詳細的研究,引入了兩個無量綱的參數(shù)x和y:
式中,ω是散射光角頻率的變化(;k)是相互作用的波矢量的絕對值,;θ為信號散射角;p為氣體壓強;η是剪切黏度.根據(jù)y值的變化,可以將瑞利-布里淵散射分為三個區(qū)域.當(dāng)y? 1時,散射發(fā)生在克努森區(qū)域,該區(qū)域內(nèi)的瑞利-布里淵散射光譜可由高斯函數(shù)來近似描述;當(dāng)y? 1時,散射發(fā)生在流體力學(xué)區(qū)域,在該區(qū)域的光譜包含一個中間的瑞利峰和兩個邊帶的布里淵峰;當(dāng)y~1時,該區(qū)域為動力學(xué)區(qū)域,該區(qū)域的光譜是克努森區(qū)域的光譜與流體力學(xué)光譜的復(fù)雜混合,此時需要解玻爾茲曼方程[14].通過對線性近似的玻爾茲曼方程求解,Tenti通過6個基本的函數(shù)實現(xiàn)了對分子碰撞積分的描述.這6個函數(shù)包含了氣體的物理參數(shù),如熱傳導(dǎo)Tr、剪切黏度η、體黏度ηB及分子的比熱容Cint.在我們的S6模型擬合程序中,令Tr=25.2×10?3W/(m·K),η=17.6×10?6kg·m?1·s?1,ηB=14.0× 10?6kg·m?1·s?1,Cint=1[15].
實驗中儀器測得的瑞利-布里淵光譜是理想的散射光譜與儀器傳遞函數(shù)卷積的結(jié)果,其可以表達為
式中I(f)為測量的瑞利-布里淵散射光譜,S(f)為理想的瑞利-布里淵散射信號,由Tenti S6模型來描述,?表示卷積運算,w(f)表示儀器的傳遞函數(shù)(Fabry-Perot(F-P)掃描干涉儀).
考慮鏡面缺陷的F-P掃描干涉儀的傳遞函數(shù)w(f)可由改進后的艾里函數(shù)給出[15]:
式中FSR為F-P掃描干涉儀的自由光譜范圍,σg為鏡面缺陷參數(shù).
由于在實驗中引入了氣溶膠,(3)式中瑞利-布里淵散射信號需考慮氣溶膠散射的影響.根據(jù)已有的文獻[16],米散射可由Dirac delta函數(shù)來描述,因此,(3)式修改為
式中Smol(f)是Tenti S6模型描述瑞利-布里淵散射的歸一化的光譜線型,Spar(f)是描述粒子散射的Dirac delta函數(shù),Imol和Ipar是Smol(f)與Spar(f)各自對應(yīng)的強度.
根據(jù)(3)和(4)式及S6模型對壓強為2,4和6 atm的理想純凈N2氣體進行仿真,圖1顯示了不同壓強下Tenti S6模型與系統(tǒng)儀器函數(shù)卷積的結(jié)果.由于瑞利-布里淵散射光譜具有對稱性,我們僅給出了頻率大于0部分的光譜.從圖1可以看出,改變氣體壓強可以增強散射信號強度,而對布里淵頻移幾乎沒有影響,這主要是由于散射信號強度與壓強成正相關(guān),而受溫度的影響較小,而布里淵頻移與溫度成正相關(guān),與壓強無關(guān).
圖1 壓強為2,4和6 atm下卷積后的S6模型的光譜仿真Fig.1.Stimulation of S6 model convolved with instrument function corresponding to 2,4,and 6 atm.
高功率532 nm激光器激發(fā)的N2氣體在90°散射方向自發(fā)瑞利-布里淵散射譜測量的實驗裝置如圖2所示.由激光器(Verdi V-10,Coherent)發(fā)出的功率為10 W、線寬為5 MHz的532 nm激光被分光鏡M1按95:5的比例分成兩束光.其中一束5%的光束作為參考光用來校準(zhǔn)F-P掃描干涉儀和一系列光具組元件,在進行光譜測量時要將該光束擋住.實驗中使用的是共焦腔式的F-P掃描干涉儀,其鏡面反射率為99.5%,自由光譜范圍是10 GHz,半高全寬為105 MHz.另外95%的光束作為激發(fā)光,其經(jīng)過反射鏡M3被焦距為20cm的凸透鏡L1聚焦到散射池的中心后(焦斑直徑約60μm)與氣池中氣體相互作用產(chǎn)生瑞利-布里淵散射.90°散射方向的信號經(jīng)過透鏡組合L2,L3準(zhǔn)直以及直徑為100μm的針孔濾波器S2的空間濾波,由透鏡組合L4,L5再次準(zhǔn)直聚焦到F-P掃描干涉儀中,在掃描控制儀的控制下,F-P掃描干涉儀實現(xiàn)光譜掃描,掃描鑒頻后的信號經(jīng)焦距為10cm的雙膠合透鏡L6聚焦到光纖耦合器后被光子探測器(SPCMAQRH-14,Perkin-Elmer)探測,并經(jīng)過光子計數(shù)卡(P7882,Fast ComTech)傳輸?shù)接嬎銠C中進行光譜采集分析.散射池被設(shè)計成帶有布儒斯特窗的梯形形狀,散射池窗片雙面都鍍有增透膜,用來減少激光能量在散射池內(nèi)的損失.為了實現(xiàn)對氣池溫度的控制,氣池上裝有溫控系統(tǒng).掃描控制儀的控制信號和光子計數(shù)卡的采集信號都可以在示波器上顯示.實驗中利用標(biāo)準(zhǔn)氣溶膠粒子發(fā)生器向散射池中加入氣溶膠,通過控制發(fā)生器產(chǎn)生氣溶膠粒子的時間來控制氣溶膠的加入含量.
利用卷積后的S6模型對N2氣體在壓強為2,4和6 atm下測量光譜的擬合結(jié)果如圖3所示.為了對測量光譜中氣溶膠含量進行量化,采用米散(射的強)度與瑞利峰強度比值的百分比來表示計算得到壓強為2,4和6 atm對應(yīng)米散射的量化值為70.58%,161.48%和67.45%.由于瑞利-布里淵散射光譜具有對稱性,圖3僅給出了頻率大于0部分的擬合光譜.從圖3可以看出,卷積后的S6模型與實驗測量的光譜能夠較好地符合,誤差較小且分布比較均勻.為了詳細地評估擬合效果,進一步引入了求和歸一化的均方根誤差,其表達式為
式中Ie(fi)和Im(fi)分別為實驗測量的瑞利-布里淵散射光譜強度和卷積后的S6模型的光譜強度,Iint為歸一化后的光譜強度的積分值,N為獨立頻率的樣本數(shù).
圖2 測量N2在90°散射方向自發(fā)瑞利-布里淵散射光譜的實驗裝置Fig.2.Setup of measuring spontaneous Rayleigh-Brillouin scattering spectra of N2at the scattering angle of 90°.
圖3 卷積后的S6模型對測量光譜的擬合Fig.3.Fitting measured spectrum with the convolved S6 model at the pressure of 2,4,and 6 atm respectively.
為了消除儀器傳遞函數(shù)對測量光譜的影響,提高光譜分辨率,實現(xiàn)測量瑞利-布里淵散射光譜與理論光譜(S6模型)的直接比較,我們利用維納濾波器對測量的散射光譜進行了反卷積處理.維納濾波器的基本思想是濾波輸出與期望輸出結(jié)果可以實現(xiàn)最佳逼近[19].在本實驗中,F-P掃描干涉儀的傳遞函數(shù)w作為維納濾波器的輸入信號,其理想的輸出函數(shù)為沖擊函數(shù)δ.根據(jù)維納濾波器的原理,濾波輸出與理想輸出之間的誤差可以表達為
式中G為維納濾波因子.
當(dāng)對(7)式求一階偏導(dǎo)為0時即可得到矩陣形式的Toeplitz方程.考慮到儀器傳遞函數(shù)的自相關(guān)、δ函數(shù)的互相關(guān)及F-P掃面干涉儀的傳遞函數(shù)w的特性時,方程(7)可以進一步簡化為
式中R+=[RHR]?1RH,W為F-P掃面干涉儀的傳遞函數(shù)w構(gòu)成的矩陣,進行變換可得到
由于受到噪聲的影響,在求解Toeplitz方程的過程中會出現(xiàn)病態(tài)問題,在濾波的過程中會將誤差非線性地放大.為了解決該問題,我們對矩陣R進行了奇異值分解,但根據(jù)方程(9)中R+的表達式[20,21],可以看出在求解濾波因子G的過程中,在奇異值為零時的點處會產(chǎn)生很大誤差,因此還需結(jié)合奇異值截斷法[22]進行求解濾波因子G,其結(jié)果為
式中un和vn表示進行奇異值分解時的左右奇異值向量元.
結(jié)合方程(3)和(10)可以得到理想的瑞利-布里淵散射光譜:
利用上式對圖1中N2在壓強為2,4和6 atm情況下的卷積儀器函數(shù)后的瑞利-布里淵散射仿真光譜進行反卷積,并與對應(yīng)壓強的理想Tenti S6模型對比,結(jié)果如圖4所示.從圖4可以看出,利用維納濾波器對不同壓強下N2氣體的瑞利-布里淵散射光譜進行反卷積獲得的光譜線型基本與S6模型一致,且誤差較小,計算得到2,4和6 atm對應(yīng)的殘差的值分別為0.8%,1.3%和4.8%.離散程度隨著壓強的增大而變大,誤差主要出現(xiàn)在信號強度接近于0以及峰值和峰谷處.結(jié)果表明在一定程度上,反卷積可以消除儀器函數(shù)的影響,實現(xiàn)信號重建,可用于對瑞利-布里淵散射光譜的分析.引起誤差的原因可能是由于在采用奇異值截斷法解Toeplitz方程時,忽略奇異值而引起的反卷積光譜的抖動,而這種抖動在高壓時更加明顯.該問題可以通過Tikhonov正則化來進一步分析解決.
為了從實驗上驗證該反卷積方法的有效性和準(zhǔn)確性,并與圖3中的卷積光譜壓強反演結(jié)果進行對比,我們對圖3中測量的各壓強下的光譜進行了反卷積處理,并用理想的S6模型對其進行了擬合,結(jié)果如圖5所示.圖5表明壓強為2,4和6 atm下的反卷積后的光譜與理想的S6模型基本一致,通過計算得到各壓強對應(yīng)的殘差的值分別為5.9%,4.6%和3.8%,但在接近光譜峰谷處,存在明顯誤差.各光譜對應(yīng)N2氣體壓強反演值分別為1.91,3.93和5.73 atm,相對誤差為4.5%,1.8%和4.5%.造成光譜擬合誤差和壓強反演誤差的主要原因包括多方面.由于反卷積是對測量的光譜進行處理,測量光譜的噪聲(氣溶膠、氣池壁對入射光的反射及背景噪聲)及測量誤差會在反卷積的過程中被非線性放大并傳遞給反卷積光譜,引起更大的誤差,同時在反卷積的過程中是將理想對稱的儀器傳遞函數(shù)作為維納濾波器的輸入函數(shù),而在實際的實驗中F-P掃描干涉儀的對稱性會受到周圍環(huán)境波動(特別是溫度變化)以及F-P中的壓電陶瓷的非線性振動的影響使其不能達到完全的理想狀態(tài),造成反卷積光譜線型發(fā)生變化.
圖4 N2在2,4和6 atm情況下的瑞利-布里淵散射反卷積光譜與Tenti S6模型的比較Fig.4.Comparison between the deconvolved spectrum and Tenti S6 model corresponding to 2,4,and 6 atm.
圖5 壓強為2,4和6 atm時理想S6模型對反卷積后的光譜擬合Fig.5.Fitting deconvolved measured spectrum with ideal S6 model at the pressure of 2,4,and 6 atm respectively.
為了方便對卷積和反卷積結(jié)果進行對比,現(xiàn)將壓強為2,4及6 atm下N2卷積和反卷積對應(yīng)的值、壓強反演值及壓強反演值與實際壓力值的相對誤差列于表1.
通過表1對卷積和反卷積兩種方法得到的光譜擬合結(jié)果的對比可知,在壓強低于2 atm時,理想S6模型對反卷積作用后的光譜的擬合及壓強反演結(jié)果要優(yōu)于卷積后的S6模型對實驗光譜擬合及壓強反演的結(jié)果,而在相對高壓(大于2 atm)時,后者光譜擬合結(jié)果及壓強反演的精度要高于前者,說明將兩種方法有效結(jié)合可以通過測量瑞利-布里淵散射提高氣體壓強的反演精度.
表1 卷積和反卷積對應(yīng)的值、壓強反演值及壓強反演值與實際壓力值的相∑對誤差Table 1.The values of,retrieved pressure and relative error between retrieved pressure and real pressure.
表1 卷積和反卷積對應(yīng)的值、壓強反演值及壓強反演值與實際壓力值的相∑對誤差Table 1.The values of,retrieved pressure and relative error between retrieved pressure and real pressure.
壓強/atm 2 4 6∑ n RMS/% 6.4 2.3 3.4卷積 反演壓強/atm 1.89 3.91 5.92相對誤差/% 5.5 2.3 1.3∑ n RMS/% 5.9 4.6 3.8反卷積 反演壓強/atm 1.91 3.93 5.73相對誤差/% 4.5 1.8 4.5
本文通過卷積和反卷積實現(xiàn)了對自發(fā)瑞利-布里淵散射光譜進行理論分析和仿真,并利用兩種方法對實驗測量光譜進行了擬合及氣體壓強反演.經(jīng)過分析計算,得到采用卷積后的S6模型對實驗光譜擬合的求和歸一化均方根誤差小于6.0%,反演壓強與實際壓強的相對誤差小于6.5%,而采用理想S6模型對經(jīng)過維納濾波器反卷積處理后的光譜擬合的求歸一化的均方根誤差小于5.0%,反演壓強與實際壓強的相對誤差小于6.0%.以上結(jié)果表明兩種方法都具有較高的準(zhǔn)確性,同時反卷積在一定程度上可以消除儀器函數(shù)的影響,實現(xiàn)信號的重建.通過詳細對比,結(jié)果表明在壓強小于2 atm的條件下,利用反卷積光譜并結(jié)合理想S6模型的光譜擬合結(jié)果和壓強反演精度要高于卷積后的S6模型對測量光譜直接處理的結(jié)果,在壓強大于2 atm時,后者的光譜擬合結(jié)果及壓強反演精度要高于前者,說明將兩種方法有效結(jié)合可以通過測量瑞利-布里淵散射提高氣體壓強的反演精度.同時通過對卷積和反卷積光譜的擬合結(jié)果和壓強反演誤差的分析得到,卷積后S6模型的擬合誤差和壓強反演誤差主要是由實驗的系統(tǒng)誤差引起的,如儀器的系統(tǒng)誤差、氣池溫度波動、散射角的不確定度等.而理想S6模型對反卷積后的光譜擬合誤差和壓強反演誤差主要是由維納濾波器在信號處理過程中采用的奇異值截斷法引起的,其可以通過Tikhonov正則化來進一步地改善,而且F-P掃描干涉儀的狀態(tài)、測量光譜的噪聲等也是影響反卷積結(jié)果的重要因素.
[1]Boley C D,Desai R C,Tenti G 1972 Can.J.Phys.50 2158
[2]Ma Y,Liang K,Lin H,Ji H 2007 Acta Opt.Sin.27 962(in Chinese)[馬泳,梁琨,林宏,冀航 2007光學(xué)學(xué)報 27 962]
[3]Gu Z,Ubachs W,Marques Jr W,van de Water W 2015 Phys.Rev.Lett.114 243902
[4]Cao C L,Xu S L,Liu E W 2013 J.Univ.Sci.Tech.China 43 510(in Chinese)[曹春麗,徐勝利,劉二偉 2013中國科學(xué)技術(shù)大學(xué)學(xué)報43 510]
[5]Gu Z Y,Ubachs W,van de Water W 2014 Opt.Lett.39 3301
[6]Meijer A S,de Wijn A S,Peters M F E,Dam N J,van de Water W 2010 J.Chem.Phys.133 164315
[7]Gerakis A,Shneider M N,Stratton B C 2016 Appl.Phys.Lett.109 031112
[8]Lock J A,Seasholtz R G,John W T 1992 Appl.Opt.31 2839
[9]Pan X G,Shneider M N,Miles R B 2005 Phys.Rew.A 71 045801
[10]Witschas B,Gu Z,Ubachs W 2014 Opt.Express 22 29655
[11]Tenti G,Boley C D,Desai R C 1974 Can.J.Phys.52 285
[12]Gu Z,Witschas B,van de Water W 2013 Appl.Opt.52 4640
[13]Vieitez M O,van Duijn E J,Ubachs W 2010 Phys.Rev.A 82 043836
[14]Witschas B,Vieitez M O,van Duijn E J,Reitebuch O,van de Water W,Ubachs W 2010 Appl.Opt.49 4217
[15]Witschas B,Lemmerz C,Reitebuch O 2012 Appl.Opt.51 6207
[16]Witschas B,Lemmerz C,Reitebuch O 2014 Opt.Lett.39 1972
[17]Mielke A F,Seasholtz R G,Elam K A 2005 Exp.Fluids 39 441
[18]Wang Y Q,Yu Y,Liang K,Marques Jr W,van de Water W,Ubachs W 2017 Chem.Phys.Lett.669 137
[19]Levinson N 1946 Stud.Appl.Math.25 261
[20]Golub G H,Reinsch C 1970 Numer.Math.14 403
[21]Henry E R,Hofrichter J 1992 Meth.Enzymol.210 129
[22]Hansen P C 1990 SIAM(Soc.Ind.Appl.Math.)J.Sci.Stat.Comput.11 503