引用本文
王琦, 刘磊, 高太长, 胡帅, 曾庆伟. 基于ARTS的傅里叶红外高光谱计算模型研究及其影响因素分析[J]. 光谱学与光谱分析, 2019,39(6): 1711-1723.
WANG Qi, LIU Lei, GAO Tai-chang, HU Shuai, ZENG Qing-wei. A Study on the Computational Model for High Spectral Infrared Sounder by Fourier Transform Technique and its Influence Factors[J]. Spectroscopy and Spectral Analysis, 2019,39(6): 1711-1723.
Doi:10.3964/j.issn.1000-0593(2019)06-1711-06  
Permissions
Copyright©2019, 《光谱学与光谱分析》期刊社
《光谱学与光谱分析》期刊社 所有
基于ARTS的傅里叶红外高光谱计算模型研究及其影响因素分析
王琦, 刘磊*, 高太长, 胡帅, 曾庆伟
国防科技大学气象海洋学院, 江苏 南京 211101
*通讯联系人 e-mail: liuleidll@gmail.com

作者简介: 王琦, 1994年生, 国防科技大学气象海洋学院硕士研究生 e-mail: 896378352@qq.com

收稿日期: 2018-04-17 接受日期: 2018-09-19
基金: 国家自然科学基金项目(41575024)资助
摘要

在基于红外高光谱辐射数据进行大气遥感方面的研究中, 准确模拟红外高光谱数据是很重要的一步。 分析了红外高光谱辐射仪的测量原理, 建立了基于Atmospheric Radiation Transfer Simulator(ARTS)的考虑仪器干涉图截断与离散化处理过程的正向模型。 在该正向模型中, 首先采用高光谱辐射传输模式ARTS模拟得到离散化理想光谱, 通过逆傅里叶变换将理想光谱转化为干涉图, 对干涉图加窗截断处理, 模拟仪器响应函数对干涉图的影响, 最后采用傅里叶变换得到仪器测量光谱。 在这一过程中, 窗口函数的选择取决于仪器的干涉图截断方式。 未经过切趾处理的仪器, 其对应的窗口函数为矩形窗口; 经过切趾函数处理, 可以减少干涉图截断造成的能量泄露现象。 逆傅里叶变换与傅里叶变换过程中必须满足Nyquist采样定律。 基于已建立的正向模型, 模拟了tmospheric Emitted Radiance Interferometer (AERI)在Southern Great Plains (SGP)站点的108组晴空辐射数据, 并与AERI的实测结果进行比较分析, 结果发现理想光谱与AERI实测光谱在吸收线上差异较大, 最大残差达到35 mW·sr-1·m-2·(cm-1)-1(简称RU)以上, 增加干涉图截断过程后, 模拟光谱与实测光谱的最大残差减小到10 RU以内。 截断过程的增加对模拟光谱的精度有明显提高, 尤其在吸收线上, 模拟光谱明显被平滑, 模拟精度显著提高。 进一步分析六种常用窗口函数截断处理的结果与AERI实测数据的残差, 结果发现, 模拟过程中选择窗口函数为矩形窗口时, 模拟光谱与AERI实测数据残差最小, 基本可以约束在5 RU以内, 确定了AERI的干涉图截断方式可以近似看作矩形截断。 另外, 在理想光谱转换为干涉图的过程中, 理想光谱分辨率的选择决定了干涉图信息的采样率以及ARTS的计算效率, 因此综合考虑模型计算精度和模型计算效率, 确定最佳的理想光谱分辨率对于提高模型计算效能是非常必要的; 基于此, 本文模拟了不同理想光谱分辨率下的仪器测量光谱, 对比分析了模拟光谱与AERI实测光谱的残差分布, 并讨论了光谱分辨率对模型计算耗时的影响。 结果表明, 对于AERI, 在对应的正向模型中设置理想光谱分辨率为0.241 1 cm-1时, 可在保证模型准确度的前提下, 最大化模型计算效率。

关键词: 辐射传输模式; 正向模型; 切趾函数; 分辨率
中图分类号:P407.6 文献标志码:A
A Study on the Computational Model for High Spectral Infrared Sounder by Fourier Transform Technique and its Influence Factors
WANG Qi, LIU Lei*, GAO Tai-chang, HU Shuai, ZENG Qing-wei
College of Meteorology and Oceanography, National University of Defense Technology, Nanjing 211101, China
*Corresponding author
Abstract

In the research of atmospheric remote sensing based on the hyper-spectral infrared radiance, it is an important step to accurately simulate the hyper-spectral infrared radiance. In this paper, the measurement principle of hyper-spectral infrared radiometer is analyzed, and a forward model is established based on ARTS by concerning about the process of interferograms-truncated and discretization. In the forward model, an ideal discretization spectrum was simulated by ARTS firstly, and then the spectrum was transformed into an interference figure through the inverse Fourier transform. After the interference figure was truncated by a specific window function, the Fourier transform was further applied to obtain the simulated spectrum (called “target spectrum” in this paper). During this process, the window function type is dependent on the method of truncation of interference figure in the instrument measurement, (for example, a rectangular function corresponds to an interference figure without a truncation process by apodization function), and both the inverse Fourier transform and the Fourier transform must satisfy the Nyquist sampling law. Based on the forward model, 108 groups of clear sky radiation data in SGP site have been simulated, and the simulation results were compared with the actual measurement results of AERI. The results showed that there was notable difference between ideal spectrum and the spectrum measured by AERI on the gases absorption line, and the maximum residual reached around 35 RU. When a truncation process was added to the simulated spectrum, the maximum residual was constrained within 10 RU, indicating that the truncation process can improve accuracy of the simulated spectrum, especially in the gases absorption lines. Furthermore, the simulated spectrum obtained from six commonly used window function was compared with the spectrum measured by AERI. The results showed that the spectrum processed with rectangular window was most close to AERI measured spectrum, which meant that the window function used in AERI can be seen as a rectangle function. Because the ideal spectral resolution determines the sample rate of the interference figure and computation efficiency of ARTS in transformation of the ideal spectrum to interference figure, it is necessary to find an appropriate ideal spectral resolution to guarantee both the modeling accuracy and efficiency. For this purpose, instrument measurement spectrum were simulated with different ideal spectral resolution, and the residuals between simulated radiations and AERI measured radiations were analyzed. Meanwhile, the influence of spectral resolution on the computational time was discussed. The results showed that when the ideal spectral resolution is set as 0.241 1 cm-1 in the forward mode, the model calculation efficiency can be maximized on the premise of modeling accuracy.

Keyword: Radiative transfer model; Forward model; Apodizing function; Resolution
引言

红外高光谱辐射测量设备具有辐射测量精度高、 光谱分辨率高和灵敏度高等特点, 在大气温湿度廓线反演[1]、 数值天气预报同化、 气候变化研究[2]以及大气痕量气体探测[3]等方面有显著优势。 目前实现高光谱分辨率红外辐射测量的技术途径主要有光栅技术和干涉技术两种方式[4]。 由于光栅技术对探测器灵敏度要求高、 光栅定标难度大等特点, 目前基于干涉技术的红外高光谱测量设备得到更加广泛的应用, 一般将此类设备称之为傅里叶变换红外高光谱仪。 该设备在测量时, 通过对原始干涉信息的复原得到光谱信息。 由于干涉仪记录的干涉信息是有限光程差范围内的信号, 也就是说存在干涉图截断[5]现象, 因此获取的光谱与真实光谱之间存在一定偏差。

为了更好地应用红外高光谱辐射数据于大气遥感研究, 一方面需要重视如何建立合理有效的算法, 以抑制反演过程中的不适定和收敛性问题; 另一方面, 需要准确建立从大气参数到仪器辐射测量值的描述模型, 尤其需要关注傅里叶变换红外高光谱仪工作原理如何对应到辐射传输模拟过程中。 现有的高光谱辐射传输模式ARTS, 虽然可以准确模拟红外辐射传输过程, 但其得到的是大气发射理想红外辐射, 并未将仪器测量原理及过程纳入到模型中, 这必然导致模型计算值与测量值存在差异, 进而对反演结果造成误差。

基于此, 研究了地基傅里叶红外高光谱仪测量原理纳入到正演模型的方法, 建立了基于ARTS辐射传输模式的地基傅里叶红外高光谱干涉仪的正向模拟系统, 以AERI为实例, 对影响模型计算精度及计算效率的因素进行了分析, 确定了在AERI数据模拟中的最佳窗口函数与模型初始理想光谱分辨率。

1 基于ARTS的红外高光谱正向模型的建立
1.1 红外高光谱辐射仪测量原理

典型的红外高光谱辐射仪有Atmospheric Infrared Sounder(AIRS)[6], Infrared Atmospheric Sounding Interferometer(IASI)[7], Cross-track Infrared Sounder(CrIS)[8], High Spectral Infrared Atmospheric Sounder(HIRAS)[9], Geostationary Interferometric Infrared Sounder(GIIRS)[10]以及AERI[11]和Atmospheric Sounder Spectrometer by Infrared Spectral Technology(ASSIST)[12]。 这些设备的测量原理基本相同, 其核心仪器迈克尔逊干涉仪通过动镜沿光轴方向的水平直线运动改变干涉光束间的光程差, 从而实现光谱信号的振幅调制, 形成光源光谱干涉图, 再由计算机经过快速傅里叶变换输出光谱图。

傅里叶变换红外光谱仪的理想干涉图函数可表示为

I(x)=∫0+∞B(σ)[I+cos(2πσx)]dσ=12I(0)+∫0+∞B(σ)cos(2πσx)dσ(1)

其中, I(x)表示光强, σ 表示波数, B(σ )为波数σ 的光谱线强度, x表示两束光的光程差。 去除红外干涉图的直流分量, 则根据傅里叶变换公式, 可得理想光谱图的表达式为

B(σ)=∫-∞+∞I(x)exp(-i2πσ)dx(2)

由于计算机只能进行离散化数据的处理, 因此, 理论上红外高光谱干涉仪得到的理想光谱图为

B(k)=∑n=0n=NI(n)exp(-j2πkn/N)　k=0, 1, 2, …, N-1(3)

值得注意的是, 由于干涉仪转镜的扫描范围是有限的, 只能采集到最大光程差内的干涉图信息, 因此会造成干涉图截断现象, 一般通过切趾函数来控制干涉图截断的形状

W(x)'=W(x)T(x)(4)

其中, W(x)为理想干涉图, T(x)为控制干涉图截断形状的切趾函数。 对于矩形截断来说,

T(x)=1|x|≤L0|x|> L(5)

其中, L为干涉仪的最大光程差。

进一步对离散化干涉图W(n)'傅立叶变换得到测量光谱B(k)'

B(k)'=f(W(n)T(n))=B(k)⊗f(T(n))(6)

从式(6)可知, 这一过程本质上等价于理想光谱B(k)和与切趾函数T(n)相关联的扫描函数f(T(n))的卷积。 其中, 与矩形窗口相对应的扫描函数为Sinc函数, 即

f(T(n))=sin(k)/k(7)

值得注意的是, Sinc函数的旁瓣幅值大, 会造成能量泄露现象, 使临近强吸收线位置的光谱出现幅值振动。 干涉图截断过程中可以加入切趾函数控制截断形状, 改善能量泄露问题。

1.2 基于ARTS的红外高光谱模拟正向模型

高光谱辐射传输模式ARTS[13]是由德国汉堡大学和瑞典查尔姆斯理工大学联合开发的大气辐射传输模型, 可对微波至热红外辐射波段范围内的辐射传输进行模拟和计算, 具有较高的光谱分辨率和计算效率, 已广泛适用于卫星遥感技术和地基观测仪器的辐射传输模拟。 因此, 本节的正向模型建立过程不再对辐射传输模式模拟理想光谱的过程中的影响因素做考虑, 只基于1.1中干涉仪的工作原理建立正向模型, 如图1所示。

图1 正向模型框图Fig.1 Flowchart of a forward model

第1步, 通过辐射传输模式ARTS模拟理想离散光谱B(k)。

第2步, 理想光谱图B(k)转化为理想干涉图W(n)。 理想光谱B(k)到干涉图W(n)的转换通过快速逆傅里叶变换IFFT进行: 具体计算方法如式(8)

W(n)=f-1[B(k)]=1N∑k=0N-1B(k)e-j2πknNn=0, 1, …, N-1(8)

第3步是干涉图W(n)的截断过程。 通过干涉图W(n)与窗函数T(n)相乘实现。 由干涉仪自身扫描限制造成的截断是矩形截断, 对应窗口函数为矩形函数; 若红外高光谱干涉仪对干涉图进行了切趾处理, 则窗口函数为切趾函数对应的窗口函数。 具体公式如下

W(n)'=W(n)T(n)(9)

第4步, 截断干涉图W(n)'转换为目标光谱B(k)'。 截断干涉图W(n)'到目标光谱B(k)'通过快速傅里叶变换FFT进行。 具体计算方法如式(10)

B(k)=f[W(n)']=∑n=0N-1W(n)'ej2πknN　k=0, 1, …, N-1(10)

对于第一步, 在ARTS中输入大气背景廓线, 给定理想光谱波段范围、 分辨率以及传感器参数, ARTS可以从Hitran库中调用气体吸收谱线, 通过逐线积分算法得到理想光谱。 值得注意的是: 理想光谱B(k)覆盖波段应大于目标波段, 以防止卷积过程中对边缘波段的采样缺失导致的模拟结果不准确; 另外, 理想光谱B(k)的分辨率应参照红外高光谱干涉仪输出光谱的分辨率与模拟波段的特性, 综合考虑光谱准确度与计算效率来得到, 这一点我们在第3节中做了实例分析。

在IFFT与FFT变换的过程中, 应满足Nyquist采样定律。 对于理想光谱到干涉图的IFFT转换与截断干涉图到目标光谱的FFT转换过程中, 应分别满足:

V1=N1dv1=1/2dxL1=N1dx=1/2dv1(11)V2=N2dv2=1/2dxL2=N2dx=1/2dv2(12)

其中, V1是理想光谱平移到以0 cm-1为起点后对应的最大波数, N1是理想光谱的采样点数, dv1是理想光谱的波数步长, dx是干涉图步长, L1为干涉仪最大光程差, V2是目标光谱平移到以0 cm-1为起点后对应的最大波数, N2是目标光谱的采样点, dv2是目标光谱的波数步长, L2是截断处理后的最大光程差; N1和N2是整数。

因此, 在已知dv1, dv2, L1, L2的情况下, 必须满足

N2N1=dv1dv2(13)

若dv1/dv2是有理数, 则取m1/m2=dv1/dv2, 使m1和m2为最小正整数, 并寻找最小整数k, 使其满足m2× 2k× dv1> V2。 m1和m2越小, 则IFFT和FFT的计算效率越高。

此时, N1=m2× 2k, N2=m1× 2k, 这可以在满足Nyquist采样定律的前提下, 最大化IFFT与FFT的计算效率。

2 增加截断过程后的模型模拟精度分析

AERI是地基红外高光谱领域发展较为成熟的仪器, 具有方便运输、 鲁棒性强的优点, 是 ARM 计划的标配仪器, 已经在ARM计划中的固定站点和分布在德国、 葡萄牙、 中国、 印度、 尼日尔、 加利福尼亚等地的移动站点进行了部署[14]。 AERI以被动的方式自动地获取520~3 020 cm-1波段的大气下行红外辐射, 光谱分辨率为0.482 1 cm-1。

本节以AERI为实例, 基于已经建立的红外高光谱模拟正向模型, 实验分析增加干涉图截断过程的正向模型模拟精度。 挑选2011年SGP站点AERI测量的108组晴空辐射, 结合其对应时间的大气背景廓线, 利用ARTS辐射传输模式分别模拟与AERI测量光谱相对应的理想光谱, 理想光谱的分辨率设为0.024 1 cm-1。 在模型模拟过程中, 干涉图截断选择未考虑切趾过程的矩形窗口以及常用的切趾函数对应的窗口函数, 包括Triangle函数、 Gauss窗、 Hamming窗、 Cos函数和Beer函数等, 具体如表1所示。

表1 常用切趾函数 Table 1 Common apodization functions

图2(a)是108组AERI的实测辐射率的平均值和与实测光谱相对应的ARTS模拟理想光谱、 模型模拟光谱, 图2(b)是理想光谱、 模型模拟光谱与AERI测量光谱的残差, 可以看出, 在800~1 400 cm-1波段, ARTS模拟的理想光谱在吸收线上与AERI实测光谱差距较大, 理想光谱的辐射率明显大于AERI测量光谱, 最大残差达到35 RU以上。 经过窗口函数处理后, 模拟光谱与AERI实测光谱明显较为一致, 模拟光谱与实测光谱的残差约束在10 RU以内。

图2 光谱与平均残差Fig.2 Spectrum data and average residuals

图3是在800~1 400 cm-1波段的实测数据与模拟数据残差的平均值。 108组实验结果一致, 都表现为理想光谱与实测光谱的残差大于经过干涉图截断过程的模拟光谱与实测光谱的残差。 进一步说明干涉图截断处理使得模拟光谱的准确度提高。

图3 108组数据的平均残差Fig.3 Average residuals of 108 sets of data

通过图2和图3可以看出, 不同窗口函数下的模拟光谱与AERI实测光谱的残差有明显差别, 在800~1 400 cm-1波段, 不同窗口函数下的模拟结果最大差值接近10 RU, 常用的几种窗口函数下的模拟光谱中, 矩形窗口的处理结果与AERI实测结果最为一致, 在800~1 400 cm-1波段内, 残差基本约束在5 RU以内。 不同窗口函数下的辐射率残差的均值差别较小, 变化在1 RU以内, 但可以明显看出, 在108组实验中, 都显示矩形窗口的模拟结果与AERI实测结果最一致。 因此, AERI的窗口函数可以近似看作矩形窗口。

3 理想光谱分辨率对模型模拟精度及计算效率的影响分析

理想的干涉图对应的是连续的理想光谱, 但目前模式只能模拟较高分辨率的理想光谱, 计算机也只能处理离散化数据。 此外, 理想光谱分辨率的选择也与ARTS辐射传输模式的计算效率相关, 光谱分辨率越高, ARTS模拟耗时越多。

基于上述事实, 本节仍以AERI的108组晴空数据为实例, 通过模拟不同理想光谱分辨率下的模拟光谱, 分析理想光谱分辨率对模型模拟光谱精度的影响, 并进一步寻找AERI数据模拟中较为合适的理想光谱分辨率, 在保证模拟光谱准确度的前提下, 最大化计算效率。

模拟中考虑正向模型的计算效率, 目标光谱为模拟的理想光谱的整数倍, 理想光谱分辨率分别模拟为0.024 1, 0.048 2, 0.072 3, 0.096 4, 0.120 5, 0.146 6, 0.168 8, 0.192 9, 0.217 0, 0.241 1和0.482 1 cm-1。 由第2节已知AERI的扫描函数为Sinc函数, 因此将各分辨率的光谱用Sinc函数卷积到0.482 1 cm-1分辨率, 比较不同分辨率卷积结果。

图4(a), (b)和(c)是波数域上108组不同分辨率下的模拟值与AERI实测值的残差, 可以看出, 在800~1 400 cm-1波段, 0.482 1 cm-1分辨率的卷积结果与AERI实测值的残差明显大于其他分辨率下的结果, 最大残差接近50 RU; 其他分辨率下的卷积结果与AERI实测值的残差较为集中, 基本可以约束在10 RU以内。 图4(d)表示的是波数域上的平均残差, 可以明显看出, 随着理想光谱分辨率的降低, 模拟结果与AERI实测结果的残差越来越大, 并在理想光谱分辨率与AERI仪器分辨率达到一致时, 产生跳变。

图4 不同分辨率下的平均残差Fig.4 Average residuals for different ideal spectral resolution spectra

图5是不同分辨率下的平均残差与ARTS模拟耗时的比对, 从平均残差的变化曲线来看, 当理想光谱分辨率高于0.241 1 cm-1时, 平均残差随分辨率的变化幅度很小, 几乎可以忽略不计, 理想光谱分辨率为0.482 1 cm-1时, 平均残差产生跳变, 这与图7的分析结果保持一致。 从ARTS模拟理想光谱的耗时曲线可以看出, 分辨率越高, ARTS模拟理想光谱耗时就越多, 在理想光谱分辨率设置为0.024 1 cm-1时, 模拟一组理想光谱的时间达到4 000 s以上, 这对需要大量模拟数据的反演应用是不可取的; 理想光谱分辨率为0.048 2和0.072 3 cm-1时, ARTS的耗时有明显降低, 耗时大约为1 400 s; 理想光谱分辨率低于0.096 4 cm-1时, ARTS的时耗随分辨率的变化减小, 尤其在分辨率低于0.217 0 cm-1时, 时耗随分辨率的变化可以忽略不计。

图5 不同分辨率下的平均残差及ARTS模拟时耗Fig.5 Average residuals and time consuming of ARTS for different ideal spectral resolution spectra

综合平均残差曲线与耗时曲线来看, 在800~1 400 cm-1波段, 目标光谱分辨率为0.482 1 cm-1时, ARTS模拟理想光谱的分辨率设置为0.241 1 cm-1, 可以达到在保证光谱准确度的前提下, 最大化计算效率。

4 结论

傅里叶变换红外光谱仪在测量红外高光谱数据时, 存在干涉图截断以及计算机的离散化处理过程。 基于ARTS辐射传输模式, 建立了考虑仪器干涉图截断与离散化处理过程的正向模型, 并以AERI数据为实例, 对增加截断过程后的模型模拟精度进行分析。 结果发现, 干涉图截断过程对模型模拟精度影响很大, 尤其在吸收线上, 经过干涉图截断处理的模拟光谱与实测光谱的残差明显减小; 干涉图截断处理使得模拟光谱的准确度明显提高。

另外, 不同窗口函数处理下的模拟光谱也有明显不同。 在AERI数据模拟中, 矩形窗口截断的模拟结果更接近AERI实测数据, 因此, AERI的窗口函数可以近似看作矩形窗口, 即AERI在干涉图的处理中, 没有加入切趾函数控制干涉图截断形状。

最后, 对正向模型中理想光谱分辨率的选择做了实验分析。 在AERI数据模拟中, 0.241 1 cm-1的分辨率可以在保证模型准确度的前提下, 最大化计算效率。 实际上, 在解决大型而复杂的计算问题时, 并行计算是一种很好的计算方式, 可以有效提高计算速度, 这是后续要深入试验分析的工作。

The authors have declared that no competing interests exist.

参考文献
[1] LIU Yang, GUAN Li(刘旸, 官莉). Meteorological Monthly(气象), 2011, 37(3): 318. [本文引用:1]
[2] QI Cheng-li, GU Ming-jian, HU Xiu-qing, et al(漆成莉, 顾明剑, 胡秀清, 等). Advances in Meteorological Science and Technology(气象科技进展), 2016, 6(1): 88. [本文引用:1]
[3] Griffith D W T, Deutscher N M, Caldow C, et al. Atmospheric Measurement Techniques Discussions, 2012, 5(3): 3717. [本文引用:1]
[4] WANG Jian-yu, LI Chun-lai, JI Hong-zhen(王建宇, 李春来, 姬弘桢, 等). Journal of Infrared and Milimeter Waves(红外与毫米波学报), 2015, 34(1): 51. [本文引用:1]
[5] Lee L, Chen H, Chen S, et al. Applied Optics, 2012, 51(20): 4622. [本文引用:1]
[6] Ruzmaikin A, Aumann H H, Manning E M. Journal of the Atmospheric Sciences, 2014, 71(7): 2516. [本文引用:1]
[7] Hilton F, Armante R, August T, et al. Bulletin of the American Meteorological Society, 2012, 93(3): 347. [本文引用:1]
[8] Marshall J L, Jung J, Lord S, et al. Bulletin of the American Meteorological Society, 2015, 88(3): 329. [本文引用:1]
[9] Qi C, Gu M, Wu C, et al. Calibration Method of High Spectral Infrared Atmospheric Sounder Onboard FY-3D Satellite, 3rd International Symposium of Space Optical Instruments and Applications, 2017. [本文引用:1]
[10] Dong C, Yang J, Zhang W, et al. Bulletin of the American Meteorological Society, 2010, 90(10): 1531. [本文引用:1]
[11] Gero P J, Turner D D. Journal of Climate, 2011, 24(18): 4831. [本文引用:1]
[12] Wang G, Zhang J. Advances in Space Research, 2014, 54(1): 49. [本文引用:1]
[13] Li Shulei, Liu L, Gao T. Journal of Atmospheric & Environmental Optics, 2016. [本文引用:1]
[14] Long C N, McFarlane S A, Genio A D, et al. Bulletin of the American Meteorological Society, 2013, 94(5): 695. [本文引用:1]