趙奉奎,徐曉美,呂立亞
(1.南京林業(yè)大學 汽車與交通工程學院,江蘇 南京 210037;2.東南大學 儀器科學與工程學院,江蘇 南京 210096)
光譜本底常使得凈峰面積的估算結果過大和峰位估算結果偏移[1]。如能量色散型X射線熒光光譜(EDXRF)分析作為一種多元素分析技術,能夠同時提供元素的定性和定量分析結果,但由于其特征峰疊加在本底之上,因此,分析時為獲取特征峰信息,必須扣除本底的影響[2]。
為了扣除本底,研究人員已經(jīng)提出了多種算法,一種方法是從硬件的角度出發(fā),通過改進光路,如設計三角形光路結構抑制本底的產(chǎn)生[3],更多的研究是從軟件的角度對已生成的光譜進行本底扣除,常采用多項式擬合[4-5]、傅立葉變換[6]、削峰法[2,7]和小波變換(WT)[8-11],以及近幾年研究較多的神經(jīng)網(wǎng)絡等[12]。如劉察等[13-14]將信號變換到小波空間,以減輕譜信號中基線的變化。大多數(shù)方法依賴于特征峰和本底在頻域分布在不同區(qū)域的特性,但當兩者在頻域有重合時,一般方法不易將其分開。小波變換具有多分辨率分析的特點,能夠對信號進行多尺度分析,非常適合區(qū)分特征峰和本底,然而直接利用小波變換進行本底扣除容易導致光譜波形畸變。Galloway等[8]提出了一種迭代小波變換(IWT)方法,相比于傳統(tǒng)小波變換能夠對信號進行更加細化的分解,可根據(jù)信號形態(tài)區(qū)分頻帶重復的信號,但未研究終值迭代條件。筆者利用信息熵作為迭代小波終值條件,改進了迭代小波變換[15]方法,但計算稍復雜。
基于此,本文提出了一種更加簡單的迭代終值準則計算方法,改進了基于迭代小波變換的本底扣除方法,旨在精確扣除本底,同時避免小波變換導致的波形畸變。
小波變換能夠在時域和頻域同時進行信號分析。在小波變換過程中,利用選定的小波基,即母小波函數(shù)ψ,對信號進行分解。在實際運用時,更多的是利用離散小波變換(DWT)。為了實現(xiàn)多分辨率分析,還需要構造尺度函數(shù)φ。經(jīng)壓縮因子2j和平移因子2jn壓縮平移的小波函數(shù)和尺度函數(shù)如式(1)~(2)所示:
(1)
(2)
根據(jù)Mallat快速小波變換,信號f分解的逼近系數(shù)aj[n]和細節(jié)系數(shù)dj[n]可根據(jù)式(3)獲得:
aj[n]=〈f,φj,n〉
dj[n]=〈f,ψj,n〉
(3)
通過與濾波器h[n]與g[n]的卷積和下采樣,根據(jù)式(4)計算逼近系數(shù)aj+1[n]和細節(jié)系數(shù)dj+1[n]:
(4)
圖1 離散小波分解過程Fig.1 Process of discrete wavelet decomposition
根據(jù)收斂定義,對正整數(shù)n,存在正整數(shù)N,當n>N時,對任意指定的ε>0,都有|fn+1-fn|max<ε,則認為f收斂。對于光譜分析,如果連續(xù)3次迭代過程中fn和fn+1之間的最大誤差均小于指定ε,則認為f收斂,具體計算過程將在下述步驟中詳述。
圖2 在不同分解層次進行小波變換所得低頻逼近圖Fig.2 Approximations obtained through wavelet transform at different levels
圖3 基于改進迭代小波變換的本底扣除流程圖Fig.3 Program flow chart of background removing method based on the improved IWT
進行迭代小波變換之前,需先確定分解層數(shù)。隨著迭代小波變換過程的進行,光譜中的特征峰成分會逐步消減。為盡量減少人工干預,所有小波變換均在相同的分解層次上進行。最優(yōu)分解層次的確定通過在不同分解層次上進行小波變換,比較其低頻逼近的方法實現(xiàn)。最優(yōu)的分解層次需使信號經(jīng)小波分解后,所有的本底成分均包含在低頻逼近中,盡管其中可能仍含有一些特征峰成分。經(jīng)過試驗比較,選取DB4小波基進行小波分解。信號在不同分解層次進行小波變換所得低頻逼近如圖2所示。圖中A7、A8、A9分別表示第7層、第8層、第9層分解的低頻逼近。當對信號進行較低分解層次的小波變換時,如第7層,所得低頻本底振蕩較為劇烈,與實際本底相差較大。當對光譜信號進行較高分解層次的小波變換時,如第8層和第9層,所得本底的變化趨勢與實際本底比較一致,但第9層所得低頻逼近過于平坦,所以選擇第8層作為最優(yōu)分解層次。在后續(xù)的迭代小波變換過程中,均在第8層對信號進行小波分解。
基于迭代小波變換的本底扣除過程如圖3所示,具體流程如下:
①fm[i]用來表示經(jīng)m-1次迭代小波變換后的光譜,其中m≥1為整數(shù),表示迭代次數(shù),i表示光譜通道,原始光譜標記為f1[i]。
②在優(yōu)化分解層次上對fm[i]進行小波分解,將得到的逼近系數(shù)am[i]作為第m次估計的本底。
③計算|am[i]-am-1[i]|max=errm,a0與f1[i]相等,表示原始光譜。將errm和用戶指定的誤差ε進行比較,比較次數(shù)由times進行索引,times初始值為0。當對所有的通道均存在errm<ε時,表示相鄰兩次所估計的本底足夠一致,將times=times+1,跳至步驟④,否則,重置times=0,并轉至步驟⑤。
④如果times≥3,則連續(xù)4次小波變換所得估計本底足夠一致,認為所得估計本底已經(jīng)收斂,最后一次小波變換所得低頻逼近am[i]作為最終估計的本底,從原光譜中減去am[i]即可實現(xiàn)本底扣除。
⑤比較所有通道值的fm[i]和am[i],將各通道最小值替換fm[i],即:
(5)
通過削去較大值的方法,迭代地去除低頻逼近中的特征峰成分,令m自增,轉至步驟②。
為驗證算法的有效性,首先利用仿真光譜對算法進行分析。已有研究表明,EDXRF光譜中的特征峰可由高斯峰表征,本底可由解析函數(shù)表征[1]。因此,本研究采用高斯峰為特征峰,多項式函數(shù)為本底,構建仿真光譜。
為驗證算法是否對噪聲敏感,在光譜上增加了信噪比為38 dB的噪聲,仿真光譜如圖4A所示。
根據(jù)“1.2”所述小波分解層次選擇方法,選取第8層為最優(yōu)分解層。相鄰2次小波分解所得errm見表1。結果顯示,隨著迭代過程的進行,errm明顯衰減。由于最大峰高7 000,選取ε=10可以認為2次所得本底足夠一致。10次迭代后,連續(xù)3次errm<ε,表明所得估計本底已近似收斂。
Iteration number(m)12345678910errm6 12060018067301711864
圖5 迭代小波變換估計的光譜與傳統(tǒng)小波變換估計的光譜對比Fig.5 Comparison of the estimated spectra with IWT and non-iterative WT
迭代過程見圖4B,圖中第1個估計本底(1st Estimated bakground)表示第1次小波變換后得到的本底。由圖可知,第一次小波變換得到的本底與實際信號本底相差很大。隨著迭代小波變換的進行,所得估計本底與仿真光譜的本底逐漸靠近。7次迭代后,所得估計本底收斂到仿真光譜本底。從原始仿真光譜中減去估計的本底即可實現(xiàn)本底扣除,去除本底后光譜見圖5(IWT spectrum)。同時也可以看到,在存在噪聲的前提下,該算法仍能夠準確估算本底。
作為對比,對該仿真光譜利用多項式擬合法和傳統(tǒng)小波分解進行本底估計。利用DB4小波對仿真光譜在第8層進行分解,不進行迭代,小波逼近系數(shù)作為估計的本底。實驗發(fā)現(xiàn),從原始光譜中去除小波變換估計的本底使光譜產(chǎn)生了較大畸變(圖5 WT spectrum)。相比之下,迭代小波變換能準確扣除本底,且未引入波形畸變。
為進一步驗證算法的有效性,利用迭代小波變換對EDXRF光譜儀所測光譜進行實驗。光譜儀采用Si-Pin探測器,高壓為45 kV。通過對光譜的觀察和實驗,選擇DB6在第8層對光譜進行分析。光譜迭代分解過程見圖6A。由圖可知,經(jīng)一次小波分解后,所得估計本底與光譜實際本底相差較大,而經(jīng)小波分解迭代6次后,所得估計本底趨于收斂。
迭代小波變換過程中的errm列于表2中,考慮到特征峰最高值為8 000,ε=10時認為所得估計本底足夠接近。經(jīng)7次迭代后,err7<ε,9次迭代后,可認為迭代過程收斂。表中數(shù)據(jù)顯示的本底變化趨勢與圖6A中一致。
表2 實驗光譜相鄰兩次估計本底誤差Table 2 Errs of two consecutive estimated backgrounds of the experimental spectrum
圖6 實驗光譜迭代小波變換過程(A)和原始實驗光譜與扣除本底光譜對比圖(B)
Fig.6 Process of iteration for the experimental spectrum(A) and comparison of the original spectrum and background-subtracted spectrum(B)
從原始光譜中減去第9次小波變換所得估計本底,即可去除本底。原始實驗光譜與本底扣除光譜對比如圖6B所示,由圖可知,本底已完全扣除,去除本底后的光譜基線與水平坐標非常接近,且波形未發(fā)生畸變,峰位未偏移,即使弱峰也保持完整,因此該算法也為低含量元素的分析奠定了基礎。
本文基于迭代小波變換的本底扣除改進算法提出了一種迭代停止標準。該算法通過對仿真光譜和實驗光譜的分析,準確扣除了光譜本底,證明了算法的有效性,為光譜凈峰面積的計算和峰位的計算奠定了基礎。與傳統(tǒng)小波變換本底扣除法相比,該算法能夠對信號進行更加精細的分解,從而分離特征峰和本底重復的部分,且不會導致波形畸變。與多項式擬合法相比,本底扣除結果更加準確。本研究迭代的思想也可由傅立葉或多項式擬合來實現(xiàn),為光譜分析提供了有力的工具。