专利名称：：一种运动模糊和散焦复合模糊的图像复原方法
技术领域：
：本发明涉及一种运动模糊和散焦复合模糊的图像复原方法，其目的在于当图像中运动模糊和散焦同时存在时，建立复合降质模型，估计模糊参数，并实施图像复原，属于图像处理领域。
背景技术：
：为了改善图像质量，目前常采用基于模型的图像恢复方法进行处理，这种方法的核心是能够准确获知模糊模型和模糊参数。运动模糊和散焦是两种常见的图像模糊，运动模糊是由拍摄图像与摄像器材间的相对运动造成的，散焦则是由拍摄器材对焦不准造成的。针对这两种模糊，目前都是分别处理，分别估计其模糊参数并建立模糊模型，常用的模糊参数估计方法有基于频域的分析方法、基于边缘的分析方法、L-曲线法、自相关法等，其中基于自相关的方法由于抗干扰能力较强、计算量较小因而成为模糊参数估计中的一个热点。通过对现有技术文献的检索发现现有的模糊参数估计方法都是根据标准运动模糊和散焦模型分别建立的。当运动模糊和散焦同时出现时，没有确定复合模型可以参考，由于模糊相互叠加，上述参数估计方法无法直接应用。本发明将建立一种运动模糊和散焦同时出现的模糊模型，并提供一种对于该模糊参数估计方法，并最终实现图像复原。
发明内容技术问题本发明提供一种运动模糊和散焦复合模糊的图像复原方法，能够对运动模糊与散焦同时出现的图像进行复原。技术方案一种运动模糊和散焦复合模糊的图像复原方法，包括如下步骤步骤1降质图像去噪禾Ij用低通滤波对降质图像g(x，y)进行去噪处理，其中χ和y分别为行坐标和列坐标，且都为大于0的整数，得到去噪后的降质图像f'(X，y)，其算法如下均值为零的高斯分布P(x，y)的函数形式Pz(il，2)可以表示为_(i\2+i22)ρ(/1/2)=——e2σ1ζ2πσ2其中下标ζ表示均值为零的高斯分布的函数形式，il和i2分别为横坐标和纵坐标，且都为实数，σ为标准差，e为自然对数函数的底数，π为圆周率，设成像系统的加性噪声n(x，y)服从均值为零的高斯分布ρ(x，y)，这里σ取0.5，建立一个3X3的高斯白噪声模板n'(x,y)"0.01130.08380.0113"n'(x,y)=0.08380.61930.08380.01130.08380.0113将降质图像g(x，y)与高斯白噪声模板η'(x,y)进行卷积，得到M1XM2大小的去噪后的降质图像f'(χ,y)f'(x,y)=g(x,y)*n'(χ,y)其中*为卷积操作，步骤2模糊方向识别利用公知的傅立叶变换来估计去噪后的降质图像f'(x，y)的主方向θ和副方向θ‘，假设此时主方向θ的角度为α，-90°彡α<90°，副方向θ‘的角度为β，0°^β<180°，旋转图像至主方向θ的角度为0°，具体方法如下设散焦与运动模糊同时发生时运动模糊的方向为主方向θ，与主模糊方向垂直的方向为副方向Θ'，计算步骤1得到的去噪后的降质图像f'(x,y)的傅立叶变换F(G)1,ω2)=JJf(x,y)e-Mxe-]^ydxdy(O1和ω2为频率变量，F(G)1,ω2)为复数，j为虚数单位，去噪后的降质图像的傅立叶变换F(G)1,ω2)可以表示为Ρ(ω17ω2)=R(Q17Q2)+J'I(ω17ω2)其中R(G)1,ω2)为实部，I(G)1,ω2)为虚部，计算待估计降质图像f'(x，y)的能量谱E(G)1,ω2)Ε(ω”ω2)=R(q17ω2)2+Ι(ω17ω2)2能量谱E(G)1,ω2)中各元素为实数，将能量谱E(G)1,ω2)显示出来，在图中心位置会出现一个类似于椭圆的白色亮斑，亮斑周围会有连续的亮环，在所有穿过亮斑的弦中，最长的弦Imax所处的方向，取0°到180°之间的角做为β的值，与此垂直的方向角为α，α=β-90°，此时主方向θ的角度为α，副方向θ‘的角度为β，以垂直于去噪后的降质图像f'(x,y)，且过去噪后的降质图像f'(x,y)中心的直线为旋转轴，如果0°^α<90°，以α角顺时针旋转去噪后的降质图像f'(x，y)，如果-90°<α<0°，以α+90°角顺时针旋转去噪后的降质图像f'(x，y)，降质图像f'(x，y)旋转后，主方向θ的角度为0°，副方向Θ'的角度为90°，得到M3XM4大小的旋转后的去噪后的降质图像"(χ,y)；步骤3图像导数计算对步骤2中得到的旋转后的去噪后的降质图像f"(x,y)，计算主方向的导数Da(x,y)和副方向的导数Db(x，y)Da(x,y)="(x，y)*La(x，y)Db(χ,y)="(χ，y)*Lb(x，y)其中下标a和b分别表示主方向和副方向，La(x，y)和Lb(x，y)分别为主方向和副方向的导数卷积模板，此时主方向θ的角度为0°，副方向θ‘的角度为90°，主方向卷积模板La(x，y)=[_1，1]，副方向卷积模板忍(xj)=。1]，此时主方向的导数Da(x，y)和副方向的导数Db(x，y)大小都SM3XM4;步骤4图像导数的自相关计算步骤4.1利用公知的自相关运算，对步骤3中得到的主方向导数Da(X，y)和副方向导数Db(x，y)，分别计算主方向导数Da(x，y)的自相*Ca(x，y)和副方向导数Db(x，y)的自相关Cb(x，y)权利要求一种运动模糊和散焦复合模糊的图像复原方法，其特征在于，包括如下步骤步骤1降质图像去噪利用低通滤波对降质图像g(x，y)进行去噪处理，其中x和y分别为行坐标和列坐标，且都为大于0的整数，得到去噪后的降质图像f′(x，y)，其算法如下均值为零的高斯分布p(x，y)的函数形式pz(i1，i2)可以表示为<mrow><msub><mi>p</mi><mi>z</mi></msub><mrow><mo>(</mo><mi>i</mi><mn>1</mn><mo>,</mo><mi>i</mi><mn>2</mn><mo>)</mo></mrow><mo>=</mo><mfrac><mn>1</mn><mrow><mn>2</mn><mi>π</mi><msup><mi>σ</mi><mn>2</mn></msup></mrow></mfrac><msup><mi>e</mi><mrow><mo>-</mo><mfrac><mrow><mo>(</mo><msup><mrow><mi>i</mi><mn>1</mn></mrow><mn>2</mn></msup><mo>+</mo><msup><mrow><mi>i</mi><mn>2</mn></mrow><mn>2</mn></msup><mo>)</mo></mrow><msup><mrow><mn>2</mn><mi>σ</mi></mrow><mn>2</mn></msup></mfrac></mrow></msup></mrow>其中下标z表示均值为零的高斯分布的函数形式，i1和i2分别为横坐标和纵坐标，且都为实数，σ为标准差，e为自然对数函数的底数，π为圆周率，设成像系统的加性噪声n(x，y)服从均值为零的高斯分布p(x，y)，这里σ取0.5，建立一个3×3的高斯白噪声模板n′(x，y)<mrow><msup><mi>n</mi><mo>′</mo></msup><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow><mo>=</mo><mfencedopen='['close=']'><mtable><mtr><mtd><mn>0.0113</mn></mtd><mtd><mn>0.0838</mn></mtd><mtd><mn>0.0113</mn></mtd></mtr><mtr><mtd><mn>0.0838</mn></mtd><mtd><mn>0.6193</mn></mtd><mtd><mn>0.0838</mn></mtd></mtr><mtr><mtd><mn>0.0113</mn></mtd><mtd><mn>0.0838</mn></mtd><mtd><mn>0.0113</mn></mtd></mtr></mtable></mfenced></mrow>将降质图像g(x，y)与高斯白噪声模板n′(x，y)进行卷积，得到M1×M2大小的去噪后的降质图像f′(x，y)f′(x，y)＝g(x，y)*n′(x，y)其中*为卷积操作，步骤2模糊方向识别利用公知的傅立叶变换来估计去噪后的降质图像f′(x，y)的主方向θ和副方向θ′，假设此时主方向θ的角度为α，90°≤α＜90°，副方向θ′的角度为β，0°≤β＜180°，旋转图像至主方向θ的角度为0°，具体方法如下设散焦与运动模糊同时发生时运动模糊的方向为主方向θ，与主模糊方向垂直的方向为副方向θ′，计算步骤1得到的去噪后的降质图像f′(x，y)的傅立叶变换F(ω1，ω2)<mrow><mi>F</mi><mrow><mo>(</mo><msub><mi>ω</mi><mn>1</mn></msub><mo>,</mo><msub><mi>ω</mi><mn>2</mn></msub><mo>)</mo></mrow><mo>=</mo><mo>∫</mo><mo>∫</mo><msup><mi>f</mi><mo>′</mo></msup><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow><msup><mi>e</mi><mrow><mo>-</mo><mi>j</mi><msub><mi>ω</mi><mn>1</mn></msub><mi>x</mi></mrow></msup><msup><mi>e</mi><mrow><mo>-</mo><mi>j</mi><msub><mi>ω</mi><mn>2</mn></msub><mi>y</mi></mrow></msup><mi>dxdy</mi></mrow>ω1和ω2为频率变量，F(ω1，ω2)为复数，j为虚数单位，去噪后的降质图像的傅立叶变换F(ω1，ω2)可以表示为F(ω1，ω2)＝R(ω1，ω2)+jI(ω1，ω2)其中R(ω1，ω2)为实部，I(ω1，ω2)为虚部，计算待估计降质图像f′(x，y)的能量谱E(ω1，ω2)E(ω1，ω2)＝R(ω1，ω2)2+I(ω1，ω2)2能量谱E(ω1，ω2)中各元素为实数，将能量谱E(ω1，ω2)显示出来，在图中心位置会出现一个类似于椭圆的白色亮斑，亮斑周围会有连续的亮环，在所有穿过亮斑的弦中，最长的弦lmax所处的方向，取0°到180°之间的角做为β的值，与此垂直的方向角为α，α＝β90°，此时主方向θ的角度为α，副方向θ′的角度为β，以垂直于去噪后的降质图像f′(x，y)，且过去噪后的降质图像f′(x，y)中心的直线为旋转轴，如果0°≤α＜90°，以α角顺时针旋转去噪后的降质图像f′(x，y)，如果90°≤α＜0°，以α+90°角顺时针旋转去噪后的降质图像f′(x，y)，降质图像f′(x，y)旋转后，主方向θ的角度为0°，副方向θ′的角度为90°，得到M3×M4大小的旋转后的去噪后的降质图像f″(x，y)；步骤3图像导数计算对步骤2中得到的旋转后的去噪后的降质图像f″(x，y)，计算主方向的导数Da(x，y)和副方向的导数Db(x，y)Da(x，y)＝f″(x，y)*La(x，y)Db(x，y)＝f″(x，y)*Lb(x，y)其中下标a和b分别表示主方向和副方向，La(x，y)和Lb(x，y)分别为主方向和副方向的导数卷积模板，此时主方向θ的角度为0°，副方向θ′的角度为90°，主方向卷积模板La(x，y)＝[1，1]，副方向卷积模板此时主方向的导数Da(x，y)和副方向的导数Db(x，y)大小都为M3×M4；步骤4图像导数的自相关计算步骤4.1利用公知的自相关运算，对步骤3中得到的主方向导数Da(x，y)和副方向导数Db(x，y)，分别计算主方向导数Da(x，y)的自相关Ca(x，y)和副方向导数Db(x，y)的自相关Cb(x，y)<mrow><msub><mi>C</mi><mi>a</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow><mo>=</mo><msub><mi>D</mi><mi>a</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow><mo>⊗</mo><msub><mi>D</mi><mi>a</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow></mrow><mrow><msub><mi>C</mi><mi>b</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow><mo>=</mo><msub><mi>D</mi><mi>b</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow><mo>⊗</mo><msub><mi>D</mi><mi>b</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow></mrow>其中下标a和b分别表示主方向和副方向，表示相关操作，此时主方向导数的自相关Ca(x，y)和副方向导数的自相关Cb(x，y)的大小都为(2×M31)×(2×M41)；步骤4.2对主方向导数的自相关Ca(x，y)按照与主方向θ垂直的方向进行累加得ACa(x，y)，对副方向导数的自相关Cb(x，y)按照与副方向θ′垂直的方向进行累加得ACb(x，y)，累加方法如下<mrow><msub><mi>AC</mi><mi>a</mi></msub><mrow><mo>(</mo><mi>m</mi><mo>,</mo><mi>n</mi><mo>)</mo></mrow><mo>=</mo><msub><mi>C</mi><mi>a</mi></msub><mrow><mo>(</mo><mn>1</mn><mo>,</mo><mi>n</mi><mo>)</mo></mrow><mo>+</mo><msub><mi>C</mi><mi>a</mi></msub><mrow><mo>(</mo><mn>2</mn><mo>,</mo><mi>n</mi><mo>)</mo></mrow><mo>+</mo><mi>L</mi><mo>+</mo><msub><mi>C</mi><mi>a</mi></msub><mrow><mo>(</mo><mn>2</mn><mo>×</mo><msub><mi>M</mi><mn>3</mn></msub><mo>-</mo><mn>1</mn><mo>,</mo><mi>n</mi><mo>)</mo></mrow><msub><mo>|</mo><mrow><mi>m</mi><mo>=</mo><mn>1,1</mn><mo>≤</mo><mi>n</mi><mo>≤</mo><mn>2</mn><mo>×</mo><msub><mi>M</mi><mn>4</mn></msub><mo>-</mo><mn>1</mn></mrow></msub></mrow><mrow><msub><mi>AC</mi><mi>b</mi></msub><mrow><mo>(</mo><mi>m</mi><mo>,</mo><mi>n</mi><mo>)</mo></mrow><mo>=</mo><msub><mi>C</mi><mi>b</mi></msub><mrow><mo>(</mo><mi>m</mi><mo>,</mo><mn>1</mn><mo>)</mo></mrow><mo>+</mo><msub><mi>C</mi><mi>b</mi></msub><mrow><mo>(</mo><mi>m</mi><mo>,</mo><mn>2</mn><mo>)</mo></mrow><mo>+</mo><mi>L</mi><mo>+</mo><msub><mi>C</mi><mi>b</mi></msub><mrow><mo>(</mo><mi>m</mi><mo>,</mo><mn>2</mn><mo>×</mo><msub><mi>M</mi><mn>4</mn></msub><mo>-</mo><mn>1</mn><mo>)</mo></mrow><msub><mo>|</mo><mrow><mn>1</mn><mo>≤</mo><mi>m</mi><mo>≤</mo><mn>2</mn><mo>×</mo><msub><mi>M</mi><mn>3</mn></msub><mo>-</mo><mn>1</mn><mo>,</mo><mi>n</mi><mo>=</mo><mn>1</mn></mrow></msub></mrow>其中m和n分别为行坐标和列坐标，且都为大于0的整数，主方向导数自相关的累加ACa(x，y)的大小为1×(2×M41)，副方向导数的自相关的累加ACb(x，y)的大小为(2×M31)×1；步骤5模糊参数识别设主方向模糊长度为d，副方向模糊半径为r，对步骤4中得到的主方向导数自相关的累加ACa(x，y)和副方向导数的自相关的累加ACb(x，y)分别绘制曲线图，在主方向导数自相关的累加ACa(x，y)曲线图中将出现一个正的最高峰kap，在正最高峰两侧，出现负的左侧最低峰kap1和负的右侧最低峰kap2，下标ap，ap1和ap2分别表示主方向导数自相关的累加ACa(x，y)曲线图中正的最高峰，负的左侧最低峰和负的右侧最低峰，两负最低峰kap1和kap2之间的距离s1的值为2d+4r，在ACa(x，y)曲线图中将出现一个正的最高峰kbp，在正最高峰两侧，出现负的左侧最低峰kbp1和负的右侧最低峰kbp2，下标bp，bp1和bp2分别表示副方向导数自相关的累加ACb(x，y)曲线图中正的最高峰，负的左侧最低峰和负的右侧最低峰，两负最低峰kbp1和kbp2之间的距离s2的值为4r，计算主方向模糊长度d和副方向模糊半径r<mrow><mfencedopen='{'close=''><mtable><mtr><mtd><msub><mi>S</mi><mn>1</mn></msub><mo>=</mo><mn>2</mn><mi>d</mi><mo>+</mo><mn>4</mn><mi>r</mi></mtd></mtr><mtr><mtd><msub><mi>S</mi><mn>2</mn></msub><mo>=</mo><mn>4</mn><mi>r</mi></mtd></mtr></mtable></mfenced><mo>⇒</mo><mfencedopen='{'close=''><mtable><mtr><mtd><mi>d</mi><mo>=</mo><mrow><mo>(</mo><msub><mi>S</mi><mn>1</mn></msub><mo>-</mo><msub><mi>S</mi><mn>2</mn></msub><mo>)</mo></mrow><mo>/</mo><mn>2</mn></mtd></mtr><mtr><mtd><mi>r</mi><mo>=</mo><msub><mi>S</mi><mn>2</mn></msub><mo>/</mo><mn>4</mn></mtd></mtr></mtable></mfenced></mrow>步骤6建立复合模糊模型根据步骤5得到的主方向模糊长度d和副方向模糊半径r，通过对标准运动模糊模型hq(x，y)和标准散焦模型hp(x，y)进行线性运算和均一化处理来建立复合模糊模型h(x，y)，下标q和p分别表示运动模糊和散焦，建立复合模糊模型h(x，y)的方法如下步骤6.1对标准运动模糊模型hq(x，y)与标准散焦模型hp(x，y)进行卷积操作，得到复合模糊模型的理论模型hv(x，y)，下标v表示复合模糊模型的理论模型，具体方法如下标准散焦模型hp(x，y)的函数形式hpt(i1，i2)为其中i1，i2为横坐标和纵坐标，且都为实数，下标pt表示标准散焦模型的数学模型，r为副方向模糊半径且以r为散焦模糊半径，r为大于零的实数，标准散焦模型hp(x，y)是标准散焦模型的函数形式hpt(i1，i2)的矩阵形式，标准散焦模型hp(x，y)中的任意元素为大于或等于零的实数，假设运动模糊方向水平，标准运动模糊模型hp(x，y)的函数形式hpt(i1，i2)为其中i1，i2为横坐标和纵坐标，且都为实数，下标qt表示标准运动模糊模型的数学模型，d为主方向模糊长度且以d为运动模糊长度，d为大于零的实数，标准散焦模型hq(x，y)是标准运动模糊模型的函数模型hqt(i1，i2)的矩阵形式，标准运动模糊模型hq(x，y)中的任意元素为大于或等于零的实数，对标准散焦模型hp(x，y)和标准运动模糊模型hq(x，y)进行卷积，得到复合模糊模型的理论模型hv(x，y)，将复合模糊模型的理论模型hv(x，y)表示为函数形式hvt(i1，i2)其中下标vt表示复合模糊模型的理论模型的函数形式，步骤6.2复合模糊模型的理论模型hv(x，y)确定后，在能量守恒的前提下，对复合模糊模型的理论模型hv(x，y)进行均一化处理，得到最终的复合模糊模型h(x，y)，方法如下步骤6.2.1计算复合模糊模型的理论模型hv(x，y)中大于零的点的个数w，步骤6.2.2将复合模糊模型的理论模型hv(x，y)中大于零的点所对应的值替换为得到最终的复合模糊模型h(x，y)，将复合模糊模型h(x，y)表示成函数形式ht(i1，i2)其中下标t表示复合模糊模型的函数形式，步骤7图像复原根据步骤6得到的复合模糊模型h(x，y)，利用公知的维纳滤波进行复原，得到得到复原图像FDA0000029710520000021.tif,FDA0000029710520000033.tif,FDA0000029710520000042.tif,FDA0000029710520000043.tif,FDA0000029710520000051.tif,FDA0000029710520000052.tif,FDA0000029710520000053.tif,FDA0000029710520000054.tif全文摘要本发明提供一种运动模糊和散焦复合模糊的图像复原方法，该方法能够对运动模糊和散焦同时出现的图像进行参数估计和图像复原，它包含如下步骤(1)建立高斯白噪声模板，将降质图像与白噪声模板卷积，达到去除噪声的目的；(2)通过图像能量谱估计图像的模糊主方向和副方向(3)计算图像的主方向导数矩阵和副方向导数矩阵；(4)分别对主方向导数矩阵和副方向导数矩阵实施自相关运算和方向累加运算；(5)根据主方向导数自相关的累加曲线和副方向导数自相关的累加曲线，估计主方向模糊长度和副方向模糊长度；(6)根据获得的主方向模糊长度和副方向模糊长度，建立复合模糊模型；(7)利用维纳滤波对降质图像进行复原。文档编号G06T5/00GK101968881SQ201010522739公开日2011年2月9日申请日期2010年10月27日优先权日2010年10月27日发明者李楠,路小波申请人:东南大学