第26卷,第2期 光谱学与光谱分析Vol 126,No 12,pp37223762006年2月 Spectroscopy and Spectral Analysis February ,2006 基于卷积型小波包变换的谱线自动提取方法刘中田1,吴福朝1,罗阿里2,赵永恒21.中国科学院自动化研究所模式识别国家重点实验室,北京 1000802.中国科学院国家天文台,北京 100012摘 要 天体光谱中的谱线包含重要的天体物理信息。
文章提出一种基于卷积型小波包变换的谱线自动提取方法。
该方法由以下主要步骤组成:(1)将观测光谱进行4层卷积型小波包变换;(2)对第四层小波包系数,采用区域相关算法以及阈值处理方法进行噪声处理;(3)选择中高频小波包系数进行谱线特征重构;(4)根据重构后的谱线特征,利用谱线搜索方法,在观测光谱中提取谱线。
作者在实验中用恒星、正常星系和活动星系光谱进行谱线提取测试,结果表明该方法具有对噪声鲁棒和谱线提取准确等特点。
用该方法提取sloan digital sky survey (SDSS )光谱中的谱线后,计算了红移并与SDSS 给出的红移进行了对比,实验结果间接验证了该方法提取谱线的有效性。
主题词 卷积型小波包;区域相关算法;天体光谱;谱线提取中图分类号:TN91117 文献标识码:A 文章编号:100020593(2006)022******* 收稿日期:2004212228,修订日期:2005205228 基金项目:国家“863”计划(2003AA133060)和国家自然科学基金(60402041)资助项目 作者简介:刘中田,1979年生,中国科学院自动化研究所博士研究生引 言 谱线提取在光谱分析中起着非常重要的作用[1],作为一种预处理手段,对基于谱线的光谱分类有着直接的影响。
以往的谱线提取方法多是以人工参与的半交互方式进行的,常用的天文软件如MIDS ,FIGARO 和IRA F 都是如此。
这些处理方法需要借助专家知识来标定谱线的位置,既费时也费力。
目前我国正在建造的大天区面积多目标光纤光谱望远镜(简称L AMOST ),建成后每个观测夜可以得到2万~4万条天体光谱,面对如此巨大的海量数据,采用自动的光谱分析方法已成为必然的选择[1]。
本文旨在为L AMOST 光谱分析和分类系统研究可靠的光谱谱线提取方法。
由于各种天体光谱的形态不同,再加上连续谱和噪声的影响,自动提取谱线是一项相当困难的工作。
自动提取谱线一般的思路是先拟合连续谱,再对光谱归一化,然后进行去噪处理。
由于连续谱的物理成因比较复杂,再加上观测设备(望远镜,仪器等)内部参数的影响,使得人们很难得到连续谱的准确模型。
以往的连续谱拟合方法如中值滤波方法、多项式插值法以及小波变换方法都是通过对光谱作了一定程度的平滑而实现的,因此拟合出来的连续谱在强谱线附近不够准确,降低了谱线提取的精度[1,2]。
赵瑞珍采用小波变换零交叉点方法[2],可以提取出比较准确的发射线或吸收线,但对于既含吸收线又含发射线的光谱,不能同时提取出两种谱线。
另外由于该方法直接对光谱的一阶小波变换的高频系数进行零交叉点搜索,这一搜索算法比较容易受噪声干扰,因此不适合低信噪比光谱的谱线提取。
本文研究了一种基于卷积型小波包变换的谱线自动提取方法。
本文的主要特点有:在小波包域内结合了区域相关算法和阈值处理方法进行噪声处理;选择相应的小波包系数进行谱线特征重构;根据重构后的谱线特征,利用谱线搜索方法,在观测光谱中提取谱线。
与小波变换零交叉点方法相比,本文方法利用小波包频带细化的特点,在小波包域内进行信噪分离,有效地减小了噪声对谱线提取的影响。
另外本文方法可以同时提取出天体光谱中的吸收线和发射线。
1 卷积型小波包变换 多分辨分析可以对信号进行有效的时频分解,但由于其尺度是按二进制变化的,所以在高频频段其频率分辨率较差。
小波包分析能够为信号提供一种更加精细的分析方法[3],它将频带进行多层次划分,对多分辨分析没有细分的高频部分进一步分解,从而可以详细分析信号与噪声在中、高频方面的特征。
一般的小波包变换算法,当分解层数较大时,各序列的数据长度就会变得很短,这对后续的数据处理很不利。
文献[4]提出的卷积型小波包变换算法较好地解决了这一问题。
关于卷积型小波包变换,分解公式为,f j+1,2n (l )=12∑k ∈Zh (k )・f j ,n (l -2j k )f j+1,2n+1(l )=12∑k ∈Zg (k )・f j ,n(l -2jk )(1)重构公式为fj ,n(l )=12∑k ∈Zh (k )・fj +1,2n(l +2j k )+12∑k ∈Zg (k )・fj+1,2n+1(l +2j k )(2)其中j 是小波包变换尺度,h (k )称为低通滤波系数,g (k )称为高通滤波系数。
由(1)式的分解公式可以看出,卷积型小波包变换在迭代运算时并没有进行隔二抽一的采样,而只是对上一层的分解结果进行2j k 位移,所以每层中的各频道序列长度始终和原始信号相同,这一特性使得对小波包系数进行进一步的处理比较方便。
2 小波包域降噪方法 对于含噪光谱信号,经过卷积型小波包分解以后,在分解的各分量中,原始信号对应的小波包系数有较强的相关性,而随机噪声却没有这种明显的相关性。
因此,可以利用小波包系数对应点处的相关性来区分系数的类别,从而进行取舍,达到去噪的效果。
为了避免小波包系数微小的偏移可能导致所求相关系数的不准确性,本文采用区域相关系数[5]。
在区域[x k -m ,x k +m ]上的小波包系数的和N j ,i (x k )=∑x k +ml =x k -mfj ,i(l )(3)称为局部区域和系数。
式中,f j ,i 为尺度j 上的第i 频段小波包系数,m 为第x k 点的区域相关点数。
第x k 点的区域相关系数为C N j ,i (x k )=N j ,i (x k )N j ,i+1(x k )(4)第x k 点的归一化区域相关系数为N j ,i (x k )=C N j ,i (x k )(P N j ,i /P C Nj ,i)1/2(5)式中P N j ,i =∑x kN2j ,i(x k ),P C Nj ,i=∑x kC 2N j ,i (x k)。
若| N j ,i (x k )|≥M ・|N j ,i (x k )|,M 为经验系数,则认为f j ,i (x k )为信号,否则认为是随机噪声。
经过区域相关算法之后,由于信号较大边缘会附带一些噪声,可能会对谱线提取造成一定的影响,因此还需要采用阈值处理。
阈值方法分为软阈值和硬阈值两种,最初由Do 2noho [6]提出。
本文采用硬阈值方法,如下所示。
取λ= σ12ln (N ),对f j ,i (x )作阈值处理f j ,i (x )=f j ,i (x ),|f j ,i (x )|≥λ0,|f j ,i (x )|≤λ(6)其中,λ为阈值, σ为中高频小波包系数f j ,i (x )均方差的均值,N 为信号离散点数。
然后利用f j ,i (x )重构信号,就可以得到去噪后的光谱信号。
3 谱线特征分析 对于含噪光谱信号,f (x )=s (x )+n (x ),s (x)为原始光谱信号,n (x )为随机干扰噪声。
巡天观测光谱的噪声来源较为复杂,可以近似为高斯白噪声。
光谱中谱线位置对应于其极值点,由于受连续谱和噪声的影响,使得通过找光谱的极值点判别谱线比较困难。
光谱主要是由连续谱、谱线和噪声构成的,在小波包变换域内一般可以较好地提取出谱线对应的特征。
谱线和噪声的能量都集中于中高频部分,分离它们比较困难,我们采用了区域相关算法和阈值处理方法,尽量将噪声能量减到最小。
实验中,我们选择了db2型的小波对应的小波包[7],采用卷积型小波包变换算法,对光谱信号进行4层小波包变换,得到第4层的小波包系数为f 4,1(x ),f 4,2(x ),…,f 4,16(x ),然后对f 4,3~f 4,16这14个小波包系数,利用区域相关算法和阈值处理方法,按照频段从高到底顺序进行处理,令f 4,1(x ),f 4,2(x )小波包系数为零,最后进行小波包重构,重构所得结果f ′(x )即为我们求得的谱线特征。
假设谱线的理想模型为一光滑函数h (x ),利用上面的变换过程重构后的谱线特征h ′(x )如图1所示。
Fig 11 Ideal model of the spectral line (a ),The original spectral line h (x );(b ),The feature of t hespectral line h ′(x ) 由上图可以看出,原始谱线h (x )的极值点a 对应于谱线特征h ′(x )上的极值点A 点。
A 点的纵坐标y A 的绝对值|y A |是|h ′(x )|的最大值,并且有y A >0 当h (x )为发射线时y A <0 当h (x )为吸收线时 定义1:设定一个阈值δA ,当|y A |>δ时,我们称极值点A 为有效极值点,它的横坐标即为一条谱线的位置。
对应于有效极值点A 的另两个极值点B ,C 我们称之为支撑点。
定义2:谱线特征h ′(x )中的D ,E 两点称为有效极值点A 的特征边界点,区间[x D ,x E ]为特征区间,可以通过支撑点B 向左搜索过零点得到D 点,通过支撑点C 向右搜索过零点得到E 点,过零点x k 由h ′(x k )・h ′(x k +1)≤0定义[2]。
如图1所示,光谱信号的谱线范围[x b ,x c ]包含于特征区间[x D ,x E ]中,在f ′(x )中检测到一条谱线后,一般要对373第2期 光谱学与光谱分析原始光谱f (x )在区间[x b ,x c ]内进行端点线性插值f (x )=f (x c )-f (x b )x c -x b(x -x b )+f (x b ),x ∈(x b ,x c )(7)并且把特征区间[x D ,x E ]内的f ′(x )值设为零,然后重新进行谱线的搜索。
当谱线比较密集的时候,特征区间之间会有重叠,增加了谱线搜索的难度,本文把有效极值点和谱线范围[x b ,x c ]的搜索分开进行,较好地解决了密集谱线搜索的问题。
4 谱线提取算法 如图1,如果我们能够找到b,c 两点,就可以准确地提取出谱线。
然而,实际中的谱线并不象h (x )这么光滑,而且谱线中还可能叠加噪声,且噪声一般很难完全滤除,因此我们选择如下的策略进行谱线范围的搜索。
(1)若y A >0,则在区间[x D ,x B ],[x C ,x E ]中对原始光谱f (x )分别搜索最小的极小值点(若无则取最小值点)得到b,c 两点;(2)若y A <0,则在区间[x D ,x B ],[x C ,x E ]中对原始光谱f (x )分别搜索最大的极大值点(若无则取最大值点)得到b,c 两点。