一种自适应层进式Savitzky-Golay光谱滤波算法及其应用
鲁一冰1,2, 刘文清1,2, 张玉钧1,2,*, 张恺1,2, 何莹1,2, 尤坤1,2, 李潇毅1, 刘国华1,2, 唐七星1,2, 范博强1,2, 余冬琪1,2, 李梦琪1,2
1. 中国科学院环境光学与技术重点实验室, 安徽光学精密机械研究所, 安徽 合肥 230031
2. 中国科学技术大学研究生院科学岛分院, 安徽 合肥 230026
*通讯联系人 e-mail: yjzhang@aiofm.ac.cn

作者简介:鲁一冰, 1991年生, 中国科学技术大学博士研究生 e-mail: yblu@aiofm.ac.cn

摘要

可调谐半导体激光吸收光谱技术(TDLAS)利用半导体激光器的可调谐和窄线宽特性, 通过选择特定气体的单条吸收线, 排除其余气体的干扰, 可以实现高精度、 高选择性的气体浓度测量, 在气体浓度检测系统中具有广泛的应用前景。 在不同的应用条件和环境下, 需要解决相应的硬件和数据处理方面的技术问题。 主要研究TDLAS技术机动车尾气CO组分浓度遥测系统中的光谱数据处理问题, 该系统利用路面漫反射回波信号遥测行驶中的机动车尾气CO组分浓度。 由于激光扫描光谱回波信号受到漫反射面情况变化、 空气环境变化、 尾气湍流影响等因素影响, 探测器收集到的信号不仅较弱同时也夹杂着多种噪声, 即测量光路信噪比较差, 故提出一种自适应层进式Savitzky-Golay(S-G)平滑滤波算法, 实现了对光谱进行滤波处理从而更加准确地反演CO浓度。 S-G滤波算法因其原理简单、 功能强大、 只需设置两个参数(窗口大小、 拟合阶数)等优点, 已广泛应用于光谱处理。 如何正确设置S-G算法参数使滤波效果在去噪不足和过度滤波之间找到平衡点, 是该滤波算法应用的一大难题。 设计的检测系统中, 测量光路光谱信号为非平稳信号, 噪声和有效信号幅度时变, 最佳窗口大小和多项式阶数随信号动态而变化, 且变化区间较大, 使用固定参数的S-G滤波器难以达到最佳效果。 提出的自适应层进式S-G平滑滤波算法, 通过逐层将测量光路光谱信号经过S-G滤波后, 与参考光路的光谱信号设置的参考段比对信号相关系数和信号一阶导相关系数的和, 以自适应得到逐层最优参数。 通过对信噪比从9.81~29.77的10组不同带噪光谱分析验证了该算法的有效性, 自适应层进式S-G算法能较好地去除噪声并还原带噪信号所携带的待测气体浓度信息, 与带噪光谱对比, 吸收光谱峰值最大误差由25.152%降至5.917%, 积分吸光度最大误差由18.1%降至3.9%。 在实现的系统中, 使用自适应层进式S-G算法对测量光路进行滤波处理, 并对不同车型、 不同排量、 燃烧不同油品的机动车在怠速和缓速通过(5 km·h-1)系统时其排放的CO浓度进行实时在线监测。

关键词: 可调谐半导体激光吸收光谱技术; 自适应层进式Savitzky-Golay平滑滤波; 机动车尾气遥测
中图分类号:TN219 文献标志码:A
An Adaptive Hierarchical Savitzky-Golay Spectral Filtering Algorithm and Its Application
LU Yi-bing1,2, LIU Wen-qing1,2, ZHANG Yu-jun1,2,*, ZHANG Kai1,2, HE Ying1,2, YOU Kun1,2, LI Xiao-yi1, LIU Guo-hua1,2, TANG Qi-xing1,2, FAN Bo-qiang1,2, YU Dong-qi1,2, LI Meng-qi1,2
1. Key Laboratory of Environmental Optics & Technology, Anhui Institute of Optics and Fine Mechanics, Chinese Academy of Sciences, Hefei 230031, China;
2. University of Science and Technology of China, Science Island Branch of Graduate School, Hefei 230026, China
*Corresponding author
Abstract

The tunable diode laser absorption spectroscopy (TDLAS) has the narrow linewidth characteristic of tunable diode lasers, so the gas concentration measurement of high precision and selectivity can be achieved by selecting the single absorption line of the specific gas to eliminate the interference of other gases, which has wide applications in the gas concentration detection. However, under different application conditions and environments, the people need to solve corresponding technical problems in hardware and data processing. In this paper, the TDLAS spectral data processing problem in the telemetry system for the vehicle exhaust CO concentration has been mainly studied. The system telemetered the exhaust CO concentration of the driving vehicle with the echo signals from the road diffuse reflection. Because the echo signal of laser scanning spectral is affected by factors such as the variation of the diffuse reflection surface, the change of the air environment, and the influence of exhaust turbulence, the signal collected by the detector is not only weak but also mixed with various noises, which means the SNR of the measured optical path is comparatively weak, so an adaptive hierarchical Savitzky-Golay (S-G) smoothing filter algorithm has been proposed in this paper, which can realize the spectral filtering processing to inversion the CO concentration more accurately. The S-G filtering algorithm has been widely used in spectral processing due to its advantages of such as simple principles, powerful functions and only two parameters setting (the window size and the fitting order). But how to set the parameters of the S-G algorithm correctly to balance the filtering effect between insufficient denoising and excessive filtering is a big problem for its application. In the designed detection system, the spectral signal of the measured optical path is non-stationary signal, and the amplitudes of the noises and effective signals are time-varying. So the optimal window size and polynomial order are changing with the signal dynamics of large range. As a result, it’s difficult to achieve the optimal filtering effect through S-G filters with fixed parameters. With the adaptive hierarchical S-G smoothing filter algorithm proposed in this paper, the sum of the signal correlation coefficient and the first derivative of the signal from measured light path spectrum signal after S-G filtering layer by layer and the reference section set by the spectrum signal of the reference light path have been compared, and then the optimal parameters of each layer can be obtained adaptively. With the analysis on 10 groups of band noise spectrum of which the signal to noise ratios(SNR) are from 9.81 to 29.77, the algorithm could effectively restore the concentration information carried by the band noise signals of the gas to be measured. Compared with the band noise spectrum, the maximum error of the absorption spectrum peak has dropped from 25.152% to 5.917%, and the maximum error of the integral absorbance has decreased from 18.1% to 3.9%. In the realized system, the adaptive algorithm has been used for the filtering processing of the measured optical path. The CO concentration emitted by motor vehicles of different models, displacement and oil product use has been monitored online in real time when they passed the system at idle and low speed(5 km·h-1).

Keyword: Tunable diode laser absorption spectroscopy; Adaptive hierarchical Savitzky-Golay algorithm; Telemetry measurement system of motor vehicle exhaust
引言

机动车发动机的燃烧过程产生的尾气排放污染物中有数十种不同的化合物, 其中, 对环境影响最大的是CO, CO2, NOx, HC和颗粒物。 对路网中机动车排放的CO进行准确的测量可用于分析机动车发动机燃烧过程、 并且有助于排放控制策略的制定[1]。 目前用于环境污染气体检测的方法主要分为非光学法和光学法, 其中光学法是利用光与气体分子相互作用的原理进行探测, 包括FTIR, NDIR[2], TDLAS[3]等方法。 TDLAS(可调谐半导体激光吸收光谱)方法因其相应快、 精度高、 无侵入、 光程远等特点, 已成为气体检测领域最常用方法之一。

提出一种基于TDLAS的垂直式机动车尾气排放CO组分遥测系统, 可无接触地实现对通过系统光路的机动车尾气中CO组分进行实时浓度遥测。 由于系统采用垂直收发结构, 相比于传统的水平式收发结构, 适用性更强, 可在实际路网中独立测量车道内车辆尾气排放组分CO浓度信息。 但相比于水平收发结构, 探测器收集到的回波信号经过路面粗糙漫反射后光强较弱, 又因为开放空间中测量CO浓度受大气湍流、 杂散光等影响, 噪声较大, 故系统信噪比较低, 需考虑在系统中使用光谱数据处理算法以保证气体反演浓度精度。

Savitzky-Golay滤波算法由于滤波效果强大、 只需设置多项式拟合阶数和窗口大小两个参数, 在光谱数据处理中, 应用已十分广泛。 针对相对平稳的信号, 可以通过离线分析滤波效果选择固定最优参数, 而在TDLAS遥测系统中, 测量光路信号、 噪声幅值时变, 使用固定参数的S-G算法难以达到理想滤波效果, 故需使用自适应S-G滤波算法, 旨在动态寻找合适的参数。 其中, Gabriel等认为S-G滤波器滤除的噪声与仪器噪声自相关最高时为最优[4], Phillip Barak等将拟合信号的残差平方和作为平滑效果的评判标准[5]。 本文提出了一种自适应层进式S-G滤波算法, 使用带噪的测量光路光谱信号通过逐层滤波并与参考光路所设置的参考段的相关系数及一阶导相关系数比对, 选取每层最优参数, 最终输出去噪后的信号。

在实验部分, 使用搭载了层进式S-G滤波算法的机动车尾气CO遥测系统, 对不同车型、 不同排量、 燃烧不同油品种类的机动车在怠速和缓速通过遥测系统时尾气排放的CO组分浓度进行了测量, 验证了系统的可行性及所提出算法的有效性。

1 系统设计及原理
1.1 垂直式机动车尾气CO组分浓度遥测系统

所设计的基于TDLAS的垂直式机动车尾气排放CO组分遥测系统结构如图1所示, 其中, 使用温度控制模块和电流控制模块对中心波长位于2 333 nm的DFB二极管激光器进行调谐, 并通过控制锯齿波信号发生模块实现对CO在该波长的吸收线扫描, 同时输出触发信号控制数据采集。 由激光器出射的激光通过分束器后, 光强较弱的一束光通过10 cm的参考光路气室入射到探测器靶面上, 其中参考光路气室中充满5%浓度的CO气体, 该分路用于数字锁线、 积分浓度定标以及自适应层进式S-G算法选参标准, 称其为参考光路。 另一束光通过长光纤接入龙门架上的遥测收发装置, 经过准直后的光束垂直入射于路面上, 经过路面漫反射, 后向散射光通过前置在望远镜系统的平凸透镜收集并聚焦于探测器靶面上, 此路称为测量光路。 探测器将光信号转换为电流信号后通过长导线进入主机中的I/V转换、 放大模块后, 进入数据采集卡进行AD转换, 最后由上位机程序对光谱信号进行滤波、 浓度反演等处理。

图1 系统结构图Fig.1 System structure

1.2 TDLAS浓度反演算法原理[6, 7]

根据Lambert-Beer定律, 当红外光通过一段气体介质时, 当激光器的发射波长和待测气体某个吸收谱线的中心波长相同时, 激光能量发生衰减且光的衰减程度与待测气体的浓度以及光程成正比关系, 激光能量衰减满足关系

It(ν)=I0(ν)exp(-kνL)(1)

式(1)中, I0(ν )为入射光强(单位: mW), It(ν )为透射光强(单位: mW), ν 为入射光的频率(单位: cm-1), L为有效吸收光程(单位: cm), kν 为光谱吸收系数(单位: cm-1), 实际中常利用吸光度来描述目标气体的吸收特性。

在TDLAS技术直接吸收方法中, 通常使用吸光度曲线积分面积反演气体浓度, 即积分吸光度A

A=aνdv=S(T)PcL(2)

其中, α ν 为吸光度, P为气体压强(atm), S(T)为温度T时的吸收线强(单位: cm-2· atm-1, 1 atm=101 325 Pa), ϕ (ν )为吸收线的归一化线型函数(单位: cm), c为吸收气体的摩尔分数。

积分吸光度在、 压强、 气体吸收线强相等的情况下, 吸光度曲线积分面积与待测气体浓度、 测量光程呈正比关系, 故

c=AL×LA×c(3)

其中, A为标气的积分吸光度, L为标气的有效测量光程, A为待测气体吸收线型积分吸光度, L为待测气体的有效测量光程。

2 自适应层进式S-G滤波

Savitzky-Golay平滑滤波算法(S-G)通过一个移动窗口, 对窗口内部元素做多项式最小二乘拟合, 得到窗口中心位置处元素平滑后的值[8]。 该方法即能够实现良好的去噪效果, 又可以较好地保留光谱信号中的有用信息。

S-G算法另一大优势在于只有两个参数需要设定, 即窗口大小(本文中所有关于窗口大小的叙述均为单侧数据点数)和多项式拟合阶数, 对于给定的信号, 这两个参数选择的正确与否会直接导致滤波效果的不同, 低阶大窗口会造成信号失真, 并且削弱吸收峰强度, 同时拉宽吸收线型, 难以保留所需要的信息; 高阶小窗口虽可以较好的保留信号信息, 但同样对噪声的滤除效果较弱。

1.1节系统设计中, 描述了激光器通过电流和温度调谐后, 输出信号经过分束器分为参考光路和测量光路。 由于参考光路和测量光路的光束是从同一激光器发射, 故由激光器带来的噪声可视为等同。 同时由于参考光路较短, 激光通过准直器直射探测器靶面, 相比于测量光路信号是通过开放空间并经由路面漫反射而来, 受到漫反射面情况变化、 空气环境变化、 尾气湍流等因素影响, 探测器收集到的信号不仅较弱同时也夹杂着多种噪声, 当两路信号经过放大后参考光路噪声影响更加明显。 自适应层进式S-G滤波器算法流程图如图2所示。

图2 自适应层进式S-G滤波器算法流程图Fig.2 Flow chart of adaptive hierarchical S-G filter algorithm

由于无法从测量光路信号中自适应的分辨噪声和有效信号, 可以使用信噪比较高的参考管光路信号作为参考, 与滤波后的信号进行比对, 观察滤波效果。 现选取参考光路光谱无吸收段作为参考段, 并设置每层S-G滤波阶数不超过5阶, 单侧数据点数不超过参考段横坐标最小值(WSmax), 遍历阶数和窗口单侧数据点数对测量光路信号进行滤波, 定义参考段相似度为S如式(4)

S=R+R1d(4)

其中, R为参考光路参考段与测量光路当前层滤波后参考段的线性相关度, R1d为参考光路参考段与测量光路当前层滤波后参考段一阶导的线性相关度。

如果当前层参考段相似度大于上一层, 则保留当前层滤波后结果, 在实际应用中, 为考虑时间成本, 可适当设置参考段参考段相似度增幅阈值STH, 即如果当前层相似度增幅没有超过STH, 则停止进入下一层滤波并输出当前层滤波结果。

现使用参考光路信号作为标准信号, 通过加入不同幅度的噪声, 观察自适应递进式S-G算法效果, 定义系统信噪比SNR

SNR=10log10PSPN=10log10n=1LS(n)2n=1LN(n)2(5)

其中, PSPN分别为标准信号和噪声信号有效功率, S(n)和N(n)分别为标准信号光谱及噪声信号, L为信号长度。

逐层滤波结果如图3所示, 其中, 添加噪声后SNR为15.56。 图3中内插图为所选用的参考段。

原始光谱信号及其一阶导如图3(a)和(b)所示, 经过加噪后, 参考段线性相关度为0.983 196, 参考段一阶导数相关度为0.023 491 7, 带噪光谱信噪比为15.56。 如图3(e)和(f)所示, 当带噪光谱信号经过第一层阶数为3, 单侧窗宽为87的S-G滤波后, 与参考段线性相关度增至0.998 105, 参考段一阶导线性相关度增至0.998 421, 此时光谱信号已经有较大优化, 信噪比提升至32.732, 但从信号参考段及信号一阶导谱线线型存在的毛刺可以判断少量噪声没有滤除。 经过第二层2阶, 单侧窗宽为62的S-G滤波后, 与标准光谱线性相关度增至0.998 958, 一阶导线性相关度增至0.999 815, 此时, 噪声已基本滤除, 信噪比为33.855, 且由图3(g)和(h)所示, 光谱信号及其一阶导平滑并且没有发生形变。

图3 信噪比为15.56的吸收光谱信号使用自适应层进式S-G滤波算法进行光谱数据处理
(a): 原始信号; (b): 原始信号一阶导; (c): 原始信号加噪; (d): 加噪信号一阶导; (e): 加噪信号经过第一层自适应层进式S-G滤波后波形; (f): 加噪信号经过第一层滤波后信号一阶导: (g): 加噪信号经过第二层自适应层进式S-G滤波后波形; (h): 加噪信号经过第二层滤波后信号一阶导
Fig.3 Absorption spectrum (SNR=15.56) processing by using adaptive hierarchical S-G filter
(a): Raw sigal; (b): Fist derivative of raw signal; (c): Signal with noise; (d): First derivative of signal with noise; (e): Signal with noise using first layer of S-G filter; (f): First derivative of signal using first layer of S-G filter; (g): Signal using second layer of S-G filter; (h): First derivative of signal using second layer of S-G filter

通过以下三个指标评价自适应递进式S-G滤波方法对不同信噪比下光谱信号去噪的效果:

①标准光谱Sstd与去噪后光谱S'的均方根误差RMSE, 如式(6)。 该指标可以表征去噪后光谱逼近标准光谱的程度。

RMSE=n=1L[Sstd(n)-S'(n)]2/L(6)

②去噪后光谱S'与标准光谱Sstd吸光度曲线经Voigt线型拟合后峰值接近率Rp。 该指标可以表征去噪后光谱吸收峰值与原峰值近似程度, Rp越接近1, 表示滤波后吸收谱谱峰退化程度越小。

③去噪后光谱S'与标准光谱Sstd吸光度曲线经Voigt线型拟合后积分吸光度比值RI。 该指标表征去噪后光谱积分吸光度与标准光谱积分吸光度的近似程度, RI越接近1, 表示自适应递进S-G滤波器结果对吸光度曲线回归效果越好。

对标准吸收光谱信号添加不同幅度的噪声, 得到了信噪比为9.81~29.77的10组不同带噪光谱, 分别使用自适应递进式S-G滤波算法, 通过上述3个指标, 分析在不同信噪比下的滤波效果。 其中, 不同信噪比光谱信号滤波前后信号RMSE, Rp, RI如图4(a), (b)和(c)所示, 伴随光谱信噪比降低, 带噪信号与标准信号的均方根误差、 吸收线型峰值比值、 积分吸光度比值均值均逐渐增大, 经过自适应递进S-G滤波后, 均方根误差RMSE虽有升高但稳定在0.1以内。 由于吸光度曲线峰值和积分吸光度可作为评判滤波器对反演待测浓度影响的标准, 由图4(b)和(c)可知, 自适应递进式S-G算法能较好的还原带噪信号所携带的待测气体浓度信息, 吸收峰值最大误差由25.152%降至5.917%, 积分吸光度最大误差由18.1%降至3.9%。

图4 不同信噪比光谱信号经过自适应层进式S-G算法前后的(a)均方根误差RMSE、 (b)吸收峰近似率Rp、 (c)积分吸光度比值RIFig.4 RMSE (a), Rp (b), RI (c) of different SRN absorption spectrums using adaptive hierarchical S-G filter

3 实验部分

第1节系统设计中, 详细介绍了基于TDLAS的垂直式机动车尾气CO检测系统组成及CO检测方法, 实测系统如图5所示, 其中龙门架内有效长宽均为3 m, 望远镜系统长0.5 m, 与龙门架衡量垂直悬挂, 测量机动车尾气CO组分有效光程为单程2.5 m。

图5 基于TDLAS的机动车尾气排放CO组分浓度遥测系统Fig.5 Remote CO concentration measurement system based on TDLAS of motor vehicle exhaust

使用搭载了自适应层进式S-G滤波算法的机动车尾气CO组分浓度遥测系统, 对表1中的3种车型分别在怠速和5 mk· h-1缓速两种工况下, 检测其排放的CO组分浓度。 其中, 为验证所搭建的机动车尾气排放CO组分浓度遥测系统的有效性, 在怠速状态下, 同时使用遥测系统和抽取式的基于NDIR原理的CO检测系统进行对比实验, 怠速状态下3车的CO排放如图6所示。

表1 测试车辆基本信息 Table 1 Basic information of test vehicles

图6 三辆机动车怠速状态下尾气CO组分浓度结果Fig.6 CO concentrations of 3 vehicles in idle status

三辆机动车怠速状态下排放CO组分浓度检测结果如图6所示。

使用所研制的机动车尾气排放遥测系统对机动车排放的CO气体浓度进行遥测时, 所检测的浓度结果与基于NDIR原理的抽取式CO检测模块检测所得浓度虽在取样方式上有较大差异, 导致检测结果不同(遥测系统所得浓度值为测量光路上包含的CO浓度)但是两种不同的系统反演所得的CO浓度趋势变化有高度一致性, 证明了所研制的系统CO模块测量结果的有效性, 少量浓度结果有偏差是因为取样方式的不同造成, 相比于抽取式, 遥测方式更容易受气象影响。 图6为三辆机动车在怠速状态下30次采样浓度检测结果, 轻型柴油车A在怠速状态下平均CO排放为59.06 mg· L-1, 较小排量的汽油车B平均CO排放为61.37 mg· L-1, 而较大排量的汽油车C尾气排放中的CO浓度较高, 平均为64.05 mg· L-1

图7为三辆机动车10次缓速通过系统测试光路时尾气排放CO浓度信息, 柴油车A10次通过测试光路平均CO浓度为58.007 mg· L-1, 机动车B平均排放CO浓度为57.143 mg· L-1, 机动车C平均排放CO浓度为65.895 mg· L-1

图7 三辆机动车怠速状态下尾气CO组分浓度结果Fig.7 CO concentrations of 3 vehicles in speed of 5 km· h-1

4 结 论

针对基于TDLAS的机动车尾气遥测CO组分遥测系统测量光路信噪比较差的问题, 提出了一种自适应层进式Savitzky-Golay平滑滤波算法, 通过分析对不同信噪比吸收光谱信号使用该方法滤波前后信号的均方根误差、 吸收峰线型比值、 积分吸光度误差, 证明了自适应层进式S-G算法的有效性, 并使用该算法搭载于所设计的机动车尾气排放CO组分浓度遥测系统, 对多种车辆在怠速和缓速通过遥测系统时CO排放进行了测量。

参考文献
[1] TANG Yuan-yuan, LIU Wen-qing, KAN Rui-feng, et al(汤媛媛, 刘文清, 阚瑞峰, ). Chinese Journal of Lasers(中国激光), 2011, 38(12): 221. [本文引用:1]
[2] Chen C, Zhang Y J, He Y. Infrared Technology, 2017, 39(6): 567. [本文引用:1]
[3] Arimoto H, Takeuchi N, Mukaihara S, et al. International Journal of Technology, 2011, 2(1): 9. [本文引用:1]
[4] Vivó-Truyols G, Schoenmakers P J. Analytical Chemistry, 2006, 78(13): 4598. [本文引用:1]
[5] Phillip Barak. Analytical Chemistry, 1995, 67(17): 2758. [本文引用:1]
[6] GAO Yan-wei, ZHANG Yu-jun, CHEN Dong, et al(高彦伟, 张玉钧, 陈东, ). Acta Optica Sinica(光学学报), 2016, 36(3): 0330001. [本文引用:1]
[7] YAO Lu, LIU Wen-qing, LIU Jian-guo, et al(姚路, 刘文清, 刘建国, ). Chinese Journal of Lasers(中国激光), 2015, 42(2): 305. [本文引用:1]
[8] Savitzky A, Golay M J E. Anal. Chem. , 1964, 36: 1627. [本文引用:1]