斜变滤波器时域加窗对图像重建的影响

2014-06-07 10:02乔志伟
计量学报 2014年3期
关键词:脉冲响应时域矩形

乔志伟

(中北大学电子与计算机科学技术学院,山西太原 030051)

斜变滤波器时域加窗对图像重建的影响

乔志伟

(中北大学电子与计算机科学技术学院,山西太原 030051)

提出了斜变滤波器的单位脉冲响应的时域加窗方法。斜变滤波器的单位脉冲响应应取奇数个点并偶对称,以保证其严格零相位。可以用来截取单位脉冲响应的窗函数有矩形窗、三角窗、汉宁窗、海明窗、布莱克曼窗等。仿真实验表明,矩形窗精度最好,三角窗因过度地压低了第一旁瓣,精度最差,其余窗对应的CT图像的精度介于两者之间。在对斜变滤波器的单位脉冲响应截断时应该使用矩形窗。

计量学;斜变滤波器;窗函数;图像重建;滤波反投影算法

1 引 言

CT(Computed Tomography)即计算机断层成像技术是最好的无损检测手段之一,已经广泛应用到医学、工业、地震学和考古学等领域[1,2]。X射线工业CT算法分为解析式算法和代数式算法2种,其中解析式算法居于优势地位[3,4]。在解析式重建算法中,又分为二维CT算法和三维CT算法,二维算法的经典代表是滤波反投影(filtered back projection,FBP)算法,而三维算法的经典代表是FDK算法[5,6]。从本质上说,由于FDK算法就是一种近似的滤波反投影算法,因此要研究影响此类算法的各种因素时,往往以FBP算法为研究对象。

在滤波反投影算法中,斜变滤波器的设计方法为【7,8】:

(1)通过对理想斜变滤波器ω加窗形成实用滤波器的频率响应H(ω)。加矩形窗的叫R-L滤波器;加sinc窗的叫S-L滤波器。

(2)由H(ω)计算得到单位冲激响应h(t)。

(3)根据冲激响应不变法,计算得到单位脉冲响应h(n)。

(4)根据投影信号和滤波投影信号均为有限长度的特点,截取h(n)为有限长度的单位脉冲响应。

在经典的FIR滤波器的窗函数法设计中,要用不同于矩形窗的更加平滑的窗去截取有限长度,以减少频谱泄露。对于上述第(4)步,各种文献均未见论及加窗的问题,即相当于均用矩形窗截断的方法。而根据信号处理的基本理论可知,加三角窗、汉宁窗、海明窗等可以有效减少频谱泄露。因此,有必要研究时域加窗对图像重建效果的影响。本文以S-L滤波器为例研究此问题。

2 斜变滤波器单位脉冲响应长度确定

根据离散时间信号卷积理论,2个有限长的信号做卷积,则卷积和信号的长度为2个信号的长度之和减1。

设离散的投影信号p(n)的取值范围为[-N,N],S-L滤波器的单位脉冲响应h(n)的取值范围为[-∞,∞],则离散的滤波投影信号g(n)的取值范围为[-∞,∞]。但在实际的CT工程实践中,滤波投影信号没有必要取无限长,而只需在[-N,N]这个范围内取值[9,10]。

根据卷积的定义,可知g(n)在2个边界点-N和N的值为

可见,h(n)的取值范围为[-2N,2N],共4N+1个点,就可以求得所需范围内的滤波投影信号g(n)。

总之,若投影信号和滤波投影信号的长度均取M点,则斜变滤波器的单位脉冲响应h(n)只需要取2M-1点。

需要注意的是,截取后的斜变滤波器h(n)必须是严格零相位的,因此,在截取2M-1个点时,必须以0为对称中心,即取值范围应该为[-(M-1),M-1]。

3 时域加窗方法

如果直接取[-(M-1),M-1]范围的h(n)作为滤波器的单位脉冲响应,那么这种方式相当于在无限长的h(n)之上加了1个以0为中心的奇数个点的矩形窗q(n)。

矩形窗是最简单的一种截取无限长单位脉冲响应为有限长单位脉冲响应的窗函数,因为是突然截断的,会引起强烈的频谱泄露。为了减小频谱泄露,就需要加一些不突然截断的窗函数。这些窗函数有三角窗、汉宁窗、海明窗、布莱克曼窗等,其性能参数见表1。-80

表1 6种窗函数的基本参数比较

显然,任何一种窗函数均不可能做到过渡带最窄,同时阻带衰减最大。过渡带宽度和阻带衰减程度是2个矛盾的参数,在工程中一般是根据对过渡带宽和阻带衰减的要求来选择窗函数。

可以认为表1中的后5种窗函数均为对矩形窗函数的优化和校正,其中汉宁窗、海明窗、布莱克曼窗其实均为余弦类的窗函数,只是各自性能有差异,而凯泽窗的波形也与海明窗近似。

在FIR滤波器的设计理论中,因为FIR滤波器的单位脉冲响应是因果序列,窗函数是加在[0,N-1]范围内的。而在对斜变滤波器的单位脉冲响应截取时,信号要偶对称,且选取奇数个点,故窗函数的表示方法与FIR滤波器设计理论中的形式不同。

为了表示方便,现假设要对h(n)在[-N,N]之间截取,则相应的窗函数为

(1)三角窗

三角窗变化比较剧烈,呈线性变化,而汉宁窗、海明窗和布莱克曼窗因均为升余弦窗,其变化较平滑。

4 仿真实验及结果分析

4.1 仿真模型及实验条件

采用滤波反投影算法对一个由4个细圆柱插到1个粗圆柱的模型的某断层进行重建。探测器分别在每个角度采集101点(采样间隔为0.1 cm),共采集180个角度的投影(角度步进间隔为1°),斜变滤波器采用S-L滤波器;采用卷积运算方式滤波,像素驱动方式为反投影。重建的图像为以旋转中心为图像中心的101×101的图像。

模型及其采集示意图见图1,模型参数见表2。

图1 实验用仿真模型及0°采集示意图

表2 模型参数

对S-L滤波器的h(n)采用5种不同的窗函数截断,研究不同窗函数对重建精度的影响。

对重建图像的质量评估采用归一化均方距离判据d和归一化平均绝对距离判据r这2个参数来评估。误差的公式分别为

4.2 实验结果

窗函数对重建精度的影响见表3,重建的图像见图2。

表3 S-L滤波器时域不同窗函数对重建精度影响

图2 时域加窗重建的图像对比

4.3 结果分析

由表3可以看出,按照精度从高到低的顺序为:矩形窗、海明窗、汉宁窗、布莱克曼窗、三角窗。

由图2可以看出,只有三角窗对应的CT图像失真明显,与其它几个区别较大,而其它几个重建图像的区别用肉眼不好分辨。三角窗的失真主要表现在中心圆变亮,同时中心圆周围的区域也明显变亮。

表3和图2反映出来的规律一致,能得出共同的结论。只是这个结论应该引起特别的注意,因为这里的加窗性能与FIR滤波器设计方法中的加窗性能不一致。

在FIR滤波器的设计中,为了减小吉布斯效应,应该使用过渡平滑的除矩形窗以外的其它滤波器。在牺牲过渡带宽不大的情况下,阻带性能从好到差的顺序为布莱克曼窗、海明窗、汉宁窗、三角窗、矩形窗。一般在工程中认为,这个顺序是一个性能由高到低的顺序。

然而,旁瓣衰减迅速(阻带衰减性能好)的窗函数应用到图像重建中,图像精度并没有提高。

采用过渡平滑的窗函数的优点是可以使得吉布斯效应减小,然而它的一个显著的缺点是窗函数的形状改变了信号时域值的真实表示。平滑可以减小吉布斯效应,从而提高精度;另一方面,加窗引起了信号时域的失真,从而降低精度,这是相互制约的一对矛盾。

当因为加窗引起的精度降低程度超过了平滑引起的精度提高的程度,总的加窗精度就要比不加窗(加矩形窗)的精度要低,加窗就没有意义了。

由表3可以看出,矩形窗的精度最高,三角窗的精度最低,其它升余弦类的窗函数的重建精度介于其间。这说明在对斜变滤波器的单位脉冲响应时域加窗时,加窗引起的精度降低程度超过了平滑过渡引起的精度提高程度。总的来说,加窗反而引起了精度的降低。

这种现象从h(n)的特征也可以解释。h(n)的第一旁瓣特别重要,它的大小对重建图像的精度影响重大。而无论是三角窗还是海明窗等,均对第一旁瓣做了较大的衰减,从而引起较大的精度损失。图3给出了加矩形窗的S-L滤波器的单位脉冲响应与加三角窗的S-L滤波器的单位脉冲响应的对比。

图3n取[-100,100]h(n)在[2,5]范围的波形

由图3可以看出,加三角窗后的单位脉冲响应绝对值变小,在[2,5]范围内三角窗的波形整体提升到了矩形窗对应的波形上面,这相当于对每一个旁瓣的幅值均做了压缩。而这种压缩必然导致较大的误差,所以三角窗形式的斜变滤波器重建出来的CT图像误差较大。而其它的余弦类的窗函数也均有不同程度地改变了旁瓣幅值,导致了精度的降低。

可见,在斜变滤波器的设计中,精度最好的加窗方式是矩形窗,没有必要加起平滑作用的窗函数。

5 结 论

在用FIR滤波器的窗函数法设计滤波器时,根据过渡带宽和阻带衰减速度的要求,可以适当选择三角窗、汉宁窗、海明窗等不是突然截断的窗函数,以减少吉布斯效应。

但是,在用窗函数截取无限长的单位脉冲响应为有限长序列时,只能用矩形窗,而不能用三角窗和升余弦类窗。其原因在于,斜变滤波器的单位脉冲响应的第一旁瓣对重建精度有显著影响,非矩形窗使得第一旁瓣的幅值被压缩,从而增大了误差。

总之,在设计滤波反投影算法的斜变滤波器的单位脉冲响应时,只能用矩形窗。

[1] 王召巴.任意角度入射的三维CT投影方法[J].计量学报,2002,23(1):1-5.

[2] 张朝宗,郭志平,张朋,等.工业CT技术和原理[M].北京:科学出版社,2009.

[3] 张顺利,张定华,赵歆波.一种射束与像素的快速遍历和求交算法[J].中国图象图形学报,2009,14(10):1961-1965.

[4] 王东明,卢虹冰,张军英.基于统计特性的小波噪声抑制在低剂量CT中的应用[J].中国图象图形学报,2008,13(05):876-881.

[5] 马建华,陈武凡.不含旋转角度微分的螺旋锥束CT重建[J].中国图象图形学报,2008,13(04):647-653.

[6] 张全红,路宏年,杨民.锥束工业CT中Feldkamp重建算法的快速实现[J].计算机工程与设计,2006,27(6):931-933.

[7] 乔志伟,魏学业,韩焱.解析法图像重建中的插值技术研究[J].计算机工程与设计,2009,30(9):2213-2216.

[8] 曾更生.医学图像重建[M].北京:高等教育出版社,2010.

[9] 李保磊,傅健,魏东波,等.工业计算机断层成像系统转台旋转中心的确定[J].航空动力学报,2009,24(07):1544-1548.

[10] 李保磊,傅健,黄巧珍,等.一种基于正弦图的工业CT系统转台旋转中心自动确定方法[J].航空学报,2009,30(07):1341-1345.

Im pact of Time Domain W indow-adding for the Ramp Filter on Image Reconstruction

QIAO Zhi-wei
(School of Electronics and Computer Science and Technology,North University of China,Taiyuan,Shanxi030051,China)

The time domain window-adding method for the unit impulse response of the ramp filter is proposed.the unit impulse response of the ramp filter should be symmetrical and its length should be odd number to ensure that it is strict zero-phase.There are many windows to truncate the unit impulse response,for example,rectangle window,triangle window,Hanning window,Hamming windows,Blackman window etc.The simulation experiment demonstrates that the rectangle window is the best and the triangle window is the worst,for the triangle window depresses the first side petal overly.The precision of the CT images reconstructed by the otherswindows is in-between.The rectanglewindow should be used when the unit impulse response of the ramp filter is been truncated.

Metrology;Ramp filter;Window function;Image reconstruction;Filtered back projection algorithm

TB973

A

1000-1158(2014)03-0268-04

10.3969/j.issn.1000-1158.2014.03.15

2012-01-16;

2013-08-30

国家自然科学基金(61071193,61171179)

乔志伟(1977-),男,山西洪洞人,中北大学副教授,博士,主要研究方向为信息获取及处理技术。zhiweiqiaook@nuc.edu.cn

猜你喜欢
脉冲响应时域矩形
基于重复脉冲响应的发电机转子绕组匝间短路检测技术的研究与应用
两矩形上的全偏差
化归矩形证直角
基于时域信号的三电平逆变器复合故障诊断
从矩形内一点说起
基于极大似然准则与滚动时域估计的自适应UKF算法
基于时域逆滤波的宽带脉冲声生成技术
脉冲响应函数下的我国货币需求变动与决定
基于有限元素法的室内脉冲响应的仿真
基于时域波形特征的输电线雷击识别