二维分数阶卡尔曼滤波及其在图像处理中的应用

1、第 32卷第 12期 电 子 与 信 息 学 报 Vol.32No.12 2010年 12月 Journal of Electronics & Information Technology Dec. 2010二维分数阶卡尔曼滤波及其在图像处理中的应用左 凯 孙同景 李振华 陶 亮(山东大学控制科学与工程学院 济南 250061摘 要:该文研究了二维分数阶卡尔曼滤波及其在图像增强与滤波中的应用问题。 首先基于分数微积分的定义, 建 立了二维线性离散系统的分数阶差分状态空间模型。 然后, 提出了一种可应用于图像信息处理的二维分数阶卡尔曼 滤波算法, 并通过实验验证了该文提出算法的有效性。仿

2、真结果证明,该算法增强了图像中的细节特征, 同时消弱 了图像中的背景噪声。关键词:图像增强;图像去噪;分数阶离散状态空间;二维分数阶卡尔曼滤波中图分类号:TN911.73 文献标识码:A 文章编号:1009-5896(201012-3027-05 DOI : 10.3724/SP.J.1146.2010.002322D Fractional Kalman Filter and Its Application to Image ProcessZuo Kai Sun Tong-jing Li Zhen-hua Tao Liang(School of Control Science and Engi

3、neering, Shandong University, Jinan 250061, China Abstract : This paper deals with the issue of 2D Fractional Kalman Filter (2DFKF and its applications to image enhancement and recognition. With the introduction of 2D fractional differential, 2DFKF recursive equation is first presented. Next, a stat

4、e space model of a image given and, based on this, the 2DFKF algorithm is proposed. Finally, an example is given to demonstrate the effectiveness of proposed algorithm and the simulation result shows that the details of the image are enhanced, while the background noise of the image is efficiently a

5、ttenuated.Key words:Image enhancement; Image denoising; Discrete fractional state-space systems; 2D Fractional Kalman Filter1引言将卡尔曼滤波器扩展到二维空间,可推广应用 于图像去噪。但是对于对比度低、边缘模糊、背景 起 伏 大 的 图 像 , 应 用 二 维 卡 尔 曼 滤 波 (Two- Dimensional Kalman Filter, 2DKF 在处理图像过 程中容易导致细节与边缘等重要信息损失。如何确 保图像滤波过程中的细节信息保护与增强,同时有 效滤除噪声一

6、直是该领域的热点研究问题。另一方 面,将分数微积分应用于图像处理,可以有效地增 强图像的纹理特征,主要方法有分数傅里叶变换 (Fractional Fourier Transform, FRT 、分数阶微分 梯度算子、分数阶微分掩模等 15 。但上述方法在 增强图像的同时,会降低图像的信噪比。 2006年文 献 6,7将分数阶微积分引入卡尔曼滤波,提出了分 数阶卡尔曼滤波算法 (Fractional Kalman Filter, FKF , 提高了信号估计精度, 但仅考虑了一维空间, 对于二维空间模型描述的图像信号系统的高性能滤 波问题尚待进一步深入研究。本文将综合应用分数 微积分和 2DKF

7、 ,针对 2DKF 在去噪过程中图像纹2010-03-11 收到, 2010-09-27改回国防科工委基础科研 (B142008.0209-08资助课题通信作者:左凯 zuokai 理细节及边缘信息模糊化的问题,研究提出了二维 分 数 阶 卡 尔 曼 滤 波 算 法 (Two-Dimensional Fractional Kalman Filter, 2DFKF , 并将该算法对 加噪图像进行去噪和增强处理, 验证算法的有效性。 2 二维线性离散分数阶差分状态空间模型 考虑如下传统的二维线性离散状态空间模型:,(, (, (, (, (, (, (, klk l Wm n m k n l m

8、n m n m n m n m n = + =+ s A s Bu wr Cs v(1 其中 (, m ns 和 (, m nr 为系统状态和观测值, (, m n u为控制输入, (, m nw 为系统噪声, (, m nv 为观测噪声。klA 为相应的系统参数矩阵, C 为转移矩阵, W 为邻域宽度。根据 Grunwald-Letnikv(G-L 的分数微积分定 义 8,可推导出系统式 (1中的 a 阶分数阶差分状态 方程为1, 1(, (1, 1 klm k n kjk ljm k n lm k j n l j = = A sA s (2 其中3028 电 子 与 信 息 学 报 第 3

9、2卷 0111, 0diag (1 , ;1diag (1 , ;01diag (1 , ;1kl kl mn i j m i j j kl n i j z i j a k l j a k l j j a k l j a j += = + = >> + < + A A A (1 (1 , 1,2, 3, ! z z z i i i a a a j j j + = " " (, , z i a R z m n mn + =为 (, i m n s 的分数阶数, j 为 数据序列长度。注意到,随着数据序列长度 j 增加,式 (2的计算量将不断增大。因此,在保证

10、精度的情 况下,通常需要将回溯计算长度 L 加以限定。 将式 (2代入系统式 (1则得二维线性离散分数 阶差分状态空间为1, 1, 00(, (1, 1 (, (, (, (, (, m k n k W j k l k l j m n m k j n l j m n m n m n m n m n = = + =+ s A s Bu w r Cs v(3 3 二维分数阶卡尔曼滤波基于系统描述式 (3, 2DFKF 算法为 (, (, (, (, (, m n m n m n m n m n =+ s s K r Cs (4a 1, 1, 0(, (1,1 (, m k n l W j k l

11、k l j m n mk j n l j m n = +sA s Bu (4bi i T T 1, (, (, (, m n m n m n m n =+K PC CP C R (4c i (, (, (, m n m n m n = P I K C P(4d i (' ' 1, 1T, 01, 1, , 0T, , (, (, Cov(1, 1 (1, 1(4em k n l W j j k l klk l j m k n lWj k l k l j j m n k l m n m k n l m k j n l j m k j n l j = = + + PA P A A

12、s s A Q式中 (, m n K 为增益矩阵, , m n R 为观测噪声方差, , m n Q 为过程噪声方差。算法证明如下 :( i ( ( ( ( 1TT 1, (, arg min(, (, (, , (, sm n m n m n m n m n m n m n = + ss s P ss r Cs R r Cs其中 ( ( ( *1, =, , , =(, |m n E m n z m n E m n z s s s s , *1Z 和 *Z 分别为前一步和当前观测序列。对上式求一阶导数,则i ( ( ( ( ( 1T 1, , , , (, , 0m n m n m n m

13、n m n m n + =P s s C R r Cs(5 ( i ( (i ( ( (i ( i ( i ( i ( i ( ( ( 11T 1, 1T 1, T T 11, T 1, , , , , , (, , (, , (, , , m n m n m n m n m n m n m n m n m n m n m n m n m n m n m n m n =+ += +s P C R CPs C R r PP C CP C R CPP s C R r (6 令 i ( i ( T T 1, (, =, (, + m n m n m n m n K PC C P C R , 由式 (

14、5和式 (6,得( ( i ( ( i ( ( ( ( ( T 1, T 1, , =, +(, , (, , (, , , , , , m n m n m n m n m n m n m n m n m n m n m n m n m n m n =+ s s P C R r K C K CP C R r s K r Cs 式 (4a和式 (4c得证。 ( ( ( ( (*11, 1, , 00*11, 1, , 0, , |(1, 1 , , (1, 1 , m k n l W j k l k l j m k n lW j k l k l j m n E m n z E m k j n l

15、 j m n m n z m k j n l j m n = = = + = + s s A s Bu w A sBu 式 (4b得证。( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( ( i ( ( ( TT T TT T T , =, , , , , , , , , (, , , , , , , , , , , , , , , =, , , m n E m n m n m n m n E m n m n m n m n m n m n m n m n m n m n m n E m n m n m n m n m n m

16、 n E m n m n m n m n m n m n =+ ×+ = × + +P s s s s s K r Cs s sK r Css I K C ss ss I K C K v v K I K C PI K C ( ( ( i ( T , , , , , m n m n m n m n m n = K R K I K C P第 12期 左 凯等:二维分数阶卡尔曼滤波及其在图像处理中的应用 3029式 (4d得证。结合文献 9提出的状态值与噪声值的 协方差计算方法以及结论且设: T T , , (, (, ,=(, (, m n m n E m n m n E m

17、n m n =R v v Q w w令 ( ( ( , , , m n m n m n = s ss,则 i ( ( ( ( ( ( ( ( T 1, 1, , 001, 1T, , 0T (, , , , , (1, 1(1, 1 , , m k n lW j k l k l j m k n lW j k l k l j m n E m n m n m n m n E m k j n l j m kj n l j E m n m n = = = = × + P s s s s A s A s w w 1, 1, , 01, 1T, , , 0T , , (1, 1Cov( 1, 1

18、 (1,1 ' ' m k n lW j k l k l j m k n lW j j k l k l k l j j m nk l m k j n l j m k j n l j m k j n l j = = + × +A P A A s s A Q式 (4e得证。至此,算法证明完毕。将 2DFKF 算法应用于图像增强和去噪,首先 要确定回溯长度 L 和分数阶参数 a ,其次将原二维 状态空间方程改写成二维线性离散分数阶差分状态 空间,然后应用算法对图像进行处理。4 仿真实验4.1 回溯长度 L 的确定给出以下系统参数矩阵,并分别令 10,20, L =50, 3

19、00, 000.1, , 0.1 0.3, 0.10.010.02 = A B C 0.71.2 = N ,其中 N 为系统分数阶矩阵。对比结果见图 1,方差 k 见表 1。由图 1和表 1可知,该系统回溯长度 50L =基 本上就可以满足精度。图 1 不同回溯长度 L 输出结果对比表 1 不同 L 下的方差对比L10 20 50 100 300 k 0.650730.13380.06759 0.06808 0.076214.2 分数阶参数 a 的确定利用 4.1节中的状态空间,将分数阶矩阵设定 为 1.2a = N , 0.2, 0.4, 0.6, 0.7, 0.9a =,考察系统误差。实验

20、结果如图 2所示。方差对比见表 2。表 2 不同 a 下的方差对比a0.2 0.4 0.6 0.7 0.9 k 0.066550.025860.00661 0.003120.00091当 a 设定过高时,若有阶跃信号时将会出现较 大误差或系统输出将会发散,详见图 3和图 4。综 上所述,选择分数阶系数 a 介于 0.5到 0.6之间。 4.3 图像的滤波与增强实验采用常见的 NSHP(Non-Symmetric Half- Plane 模型用于描述图像, 邻域宽度确定为 3W =, 状态空间模型如下:( ( ( ( ( ( ( ( 0,01,12,21,02,02,10,10,21,2(, 1

21、, 12, 2 3, 32, 3, 3, 1 , 2, 3 1, 3(, (, (, (, m n m n m n m n m n m n m n m n m n m n m n m n m n m n = + + + + + + + + +=+s A s A s A s A s A s A s A s A s A s w r Hs v 图 2 不同分数阶 a 输出结果对比 图 3 0.9a =输出结果 图 4 1.1a =输出结果3030 电 子 与 信 息 学 报 第 32卷基于该模型,其二维线性离散分数阶差分状态 空间可以直接改写为( ( ( ( 100,00101,10102,2101

22、,0010102,02,10,1(, 1, 1 2, 2 3, 3 2, 3, ( 3, 1 jj jj jj jj j jj j jm n m j n j m j n j m j n j m j n mj n m j n = + + + + + + s A s A s A s A s A s A s A s ( ( 1010100,21,20, 2 , 3(1, 3 (, j j jj j m n j m n j m n j m n = + + + A s A s w(, (, (, m n m n m n =+r Hs v其中 kl A 通过递归最小二乘法 (RLS估计, H 为适当 维数

23、的单位矩阵。设定回溯长度 10L =,分数阶 0.6a =。使用直方图均衡化和分数阶微分掩模方法对图 像增强, 效果见图 5和图 6。 利用本文提出的 2DFKF 算法对加有高斯白噪声 (均值为零,方差为 0.1 的图 像进行滤波增强,图 7显示了几种滤波方法的图像 滤波效果。表 3为上述滤波方法的信噪比。实验结果表明,二维分数阶卡尔曼滤波算法在 图像去噪能力上明显优于其它传统的滤波方法,且 具有提高图像细节和增强图像质量的能力。表 3 不同滤波方法处理后的信噪比信噪比噪声图像中值 滤波维纳 滤波2DKF 2DFKF SNR(dB 3.5029 4.3210 4.4699 6.3854 6.6

24、512 PSNR(dB27.568328.3864 28.5353 30.450830.71665 结束语本文基于分数阶微分和卡尔曼滤波原理,提出 了离散状态空间的二维分数阶卡尔曼滤波算法。利 用二维分数阶卡尔曼滤波算法对含有噪声的图像进 行了处理。 该算法有效地减弱了图像中的背景噪声, 增强了图像中的细节特征,验证了算法的有效性。 图 5 无噪声条件图像增强对比 图 6 有噪声条件图像增强对比第 12期 左 凯等:二维分数阶卡尔曼滤波及其在图像处理中的应用 3031 图 7 不同滤波方法处理后的图像参 考 文 献1朱呈祥 , 邹云 . 分数阶控制研究综述 J. 控制与决策 , 2009, 2

25、4(2: 161-169.Zhu Cheng-xiang and Zou Yun. Summary of research on fractional-order controlJ. Control and Decision, 2009, 24(2: 161-169.2Bai Jian and Feng Xiang-chu. Fractional-order anisotropic diffusion for image denoisingJ. Image Processing, 2007, 16(10: 2492-2502.3Pu Yi-fei, Zhou Ji-liu, and Yuan

26、Xiao. Fractional differential mask: a fractional differential-based approach for multiscale texture enhancementJ. Image Processing, 2010, 19(2: 491-551.4艾必刚 , 罗以宁 , 蒋涛等 . 分数阶微分梯度算子在图像增强中 的应用 J. 四川大学学报 (自然科学版 , 2009, 46(2: 343-347. Ai Bi-gang, Luo Yi-ning, and Jiang Tao, et al. Application of fractional-order differential gradient operator in image enhancementJ. Journal of Sichuan University(Natual Science Editon, 2009, 46(2: 343-347.5张军 , 韦志辉 . 分数阶图像去噪变分模型及投影算法 J. 计算 机工程与应用 , 2009, 45(5: 1-6.Zhang Jun and Wei Zhi-hui. Fractional variational model an