张 静, 李 钰, 任舜文

(华东理工大学信息科学与工程学院, 上海 200237)

研究论文

基于有限脉冲响应滤波器的实时小波算法及其在色谱信号解析中的应用

张 静, 李 钰*, 任舜文

(华东理工大学信息科学与工程学院, 上海 200237)

对色谱信号进行正确、实时地解析是利用色谱仪器系统进行过程检测的关键技术。小波分析方法能对色谱信号进行有效解析,但实时性一直是制约该方法在色谱信号处理中广泛应用的一个瓶颈。多分辨率分析与重构算法(MALLAT)是小波分析的一种快速算法,但是该算法的实时性仍然需要进一步提高。针对这一问题,以MALLAT算法原理为基础,利用在分解重组过程中两组滤波器系数之间的关系提出了一种基于有限脉冲响应滤波器的实时小波分析算法。仿真结果表明,与经典的MALLAT小波快速算法相比,该算法在保持信号有效分解与重构的基础上,运行耗时明显缩短,实时性得到较大提高。

色谱信号;在线模式;小波分析;实时性

色谱信号的处理系统是色谱仪器系统的重要组成部分。色谱法具有测量范围广、分析速度快、选择性高等特点,在许多领域都有着广泛的应用。目前,在色谱信号处理中,小波变换已经成为一种重要的方法,在基线校正、平滑滤噪、谱图压缩、重叠色谱峰分析、特征提取等方面都有成功的应用[1-8]。

传统的分析方法属于离线模式的实验分析,对算法的实时性要求不高;然而,在色谱的在线分析中采用实时性高的算法却十分必要,可以根据数据分析及时调节过程参数,有利于生产过程的优化。多分辨率分析与重构算法(MALLAT)是一种针对离散正交小波变换的快速算法[9-12]。在对色谱信号进行分解重组的过程中,MALLAT算法自底到顶按不同尺度进行小波分解并进行逆向重构,该算法的不足之处在于,需在信号采集完成之后进行处理,不能对信号进行实时分解与重构[13]。近年来,一些针对实时性要求的小波快速算法相继出现。Qian等[14]提出了基于滑动时窗的实时小波算法,通过预先存贮一部分小波系数避免重复计算,以较小的存储空间为代价提高计算效率。但是,利用滑动时窗来进行数据处理时,总是需要剪取实时数据的最新一段,这对数据的处理效率有一定影响。Rao等[15]提出了垒墙式算法,采取隔一抽取的方法每进行两次下层卷积之后计算一次上层卷积,像垒墙一样,整体上纵向发展。垒墙式算法虽然在一定程度上达到了实时性的要求,但只是对数据进行逐段的处理,并没有达到最大化的实时性要求。Wang等[16]提出的基于实时小波算法的暂态扰动检测FPGA(field-programmable gate array)实现算法,利用了FPGA高集成性、高速度、延迟时间一般为ns级等特点,与经典MALLAT算法相结合,降低了数据处理时间,利用高性能的硬件特性达到实时性的要求,可是该算法并未从软件上进行改进,因此实现成本有所增加。Liu等[17]提出的压电振动陀螺仪实时小波滤波方法,由无限冲击响应低通滤波器先对信号低频部分进行截取,之后采用滑动时窗算法对数据进行实时小波处理。因为滑动窗口的宽度要保持为包含最新2J个采样数据,所以对数据的处理效率也有一定影响。

为进一步提高小波变换在色谱信号解析中的实时性,利用MALLAT算法处理过程中各组滤波器系数之间存在翻转、移位等关系,本文提出了基于有限脉冲响应滤波器的实时小波分析算法。该算法的特点是在采集数据的同时可以对数据进行更及时的处理,克服了垒墙式算法对数据进行下层卷积处理之后再进行上层卷积处理的局限性,对采集的新数据进行实时处理,实现对色谱信号的实时分解与重构。

本文首先对MALLAT小波快速算法的原理进行了分析;然后提出了基于有限脉冲响应滤波器的实时小波算法。仿真实验表明了本文提出的算法的准确性及有效性,该算法在色谱信号的解析中具有应用价值。

1 MALLAT小波快速算法的基本原理

由小波多分辨率分析的定义,假定ψ(t)和φ(t)分别为尺度空间V0和小波空间W0的一个标准正交基函数,h0(n)为低通滤波器系数,h1(n)为相对应的高通滤波器系数。由于V0∈V-1、W0∈V-1,ψ(t)、φ(t)可用V-1空间的正交基ψ-1,n(t)线性展开为:

(1)

其中,展开系数h0(n)=<ψ, ψ-1,n>, h1(n)=<φ, φ-1,n>。因此:

(2)

将此式对时间进行伸缩和平移变换,并令m=2k+n(其中k、n均为序列的偏移量，均从0开始递增，m为总的偏移量),则:

(3)

f(t)=∑kcj-1,k2(-j+1)/2ψ(2-j+1t-k)

(4)

将f(t)分别投影到Vj、Wj空间,则有:

f(t)=∑kcj,k2-j/2ψ(2-jt-k)+

∑kdj,k2-j/2φ(t)(2-jt-k)

(5)

上式中,cj,k为剩余系数或尺度系数,dj,k为小波系数:

将式(3)代入式(4),得:

cj,k=∑mh0(m-2k)∫Rf(t)2(-j+1)/2×

(6)

(7a)

(7b)

将Vj空间尺度系数cj,k进一步分解下去,可分别得到Vj+1和Wj+1空间的剩余系数cj+1,k和dj+1,k,即著名的MALLAT塔式分解算法,算法原理见图1。

图 1 MALLAT分解算法原理Fig. 1 Decomposition algorithm based on the principle of MALLAT

设函数f(t)∈Vj-1,将二尺度方程式代入式(5),得:

f(t)=∑kcj,k∑nh0(n)2-(j+1)/2ψ(2-j+1t-2k-n)+∑kdj,k∑nh1(n)2-(j+1)/2ψ(2-j+1t-2k-n)

上式两边同时对ψj-1,m(t)进行内积,并且由尺度函数的正交性得:

cj-1,m=∑kcj,kh0(m-2k)+∑kdj,kh1(m-2k)

(8)

上式即为小波变换的系数重建公式。重建算法原理见图2。

图 2 MALLAT重构算法原理Fig. 2 Reconstruction algorithm based on the principle of MALLAT

2 基于有限脉冲响应滤波的实时小波算法

利用小波系数、尺度系数与滤波器系数的关系,本文对式(7)进行变形,在分解与重构过程中,利用低通、高通滤波器具有的翻转移位等特点,可以用低通滤波器来代替高通滤波器,降低在分解与重组过程中因卷积而相乘、相加的次数,达到快速分解重组的目的。

式(7)的物理意义表明,cj-1,m(n)分别经过脉冲响应为h0(n)、h1(n)的数字滤波器的作用之后,再分别对其进行偶数抽取得到cj,m(n)、dj,m(n)。此处,h0(n)和h1(n)分别为低通滤波器H和高通滤波器G的冲击响应序列,并且g(k)=(-1)1-kh(k),h(k)=h(-k)。其中，g(k)和h(k)分别为高、低通滤波器系数序列。

设原始信号为cm(n),低通滤波器的长度为N,令t=k-2n,则k=2n+t;由于h的长度为N,所以式(7a)可以写成:

(9)

又由滤波器的性质h(k)=h(-k)可得:

cm(2n+1+t)]-h(0)cm(2n)

(10)

同理,式(7b)可以写成:

(11)

亦可得:

cm(2n+1+t)]-h(0)cm(2n-1)

(12)

对于重构算法,将式(8)写为:

cm(n)=2pm(n)+2qm(n)

(13)

a)在pm(n)中,令n-2k=t,则:

(14)

上式可以写为:

(15a)

(15b)

b)在qm(n)中,由滤波器的性质h(k)=h(-k)可写成如下形式:

(16)

令1-n+2k=t,则:

(17)

由滤波器性质h(t)=h(-t)可得:

(18a)

(18b)

因此,基于有限脉冲响应的实时小波信号重构方法为:

由式(19)可得,该算法的输入量是随着一定量数据的输入开始计算的,伴随着数据的不断输入对数据进行不断的卷积,并得出相应的小波系数和尺度系数。

用经典MALLAT小波快速算法对色谱信号数据进行处理时,在分解阶段中,由于需对色谱信号原始数据进行2∶1抽样,故其分解过程可以隔一个数据进行。因此,在分解过程中,每计算一个低频系数cj,k与一个高频系数dj,k,需要作2N次浮点乘法运算和2(N-1)次浮点加法运算。对M(M=2m,为偶数)个原始离散信号样本进行解析,一级小波分解需要进行(M/2)×M×N次浮点乘法运算与(M/2)×2×N次浮点加法运算;而在一级小波合成运算中恢复一个点需要作N次浮点乘法运算与(N-1)次浮点加法运算,则恢复M个点需要进行M×N次浮点乘法运算与M×(N-1)次浮点加法运算。由上述分析可得,用经典MALLAT算法对色谱数据进行一次分解与重组需要进行2×M×N次浮点乘法运算与2×M×(N-1)次浮点加法运算。用本文提出的基于有限脉冲响应滤波器的实时小波算法对色谱信号进行处理时,在分解阶段需要进行(M/2)×2×N次浮点乘法运算与(M/2)×2×(N-1)次浮点加法运算。由式(19)可知,分解过程中每计算一次得到一个低频系数与一个高频系数,合成过程中利用滤波器系数存在的特点可以在计算过程中一次恢复相邻的两个数值,合成阶段需要进行2×(N-1)次浮点乘法运算,恢复M个点则需要进行(M/2)×2×(N-1)次浮点乘法运算,由此可得利用改进的算法对色谱数据进行一次分解与重组处理时需要进行M×N+M×(N-1)次浮点数乘法运算与2×M×(N-1)次浮点加法运算。由上述分析可知本文提出的算法与经典算法相比可以减少M次浮点数乘法运算,随着样本量越来越大,改进算法进行浮点数运算的计算量会远远低于经典算法,其增益十分显著。

综上所述可得,虽然本文提出的算法和经典算法都可以对数据进行有效的处理,但是在对数据进行解析时,两种算法处理的数据输入量和时间复杂度不同。经典算法须在数据全部接收之后才能对数据进行解析,且时间复杂度较大,而改进算法在数据进行采集时就可以对已采集的数据进行处理,且随着色谱信号数据样本量的增加,本文提出的算法相较于经典算法处理速度有明显的提升,能更好地达到对数据进行实时性处理的目的。

3 实验仿真

3.1 色谱峰信号的建模

色谱峰信号通常可以通过Gaussion峰信号和Tsallis峰信号进行建模[18,19]。其中,Tsallis峰的参数q取值不同可以得到Gaussion峰信号波形和Lorentzian峰信号波形,Tsallis峰信号的数学表达式为:

(20)

式中,A表示峰幅值的大小,μ表示峰顶点的位置坐标,σ表示峰宽大小。显然,通过调节参数q可以得到不同的峰形。

在本文构造的色谱峰模型中,取q=1.5,A分别为2、1.5、1.5、0.85,μ分别为20、22、28、40,σ分别为2、1.5、1、1。构造的色谱峰信号见图3a。

3.2 算法仿真

本文将该改进算法用于色谱峰的分解与重构,首先对色谱峰进行建模,并对信号进行随机高斯噪声加噪处理,近似模拟出色谱峰的形态,见图3b(其中横坐标表示采样点数,纵坐标表示信号的相对强度)。

在仿真实验中,本文选用的滤波器见图4,滤波器长度为16。在小波理论中,分解与重构所用低通与高通滤波器系数都存在对称、翻转、移位等关系。因此,其他滤波器谱图都可由该低通滤波器通过对称、翻转、移位等操作得到。

使用经典MALLAT算法对色谱信号进行一级分解,得到的分解近似信号与细节信号见图5。

图 4 低通滤波器频谱Fig. 4 Spectrum of low pass filter

图 3 色谱(a)原信号及(b)其含噪信号的波形图Fig. 3 Chromatographic (a) original signal and (b) its noisy signal

图 5 两种算法分别对色谱信号进行一级分解的频谱图Fig. 5 Chromatographic signal spectrum decomposed at first stage with two different algorithms

在多分辨率分析的定义中,虽然研究了伸缩与平移的特性,但是没有研究尺度函数支撑长度的问题。在实际工程应用中,考虑到计算的有效性,即尺度函数和小波函数都是有限支撑的,其所能卷积的数据长度也是有限的。本文提出的基于有限长脉冲响应滤波器的实时小波算法,由滤波器长度来决定初次卷积的数据长度,同时把新采集的数据及时加入可卷积的数据长度之中,在采集数据的同时对采集到的与滤波器长度匹配的数据进行分解与重构。用该算法对信号进行第一次分解的频谱图见图5。

重构信号时,使用经典算法对信号进行一级重构,以及使用本文提出的实时小波算法处理结果见图6。

在仿真实验中,本文用两种不同方法在一台双核CPU(1.70 GHz, 2.40 GHz)、安装Windows 8.1系统的电脑上进行色谱信号解析,并在Matlab R2009a下进行仿真测试。分别用经典算法与本文提出的实时小波算法对本文建模的含噪色谱信号(信噪比为33.688 3 dB,均方误差为0.039 8)进行处理,得到的重构信号见图7。

经两种算法处理后,信号的信噪比变为42.651 4 dB,均方误差变为0.014 2。实验中对分解与合成所耗费的时间进行了多次统计计算并求平均值,结果见表1。

图 6 两种算法分别对色谱信号进行一级重构的频谱图Fig. 6 Chromatographic signal spectrum reconstructed at first stage with two different algorithms

ItemAveragetime-consuming(classicalalgorithm)/sAveragetime-consuming(proposedalgorithm)/sLowfrequencysignal(decomposition)4.87865×10-42.4437×10-6Highfrequencysignal(decomposition)4.17291×10-42.8103×10-6Lowfrequencysignal(reconstruction)4.43476×10-42.8653×10-5Highfrequencysignal(reconstruction)3.84254×10-43.0425×10-5Reconstructionoriginalsignal1.73227×10-51.5469×10-5

图 7 两种算法分别对色谱信号进行重构的频谱图Fig. 7 Chromatographic signal spectrum reconstructed with two different algorithms

仿真实验中,本文用两种算法对含不同信噪比(signal to noise ratio, SNR)的信号进行了处理,结果见表2。

表 2 两种算法对不同信噪比信号进行处理的实验结果对比分析

由图7、表1、表2可得,在同等条件下,改进算法可与经典算法获得相同的去噪处理效果。但改进算法对色谱信号进行解析比经典算法所用时间明显减少。改进算法降低了在卷积过程中浮点数运算的次数,随着原始色谱采样数据量的增加,改进算法的时间增益较经典算法会越来越明显,实时性的优势也会更加显著。

用经典算法与实时小波算法分别对一段自高压液相色谱仪中采集的信号进行解析,实验结果如图8所示。自色谱仪采集的含噪信号信噪比为48.858 0 dB,经处理过后信噪比为63.918 9 dB,经典算法耗时2.003 5×10-5s,实时小波算法耗时1.240 2×10-5s。

图 8 (a)含噪信号及(b)其重构信号的色谱图Fig. 8 Chromatographic (a) noisy signal and (b) its reconstruction signal

4 结论

本文提出了一种基于有限脉冲响应滤波器的实时小波算法,并将其应用于对色谱信号的去噪及实时处理的仿真实验中。与经典MALLAT算法相比,本算法对信号的分解与重组过程进行了更为清晰的解析,而且克服了经典算法在对信号进行实时处理方面的局限性,明显提高了实时性,取得了很好的解析与重构效果。

[1] Liu Z C, Chen X G, Li Y F, et al. Journal of Beijing University of Chemical Technology, 2007, 34(1): 71

刘志成, 陈祥光, 李宇峰, 等. 北京化工大学学报, 2007, 34(1): 71

[2] Wu Q-Q, Lee J-P, Park M-H, et al. Procedia Eng, 2014, 97:

819

[3] Li Y, Wang S W, Lin Z P. Journal of East China University of Science and Technology (Natural Science Edition), 2014, 40(6): 752

李钰, 王圣伟, 林兆培. 华东理工大学学报(自然科学版), 2014, 40(6): 752

[4] Wang S C, Wang G T, Wang S J, et al. Low Voltage Apparatus, 2009(5): 12

王世成, 王国涛, 王淑娟, 等. 低压电器, 2009(5): 12

[5] Duan Y Q, Zeng H, Zhuo H. Mechanical Engineering & Automation, 2016(3): 31

段益群, 曾辉, 卓辉. 机械工程与自动化, 2016(3): 31

[6] Wei Z, Sun Z L, Xie Y J, et al//International Symposium on Parallel Architectures, Algorithms and Programming, IEEE, 2010: 138

[7] Zhang B, Zhao D E, Liu J, et al. Journal of Gun Launch & Control, 2015, 36(3): 68

张斌, 赵冬娥, 刘吉, 等. 火炮发射与控制学报, 2015, 36(3): 68

[8] Li L Y, Guo H J, Liu C, et al, Electronics Quality, 2014(2): 17

李丽颖, 郭红静, 刘灿, 等. 电子质量, 2014(2): 17

[9] Suma M N, Narasimhan S V, Kanmani B. Sadhana-Acad P Eng S, 2015, 40: 173

[10] Wang S, Wu X, Liu T, et al. Journal of Computer-Aided Design & Computer Graphics, 2016, 28(1): 106

王森, 伍星, 刘韬, 等. 计算机辅助设计与图形学学报, 2016, 28(1): 106

[11] Liu L, Zhang Z J, Liu L, et al. Journal of Jilin University (Science Edition), 2015, 53(2): 267

刘磊, 张子佳, 刘雷, 等. 吉林大学学报(理学版), 2015, 53(2): 267

[12] Luo C H, Du Y, Huang J S. Science & Technology Association Forum, 2012 (12): 91

罗朝辉, 杜宇, 黄激珊. 科协论坛, 2012(12): 91

[13] Xian G M, Wang Z Y. South China University of Technology, 2005: 653

[14] Qian F, Hu G M. Signal Processing, 2007, 23(3): 361

钱峰, 胡光珉. 信号处理, 2007, 23(3): 361

[15] Rao G A, Kang Y H, Chen L J, et al. Journal of Instrument, 2005, 26(2): 181

饶贵安, 康宜华, 陈龙驹, 等. 仪器仪表学报, 2005, 26(2): 181

[16] Wang K, Wang M, Wang J M. Application of Electronic Technique, 2012, 38(10): 6

王康, 王蒙, 王健民. 电子技术应用, 2012, 38(10): 6

[17] Liu Y, Duan Y Y, Zeng L L, et al. Piezoelectric & Acoustooptics, 2012, 34(1): 33

刘宇, 段耀宇, 曾燎燎, 等. 压电与声光, 2012, 34(1): 33

[18] Duan Z H. Guangdong Chemical Industy, 2015, 42(10): 56

段振华. 广东化工, 2015, 42(10): 56

[19] Lin Z P, Li Y, Wu H W. Journal of East China University of Science and Technology, 2014, 40(1): 91

林兆培, 李钰, 吴慧文. 华东理工大学学报, 2014, 40(1): 91

A real-time wavelet algorithm and its application in chromatographic signal analysis based on finite impulse response filter

ZHANG Jing, LI Yu*, REN Shunwen

(SchoolofInformationScienceandEngineering,EastChinaUniversityofScienceandTechnology,Shanghai200237,China)

The precise real-time analysis of the chromatographic signal is the key technology in process detection using chromatographic instrument system. Chromatographic signal can be effectively analyzed with the wavelet transform. However, the performance of real-time is not satisfactory and become a bottleneck of its wide applications. As a fast algorithm of wavelet analysis, the multi-resolution analysis (MALLAT) algorithm can still not satisfy the requirements of real-time performance. Considering the MALLAT algorithm and the relationship between filter coefficients of the decomposition and reconstruction, a novel real-time wavelet analysis method based on finite impulse response filter was proposed in this paper. Simulation results showed that the real-time performance of proposed algorithm is greatly improved on the basis of same effects compared with the classical MALLAT wavelet algorithm.

chromatographic signal; online mode; wavelet analysis; real-time performance

10.3724/SP.J.1123.2016.09044

2016-09-27

O658

A

1000-8713(2017)04-0368-07

* 通讯联系人.E-mail:liyu@ecust.edu.cn.