CN104021529B - 一种模糊图像非盲复原方法 - Google Patents

Abstract

本发明公开了一种基于泊松概率模型和高斯尺度混合型马尔科夫专家场的模糊图像非盲复原方法，其特征在于其实施步骤如下：1)用泊松概率模型对噪声进行建模；2)用高斯尺度混合型马尔科夫专家场对复原图像进行建模；3)将上述两个模型相乘，得到模糊图像非盲复原的贝叶斯后验概率模型，经过负自然对数运算，将其转化为最大后验估计问题；4)采用嵌套式二次惩罚函数法对步骤3)中的最大后验估计问题求解。本发明能够对模糊图像进行有效复原，获得高质量的复原图像，优于传统方法。

Description

一种模糊图像非盲复原方法

技术领域

本发明涉及计算机图像处理技术，特别涉及一种基于泊松概率模型和高斯尺度混合型马尔科夫专家场的图像去模糊方法。

背景技术

在日常摄影中，经常会遇到拍摄出的图像模糊不清，细节难以分辨的情况。这主要是由于在拍摄过程中照明条件较差，因此需要提高相机的ISO或延长曝光时间以使成像器件获得充足的曝光量。然而提高相机的ISO会使所得图像中包含大量噪声，影响图像的感观效果，而延长曝光时间，则容易造成相机和所拍摄景物之间的相对运动，产生图像模糊。

解决图像模糊的方法多种多样，最常见的是采用稳像设备，如三脚架或镜头稳像器等，但是稳像设备通常较笨重或价格昂贵。除此之外，图像复原的另一条途径是采用数学方法，又称为图像反卷积方法，它能够避免稳像设备的诸多缺点，具有很高的应用价值。

理论上，图像模糊可以表述为清晰图像与模糊核的卷积，同时由于电子器件等外界因素的影响，所得模糊图像中会或多或少的包含噪声。图像复原是图像模糊的逆过程，从模糊图像的形成原理可知，模糊核是图像复原的一个关键因素，在某些情况下，可通过刃边法等测量手段获得近似的模糊核，此时的图像复原称为非盲复原，其数学形式较为简单，但它却是一个典型的病态问题，若采用直接复原方法，即使很少量的噪声，也会使复原结果中出现大量的噪声和振铃，大幅偏离真实情况。解决病态问题的方法称为正则化方法，其原理是对原问题进行修正，增加一些约束条件，使修正后的问题的解尽可能的接近真实情况。

实施正则化方法的途径有很多种，通常的方法是在贝叶斯后验概率框架下对图像复原问题进行重新建模，并将其转化为最大后验估计问题。问题修正过程中，需要对噪声发生的概率分布和复原图像的先验概率分布分别进行建模，常用的噪声概率模型有高斯概率模型和泊松概率模型，由于高斯噪声是二次的，便于优化求解，因此其应用范围远超过泊松概率模型，而本专利则采用了泊松概率模型对噪声进行建模，目的是证明采用该模型，通过适当的建模和优化过程，同样能够获得高质量的复原效果。复原图像的先验概率模型则相当于正则项，即约束条件，其选择的优劣直接关系到复原效果的好坏，常用的先验概率模型有高斯模型、稀疏模型等等，本专利选择了一种高斯尺度混合型马尔科夫专家场，该模型是采用自然图像数据库，利用特定的优化算法训练得到的，较传统模型相比，能够更加准确的描述图像的概率分布特征。

最大后验问题的优化方法同样对图像复原效果有着至关重要的影响。目前，基于噪声的高斯概率模型的图像复原问题的求解方法已经较为成熟，而基于泊松概率模型的图像复原问题的求解还有待研究，虽然常用的优化迭代方法如RL算法、One-Step-Late方法可以用于求解这类问题，但优化过程的收敛性难以保证，特别是对于复杂的正则项，其鲁棒性较差，很难得到令人满意的复原效果。

发明内容

本发明解决的技术问题是提供一种能够在已知模糊图像和模糊核的情况下，去除图像模糊，提高图像对比度和清晰度，同时增加图像细节，复原效果好的模糊图像复原算法。

为了解决上述技术问题，本发明公开了一种基于泊松概率模型和高斯尺度混合型马尔科夫专家场的图像去模糊方法，其实施步骤如下：

1)用泊松概率模型对模糊图像中由成像器件引入的噪声进行建模；

2)用高斯尺度混合型马尔科夫专家场对复原图像进行建模；

3)将上述两个模型相乘，得到模糊图像非盲复原的贝叶斯后验概率模型，经过负自然对数运算，将其转化为最大后验估计问题；

4)采用嵌套式二次惩罚函数法对步骤3)中的最大后验估计问题求解。

所述步骤1)中的泊松概率模型的表达式为：

$P\left(g\right|f)=\frac{\exp [-{\mathrm{\Σ}}_j{(h\⊗f)}_j]{\mathrm{\Π}}_j{(h\⊗f)}_j^{g_j}}{{\mathrm{\Π}}_j\left(g_j\right)!}$

其中，h表示已知模糊核，表示卷积运算，！表示阶乘运算，表示h和f卷积所得的图像，j表示像素索引，表示h和f卷积所得的图像的第j个像素，g_j表示g的第j个像素，表示一维逐项求和运算，表示一维逐项求积运算，表示的g_j次幂，(g_j)！表示g_j的阶乘。

所述步骤2)中的高斯尺度混合型马尔科夫专家场的表达式为：

$P\left(f\right)\∝\underset{i,j}{\mathrm{\Π}}\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right],$

其中，表示二维逐项求积运算，ψ表示构成P(f)的函数，

$\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right]\∝\underset{m}{\mathrm{\Σ}}\frac{{\mathrm{\π}}_m}{{\mathrm{\σ}}_m}\exp \left(\frac{-{{(w_i\⊗f)}_j}^2}{{2\mathrm{\σ}}_m^2}\right),$

m表示构成函数ψ的指数函数的索引，构成ψ的指数函数的总数M的值为8，每一组加权系数π_m和标准差σ_m均为固定常数，并且

π_1～8＝{0.1940,0.0906,0.4631,0.0423,0.1021,0.0816,0.0028,0.0234}，

σ_1～8＝{0.0032,0.0147,0.0215,0.0316,0.0464,0.0681,0.100,0.1468}；

每一个w_i表示一个高通滤波器，i表示构成P(f)的高通滤波器索引，高通滤波器总数是25。

所述步骤3)中的模糊图像非盲复原的贝叶斯后验概率模型的表达式为：

P(f|g)∝P(g|f)P(f)，

其中，P(f|g)表示在g发生的情况下f发生的概率，P(g|f)表示噪声发生的概率，P(f)表示复原图像发生的先验概率。

所述步骤3)中的最大后验估计问题的表达式为：

$f=\arg {\min }_f\left\{\mathrm{\λ}\underset{j}{\mathrm{\Σ}}\right[{(h\⊗f)}_j-g_j\ln {(h\⊗f)}_j]-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right]\left\}\right\}$

其中，λ为正则化系数，表示二维逐项求积运算，表示的自然对数，表示的自然对数。

所述步骤4)中的嵌套式二次惩罚函数法的详细步骤如下：

a)引入辅助变量u和惩罚系数β₁，将u赋值为g，将β₁初始化为1，将所述步骤3)中的问题转化为：

$f=\arg \mathrm{mi}{n}_f\left\{\mathrm{\λ}\underset{j}{\mathrm{\Σ}}\right[u_j-g_j\ln u_j]+\frac{{\mathrm{\β}}_1}{2}\underset{j}{\mathrm{\Σ}}{[u_j-{(h\⊗f)}_j]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right]\left\}\right\}$

b)固定f，对u进行优化，u的第j个像素u_j的优化结果是下式所示的二次方程的非负根：

${\mathrm{\β}}_1{u}_j^2+[\mathrm{\λ}-{\mathrm{\β}}_1{(h\⊗f)}_j]u_j-{\mathrm{\λg}}_j=0$

对u进行逐像素优化，得到其优化结果。

c)当得到所有u_j的估计值后，转而求解f，此时问题转化为：

$f=\arg \mathrm{mi}{n}_f\{\frac{{\mathrm{\β}}_1}{2}\underset{j}{\mathrm{\Σ}}{[u_j-{(h\⊗f)}_j]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right]\left\}\right\}$

求解该问题的步骤如下：

c-1)引入与P(f)表达式中滤波器相同数量并一一对应的辅助变量v_i，i＝1,2,...,25和一个惩罚系数β₂，将β₂初始化为1，将步骤c)中的问题转化为：

$f=\arg {\min }_f\{\frac{\mathrm{\λ}}{2}\underset{j}{\mathrm{\Σ}}{[u_j-{(h\⊗f)}_j]}^2+\frac{{\mathrm{\β}}_2}{2}\underset{i,j}{\mathrm{\Σ}}{[{\left(v_i\right)}_j-{(w_i\⊗f)}_j]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{\left(v_i\right)}_j\right]\left\}\right\}$

c-2)固定f，并采用牛顿—拉弗森迭代方法对所有变量v_i进行逐一优化求解；(参见：韩丹夫，吴庆标，数值计算方法，浙江大学出版社，2006年第一版，p.145-151)

c-3)当得到所有v_i的估计值后，转而求解f，采用在频域中得到f的解析解，经傅里叶逆变换得到f，然后为β₂乘以放大因子R，R＝2；

c-4)判断β₂是否小于β_max，β_max＝2²⁰，若β₂＜β_max，则转而执行步骤c-2)，否则执行步骤c-5)；

c-5)为β₁乘以放大因子R。

d)判断β₁是否小于β_max，若β₁＜β_max，则转而执行步骤b)，否则执行步骤e)；

e)输出f作为最终复原结果。

本发明具有以下优点：本发明基于模糊图像非盲复原的贝叶斯后验概率框架，分别采用了泊松概率模型和一种高斯尺度混合型马尔科夫专家场对噪声和复原图像进行建模，并通过负自然对数运算将其转化为最大后验估计问题。在对所得的最大后验估计问题求解的过程中，采用了一种嵌套式二次惩罚函数法求解，能够获得清晰度高、细节丰富的复原图像，能够保证收敛于理想的复原结果。

附图说明

下面结合附图和具体实施方式对本发明做更进一步的具体说明，本发明的上述和/或其他方面的优点将会变得更加清楚。

图1为本发明实施例的总流程图。

图2为本发明实施例步骤4)中的嵌套式二次惩罚函数优化方法的流程图。

图3为本发明实施例的模糊图像。

图4为本发明实施例的已知模糊核。

图5为本发明实施例的利用经典的RL算法所得的复原效果。

图6为本发明实施例的复原效果图。

具体实施方式

如图1所示，本实施例基于泊松概率模型和高斯尺度混合型马尔科夫专家场的模糊图像非盲复原方法的实施步骤如下：

1)用泊松概率模型对噪声进行建模；

2)用高斯尺度混合型马尔科夫专家场对复原图像进行建模；

3)将上述两个模型相乘，得到模糊图像非盲复原的贝叶斯后验概率模型，经过负自然对数运算，将其转化为最大后验估计问题；

4)采用嵌套式二次惩罚函数法对步骤3)中的最大后验估计问题求解。

如图3所示，是一幅典型的模糊图像，其形成过程可以表述为清晰图像与模糊核的卷积，同时包含了大量的泊松噪声，其表达式如下：

$g=h\⊗f+n$

其中，g表示已知模糊图像，h表示已知模糊核，f表示清晰图像，表示理想的模糊图像，n表示模糊图像中由成像器件引入的噪声，在本实施例中假设噪声符合泊松概率模型。

非盲图像复原是图像模糊的逆过程，它是一个病态问题，由于噪声n的影响，经典的直接复原方法，如Richardson-Lucy(RL)算法会在复原的结果中引入大量的噪声和振铃等负面效应。因此，需要对采取一些修正措施，使问题的解趋于正确。

解决病态问题的方法称为正则化，其原理是为原问题引入一定的限制条件，对原问题进行修正，使新问题的解尽可能的接近真实解。经典的方法是在贝叶斯后验概率框架下对非盲复原问题进行重新建模，其表示式如步骤1)和步骤2)中所示：

P(f|g)∝P(g|f)P(f)

其中，f和g分别表示复原图像和已知模糊图像，P(f|g)表示在g发生的情况下f发生的概率，P(g|f)表示噪声发生的概率，P(f)表示复原图像发生的先验概率。

在对P(g|f)所表示的噪声进行建模时，本实施例选用了如步骤1)所述的泊松概率模型，表达式如下：

$P\left(g\right|f)=\frac{\exp [-{\mathrm{\Σ}}_j{(h\⊗f)}_j]{\mathrm{\Π}}_j{(h\⊗f)}_j^{g_j}}{{\mathrm{\Π}}_j\left(g_j\right)!}$

其中，h表示已知模糊核，表示卷积运算，！表示阶乘运算，表示h和f卷积所得的图像，j表示像素索引，表示h和f卷积所得的图像的第j个像素，g_j表示g的第j个像素，表示一维逐项求和运算，表示一维逐项求积运算，表示的g_j次幂，(g_j)！表示g_j的阶乘。

在对P(f)进行建模时，本实施例选用了如步骤2)所述的高斯尺度混合型马尔科夫专家场，表达式如下：

$P\left(f\right)\∝\underset{i,j}{\mathrm{\Π}}\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right],$

其中，表示二维逐项求积运算，ψ表示构成P(f)的函数，

$\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right]\∝\underset{m}{\mathrm{\Σ}}\frac{{\mathrm{\π}}_m}{{\mathrm{\σ}}_m}\exp \left(\frac{-{{(w_i\⊗f)}_j}^2}{{2\mathrm{\σ}}_m^2}\right),$

m表示构成函数ψ的指数函数的索引，构成ψ的指数函数的总数M的值为8，每一组加权系数π_m和标准差σ_m均为固定常数，并且

π_1～8＝{0.1940,0.0906,0.4631,0.0423,0.1021,0.0816,0.0028,0.0234}，

σ_1～8＝{0.0032,0.0147,0.0215,0.0316,0.0464,0.0681,0.100,0.1468}；

每一个w_i表示一个高通滤波器，i表示构成P(f)的高通滤波器索引，高通滤波器总数是25。

如本实施例步骤3)所述，将上述两个模型相乘后并作负自然对数运算，得到所对应的最大后验估计问题，如下式所示：

$f=\arg {\min }_f\left\{\mathrm{\λ}\underset{j}{\mathrm{\Σ}}\right[{(h\⊗f)}_j-g_j\ln {(h\⊗f)}_j]-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right]\left\}\right\}$

其中，λ为正则化系数，表示二维逐项求积运算。

这是一个非线性优化问题，传统方法很难得到令人满意的结果，本实施例采用了一种嵌套式二次惩罚函数法对其求解，图2是其流程图，详细求解步骤如下：

a)引入辅助变量u和惩罚系数β₁，将u赋值为g，将β₁初始化为1，将所述步骤3)中的问题转化为：

$f=\arg \mathrm{mi}{n}_f\left\{\mathrm{\λ}\underset{j}{\mathrm{\Σ}}\right[u_j-g_j\ln u_j]+\frac{{\mathrm{\β}}_1}{2}\underset{j}{\mathrm{\Σ}}{[u_j-{(h\⊗f)}_j]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right]\left\}\right\}$

b)固定f，对u进行优化，u的第j个像素u_j的优化结果是下式所示的二次方程的非负根：

$u=\arg {\min }_u\left\{\mathrm{\λ}\underset{j}{\mathrm{\Σ}}\right[u_j-g_j\ln u_j]+\frac{{\mathrm{\β}}_1}{2}\underset{j}{\mathrm{\Σ}}{[u_j-{(h\⊗f)}_j]}^2\}$

这是一个凸优化问题，对其求倒数并令结果为零，可知其解是下式所示的二次方程的非负根：

${\mathrm{\β}}_1{u}_j^2+[\mathrm{\λ}-{\mathrm{\β}}_1{(h\⊗f)}_j]u_j-{\mathrm{\λg}}_j=0$

对u进行逐像素优化，得到其优化结果。

c)固定优化所得的u，转而求解f，此时问题转化为：

$f=\arg \mathrm{mi}{n}_f\{\frac{{\mathrm{\β}}_1}{2}\underset{j}{\mathrm{\Σ}}{[u_j-{(h\⊗f)}_j]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{(w_i\⊗f)}_j\right]\left\}\right\}$

求解该问题的步骤如下：

c-1)引入与P(f)表达式中滤波器相同数量并一一对应的辅助变量v_i，i＝1,2,...,25和一个惩罚系数β₂，将β₂初始化为1，将步骤c)中的问题转化为：

$f=\arg {\min }_f\{\frac{\mathrm{\λ}}{2}\underset{j}{\mathrm{\Σ}}{[u_j-{(h\⊗f)}_j]}^2+\frac{{\mathrm{\β}}_2}{2}\underset{i,j}{\mathrm{\Σ}}{[{\left(v_i\right)}_j-{(w_i\⊗f)}_j]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{\left(v_i\right)}_j\right]\left\}\right\}$

c-2)固定f，并对所有变量v_i进行逐一优化求解，此时的优化问题为：

$v_i=\arg {\min }_{v_i}\{\frac{{\mathrm{\β}}_2}{2}\underset{i,j}{\mathrm{\Σ}}{[{\left(v_i\right)}_j-{(w_i\⊗f)}_j]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\left[{\left(v_i\right)}_j\right]\left\}\right\}$

可采用牛顿—拉弗森迭代方法对v_i求解。

c-3)当得到所有v_i的估计值后，转而求解f，采用在频域中得到f的解析解，如下式所示：

$F\left(u\right)=\frac{{\mathrm{\λH}}^{*}\left(u\right)G\left(u\right)+\mathrm{\β}{\mathrm{\Σ}}_i{W}_i\left(u\right)V_i\left(u\right)}{{\mathrm{\λH}}^{*}\left(u\right)H\left(u\right)+\mathrm{\β}{\mathrm{\Σ}}_i{W}_i^{*}\left(u\right)V_i\left(u\right)}$

其中，大写字母表示对应变量的傅立叶变换，u表示某一频率成份，上标*表示复共轭。F(u)经傅里叶逆变换得到f，然后为β₂乘以放大因子R，R＝2；

c-4)判断β₂是否小于β_max，β_max＝2²⁰，若β₂＜β_max，则转而执行步骤c-2)，否则执行步骤c-5)；

c-5)为β₁乘以放大因子R。

d)判断β₁是否小于β_max，若β₁＜β_max，则转而执行步骤b)，否则执行步骤e)；

e)输出f作为最终复原结果。

如图4所示为本实施例的已知模糊核，图5为本实施例采用经典的RL算法复原的效果，可见其中包含了大量的噪声，复原效果严重偏离了真实情况。图6为本实施例采用嵌套式二次惩罚函数方法的求解结果，可见复原效果较清晰，含有丰富的细节，噪声和振铃等负面效应得到有效抑制，有效提高了图像质量，明显优于图5。

本发明提供了一种模糊图像非盲复原方法，具体实现该技术方案的方法和途径很多，以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。本实施例中未明确的各组成部分均可用现有技术加以实现。

Claims (4)

1.一种基于泊松概率模型和高斯尺度混合型马尔科夫专家场的模糊图像非盲复原方法，其特征在于，包括如下步骤：

1)用泊松概率模型对模糊图像中的噪声进行建模；

2)用高斯尺度混合型马尔科夫专家场对复原图像进行建模；

3)将上述步骤1)和步骤2)两个模型相乘，得到模糊图像非盲复原的贝叶斯后验概率模型，经过负自然对数运算，将其转化为最大后验估计问题；

4)采用嵌套式二次惩罚函数法对步骤3)中的最大后验估计问题求解；

步骤1)中所述噪声为成像器件引入的噪声；

所述步骤1)中的泊松概率模型对噪声进行建模的表达式为：

$P\left(g\right|f)=\frac{\exp \[-{\mathrm{\Σ}}_j{(h\⊗f)}_j\]{\mathrm{\Π}}_j{(h\⊗f)}_j^{g_j}}{{\mathrm{\Π}}_j\left(g_j\right)!},$

其中，P(f|g)表示在g发生的情况下f发生的概率，h表示已知模糊核，f和g分别表示复原图像和已知模糊图像，表示卷积运算，！表示阶乘运算，表示h和f卷积所得的图像，j表示像素索引，表示一维逐项求和运算，表示一维逐项求积运算，表示h和f卷积所得的图像的第j个像素，g_j表示g的第j个像素，表示的g_j次幂，(g_j)！表示g_j的阶乘；

所述步骤2)中的高斯尺度混合型马尔科夫专家场对复原图像进行建模的表达式为：

$P\left(f\right)\∝\underset{i,j}{\Π}\mathrm{\ψ}\[{(w_i\⊗f)}_j\],$

其中，表示二维逐项求积运算，ψ表示构成P(f)的函数，

$\mathrm{\ψ}\[{(w_i\⊗f)}_j\]\∝\underset{m}{\Σ}\frac{{\mathrm{\π}}_m}{{\mathrm{\σ}}_m}\exp \left(\frac{-{{(w_i\⊗f)}_j}^2}{2{\mathrm{\σ}}_m^2}\right),$

m表示构成函数ψ的指数函数的索引，m取值范围1～8，构成ψ的指数函数的总数M的值为8，每一组加权系数π_m和标准差σ_m均为固定常数，并且

π_1～8＝{0.1940,0.0906,0.4631,0.0423,0.1021,0.0816,0.0028,0.0234}，

σ_1～8＝{0.0032,0.0147,0.0215,0.0316,0.0464,0.0681,0.100,0.1468}；

每一个w_i表示一个高通滤波器，i表示构成P(f)的高通滤波器索引，i取值范围1～25，高通滤波器总数是25。

2.根据权利要求1所述的模糊图像非盲复原方法，其特征在于，所述步骤3)中的模糊图像非盲复原的贝叶斯后验概率模型的表达式为：

P(f|g)∝P(g|f)P(f)，

其中，P(f|g)表示在g发生的情况下f发生的概率，P(g|f)表示噪声发生的概率，P(f)表示复原图像发生的先验概率。

3.根据权利要求2所述的模糊图像非盲复原方法，其特征在于，所述步骤3)中的最大后验估计问题的表达式为：

$f=\arg {\min }_f\{\mathrm{\λ}\underset{j}{\mathrm{\Σ}}\[{(h\⊗f)}_j-g_j\ln {(h\⊗f)}_j\]-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\[{(w_i\⊗f)}_j\]\left\}\right\},$

其中，λ为正则化系数，1000≤λ≤10000，表示二维逐项求积运算，表示的自然对数，表示的自然对数。

4.根据权利要求2所述的模糊图像非盲复原方法，其特征在于，所述步骤4)中的嵌套式二次惩罚函数法的包括如下步骤：

a)引入辅助变量u和惩罚系数β₁，将辅助变量u赋值为g，将惩罚系数β₁初始化为1，将步骤3)中的问题转化为：

$f=\arg {\min }_f\{\mathrm{\λ}\underset{j}{\mathrm{\Σ}}\[u_j-g_j\ln u_j\]+\frac{{\mathrm{\β}}_1}{2}\underset{j}{\mathrm{\Σ}}{\[u_j-{(h\⊗f)}_j\]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\[{(w_i\⊗f)}_j\]\left\}\right\},$

b)固定f，对u进行优化，u的第j个像素u_j的优化结果是下式所示的二次方程的非负根：

${\mathrm{\β}}_1{u}_j^2+\[\mathrm{\λ}-{\mathrm{\β}}_1{(h\⊗f)}_j\]u_j-{\mathrm{\λg}}_j=0,$

对u进行逐像素优化，得到其优化结果；

c)固定所得的u的优化结果，对f进行优化，此时问题转化为：

$f=\arg {\min }_f\{\frac{{\mathrm{\β}}_1}{2}\underset{j}{\Σ}{\[u_j-{(h\⊗f)}_j\]}^2-\underset{i,j}{\Σ}ln\{\mathrm{\ψ}\[{(w_i\⊗f)}_j\]\left\}\right\},$

求解该问题的步骤如下：

c-1)引入与P(f)表达式中滤波器相同数量并一一对应的辅助变量v_i，i＝1,2,…,25和一个惩罚系数β₂，将β₂初始化为1，将步骤c)中的问题转化为：

$f=\arg {\min }_f\{\frac{\mathrm{\λ}}{2}\underset{j}{\mathrm{\Σ}}{\[u_j-{(h\⊗f)}_j\]}^2+\frac{{\mathrm{\β}}_2}{2}\underset{i,j}{\mathrm{\Σ}}{\[{\left(v_i\right)}_j-{(w_i\⊗f)}_j\]}^2-\underset{i,j}{\mathrm{\Σ}}\ln \{\mathrm{\ψ}\[{\left(v_i\right)}_j\]\left\}\right\},$

c-2)固定f，并采用牛顿—拉弗森迭代方法对所有变量v_i进行逐一优化求解；

c-3)当得到所有v_i的优化结果后，求解f，采用在频域中得到f的解析解，经傅里叶逆变换得到f，然后将β₂乘以放大因子R，R＝2；

c-4)判断β₂是否小于β_max，β_max＝2²⁰，若β₂＜β_max，则转而执行步骤c-2)，否则执行步骤c-5)；

c-5)将β₁乘以放大因子R；

d)判断β₁是否小于β_max，若β₁＜β_max，则转而执行步骤b)，否则执行步骤e)；

e)输出f作为最终复原结果。