成像差分吸收光谱技术的软件研发与数据反演

曹子昊,曾议,鲁晓峰,廖捷,杨东上,常振,司福祺,奚亮

(1合肥学院生物食品与环境学院,安徽 合肥 230601;2中国科学院合肥物质科学研究院安徽光学精密机械研究所,中国科学院环境光学与技术重点实验室,安徽 合肥 230031;3中国科学技术大学,安徽 合肥 230026)

0 引 言

近年来,随着我国工业化程度的提升及城市化进程的加快,大气污染物的排放量也在不断地增多,这种人为排放源加剧了环境污染,极大地影响了人们的正常生活。因此需要重视大气环境现状并分析影响大气质量的各种因素,实现对大气污染物的监测,以便有效地控制及治理。

污染气体的浓度分布及污染源排放情况是近些年来大气环境监测研究的重点。差分吸收光谱(DOAS)技术通过遥感的方式得到目标区域的光谱信息,利用大气中的气体分子对光在不同波段的特征吸收以实现气体定性及定量测量[1],因其具有易于操作、价格低廉、灵敏度高等优点,越来越广泛地应用到大气环境监测中。成像差分吸收光谱(IDOAS)技术是DOAS技术与成像光谱技术结合产生的,利用成像技术可以直接获取目标区域的空间分布信息[2-4],进而分析得到目标区域的分布特征及轮廓形状。此外,将目标区域的辐照按照波长分解可以得到目标区域任一点的光谱数据,结合DOAS算法进行反演即可获得痕量气体的浓度信息,进一步得到痕量气体浓度的二维分布图。

在国外,德国海德堡大学使用IDOAS技术对墨西哥火山烟羽中的SO2、BrO和OClO进行了二维成像观测[5];美国加利福尼亚大学使用IDOAS技术对休斯顿加尔维斯顿区污染排放产生的HCHO和SO2进行了观测研究[6];德国不来梅大学环境物理研究所开展机载IDOAS实验并获得区域NO2的分布情况[7]。在国内,IDOAS技术的研究也取得了一定进展,刘进等[8]利用IDOAS技术在机载和地基平台都对污染气体的二维分布进行了相关的研究;Cheng等[9]用IDOAS技术监测了船舶污染源的排放情况,得到了NO2二维分布情况。

本文首先介绍了地基IDOAS工作原理及系统构成,着重研究了数据反演中所使用的软件工具(QDoas),解析其代码构成及框架功能,为实现更高效更灵活的气体探测,完成了反演软件的改进工作,以满足个性化的数据处理要求。该软件成功在Windows系统上编译运行,为数据反演提供了更快速和精准的工具。研究地基IDOAS数据反演方法,使用新研发的软件在铜陵富鑫钢铁厂开展实验,成功反演得到烟囱出口附近的SO2及NO2的二维浓度分布图。

1 测量原理

1.1 IDOAS工作原理

图1为IDOAS的原理示意图。该技术采用面阵CCD记录目标区域每个像元的高分辨率光谱信息,包括光谱维和空间维信息。工作方式为帧转移,其优点是感光面的利用率和分辨率高,广泛地应用于成像光谱仪。将目标区域按列划分,以每列为单位对目标区域扫描,通过角度的转动实现目标区域的整体扫描[10]。对收集到的光谱信息进行反演得到多组气体浓度值,将得到的气体浓度值与面阵CCD每个像元利用空间分布信息进行精准匹配,实现目标气体的二维分布成像[11,12]。该技术可同时测量多种痕量气体。

图1 IDOAS原理示意图Fig.1 Schematic diagram of IDOAS principle

1.2 仪器与扫描测量

基于以上原理,成功搭建地基IDOAS扫描系统,该系统主要包括二维转台、紫外镜头、高光谱成像光谱仪、控制与采集系统等[3],如图2所示。太阳散射光经紫外镜头汇集进入成像光谱仪的入射狭缝,经过反射后成像到面阵CCD上。面阵CCD的空间维像素为512 pixel,光谱维像素为2048 pixel,即每次测量可以得到A1-A2纵列方向上的512条具有2048 pixel的光谱。通过转动二维转台电机,以推扫的方式从图2中A1-A2至B1-B2的方向进行扫描,得到目标区域的光谱信息及空间信息,可以实现方位角从0°～360°,俯仰角从0°～90°的扫描[13]。测量过程中选择目标区域中视场范围内没有污染物且能见度高的区域的测量光谱作为参考光谱,用以后续的数据反演。图3为地基IDOAS扫描系统实物图。

图2 地基IDOAS扫描系统工作原理图Fig.2 Working principle diagram of ground-based IDOAS scanning system

图3 地基IDOAS扫描系统实物图Fig.3 Physical map of the ground-based IDOAS scanning system

2 数据反演及软件研发

2.1 数据预处理

使用面阵CCD在每次成像时可以同时获得空间维和光谱维的信息,针对空间维上的每个点使用QDoas软件反演,得到该点的污染气体浓度信息。推扫成像后,可获取扫描区域内的污染气体浓度信息。

因此,为适应QDoas反演,提取空间维每个点在时间维上的所有光谱组成一个文件,然后利用该文件进行反演,得到该空间维点在时间维上的变化信息。IDOAS每次成像获取的文件为像素值数据矩阵,如图4所示。

图4 IDOAS数据立方体示意图Fig.4 Schematic diagram of IDOAS data cube

对数据立方体进行坐标转换,图4所示空间维每个点代表一个文件,此文件包含数据的光谱信息,不同行所对应的光谱信息为不同时间段的数据。进行空间维的像元合并,能够提升信噪比和反演精度。

2.2 数据反演及二维成像

图5为DOAS拟合的流程图。在外场实验中,仪器配有光谱文件,该文件标识了每行光谱信息及对应的波长值。进行DOAS拟合时,首先添加太阳光谱并设置所需各种气体的吸收截面及狭缝函数太阳参考谱,然后提取原始文件的遥感影像像元亮度(DN)值导入QDoas软件中。使用选择好的光谱仪光谱校正文件对数据对应的光谱波长进行校正,经过多次迭代,最后生成的校正(shift)图的值保持在±0.1之间,得到校正文件及Ring截面文件(用于消除大气Ring效应的影响)。

图5 DOAS拟合流程图Fig.5 DOAS fitting flowchart

实验过程中需要使用最小二乘拟合,以确保参考截面与测量光谱的分辨率相同。通过Qdoas中的卷积工具实现对参考截面的卷积,通过软件QDoas中所提供的Ring工具对上述保存的Ring截面文件与生成的迭代文件及太阳光谱截面文件进行卷积,得到所需的Ring截面文件,该过程视为一种吸收成分参与光谱拟合。对于污染气体NO2,根据其光谱结构选择337～370 nm的波段进行反演,参与反演的气体包括NO2、HCHO、O4、O3。对于污染气体SO2,根据其光谱结构选择308～330 nm的波段进行反演,参与反演的气体包括SO2、O3、NO2。使用QDoas反演软件对测量光谱逐条反演,根据得到的多组污染气体斜柱浓度信息,经过二维拼图得到污染气体的分布信息图。

2.3 软件优化研究

DOAS监测技术的应用已经较为成熟,但在使用QDoas软件进行浓度反演时,由于其软件固有的操作模式和功能,存在数据处理繁琐、反演耗时较长等问题。为更加高效精确的进行污染气体浓度反演,对软件的关键代码部分做进一步的研究和开发。基于C++语言使用QT软件对QDoas的源代码进行了解析、重整,提取和优化原有功能,编写定制了适配Windows系统的新反演软件。

图6为在QT软件上的源码图,其中all.pro为all中各模块的公共配置文件,结合软件功能对代码进行分析。Convolution、qdoas、ring和usamp为软件的界面设计模块,qdoas为主体,通过按钮调用其它功能的界面模块。Cmdline为命令行操作模块,通过命令实现程序的功能,各参数的设置则通过读取xml文件实现,所以也包括了xml文件的解析工作。Commom具有消息弹出框的编写及绘制图像的功能。Engine为软件的关键部分,主要包括算法方面的编写以及各种格式文件的读取。Mediator可以解释为一个中介系统,通过它链接各模块之间的功能。

图6 软件源码图Fig.6 Diagram of software source code

程序运行需要用到以下库文件:qwt库(生成图表)、gsl库(科学计算)、coda库(读取GOME-2数据)、HDF-EOS2库及HDF5库。qwt库在Windows系统中能够直接找到并使用,而其它的几个库在Windows系统中没有直接提供,需要通过cmake+mingw编译所提供文件的方式才能得到。若要实现程序的基本功能,所需要研究的模块为cmdline、engine、mediator这三个模块,若为后续的图形界面开发考虑,也可以将common部分保留,外部库仅调用gsl库文件就可以满足基本需求,研究编写相关代码,实现需求功能。在编译过程中生成libcommon.a、libengine.a和libmediator.a三个静态库供cmdline模块调用,编译完成后生成可执行程序,能够实现污染气体的反演。

原先处理多个文件的项目时,需要多次运行程序,为方便处理,编写脚本重复调用运行程序,但需要等待一个程序结束后才能运行下一个程序,比较耗时。而将所需源代码成功编译后,优化后的代码使用多线程的方式同时处理多个文件,大大减少反演时间,并且在出现错误或实现某些定制的需求功能时,能够很好地发现并改正错误以及添加想要的定制功能。

3 结果与讨论

2019年11月6日,采用地基IDOAS系统对铜陵富鑫钢铁厂开展了现场观测实验。图7为富鑫钢铁的卫星地图,标识IDOAS处即为测量系统的摆放位置。

图7 富鑫钢铁卫星图Fig.7 Satellite picture of Fuxin Steel

图8为获取的目标区域的灰度图像,于上午10:25:32、10:34:53及10:44:17对目标区域进行了三次测量,工作过程中将二维转台设置为在0°～90°范围内每2秒旋转1°,使用研发的软件对所有测得的光谱进行反演,得到不同时段污染气体NO2(图9)和SO2(图10)的二维浓度分布图。

图8 观测区域示意图Fig.8 Schematic diagram of observation area

图9 2019年11月6日不同时间点NO2的二维浓度分布图。(a)10:25:32;(b)10:34:53;(c)10:44:17Fig.9 Two dimensional concentration profile of NO2at different time points on November 6,2019.(a)10:25:32;(b)10:34:53;(c)10:44:17

图10 2019年11月6日不同时间点SO2的二维浓度分布图。(a)10:25:32;(b)10:34:53;(c)10:44:17Fig.10 Two dimensional concentration profile of SO2at different time points on November 6,2019.(a)10:25:32;(b)10:34:53;(c)10:44:17

与目标区域的观测示意图进行对比,比较直观且清晰地看出污染气体NO2和SO2的斜柱浓度(SCD)分布信息。NO2主要集中于工厂建筑附近区域,浓度呈向上逐步递减的趋势,且烟囱左边区域NO2浓度明显高于右边区域NO2浓度,左边的NO2斜柱浓度均保持在1×1017molecules·cm-2左右,右边区域NO2斜柱浓度均保持在7×1016molecules·cm-2左右。而SO2浓度最高值主要集中在烟囱附近,其斜柱浓度保持在9×1016～1×1017molecules·cm-2,存在随风扩散的情况,与目标区域示意图中所示结果比较一致。

对铜陵富鑫钢铁厂烟囱排放口的观测实验成功观测到SO2与NO2污染气体的实时分布信息,进一步也可对其他形式排放(如汽油、柴油燃烧)且存在于大气边界层的污染气体进行扫描观测。成像系统能够很快得到目标区域的污染气体浓度分布图,具有很强的直观性和时效性,对分析城市边界层污染气体分布和传输提供了强有力的支撑,有利于相关研究人员对地区的污染情况做更近一步的研究及分析,为大气环境监测提供了可靠的支持,为应急管理提供了科学的数据支撑。

4 结 论

介绍了地基IDOAS扫描系统的基本结构和工作原理。在使用QDoas光谱反演软件的过程中,发现数据处理繁琐、反演耗时较长及出现错误难以分析原因等问题,为更高效准确地使用DOAS系统测量各种大气痕量气体,定制研发了新的反演软件,基于C++语言使用QT软件对代码做了进一步的改进和优化,并成功移植到Windows平台,满足数据反演的需求,为更加深入的算法反演研究提供了软件支持。对铜陵富鑫钢铁厂开展了现场观测实验,并通过研发的定制软件反演得到污染气体NO2和SO2的斜柱浓度值,对图像像元精准匹配得到污染气体二维分布信息图,为大气环境监测提供了更直观和实时的测量手段。