CN108564544A - 基于边缘感知的图像盲去模糊组合稀疏优化方法 - Google Patents

Abstract

本发明公开了一种基于边缘感知的图像盲去模糊组合稀疏优化方法，在图像L0稀疏先验的基础上引入相对总变分正则项，对自然图像进行盲去模糊，按照如下步骤进行：输入模糊图像y，其模糊核为k，待求解的清晰图像为x；初始化待求解的图像x为模糊图像y，初始化参数λ，σ等参数；利用模糊核和图像交叉估计的方法由粗到细求解模糊核k；根据步骤3)最后计算出来的模糊核k对模糊图像y进行非盲去模糊，求出清晰的图像x^(L)；对清晰的图像做最后的处理，得到最终的清晰图像。

Description

基于边缘感知的图像盲去模糊组合稀疏优化方法

技术领域：

本发明涉及一种基于边缘感知的图像盲去模糊组合稀疏优化方法，属于图像处理技术领域。

背景技术：

现代社会科学技术发展迅速，特别是计算机和多媒体技术的发展，人们认识到信息处理技术的重要性。图像作为信息传递的载体，随着摄像机、照相机在普通家庭的普及，能够普遍轻松地获得，因此如何简单高效地确保获得的图像高清变得越来越重要。然而在日常的拍摄过程中，我们会发现在获取、记录、处理、传输过程中由于成像设备、天气因素、处理和传输方式等原因，往往导致获得的图像是模糊的。直接利用降质图像信息可以估计清晰图像，但是由于成像时很多参数都是未知的，所以我们需要同时估计图像的模糊核和清晰的数字图像。这种在不知道模糊核信息的情况下估计清晰图像就是图像的盲恢复。在许多实际中，模糊核往往都是未知的，这时对退化的模糊图像进行有效的恢复是个具有挑战性的问题。

图像去模糊技术的理论依据是：根据图像退化的原因进行对图像进行模型化处理，将降质图像恢复或重建为原始图像，具体可以表示成：y＝k*x+n，其中这里y表示观测到的模糊图像(退化图像)，k代表模糊核(也可以称为模糊函数、点扩算函数或退化因子等)，*代表卷积运算，x则是需要复原的清晰图像(原始图像)，n为加性噪声(通常假设为高斯噪声)，所以图像去模糊问题从数学角度来说就是如何从模糊图像y得到出清晰图像x，这需要进行解卷积。盲去卷积问题是更具挑战性的病态逆问题。因为盲去模糊问题中模糊核是未知的，所以其算法一般采用对模糊核和清晰图像分开交替估计的步骤：1)首先对图像模糊核k进行估计，即利用初始恢复的清晰图像估计得到图像的模糊核；2)再对清晰图像x进行估计，利用先前估计得到的模糊核对图像进行非盲估计，得到去模糊的图像。这两个过程依次交替，最终得到清晰的图像。

纵观前期的去模糊工作，图像盲去模糊的成功都是尽可能高效地利用图像的先验知识和用于核估计的边缘检测，提出正则化项约束解空间，迭代求解清晰图像。传统的方法一般通过图像的梯度的L0范数作为先验知识实现盲去模糊，这种方法对于大多数的自然图像效果都是不错。最近，Pan等人通过结合图像亮度和图像梯度的L0范数先验实现文本图像的盲去模糊，图像先验源于观察文本图像的不同属性，基于这个先验，产生核估计的可靠中间结果。不需要检测突出边缘。在最后的图像恢复步骤，去除artifacts(伪像，人工，噪声)并且去模糊，进一步优化了自然图像去模糊效果。该方法在文本图像盲恢复问题上基本达到了目前已知的最优结果，实际上该方法对于自然图像、低照度场景的图像也可以获得了很好的去模糊效果。根据文本图像的统计结果可知，文本图像主要是双色的(梯度值有两个)，不符合自然图像的重尾分布的梯度统计，如果图像产生模糊，则会出现其他的梯度值，Pan的方法正是利用文本图像的这一先验知识进行正则化约束，得到清晰的图像结果。但是这一统计特性在自然图像上又不存在了，因为自然图像的亮度值丰富多样，无法利用文本图像的统计特性(只有两个值)作为约束条件。当然图像盲去模糊还有很多其他的方法，不过相对来说计算的复杂性更高。最近Pan的另一项工作引起了很大的关注：结合图像梯度L0范数和图像暗通道L0范数实现了鲁棒性更强、精确性更高的盲去模糊效果，这个工作给研究盲去模糊问题带来了一个新的思考方向：解决图像去模糊的问题是否可以借鉴到处理其他图像问题的方法，以及是否还存在更好的图像先验表示方法。

公开于该背景技术部分的信息仅仅旨在增加对本发明的总体背景的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域一般技术人员所公知的现有技术。

发明内容：

本发明的目的在于提供一种基于边缘感知的图像盲去模糊组合稀疏优化方法，提高算法对模糊核和清晰图像估计的结果,从而克服上述现有技术中的缺陷。

为实现上述目的，本发明提供的方案为：在图像L0稀疏先验的基础上引入相对总变分正则项，对自然图像进行盲去模糊，这个方法更好地保留了图像的边缘信息，从而获得了不错的结果。根据相对总变分的定义，在模糊图像中相对总变分项的值大于清晰图像，以这个作为盲去模糊时的先验对图像进行约束，具体是按照如下步骤进行：

S1.输入模糊图像y，其模糊核为k，待求解的清晰图像为x；

S2.初始化待求解的图像x为模糊图像y，初始化参数λ，σ等参数；

S3.利用模糊核和图像交叉估计的方法由粗到细求解模糊核k；

S4.根据步骤3)最后计算出来的模糊核k对模糊图像y进行非盲去模糊，求出清晰的图像x^(L)；

S5.对清晰的图像做最后的处理，得到最终的清晰图像。

本方面进一步限定的技术方案为：

优选地，上述技术方案中，步骤3具体为：

①对模糊核粗略估计为k⁽⁰⁾；初始化进行非盲去模糊的次数为l＝0；

②根据上面粗略估计的模糊核k,利用下面的公式1求解清晰的中间图像：

其中x^(l)表示的是第lλ次迭代求解后获得的待求解的清晰图像，表示目标函数取最小值x、u、g所取的值，u、g是人工引入的辅助变量，分别对应结构图像S和图像梯度x^(l-1)表示第l-1次迭代求解后的待求解的清晰图像，k^(l-1)表示第l-1次迭代求解后的模糊核，y表示模糊图像，*表示卷积运算，表示2范数的平方，表示图像x^(l-1)的梯度，||·||₀表示0范数，S^(l-1)表示x^(l-1)的结构图像，λ、σ是正则项系数，μ、β是引入变量的系数，初始化β＝2λσ、μ＝2λ，表示的是图像x^(l-1)的相对总变分正则项，其中D_h(p)、D_v(p)是窗口全变分，L_h(p)、L_v(p)是窗口固有变分，p是图像上的一个像素点，h、v分别代表的是水平和垂直方向，ε是一个非常小的正数，避免分母L_h(p)和L_v(p)为0；

③根据②估计出来的清晰图像x^(l)，可以直接用下面的公式2求解模糊核：

k^(l)表示的是第l次求解后获得的模糊核，表示目标函数取最小值k所取的值，x^(l)表示的第l次求解后的清晰图像，k^(l-1)表示第l-1次求解后的模糊核，y表示模糊图像，*表示卷积运算，表示2范数的平方，γ是正则项参数；

④更新λ的值，判断l是否小于最大的循环次数L，更新l＝l+1，重复②、③，得到最后的模糊核k＝k^(L)。

优选地，上述技术方案中，步骤3中构建的L0+RTV正则化模型；

基于图像梯度L0范数的显著边缘先验和相对总变分的正则化模型如下：

公式3，

令u＝S，当μ和β趋于无穷时，公式3变化为公式4：

令u，g的值趋近于零，则计算公式5：

通过下述公式6即可获得清晰图像x：

式中：F(·)和F^-1(·)代表的是傅里叶变换和傅里叶反变换，是傅里叶变换的复共轭形式， 分别代表水平和垂直方向的微分算子；

设定清晰图像x，分别通过公式7、8、9计算u、g的值：

则：

由上面的公式计算u和g，结合公式6即可得到第l次迭代时的清晰中间图像x^(l)。

优选地，上述技术方案中，加入的RTV正则项，计算结构图像的公式10：

其中，是目标函数取最小值S的取值，S是结构图像，I是输入的图像，是对输出和输入图像像素点差值求和，保证输出的结构图像准确性，η是正则化参数，正则项是相对总变分，D_x(p)、D_y(p)是窗口全变分，具体表示L_x(p)、L_y(p)是窗口固有变分，其与梯度方向没有关系，具体表示 g_p,q是根据空间相关性定义的权重函数，其中α控制窗口尺度，p是变分区域R(p)的中心像素点，q是变分区域R(p)中任意一点，ε、ε_S是极小的正数。

关于相对总变分项的具体计算如下公式11：

其中公式12和13分别为：

U_xp表示每个像素点组合了该点领域内的梯度信息，W_xp表示仅和该点的梯度有关；类似的，y方向上的表达形式如下公式14和15：

通过对相对总变分项的拆分，10式可以写成下面的矩阵形式16：

其中，C_x、C_y是利用向前查分来逼近离散梯度算子所得的托普利兹矩阵，v_S、v_I分别是S和I的向量表示法，U_x、U_y、W_x、W_y是对角矩阵，对角线上的值为U_x[i,i]＝u_xi，U_y[i,i]＝u_yi，W_x[i,i]＝w_xi，W_y[i,i]＝w_yi；利用欧拉—拉格朗日方程计算12式，并将极小化问题转化成线性问题公式17：

其中E是单位矩阵，是根据结构矢量计算的权重矩阵，(E+ηL^t)是个非负对称的拉普拉斯矩阵；

通过多次迭代计算12，13，14，15，17可以得到结构图像S。

与现有技术相比，本发明具有如下有益效果：

第一，本发明基于边缘感知，组合L0稀疏和相对总变分，更好地保留了图像的边缘细节信息，对模糊核估计更加精准，从而提高了估计的清晰图像效果。

第二，通过大量的实验发现，跟原方法比较，本发明在自然图像上去模糊效果不错，在Levin数据集的效果基本高于目前的已知的盲去模糊方法。

附图说明：

图1为本发明实现的流程图；

图2为待处理的模糊图像；

图3a是用Li Xu,Cewu Lu等人在文章“mage Smoothing via L0GradientMinimization”中公开的去模糊方法对图2的去模糊结果；

图3b是用Jinshan Pany,Zhe Hu等人在文章“Deblurring Text Images via L0-Regularized Intensity and Gradient Prior”中公开的去模糊方法对图2的去模糊结果；

图3c是用本发明方法对图2的去模糊结果；

图4a模糊的文本图像；

图4b是用Li Xu,Cewu Lu等人的方法对图4a的去模糊结果；

图4c是用Jinshan Pan,Zhe Hu等人的方法对图4a的去模糊结果；

图4d是用本发明方法对图4a的去模糊结果。

具体实施方式：

下面对本发明的具体实施方式进行详细描述，但应当理解本发明的保护范围并不受具体实施方式的限制。

除非另有其它明确表示，否则在整个说明书和权利要求书中，术语“包括”或其变换如“包含”或“包括有”等等将被理解为包括所陈述的元件或组成部分，而并未排除其它元件或其它组成部分。

具体实施方式

下面将结合附图对实验进行进一步详细说明。

参考图1，本发明的具体步骤如下：

步骤1：输入模糊图像：

输入模糊图像y，其模糊核为k，待求解的清晰图像为x；

步骤2：初始化参数：

初始化待求解的图像x为原始模糊图像y，初始化λ，σ等参数；

步骤3：计算模糊核：

利用模糊核和图像交叉估计的方法由粗到细求解模糊核：

①初始化模糊核，粗略估计为k⁽⁰⁾，初始化非盲去模糊地次数为l＝0；

②初始化β＝2λσ，初始化μ_max、β_max值，根据上面粗略估计出来的模糊核k，利用公式(8)计算辅助变量u的值。

其中结构图像S，具体的计算可以用下面的代码实现：

Step1：输入图像I，参数η、α；

Step2：初始化t＝0，S⁰＝I；

Step3：利用公式(12)，(13)，(14)，(15)计算权重U_x、U_y、W_x、W_y；

Step4：解线性方程(17)；

Step5：按要求多次迭代step3，step4；

Step6：输出结构图像S。

③利用公式(9)计算辅助变量g的值

④利用公式(6)计算清晰图像

x^(l)表示的是第lλ次循环求解后获得的待求解的清晰图像，u、g是辅助变量，分别对应S和k^(l-1)表示第l-1次循环求解后的模糊核，y表示模糊图像，F(·)和F^-1(·)代表的是傅里叶变换和傅里叶反变换，是傅里叶变换的复共轭形式， 分别代表水平和垂直方向的微分算子。

⑤根据交叉估计的思想，利用③估计出来的清晰图像x^(l)，求解模糊核：

k^(l)表示的是第l次求解后获得的模糊核，表示目标函数取最小值k所取的值，x^(l)表示的第l次求解后的清晰图像，k^(l-1)表示第l-1次求解后的模糊核，y表示模糊图像，*表示卷积运算，表示2范数的平方，γ是正则项参数；

⑥更新λ的值，判断l是否大于最大的循环次数L，更新l＝l+1，更新图像的结构S，重复上面交叉估计的步骤，得到最后的模糊核k＝k^(L)；

步骤4：计算清晰图像

根据步骤3)中最后迭代计算出来的模糊核k，对原始模糊图像y进行非盲去模糊，计算清晰的图像x^(L)；

步骤5：得到清晰图像

对步骤4得到的图像做一些必要的处理，例如去除人工伪像等工作，得到

最后的清晰图像。

下面结合实验条件对本发明的去模糊效果进行进一步的描述。

1.实验条件

本发明的实验运行***是Inter(R)Core(TM)[email protected],64位Windows操作***，使用的仿真软件为MATLAB(r2014a)。实验使用的模糊图像来源是A.Levin,Y.Weiss,F.Durand,W.T.Freeman等人发表的文章“Understanding and evaluating blinddeconvolution algorithms”中模糊的自然图像数据库，如图2所示。

2.实验内容

用本发明和现有的三种去模糊方法对图2进行盲去模糊处理，结果如图3所示。

图3a是用Li Xu,Cewu Lu等人在文章“mage Smoothing via L0 GradientMinimization”中公开的去模糊方法对图2的去模糊结果；

图3b是用Jinshan Pany,Zhe Hu等人在文章“Deblurring Text Images via L0-Regularized Intensity and Gradient Prior”中公开的去模糊方法对图2的去模糊结果；

图3c是用本发明方法对图2的去模糊结果；

3.仿真结果分析

仿真实验中，采用峰值信噪比PSNR指标来评价实验结果，PSNR分别定义为：

其中：x代表原始清晰图像，Z代表恢复出来的去模糊图像。峰值信噪比越大说明去模糊性能越好。

图3a的PSNR值为32.13dB；

图3b的PSNR值为31.66dB；

图3c的PSNR值为33.08dB。该部分数据是三种方法对图2进行实验的结果，图2是Levin数据集里面的具体的一张模糊图片，Levin数据集里面包含八个不同的模糊核，表格的实验结果是分别对原清晰图像对应八个模糊核的模糊图片，进行盲去模糊后的结果。

为了更好地比较这三种方法的去模糊性能，我们对Levin数据集上一组标准模糊图片进行实验，结果如表1所示。

表1为模糊图片实验结果：

由所得的数据可以知道，本发明方法的去模糊结果的PSNR值高于其他方法的去模糊结果的PSNR值，即本发明比现有技术具有更好地去模糊效果。

在实验过程中，我们发现本发明方法对文本图像去模糊也能取得很好的结果，结果如图4d所示。模糊文本图像来源是Jinshan Pan,Zhe Hu等人发表的文章“DeblurringText Images via L0-Regularized Intensity and Gradient Prior”中的数据库。

图4a模糊的文本图像；

图4b是用Li Xu,Cewu Lu等人的方法对图4(a)的去模糊结果；

图4c是用Jinshan Pan,Zhe Hu等人的方法对图4(a)的去模糊结果；

图4d是用本发明方法对图4a的去模糊结果。

前述对本发明的具体示例性实施方案的描述是为了说明和例证的目的。这些描述并非想将本发明限定为所公开的精确形式，并且很显然，根据上述教导，可以进行很多改变和变化。对示例性实施例进行选择和描述的目的在于解释本发明的特定原理及其实际应用，从而使得本领域的技术人员能够实现并利用本发明的各种不同的示例性实施方案以及各种不同的选择和改变。本发明的范围意在由权利要求书及其等同形式所限定。

Claims (4)

1.基于边缘感知的图像盲去模糊组合稀疏优化方法，其特征在于，在图像L0稀疏先验的基础上引入相对总变分正则项，对自然图像进行盲去模糊，按照如下步骤进行：

S1.输入模糊图像y，其模糊核为k，待求解的清晰图像为x；

S2.初始化待求解的图像x为模糊图像y，初始化参数λ，σ等参数；

S3.利用模糊核和图像交叉估计的方法由粗到细求解模糊核k；

S4.根据步骤3最后计算出来的模糊核k对模糊图像y进行非盲去模糊，求出清晰的图像x^(L)；

S5.对清晰的图像做最后的处理，得到最终的清晰图像。

2.根据权利要求1所述的基于边缘感知的图像盲去模糊组合稀疏优化方法，其特征在于，步骤3具体为：

①对模糊核粗略估计为k⁽⁰⁾；初始化进行非盲去模糊的次数为l＝0；

②根据上面粗略估计的模糊核k,利用下面的公式1求解清晰的中间图像：

其中x^(l)表示的是第λ次迭代求解后获得的待求解的清晰图像，表示目标函数取最小值x、u、g所取的值，u、g是人工引入的辅助变量，分别对应结构图像S和图像梯度x^(l-1)表示第l-1次迭代求解后的待求解的清晰图像，k^(l-1)表示第l-1次迭代求解后的模糊核，y表示模糊图像，*表示卷积运算，表示2范数的平方，表示图像x^(l-1)的梯度，||.||₀表示0范数，S^(l-1)表示x^(l-1)的结构图像，λ、σ是正则项系数，μ、β是引入变量的系数，初始化β＝2λσ、μ＝2λ，表示的是图像x^(l-1)的相对总变分正则项，其中D_h(p)、D_v(p)是窗口全变分，L_h(p)、L_v(p)是窗口固有变分，p是图像上的一个像素点，h、v分别代表的是水平和垂直方向，ε是一个非常小的正数，避免分母L_h(p)和L_v(p)为0；

③根据②估计出来的清晰图像x^(l)，可以直接用下面的公式2求解模糊核：

k^(l)表示的是第l次求解后获得的模糊核，表示目标函数取最小值k所取的值，x^(l)表示的第l次求解后的清晰图像，k^(l-1)表示第l-1次求解后的模糊核，y表示模糊图像，*表示卷积运算，表示2范数的平方，γ是正则项参数；

④更新λ的值，判断l是否小于最大的循环次数L，更新l＝l+1，重复②、③，得到最后的模糊核k＝k_(L)。l＝0。

3.根据权利要求1所述的基于边缘感知的图像盲去模糊组合稀疏优化方法，其特征在于，步骤3中构建的L0+RTV正则化模型；

基于图像梯度L0范数的显著边缘先验和相对总变分的正则化模型如下：

公式3，

令u＝S，当μ和β趋于无穷时，公式3变化为公式4：

令u，g的值趋近于零，则计算公式5：

通过下述公式6即可获得清晰图像x：

式中：F(·)和F^-1(·)代表的是傅里叶变换和傅里叶反变换，是傅里叶变换的复共轭形式， 分别代表水平和垂直方向的微分算子；

设定清晰图像x，分别通过公式7、8、9计算u、g的值：

则：

由上面的公式计算u和g，结合公式6即可得到第l次迭代时的清晰中间图像x^(l)。

4.根据权利要求1所述的基于边缘感知的图像盲去模糊组合稀疏优化方法，其特征在于，加入的RTV正则项，计算结构图像的公式10：

其中，是目标函数取最小值S的取值，S是结构图像，I是输入的图像，是对输出和输入图像像素点差值求和，保证输出的结构图像准确性，η是正则化参数，正则项是相对总变分，D_x(p)、D_y(p)是窗口全变分，具体表示L_x(p)、L_y(p)是窗口固有变分，其与梯度方向没有关系，具体表示 g_p,q是根据空间相关性定义的权重函数，其中α控制窗口尺度，p是变分区域R(p)的中心像素点，q是变分区域R(p)中任意一点，ε、ε_S是极小的正数。

关于相对总变分项的具体计算如下公式11：

其中公式12和13分别为：

U_xp表示每个像素点组合了该点领域内的梯度信息，W_xp表示仅和该点的梯度有关；类似的，y方向上的表达形式如下公式14和15：

通过对相对总变分项的拆分，(10)式可以写成下面的矩阵形式16：

其中，C_x、C_y是利用向前查分来逼近离散梯度算子所得的托普利兹矩阵，v_S、v_I分别是S和I的向量表示法，U_x、U_y、W_x、W_y是对角矩阵，对角线上的值为U_x[i,i]＝u_xi，U_y[i,i]＝u_yi，W_x[i,i]＝w_xi，W_y[i,i]＝w_yi；利用欧拉—拉格朗日方程计算12式，并将极小化问题转化成线性问题公式17：

其中E是单位矩阵，是根据结构矢量计算的权重矩阵，(E+ηL^t)是个非负对称的拉普拉斯矩阵；

通过多次迭代计算12，13，14，15，17可以得到结构图像S。