許皓文，陸 浩，王振占

（1. 中國科學(xué)院國家空間科學(xué)中心微波遙感技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京 100190；2. 中國科學(xué)院大學(xué)，北京 100049）

1 引言

微波臨邊探測是一種新型的觀測地球大氣的方式，主要通過天線從天底往上垂直掃描大氣臨邊切點(diǎn)時測量大氣毫米波、亞毫米波的熱發(fā)射，得到各種大氣變量的遙測值［1］. 在微波臨邊探測中，高光譜微波輻射計作為一種被動微波遙感器件，主要用來探測大氣的痕量氣體成分［2］，而數(shù)字譜儀是高光譜微波輻射計的核心部件，主要承擔(dān)了后端的譜分析功能［3，4］.

在目前已經(jīng)發(fā)射成功并獲取到觀測數(shù)據(jù)的臨邊探測衛(wèi)星中，最早搭載微波臨邊探測器（Microwave Limb Sounder，MLS）的是美國于1991年發(fā)射的高層大氣研究衛(wèi)星（Upper Atmosphere Research Satellite，UARS），其搭載的是濾波器組譜儀（Filter Banks Spectrometer，F(xiàn)BS）［5，6］.2001 年，瑞典航天公司聯(lián)合其他幾個國家發(fā)射了高層大氣物理研究和天文學(xué)聯(lián)合任務(wù)衛(wèi)星（Odin），其搭載的是數(shù)字自相關(guān)譜儀（Digital Auto-Correlator Spectrometer，DACS）和聲光譜儀（Acousto-Optical Spectrometer，AOS）［7，8］. 美國對地觀測系統(tǒng)（Earth Observing System，EOS）中的Aura 衛(wèi)星于2004 年發(fā)射，是UARS的繼任者，其搭載的是FBS和DACS［9，10］.日本在2009年底發(fā)射了超導(dǎo)亞毫米波臨邊輻射探測儀（JEM/SMILES），其后端是由AOS組成［11，12］.

隨著數(shù)字信號處理技術(shù)和半導(dǎo)體技術(shù)的快速發(fā)展，數(shù)字譜儀也逐漸成熟并展示出了其特有的優(yōu)勢. 數(shù)字譜儀主要包括了DACS 和快速傅里葉變換譜儀（Fast Fourier Transform Spectrometer，F(xiàn)FTS），主要特點(diǎn)是利用超大規(guī)模集成電路實(shí)現(xiàn)寬頻帶信號的并行處理，而且還具有抗干擾能力強(qiáng)、穩(wěn)定性高、體積小和更加靈活等優(yōu)點(diǎn)［13］.

考慮到衛(wèi)星上搭載時的限制因素，以及系統(tǒng)對輸入帶寬和譜分辨率的要求等，本文提出了一種新型寬帶、實(shí)時快速傅里葉變換類型的數(shù)字譜儀設(shè)計方案. 與DACS相比，F(xiàn)FTS在實(shí)現(xiàn)方式上具有更大的優(yōu)勢，可以降低系統(tǒng)復(fù)雜度，節(jié)約硬件資源. 本文所設(shè)計的FFTS分析帶寬2 GHz，通道數(shù)為1 024（1k），與各個衛(wèi)星上搭載譜儀的參數(shù)指標(biāo)對比如表1所示. 由于各類譜儀體制不同，AOS以像素點(diǎn)進(jìn)行衡量，其余類型均以通道數(shù)進(jìn)行衡量.

從表1可以看出，本文所設(shè)計的數(shù)字譜儀在分析帶寬和通道數(shù)上均有所提高. 在實(shí)際設(shè)計中一般要求分析帶寬盡可能大，通道數(shù)盡可能多，但是從數(shù)字譜儀研制的角度來看，高帶寬和多通道是很難同時兼顧的，所以在微波臨邊探測中，數(shù)字譜儀要充分平衡分析帶寬和通道數(shù)之間的關(guān)系，針對不同的探測需求對譜分辨率進(jìn)行動態(tài)的調(diào)整. 例如在大氣的中高層區(qū)域，降低譜分辨率會在一定程度上降低反演的精度，過低的譜分辨率無法辨別大氣上層的氣體吸收峰，所以在保持足夠靈敏度的前提下，2 MHz的譜分辨率有助于提高中高層大氣的反演精度［14］.

表1 譜儀參數(shù)指標(biāo)對比

2 系統(tǒng)原理及核心算法設(shè)計

2.1 高光譜微波輻射計系統(tǒng)原理

高光譜微波輻射計是一種能夠獲取大量窄帶連續(xù)頻譜通道的微波遙感系統(tǒng)，其主要原理是測量天線的輻射溫度，之后對射頻前端接收機(jī)產(chǎn)生的中頻模擬信號進(jìn)行數(shù)字化處理，通過后端數(shù)字譜儀實(shí)現(xiàn)輸入信號譜分析的功能.

高光譜微波輻射計系統(tǒng)框圖如圖1所示，系統(tǒng)中主要集成了天線、射頻前端接收機(jī)以及后端數(shù)字譜儀. 接收機(jī)又包括了射頻放大器、混頻器和中頻放大器等，其主要作用是對中心頻率包含在帶寬內(nèi)的輸入信號進(jìn)行放大及帶外抑制，輸出同樣帶寬的中頻信號［15］. 數(shù)字譜儀主要是利用集成的模擬數(shù)字轉(zhuǎn)換器（Analog to Digital Converter，ADC）和現(xiàn)場可編程門陣列（Field Program?mable Gate Array，F(xiàn)PGA）等芯片，通過模數(shù)轉(zhuǎn)換，快速傅里葉變換（Fast Fourier Transform，F(xiàn)FT），從而得到功率譜密度（Power Spectral Density，PSD）的結(jié)果. 通過對數(shù)字譜儀獲取的大氣中痕量氣體的功率譜進(jìn)行定標(biāo)分析，可以得到大氣成分的輻射亮溫譜分布，從而反演大氣參數(shù)和分布變化，最終得到大氣成分的分布規(guī)律，以便對其進(jìn)行合理的科學(xué)研究.

圖1 高光譜微波輻射計系統(tǒng)結(jié)構(gòu)

2.2 核心算法設(shè)計

基16實(shí)時復(fù)數(shù)FFT 算法是實(shí)現(xiàn)數(shù)字譜儀功能的核心算法，其是在離散傅里葉變換（Discrete Fourier Trans?form，DFT）的基礎(chǔ)上經(jīng)過并行處理［16］和復(fù)數(shù)處理［17］的改進(jìn)和結(jié)合而得到的. 并行處理可以使FFT 的計算更加快速，同時降低計算復(fù)雜度，提高運(yùn)算效率；而復(fù)數(shù)處理更加適用于硬件的實(shí)現(xiàn)過程，可以節(jié)約硬件資源，提高運(yùn)算速度.

2.2.1 并行處理

對于輸入的離散序列，可以利用DFT 的快速運(yùn)算將原本的串行計算轉(zhuǎn)化為并行計算，通過數(shù)據(jù)分解、矩陣映射、一次計算、對應(yīng)相位因子乘積和二次計算就可以完成并行FFT的處理.

假定離散序列x（n）包含N點(diǎn)數(shù)據(jù)，如果N可以分解為2 個整數(shù)R和T的乘積，如式（1）所示，且R和T均為2的整數(shù)次冪的形式，那么對x（n）做DFT可以得到長度為N的復(fù)序列X（k），如式（2）所示.

其中

如果將x（n）映射為矩形數(shù)組x（r，t），將DFT的結(jié)果X（k）映射為相應(yīng)的矩形數(shù)組X（p，q），在此過程中，選擇對x（n）采用式（4）的按列映射的方式，而對X（k）采用式（5）的按行映射的方式.

因此可以得到

其中

化簡后可得

可以將式（8）分解為以下3 步，如式（9）～（11）所示.

（1）計算T點(diǎn)的DFT，即

對矩形數(shù)組x（r，t）按行計算有r=0，1，…，R-1.

（2）定義新的矩形數(shù)組G（r，q），即

（3）計算R點(diǎn)的DFT，即

對矩形數(shù)組G（r，q）按列計算有q=0，1，…，T-1.

2.2.2 復(fù)數(shù)處理

對于輸入的實(shí)數(shù)序列，可以利用DFT 的對稱性將2N點(diǎn)實(shí)數(shù)序列的計算轉(zhuǎn)化為N點(diǎn)復(fù)數(shù)序列的計算，之后將輸出的N點(diǎn)復(fù)數(shù)序列進(jìn)行適當(dāng)?shù)呐判?，就可以得到原始?N點(diǎn)實(shí)數(shù)序列FFT的結(jié)果.

假定x（n）為一個2N點(diǎn)的實(shí)數(shù)序列，令x1（n）代表該序列中的偶數(shù)分量，x2（n）代表該序列中的奇數(shù)分量，則復(fù)序列y（n）可被定義為

由于DFT是線性的，所以y（n）的DFT可以表示為

其中，YRe（k）和YIm（k）分別代表了Y（k）的實(shí)部和虛部，序列x1（n）和x2（n）可以用y（n）表示為

其中，y*為y的共軛對稱，由x1（n）和x2（n）的DFT結(jié)果可得

將式（13）分別代入式（16）和式（17），可得

根據(jù)式（18）和式（19），可以定義4 組中間因子YRe1（k），YRe2（k），YIm1（k），YIm2（k），如式（20）所示：

對于2N點(diǎn)的實(shí)數(shù)序列x（n），對其進(jìn)行FFT 還可以表示為

其中

綜合以上公式，可以得到

根據(jù)式（23）可以得到2N點(diǎn)實(shí)數(shù)序列經(jīng)過FFT運(yùn)算后前一半的結(jié)果，后一半與前一半對稱分布.

3 數(shù)字譜儀設(shè)計及實(shí)現(xiàn)

3.1 整體設(shè)計流程

數(shù)字譜儀主要通過ADC 和FPGA 等芯片，利用核心的基16 實(shí)時復(fù)數(shù)FFT 算法完成相應(yīng)的寬帶信號處理，累加后從而得到功率譜的結(jié)果.2 GHz 帶寬數(shù)字譜儀整體數(shù)據(jù)流程如圖2 所示，4.8 GHz 高速數(shù)據(jù)在ADC內(nèi)部經(jīng)過1∶8的多路復(fù)用（DMUX）之后［18］，輸出的數(shù)據(jù)和時鐘為低電壓差分信號（Low Voltage Differential Sig?nal，LVDS），輸入至FPGA 的LVDS 數(shù)據(jù)和時鐘首先通過Buffer 轉(zhuǎn)換成單端數(shù)據(jù)和時鐘，隨后數(shù)據(jù)進(jìn)入FPGA的1∶2 DMUX 模塊，之后通過核心算法進(jìn)行處理，經(jīng)過一系列的拆分、整合和計算后可以得到每一包數(shù)據(jù)的結(jié)果，累加后通過串口上傳至上位機(jī).

圖2 數(shù)字譜儀整體數(shù)據(jù)流程

3.2 核心算法實(shí)現(xiàn)

在本文數(shù)字譜儀核心算法的實(shí)現(xiàn)過程中，根據(jù)式（9）～（11），將T和R的值分別設(shè)定為64和16. 數(shù)字譜儀核心算法架構(gòu)如圖3 所示，實(shí)際架構(gòu)包括了16 路并行數(shù)據(jù)的處理，這里僅以其中2 路作為代表進(jìn)行說明.核心算法架構(gòu)主要包括了64 點(diǎn)FFT 模塊、旋轉(zhuǎn)因子存儲模塊、數(shù)據(jù)相乘模塊、基16-FFT模塊、數(shù)據(jù)排序模塊、中間因子模塊、數(shù)據(jù)計算模塊、數(shù)據(jù)累加模塊和并轉(zhuǎn)串模塊［17］.

圖3 數(shù)字譜儀核心算法架構(gòu)

64 點(diǎn)FFT 模塊首先利用FFT IP 核同時完成16 路，每一路64 點(diǎn)的FFT，IP 核設(shè)置為流水線的方式可以保證數(shù)據(jù)的實(shí)時性. 旋轉(zhuǎn)因子存儲模塊是通過FPGA ROM 存儲復(fù)數(shù)旋轉(zhuǎn)因子，為了節(jié)約資源，本設(shè)計中將復(fù)數(shù)旋轉(zhuǎn)因子的實(shí)部和虛部合并為一個數(shù)據(jù)進(jìn)行存儲.數(shù)據(jù)相乘模塊對前面的結(jié)果進(jìn)行復(fù)數(shù)的乘法運(yùn)算. 中間因子模塊和數(shù)據(jù)計算模塊分別按照式（20）和式（23）進(jìn)行相應(yīng)的計算，累加后通過并轉(zhuǎn)串最終可以得到1 024個復(fù)數(shù)點(diǎn)的結(jié)果.

本文在算法實(shí)現(xiàn)中為了達(dá)到實(shí)時處理的目的，特別對基16-FFT 模塊和數(shù)據(jù)排序模塊進(jìn)行了改進(jìn)和優(yōu)化. 在基16-FFT 模塊，根據(jù)式（11），采取了最直接的矩陣相乘的運(yùn)算方式. 在數(shù)據(jù)排序模塊，采取了“乒乓”RAM 的數(shù)據(jù)處理方式，從而保證了整個設(shè)計過程中實(shí)時計算的持續(xù)性.

4 仿真結(jié)果分析

通過模擬正弦信號以及真實(shí)亮溫信號可以對核心算法完成前期的設(shè)計仿真以及驗(yàn)證，同時根據(jù)真實(shí)亮溫信號仿真得到的功率譜和亮溫譜還可以對數(shù)字譜儀的性能進(jìn)行評估.

利用MATLAB 產(chǎn)生2 GHz 正弦信號，采樣率設(shè)定為4.8 GHz，將產(chǎn)生的信號通過核心算法進(jìn)行運(yùn)算，經(jīng)過一定次數(shù)的累加后輸出，并和原始信號直接進(jìn)行FFT 運(yùn)算，并經(jīng)同樣累加次數(shù)的結(jié)果進(jìn)行對比. 圖4（a）給出了同樣累加次數(shù)的情況下，2 GHz 正弦信號直接進(jìn)行FFT 運(yùn)算和通過核心算法運(yùn)算后的功率譜結(jié)果和峰值通道對比，圖4（b）為2 種算法對應(yīng)通道的功率結(jié)果差值.

從圖4中可以看出，對于輸入的正弦信號，2種算法的功率譜結(jié)果高度重合，且峰值通道完全一致，每個通道對應(yīng)差值的量級為10-3，在峰值處因?yàn)樾盘柲芰枯^大，所以該處附近差值的散點(diǎn)相對聚集，而其他地方差值的散點(diǎn)相對分散.

圖4 2GHz正弦信號仿真結(jié)果

除了模擬正弦信號外，還可以利用大氣輻射傳輸模擬器（Atmospheric Radiative Transfer Simulator，ARTS）模擬臨邊探測時的真實(shí)大氣亮溫. 將由ARTS 產(chǎn)生的亮溫值轉(zhuǎn)化為進(jìn)入輻射計的熱噪聲，從而可以根據(jù)譜輻射計的仿真模型得到功率譜的結(jié)果，對功率譜進(jìn)行定標(biāo)可以得到真實(shí)的亮溫譜分布［19］，這樣就可以從實(shí)際的目標(biāo)場景角度對寬帶數(shù)字譜儀進(jìn)行仿真驗(yàn)證. 這里以240 GHz 輻射計的目標(biāo)亮溫譜為例進(jìn)行說明，仿真條件設(shè)置為輻射計噪聲溫度1 000 K 以及積分時間50 ms.

圖5 為240 GHz 輻射計中利用ARTS 模擬的天底角大氣亮溫，其中包含了所探測目標(biāo)O3和HNO3的吸收峰. 輻射計接收的是目標(biāo)的熱噪聲，根據(jù)譜輻射計的仿真模型，經(jīng)過熱輻射噪聲信號模型、射頻前端模型、數(shù)字后端模型和兩點(diǎn)定標(biāo)模型可以得到功率譜結(jié)果和亮溫譜結(jié)果.

圖5 240 GHz輻射計中模擬的天底角大氣亮溫

因?yàn)橹苯覨FT 算法和核心算法對應(yīng)的功率譜和亮溫譜都高度重合，所以這里僅列出2 種算法對應(yīng)通道的功率結(jié)果差值和對應(yīng)通道的亮溫結(jié)果差值，分別如圖6、圖7所示.

圖6 臨邊探測中真實(shí)場景目標(biāo)的功率偏差仿真結(jié)果

從圖6 和圖7 可以看出，2 種算法對應(yīng)通道的功率偏差量級都在10-5左右，2 種算法對應(yīng)通道的亮溫偏差最大值為0.02 K 左右，通過模擬不同的信號從多個維度驗(yàn)證了核心算法的合理性和正確性.

圖7 臨邊探測中真實(shí)場景目標(biāo)的亮溫偏差仿真結(jié)果

靈敏度是輻射計的一項重要指標(biāo)，可以將它定義為接收機(jī)在頻域上最小可分辨的亮溫，通過標(biāo)準(zhǔn)差的分析方法可以獲得輻射計的靈敏度. 此外，靈敏度的計算公式可以簡化為

其中，Tsys表示輻射計的系統(tǒng)溫度，β表示等效帶寬，dτ表示單次測量的積分時間. 將2 GHz 全帶內(nèi)以自相關(guān)方式計算得到的靈敏度認(rèn)為是理想全帶靈敏度，在系統(tǒng)溫度、積分時間、通道數(shù)等條件保持不變的情況下，理想子帶靈敏度和理想全帶靈敏度應(yīng)該存在C之間的關(guān)系，其中C表示通道數(shù).

通過對功率譜的結(jié)果進(jìn)行定標(biāo)分析，建立定標(biāo)方程后多次仿真得到的理想全帶靈敏度的均值為0.07 K左右，經(jīng)過核心算法后子帶靈敏度的均值為2.62 K 左右. 多次仿真結(jié)果表明，核心算法子帶靈敏度與理想子帶靈敏度相比惡化程度低于20%. 圖8為其中一次仿真后核心算法子帶靈敏度的結(jié)果.

圖8 核心算法子帶靈敏度仿真結(jié)果

5 測試結(jié)果分析

數(shù)字譜儀的譜分辨率Δf、信噪比（Signal to Noise Ratio，SNR）、信納比（Signal to Noise and Distortion Ra?tio，SINAD）和有效位數(shù)（Effective Number Of Bits，ENOB）是衡量數(shù)字譜儀工作性能的重要指標(biāo). Δf，SNR，SINAD 和ENOB 的定義分別如式（25）～（28）所示［20，21］. 式（25）中fs和K分別代表采樣頻率和采樣點(diǎn)數(shù). 式（26）中PS和PN分別代表了信號和噪聲的功率.式（27）中PD代表失真功率，式（28）中根據(jù)SINAD 的結(jié)果可以計算ENOB.

因?yàn)楦吖庾V微波輻射計本身測量的就是噪聲信號，但是在實(shí)際測試中直接產(chǎn)生寬帶的噪聲信號是較為困難的，所以同前期的仿真保持一致，實(shí)測中同樣產(chǎn)生2 GHz的正弦信號來進(jìn)行數(shù)字譜儀的性能分析.2 GHz正弦信號實(shí)測譜線的頻率-功率圖如圖9所示，圖9的子圖中是對應(yīng)的通道-功率圖.

圖9 2GHz正弦信號實(shí)測的功率譜結(jié)果

從圖9 的通道-功率圖中可以看出，信號實(shí)測功率譜和仿真功率譜的峰值通道保持一致，而且在該處附近能量較大. 因?yàn)槭艿紸DC 量化噪聲、內(nèi)核相位偏差以及諧波的影響，在信號之外會產(chǎn)生一些小的峰值，由于其功率都在-40 dB以下，所以對數(shù)字譜儀ENOB的影響也很小. 從圖9 頻率-功率圖中可以看出，對于2 GHz點(diǎn)頻信號，實(shí)測功率譜的峰值頻率為1 999 MHz，和輸入的正弦信號相比偏差在1 MHz 左右. 而實(shí)際實(shí)際情況中，氣體分子的每條吸收譜線，并不是只在一個頻率上存在吸收的直線形式，一些效應(yīng)會使吸收譜線發(fā)生譜線展寬，從而具有一定的寬度. 所以在譜分辨率恒定的情況下，1 MHz 的頻率偏差可以滿足本文中2.3 MHz 譜分辨率的設(shè)計要求.

根據(jù)得到的譜線結(jié)果，按照式（26）～（28），可以分別計算SRN，SINAD 和ENOB. 本設(shè)計方案中，結(jié)合不同頻率的測試結(jié)果，實(shí)測的ENOB 至少為3 bit，同時SNR優(yōu)于20 dB，同等條件下和仿真結(jié)果相比SNR 降低了2%左右.

6 結(jié)論

高光譜微波輻射計主要是利用后端的數(shù)字譜儀來實(shí)現(xiàn)模擬信號的譜分析，從而得到各種遙測大氣成分的功率譜結(jié)果. 本文提出了一種新型的基于高速ADC和FPGA-Kintex7 的寬帶、實(shí)時數(shù)字譜儀的設(shè)計和實(shí)現(xiàn)方案，并通過改進(jìn)、結(jié)合FFT 的并行處理和復(fù)數(shù)處理，形成了核心的基16 實(shí)時復(fù)數(shù)FFT 算法，可以完成并行信號的實(shí)時處理. 對正弦信號和ARTS 模擬的真實(shí)亮溫信號進(jìn)行仿真，結(jié)果表明，核心算法和直接FFT算法后對應(yīng)通道的亮溫偏差最大值為0.02 K 左右，而且與理想子帶靈敏度相比，核心算法的子帶靈敏度惡化程度低于20%，從多個維度驗(yàn)證了核心算法的合理性和正確性. 本文所設(shè)計的2 GHz 帶寬、1k 通道的數(shù)字譜儀譜分辨率恒定為2.3 MHz，實(shí)測的ENOB 至少為3 bit，同時SNR 優(yōu)于20 dB，可以滿足中高層大氣譜線探測的需求，同時其硬件設(shè)計與核心算法設(shè)計的方案為實(shí)現(xiàn)更高帶寬的數(shù)字譜分析系統(tǒng)提供了參考，對于臨邊探測中高層大氣的痕量氣體觀測具有重要的意義.