基于氧氣A波段發射譜線臨近空間大氣溫度的反演及分析

楊曉君，王后茂, 李葉飛，王詠梅，胡秀清

1. 中國科學院國家空間科學中心, 空間環境探測研究室，北京 100190 2. 天基空間環境探測北京市重點實驗室，北京 100190 3. 中國科學院空間環境態勢感知技術重點實驗室，北京 100190 4. 中國科學院大學天文與空間科學學院，北京 100049 5. 上海衛星工程研究所，上海 200240 6. 中國氣象局國家衛星氣象中心，北京 100081

引 言

地球大氣臨近空間區域的物理和化學特性研究需要準確的了解該區域的溫度分布，在臨近空間區域發生的許多過程都具有明顯的溫度依賴性。此外，臨近空間大氣參數很大程度上也影響著空間天氣預報以及航天器發射與再入軌過程的安全性預估等[1]，因此對臨近空間大氣的觀測與研究具有重要的科學意義和應用價值。雖然火箭和地面測量已經可以獲得準確的溫度廓線，但很難用它們建立一個全球的溫度分布數據庫，星載探測可以在較短時間內進行全球或者區域的大范圍觀測。美國空間物理研究實驗室基于搭載在URAS衛星上的高分辨率多普勒成像儀HRDI(The high-resolution Doppler imager)通過觀測氧氣A波段的發射譜線進行溫度的測量[2], 采用基于擾動理論的最優估計法對80～100 km高度的大氣進行反演，利用三條譜線協同反演得到誤差約為7 K的溫度廓線，首次提供了該高度范圍全球溫度測量結果。之后, 2001年搭載在TIMED衛星上的TIDI(TIMED Doppler Interferometer)多普勒成像儀基于氧氣A波段發射譜線對60～300 km高度的臨近空間大氣溫度進行了遙感觀測，基于譜線的觀測值采用約束非線性最小二乘擬合得到溫度廓線[3]。2001年2月20日，搭載在Odin衛星上的光譜儀和紅外成像系統OSIRIS(optical spectrograph and infraRed imaging system)對氧氣A帶的發射光譜進行了探測，將90～110 km內的氧氣A波段光譜與模擬光譜進行逐像素比較，用最小二乘法使波段譜差的平方和最小，從而得到中層-低熱層區域的溫度，在90 km附近的溫度反演精度約為±2 K，在較高高度的精度為±6 K[4]。先進地球觀測衛星ADEOS上搭載的大氣光譜儀ILAS測量了包括以762 nm為中心的分子氧A波段在內的753～784 nm波段的大氣吸收光譜，通過衛星太陽掩星對氧氣A波段吸收光譜的測量，采用層析法獲得了整個平流層的溫度分布[5]，與儀器功能相關的誤差是系統不確定性中最大的誤差源, 與所考慮的系統不確定性相關的溫度的均方根誤差估計為4 K。全球高分辨率熱層成像邁克爾遜干涉儀MIGHTI(Michelson Interferometer for Global High-resolution Thermospheric Imaging)搭載在NASA的電離層連接探測器ICON上于2019年10月11日發射[6]，MIGHTI被用來測量地球熱層的中性溫度[7]，反演方法是使用離散波長測量A波段的光譜形狀，從波段的旋轉包絡線推斷出環境溫度，在白天90～105 km之間的測量精度為1K，而在夜間，測量精度在90 km時為1 K，在105 km時增加到了3 K[8]。

以上儀器選取的譜線都不同，而影響反演結果非常重要的一點就是譜線的選擇。我們在氧氣A帶所有譜線中對比不同譜線特性，尋找譜線選擇的判斷條件，然后選取合適的譜線基于最優估計的方法進行溫度的反演，并進一步探討了譜線線強對反演結果的影響。

1 正演模型

產生氧氣A波段的激發態氧氣分子是由一系列復雜的光化學反應產生的，當激發態氧氣分子自發輻射產生氧氣A波段時會發出光子，波數ν處的特定轉動譜線發出的光子數可以表示為fDη，D是多普勒線型，η是體發射率，f是譜線強度，是溫度的函數

(1)

式(1)中，fref是在參考溫度Tref=296°K時的譜線強度，E′是轉動譜線的上能級能量。式(2)是以vc為譜線中心的多普勒線型

(2)

(3)

其中，k是玻爾茲曼常數，c是光速，m是氧氣分子量。

綜上所述，衛星觀測到的光譜亮度可以從切點高度Zt沿視線進行路徑積分得到

(4)

其中，s是視線路徑，n為氧氣數密度，σ為吸收截面。由式(4)可得，光譜輻射主要由氧氣體發射率和溫度決定，而光譜輻射對溫度的依賴主要是通過自吸收截面等的作用[9]。

為了更真實的模擬星載測量，引入儀器函數，因此模擬的儀器計數值用式(5)獲得

(5)

式(5)中，I為卷積的歸一化儀器函數，是半高寬(FWHM)為20 cm-1的高斯函數。

2 反演模型

2.1 擾動理論

衛星的觀測結果可表示為

C=KT+ε

(6)

式(6)中，C是衛星觀測的輻射強度向量；T是待反演的溫度向量；K是描述大氣物理過程正演模型的權函數矩陣；ε是測量誤差向量。

首先，分別模擬臨邊輻射強度測量值和用體發射率廓線、溫度廓線初值計算得到的正演模型值。模擬過程中的譜線線強、上能級能量等參數均來自2004年的HITRAN光譜數據庫[10]，而大氣分子濃度和溫度廓線數據均來自MISIS-00模型[11]，體發射率廓線是由光化學模型利用溫度各組分數密度廓線計算得到[12]。測量值Cmeans和用初始溫度、體發射率廓線計算得到的正演模型值C00之間的差值可以用線性積分算子表示

本文筆者將微課應用到成人繼續教育中，將信息化教學手段融入了傳統教學中，從而形成了線上線下的混合式教學模式。這種混合式教學模式是一種優秀的教學模式，它綜合了MOOC的優勢，彌補了MOOC缺乏管理機制的缺陷，利用現有大量的MOOC資源，降低了微課制作的工作量，突現了“互聯網+”時代的優勢。

(7)

其中z為切點高度(km)； 矩陣中C和T下標為視線切點高度(km)； 矩陣K下標為所在行視線積分的高度和依次改變溫度值所對應的高度； ΔT=Ttrue-T00為真實溫度值與初始溫度值的差值；Cmeans是輸入溫度為T00+ΔT的模擬值；C00是輸入溫度為T00的模擬值。

以1 km為格點對溫度廓線加入擾動，每改變一個格點就產生一條新的擾動廓線ΔT/T00，51條廓線就可以用ΔT/T00組成的矢量Tf來表示。為了得到更高的精確度，式(7)中的積分與正演模型中的一樣，積分步長都為1 km。將模擬臨邊掃描的所有測量計數值作為C矢量，這樣就得到了一個矩陣方程ΔC=KTf用來求解未知量Tf。

2.2 反演算法

理論上，溫度可以直接利用氧氣A波段的觀測輻射光譜通過求解輻射傳輸方程來獲得，然而事實上求解這一方程是不可能的，這主要是因為其解是不唯一且不穩定的； 同時由于觀測誤差的存在，使得簡單的單步反演法在實際應用時容易導致較大的反演誤差，抗干擾能力比較低?；谪惾~斯理論的最優估計法是一種求解有噪聲反演問題的非常好的手段，用實時探測數據替代非待反演參數，并利用先驗信息對反演結果進行約束和修正，得到待反演參數最合理的解，提高了反演精度和抗干擾能力。

在2.1擾動理論的基礎上，選擇了基于貝葉斯理論的最優估計法進行溫度的反演。假定測量誤差和先驗信息誤差服從高斯分布，將求解輻射傳輸方程轉化為使概率函數最大化的最優化數學問題，求解這個概率問題就可以獲得溫度的最優解。

將T看做一個多維隨機變量，其概率密度函數為

(8)

條件概率密度函數可表示為

(9)

其中，Sa是描述T0不確定性的先驗信息誤差協方差矩陣； 對角矩陣Sε是從觀測值得到的測量噪聲協方差矩陣，是描述C和正演模型不確定性的觀測誤差協方差矩陣。假設測量誤差的總體方差和探測器上的計數值一樣服從泊松分布。先驗估計協方差矩陣Sa可以模擬為

Sa=αTJ

(10)

式(10)中，z是反演的網格高度，αT=103K2。先驗估計協方差矩陣表示了一個預期的均方根變化約為5%(～10 K)的先驗溫度廓線。

根據貝葉斯定理，最有可能的大氣狀態剖面是使式(11)最大化[13]

P(Tf|C)=P(Tf)P(C|Tf)/P(C)

(11)

其中，P(Tf)是溫度參數分數變化的先驗概率分布；P(C)觀測數據的先驗概率分布，其為常數；P(C|Tf)給定溫度參數分數變化條件下具體觀測數據的似然概率；P(Tf|C)是組合先驗信息和似然概率得到的溫度參數分數變化后驗概率。

將式(8)和式(9)相乘，最有可能的T值是質式(12)最小化

(12)

式(12)最小化時的T的估計值為

(13)

其中，T0=0, 所以方程可化為

(14)

再由Test求得溫度廓線的最優估計值Ttrue。

3 反演結果及討論

3.1 正演模型

為了驗證上述反演方法的可行性和精確度，利用構建的正演模型來模擬實測值進行溫度的反演研究。首先利用一組溫度輸入值T00作為初始值來獲得一組正演模擬值，然后在T00的基礎上加上一組隨機溫度變量得到新的溫度廓線Ttrue=T00+ΔT，再由新的溫度廓線模擬一組實測輻射值，最后由該兩組輻射模擬測量值組合進行溫度反演，將反演得到的結果與初始溫度輸入值進行比較。圖1就是正演模型所用的溫度廓線初始值(T00)和模擬的實際測量值所用的溫度廓線(即溫度真實值：Ttrue)。

圖1 溫度廓線初始值和實際測量值對比

3.2 無噪聲時的反演結果分析

基于HITRAN數據庫中氧氣A帶(759～767 nm)內包含的幾十條譜線參數，我們將所有的譜線進行了輻射模擬與溫度反演分析，結果表明不同譜線的反演結果精度不同，精度較好的譜線主要集中在759～760和764～766 nm這兩個波段范圍，圖2為氧氣A帶中無噪聲情況下60～110 km范圍內的溫度反演精度較高的譜線。如圖2所示，ΔT=Ttrue-T00為真實溫度值與初始溫度值的差值，其他折線為利用不同譜線反演ΔT的結果，17條氧氣A波段的溫度反演結果在60～110 km高度范圍內效果較好，平均反演偏差為4.1 K。

圖2 無噪聲條件下譜線在60～110 km范圍內的 溫度反演結果對比

但有些譜線的反演誤差在80 km以上較小，80 km以下卻較大，如圖3所示，764.28和764.17 nm為圖2中反演精度較高的譜線，另外兩條為反演精度較差的譜線(761.72和761.25 nm)。將他們的反演結果進行進一步的對比分析，如圖3所示，80 km以上四條譜線顯示了較好的反演精度，平均反演誤差<5 K； 而在80 km以下，761.72和761.25 nm的反演結果與真實值存在較大的偏差，平均偏差達到了34.9 K。原因如下： 由式(4)可以看出，溫度通過影響線強和自吸收兩部分來影響輻射強度，且溫度對它們的影響正好相反，權函數就是用來表示溫度對輻射強度影響大小的，而反演結果的差距可從其權函數中得到規律。圖4為兩條譜線761.25和764.28 nm的權函數，為方便對比各選取一條代表線標記為紅色(71和80 km)，圖中每一點代表的都是對應高度處的溫度變動對衛星輻射強度觀測的影響與溫度的擾動量的比值，即ΔC/ΔT，溫度擾動ΔT設置為5 K，由圖4(a)可以看出，當溫度對自吸收的影響所占比重大于對線強的影響時，ΔC/ΔT會發生正負翻轉，反演精度會變差； 從圖4(b)可以看出，影響因素占主導地位的是線強，溫度主要通過對線強的影響來改變輻射強度，其對自吸收的影響比重較小，ΔC/ΔT未發生正負翻轉，反演效果比較好。同時，我們將其他譜線的反演結果進行比較得出相同的結論，因此可以根據這個規律，從氧氣A波段眾多譜線中篩選適合用于溫度反演的譜線，即溫度對自吸收影響相對較小的譜線。

圖3 無噪聲條件下四條譜線在60～110 km范圍內的 溫度反演結果對比

圖4 權函數矩陣Fig.4 Weight function matrix

3.3 添加噪聲后的反演結果分析

在實際的觀測中，除了上述考慮的因素外，還需要考慮信噪比的影響。因此，本文將圖3中的四條譜線模擬值加入隨機噪聲并進行了反演及結果分析，結果如圖5所示。764.28和764.17 nm兩條譜線的反演結果與真實值基本吻合，說明噪聲對這兩條譜線的影響非常小，該譜線的輻射強度能夠滿足溫度反演的要求。而761.72和761.25 nm這兩條譜線在80 km以下的反演偏差變得更大，最大偏差達到幾百K，原因是這兩條譜線本身在80 km以下反演誤差就非常大，微小擾動就會造成較大的誤差，所以加入噪聲后使得反演偏差變得更大。

圖5 加入噪聲情況下四條譜線在60～110 km范圍內的溫度反演誤差

為了進一步分析噪聲對反演結果的影響，本文選擇了四條較強譜線和一條弱譜線，如表1所示，其中762.2 nm為弱線，反演結果如圖6所示。圖6是762.2 nm無噪聲情況和加入噪聲時的反演結果比較，藍線為無噪聲情況下的結果，紅線為加入噪聲的結果。無噪聲情況下，弱線762.2 nm與強線764.28和764.17 nm的反演結果相對一致且都比較好； 加入噪聲后，由于本身輻射強度弱，所以信噪比低，反演效果變得比較差。原因是弱線容易受噪聲的干擾，抗干擾能力低。因此，譜線選擇時，線強也是另一個重要的判斷依據，強線更有利于提高反演精度。

表1 模型中用到的譜線參數Table 1 Constants of rotational lines of model

圖6 譜線762.2 nm在無噪聲情況和加入噪聲時60～110 km范圍內的溫度反演誤差

為了進一步確定線強對反演結果的影響，我們選擇通過增加線強的方式進行反演精度分析，圖7為分別將線強增加1～5個量級后得到的反演結果。由圖可得，80 km以上，隨著線強量級的增大，反演精度增加，反演結果更接近初始輸入值； 線強增加三個量級時，反演精度可以達到5 K。由此得出： 弱線762.2 nm在信噪比足夠大，也就是線強達到10-26的量級時，也可用來反演80 km以上的溫度并獲得精確的反演結果。

圖7 增加譜線762.2 nm線強量級對反演結果的影響

4 結 論

在臨近空間大氣溫度反演中，不同的譜線反演精度不同，基于貝葉斯反演理論選取了不同的譜線進行一系列反演試驗及結果分析，并基于權函數和譜線線強與反演精度差異的關系進行了分析。

具體結果如下：

(1)當溫度對自吸收的影響所占比重大于對線強的影響時，權函數會發生正負翻轉，反演精度變差，原因是自吸收降低了輻射強度對溫度的靈敏度，因此，權函數變化規律可以作為譜線選擇的判斷依據。

(2)在有噪聲的情況下，強線比弱線的抗干擾能力更強，反演精度較高，更適合用于溫度的反演，所以線強也是譜線選擇的另一個重要的依據。根據分析結果，當譜線線強達到10-26的量級時，譜線的輻射強度較強，此時譜線可用于80 km以上的溫度反演。