有约束最小二乘图像复原算法设计与实现

本 科 毕 业 论 文（设计）
课题名称
有约束最小二乘图像复原算法设计与实现 学 院
机械与电气工程学院 专 业
电气工程及其自动化 班级名称
XXXX 学生姓名
XXXX 学 号
XXXX 指导教师
XXXX
完成日期 XXXX
XXXXXXXX
有约束最小二乘图像复原算法设计与实现
摘要正则化图像复原方法是通过引入图像先验知识相关的正则项，将不适定问题转化为近似的适定问题，从而获得稳定的近似解的过程。

有约束最小二乘法是解决图像去模糊及加性噪声去除的通用方法。

本文针对常见的图像模糊和高斯噪声两类退化因素进行了图像复原算法的设计。

首先，实现了基于稀疏性约束的梯度投影算法，算法通过选择不同的搜索方式分为GPSR-Basic算法和GPSR-BB算法，实验结果表明GPSR-BB算法具有更优的细节保持能力及更快的求解速度；其次，实现了一种动量梯度投影算法，该算法在梯度下降方向上增加一个动量项，提高了算法的快速性和鲁棒性；最后，设计了一种自适应参数的动量梯度投影算法（A-Momentum），采用Barzila Borwein 的近似方式来获得参数的自适应设置，从而减少迭代次数。

实验结果表明在相同的截止条件下，A-Momentum算法仅需22次迭代即可到达稳定解，而GPSR-Basic算法和GPSR-BB算法则分别需要48次和47次。

关键词图像复原；梯度投影算法；动量梯度下降；自适应
ABSTRACT The regularization image restoration methods are to convert the ill-posed problem into an approximate problem by introducing the regular term related to the prior knowledge of the image, so as to obtain a stable approximate solution. Constrained least squares is a general method to solve image de-blur and additive noise removal. In this paper, we design the image restoration algorithm for the common image blur and Gaussian noise degradation factors. Firstly, the gradient projection algorithm based on sparsity constraint is realized. The algorithm is divided into GPSR-Basic algorithm and GPSR-BB algorithm by selecting different search methods. The experimental results show that the GPSR-BB algorithm has better detail retention ability and faster Secondly, a momentum gradient projection algorithm is implemented. The algorithm adds a momentum term to the descending direction of the gradient, which improves the fastness and robustness of the algorithm. Finally, a momentum gradient
projection of adaptive parameters is designed. The algorithm (A-Momentum) uses Barzila Borwein's approximation to obtain the adaptive setting of the parameters, thereby reducing the number of iterations. The experimental results show that the A-Momentum algorithm only needs 22 iterations to reach the stable solution under the same cut-off condition, while the GPSR-Basic algorithm and the GPSR-BB algorithm need 48 times and 47 times respectively.
KEY WORDS Image restoration, gradient projection algorithm, momentum gradient descent, adaptive
1 前言 (1)
1.1 研究背景及意义 (1)
1.2 图像复原方法与研究现状 (1)
1.3 本文主要内容 (3)
2 有约束最小二乘法基础 (4)
2.1 图像退化模型 (4)
2.2 有约束最小二乘法 (5)
2.2.1 图像复原的正则化模型 (5)
2.2.2 基于稀疏特性的正则化模型 (6)
2.2.3 TV正则化模型 (8)
2.2.4 正则化方法的优点 (8)
3 梯度投影算法 (10)
3.1 梯度投影（GPSR）算法 (10)
3.2 GPSR Basic算法 (11)
3.3 GPSR BB算法 (13)
3.4 实验结果分析 (15)
4 动量梯度投影法 (18)
4.1 动量梯度下降方向 (18)
4.2 动量梯度投影算法 (19)
4.3 自适应动量梯度投影算法 (21)
4.5 实验结果分析 (24)
5 结论 (28)
致谢 (29)
参考文献 (30)
附录 (31)
1.1 研究背景及意义
人类认识世界大都通过视觉、听觉、触觉等来获得信息，据统计，人类的外界信息大约有百分之七十五来自视觉系统，而视觉信息的来源是图像，图像是通过各种设备和系统采集得到的。

在图像的获取、传输以及保存过程中，因为成像设备的物理局限性、外界环境条件的限制，大气的湍流效应、摄像设备中光学系统的衍射、传感器特性的非线性、光学系统的像差等等原因，使得获取的图像会引入不同程度的模糊和噪声，造成图像质量的下降或者畸变和失真，对后续特征提取、图像分析、缺陷检测等带来影响[1]。

图像退化有很多因素，其中分为两大类：其一，噪声，在图像获得、采集、存储的过程中，由于热噪声、电子噪声、显示器件的结构噪声等等会导致图像退化。

其中噪声可以按照是否依赖于信号本身分类与按照统计特性分类，前者有加性噪声（不依赖于信号本身）、乘性噪声（依赖于信号本身），后者有高斯噪声（加性）、泊松噪声（乘性）、Gamma噪声等等。

其二，模糊，在成像设备远动模糊、聚焦不良散焦模糊、大气粒子（水下微生物粒子）影响、光反射、光折射等会造成图像模糊。

用数学来描述噪声可用点扩散函数（卷积）、高斯模糊（常考虑）、均匀模糊等等。

图像复原是数字图像复原的重要分支之一，它是利用图像退化的背景知识和图像的先验知识，建立图像退化的数学模型，对退化图像进行处理，从而复原出原始“干净”图像的过程。

图像复原技术在近几十年里，进步飞速，各种方法层出不穷。

消除噪声、改善图像的清晰度，从而提高图像处理应用的精确性、有效性是其目标。

它的应用范围已经深入到了天文观测、军事观测、医疗成像、矿产检测、数字城市监控、日常生活等等。

因此，图像复原研究具有重要的理论意义和实际应用价值。

1.2 图像复原方法与研究现状
随时代的发展，图像复原方法越来越受到人们关注，整体来看，可分为两类——确定性方法和正则化方法，确定性方法具有明确的滤波函数，如非邻域滤波法，最近邻域滤波法、维纳滤波、最小二乘滤波等，其中维纳滤波、最小二乘滤波的图像复原效果最好。

正则化方法是通过引入图像先验知识相关的正则项，将不适定问题转化为近似的适定问题，从而获得稳定的近似解[2]。

正则化方法起始于最小二乘条件下的逆滤波复原方法，其物理意义明确，数学表达式简单明了。

但逆滤波方法存在求逆的运算过程，当噪声不为零时，退化算子很小或者为零的时候，噪声会被放大，使得复原的图与原图有较大的差别。

因此，研究者们提出了有约束的最小二乘法，采用正则化约束项的设计来缩小最小二乘法的求解空间，从而获得更为精确的解。

Tikhonov正则化方法与TSVD(截断奇异值分解)方法是最早将正则化理论引入图像复原，其中TSVD方法虽然在图像复原的病态特性进行数学上的消除，但是本质上没有引入原图的先验信息。

Tikhonov方法引入了原图的二次函数的最小化问题，算法上特别好，缺点是会出现边缘模糊、寄生波纹。

压缩感知是由 E. J. Candes、J. Romberg、T. Tao 和 D. L. Donoho 等科学家于2004 年提出的，它作为一个新的采样理论，通过开发信号的稀疏特性，在远小于Nyquist 采样率的条件下，用随机采样获取信号的离散样本，然后通过非线性重建算法完美的重建信号。

压缩感知理论一经提出，就引起学术界和工业界的广泛关注。

2006年 Mario A. T. Figueiredo,，Robert D. Nowak,，Stephen J. Wright 将压缩感知理论应用到图像复原领域，提出了基于梯度投影的稀疏重建及图像复原算法。

该算法采用l1范数作为图像稀疏度的约束，刻画图像细节信息，取得了良好的复原效果。

这也是本文重点实现的算法之一。

随后，Chambolle等人提出了ISTA 算法,该算法的优点是迭代格式简单，但收敛速度较慢。

为了解决ISTA算法[11]的收敛速度问题，快速迭代收缩阈值算法(FISTA)[3]应运而生,快速迭代收缩阈值算法(FISTA)具有简单的算法的迭代格式，且其的收敛速率、数值试验的效果都比ISTA 算法要出彩的多[4]。

“外推系数带参数的邻近梯度算法”是刘紫娟等人在FISTA算法的根基上研究出的一种新的图像复原算法，“外推系数带参数的邻近梯度算法”在相应的条件下的全局收敛速率与“快速迭代收缩阈值算法(FISTA)”的全局收敛速率相同，在这个算法中，算法的外推系数具有一般的形态,同时该算法通过数值试验取得了比“快速迭代收缩阈值算法(FISTA)”更加明了的图像。

近年来,图像复原领域的学者们探讨了将图像不光滑的部分转化为光滑化部分,而且取得了很多有效且简单的图像复原算法，就如，GPSR_Basic算法与GPSR_BB算法，其中GPSR_Basic根基选取最快速下降法的步长选取相应的步长而GPSR_BB的步长要通过结合投影技术老获得。

还有，加速动量梯度投影法(MGPSR)，加速动量梯度投影法(MGPSR)是一种有效且简单的图像复原算法，加速动量梯度投影法(MGPSR)是
MaG等人[11]通过添加一个动量项到负梯度方向上而研究出的图像复原算法。

加速动量梯度投影法(MGPSR)的优缺点为：从数学领域来看，通过此方法能够在数字取得比较的结论，但其下降参数为固定值，取决于算法设计者的经验，从学习速率、动量参数来看，该算法不够灵敏，因此算法在鲁棒性和快速性方面还需要进一步提升。

本文在自适应参数设计方面对动量梯度投影算法进行了改进。

1.3 本文主要内容
本文针对常见的图像模糊和高斯噪声两类退化因素进行了图像复原算法的设计。

第二章讲关于图像复原的基础知识，包括讲图像退化模型，有约束最小二乘法的正则化模型。

第三章介绍了梯度投影法的原理，并通过实验，实现GPSR Basic算法与GPSR BB 算法。

第四章实现了动量梯度投影法[7]，并在这个算法的基础上，我改进算法，实现了自适应动量梯度投影法。

第五章进行了总结。

2 有约束最小二乘法基础
2.1 图像退化模型
导致图像退化的条件很多，如成像系统高速运动造成的模型、成像过程中的电子噪声、光电转换噪声等，如何依据实际需要，找出主要的退化因素并建立反映该退化因素的图像退化模型是至关重要的。

一般情况下，图像的退化因素可类比于线性系统的退化，从而建立系统退化模型来近似描述退化的原因。

思路为：
图2-1 建立系统退化模型的思路
依据图像信号和噪声之间的关系可将自然光成像过程中的退化因素归结为：加性热噪声引起的图像退化以及由点扩散函数造成的图像模糊退化。

在以上两类退化因素的影响下，图像退化过程见图2-2，其输入与输出关系可以用公式2-1来表示。

图2-2 退化模型 g(.)(,)*(,)(,)
x y H x y f x y n x y =+ （2-1）
其中，N M R y x ⨯∈),(f 为输入的原始干净图像，),0(~),x n 2I N y α（)为标准差为α均值为0的高斯噪声，N R y x g ∈),(为已知看到的退化图像，退化函数H 表示为PSF(点扩函数)或者降晰函数。

由上图的数学退化模型看出，在已知(x,)H y 、),(y x n 的情况下，进行反演运算退化图像，得到一个()y x ,f 的最佳估计()y x ,g ，就是图像复原。

而对图像复原问题来说，它是从退化的看到的图像估计出原始的“干净”图像的过程，是数学上典型的逆问题（Inverse Problem ）求解的研究范围。

从数学意义上来看，可以将式（2-1）所示的图像退化过程抽象为映射的形式，就是将“干净”图像映射到观测图像的过
程，其表达式[11]为：
)),((),(g y x f T y x = （2-2）
其中，M N R R T →•••
)(是成像系统特征相关的映射算子。

以就是说，其实求解一个Fredholm 算子方程[2-3]的问题就是求解图像复原问题，所有就是要发现一个逆变换，使得
)),((),(f 1y x g T y x -= （2-3）
在实际生活中，由于噪声大气湍流效应、光学系统的绕射、光学系统的像差、成像设备与物体的相对运动、传感器特性的非线性、感光胶卷的非线性和胶片颗粒噪声、摄像扫描所引起的几何失真等等干扰因素，使得在解决图像复原问题时，图像复原算法的求解空间变大，所以式（2-3）是一个病态性问题，直接求解十分困难、具有不适定性。

有约束正则化方法为解决此类不适定性问题提供了有效的途径。

2.2 有约束最小二乘法
正则化方法是通过引入图像先验知识相关的正则项，将不适定问题转化为近似的适定问题，从而获得稳定的近似解。

正则化方法的保真项为最小二乘法的二次函数项，约束项则由最初的Tikhonov [2]正则化发展至基于稀疏约束的正则化，因此正则化方法也叫做有约束的最小二乘法。

针对退化模型（2-1）的有约束最小二乘法一般可表示为如下形式[5]： 221min ()()2
x f x y Ax x τ=-+Φ （2-4）
其中22y Ax -称为数据保真项，用来限制观测数据和真实数据之间的误差大小；
:N R R Φ→为凸正则化函数，是对解空间的约束；(0,)τ∈+∞为正则化参数，用来平衡保真项和正则项在复原过程中的作用。

2.2.1 图像复原的正则化模型
图像复原早期的时候用的Tikhonov 正则化方法和截断奇异值分解（TSVD ）方法就是基于正则化理论。

利用源高质量图像是“平滑的”这一假定，并将这一假定作为先验信息引入到图像复原算法中，它的作用就是抑制近似解中的高频振荡噪声，所以得到如下模型：
22222
1)(min Lx Ax y x f x τ+-= （2-5）其中，τ是非负正则化参数，L 是一阶、二阶微分算子或者密度函数。

其中，τ的
作用是平衡保真项22Ax y -和解的平滑性度量22Lx ，如何选择τ的值，对实际应用
中是非常重要的，同时，这也是一个值得研究的问题。

如果从集合交集的角度求解式(2-4)的相关优化问题，则可以把式（2-5）的解认为是在椭圆y Q ,x 和椭圆x Q 的交集中，只要证明了当2）（E ε
τ=时，可直接进行求解就可以。

其中椭圆y Q ,x 和椭圆x Q 定义[10]为：
} |{x 22
2,x ε≤-=Ax y Q y ，}|{22
2E Lx x Q x ≤= （2-6)
其中，ε、E 为事先拟定的边界值。

在图像处理的过程中，很多图像都会呈现出不连续性，大部分图像都具有许多的边缘纹理细节，对于式（2-5）来说，对于二次正则项22Lx ，其主要功能为惩罚图像的不连续性，但过于平滑导致噪声和边缘细节信息
都同时被平滑掉，无法达到“保持边缘”和“平滑噪声”的平衡。

对于Tikhonov 正则化方法的缺点，因为基于总变分正则化模型、稀疏特性一范数模型为象征的非二次正则化模型不断出现，且因其良好细节保持特性越来越具有研究价值。

2.2.2 基于稀疏特性的正则化模型
在数据的传输过程中，存在着大量的数据冗余，这对采集数据的存储设备、硬件设备带来巨大的麻烦，且数据量太多导致难以用数学工具来进行分析计算，所以传统的图像复原方法难以处理越来越复杂的图像复原问题。

为了解决传统图像复原的问题，所以稀疏正则化方法应需求被研究出来，从数学领域来看，它是将图像的保真项、图像的稀疏约束相结合起来的优化模型，从而获得稀疏解；从模型构造来看，可将稀疏正则化模型具体的分为基于图像本身稀疏性约束、变换域稀疏性约束。

1) 基于图像本身稀疏性约束
0l 范数最小化问题是最早的稀疏约束问题，因为0l 范数的强非凸特性，导致该问题直接求解不具可行性。

直到2005年，论证了当随机测量矩阵A 满足有限等距性质(RIP)时，0l 范数模型可由1l 范数模型替换，因而0l 范数模型可转化为凸优化模型 （2-7）来求解[8]。

Ax y t s x x
=⋅⋅1min （2-7）
其中，N R ∈x 代表待求解图像，M R y ∈代表观测到的图像，N M R A ⨯∈代表测量矩阵。

∑==N i i x x 1
1是一范数，代表着图像稀疏性度量，他的作用为描画图像的边缘详细的
信息。

相比Tikhonov 正则化方法的二范数限制，一范数对图像中尖锐的边缘不灵敏，即对异常值不锐敏，所以在图像及信号处理领域能够被广泛应用。

因为Basis Pursuit De-noising (BPDN)算法将模型（2-7）扩展到了图像去噪范畴，而且BPDN 算法模型结合了“将1l 范数当作图像稀疏度的度量”与“2l 范数保真项”，创建了非常有名的12l l -模型[7]
122x 2
1min τ+-Ax y x （2-8）
其中，N R ∈x 代表待求解图像，M R y ∈代表观测到的图像，N M R A ⨯∈代表测量矩阵。

τ是非负常数。

另外，12l l -模型的两个等价有约束优化条件如下： 212min ..x t x
y Ax s t -≤ （2-9） ε≤-221..x m in Ax y t s x （2-10）
其中，ε、t 都是非负常数。

显而易见，分析于模型（2-9）和模型（2-10）可知，前者为二次约束的线性规划（QCLP ）问题，而后者为二次规划（QP ）问题。

基于“1l 范数极小化也就是说图像x 的多个分量的像素值靠近于0，蜕变出一个稀疏解”的特征，在将近的三十多年的时间内12l l -模型多次用于多个考虑稀疏特性的信号处理问题中。

2) 变换域稀疏性约束
变换域的稀疏正则化方法，从本质上来说就是对图像进行某种变换，通常变换方法是字典稀疏变换、小波变换等；从数学本质上来分析就是将稀疏约束、保真项
相结合来处理图像，其中的稀疏约束是使变换系数在统计数据上呈现稀疏分布。

12l l -模型在基于小波图像复原方法中的表达式[8]为：
1222
1)(m Wx RWx y x f in x τ+-= （2-11）其中，W 是包含了小波基的矩阵，RW 相当与模型（2-7）中的A ，x W 为未知图像的稀疏系数，R 是模糊矩阵，代表的是成像设备自身引起的退化，。

变换域的稀疏正则化方法常用于压缩感知（CS ）理论框架下的图像复原，它的过程是：首先，确立一个用来分析图像的字典库；然后，犹如模型（2-10）中的小波基一样对图像进行稀疏表示；最后，采取模型（2-8）对图像进行重建复原。

压缩感知理论广泛地应用于核磁共振成像、图像稀疏重建、医学CT 图像重建等领域，是一种在采样的同时实现数据压缩的理论框架。

2.2.3 TV 正则化模型
ROF 模型，即全变分去噪模型是1992年Fatemi 、Rudin 、Osher 提出，它是TV 正则化初始模型[6]，它的模型表达式为：
2..)(min σ=-y x t s x TV x （2-12）其中，1)(x X TV ∇=为全变分的一范数，2σ表示噪声方差。

将模型（2-12）转化为
无约束优化问题，得到如下的表达式：
)(2
1)(min x TV y x x f x τ+-= （2-13）其中，τ表示正则化参数。

该模型广泛应用于多个领域，因为其解空间为有界变差空间，能够在图像处理过程中保持着良好的图像边缘。

V ogel [7]等采取图像复原模（2-14），实现了相应的线性化快速算法，将TV 正则化模型扩展到了复原模糊噪声图像，表达式为：
)(2
1)(min x TV Ax y x f x τ+-= (2-14) 因为全变分去噪模型不完是适合图像复原处理的形态学原则，比如：获得稳态解中易形成“阶跃”现象、“分片常数”效应等等。

基于此，在模型表达式（2-14）表达的高阶全变分模型上引入高阶梯度约束，确保了平滑区域的光滑性、保证了图像复原的快速收敛，表达式为：
12y 2
1)(min x Ax x f x ∇+-=α （2-15） 随着研究的深入，TV 正则化方法慢慢发展成了图像复原处理领域中各类逆问题的多数求解方法。

2.2.4 正则化方法的优点
正则化方法灵活多变、求解简单的优点，渐渐成为了主流的图像复原方法。

其缘由为：首先，正则化方法能够限定复原解空间，因为图像复原保真项被施以正则化约束；其次，正则化方法能够减弱较小奇异值对复原效果的灾难性影响，因为原问题的解可以用一族与原问题近似的适定性问题的解来逼近；第三，随着多年研究者对正则化方法的改进，已经从不同视角对图像先验模型从事了刻画与研究，比如图像的边缘结构可以用全变分正则化模型来描画等等，所以可以根据具体的要求选择最恰当的方法；第四，因为正则化方法处理细节保持、噪声去除的问题有很好的
效果，而图像复原的主要问题是如何处理同时好平滑噪声、保持边缘和纹理细节，其中噪声去除、边缘保持的复原过程和平滑噪声要处理好对应着高频部分是相互矛盾，所以使用正则化方法是很有必要的；最后，在如今的数学领域，凸函数的求解方法很多且计算方便，如：共轭梯度法、梯度下降法等，而正则化方法的目标函数多为凸函数，且模型多为有约束的优化问题，所以正则化方法的扩展性好、理论系统完善。

3 梯度投影算法
自然光成像过程中的退化因素多归结为高斯加性噪声以及运动模糊，结合（2-1）所述的图像退化模型[8]，可归结为：
n Ax +=y （3-1） 其中，y 为观测到的退化图像，x 为待复原的原始干净图像，A 为点扩散函数得到的模糊矩阵，n 为高斯白噪声。

图像复原问题为由已知的观测图像y 来求解原始干净图像x 的过程。

1l 范数是图像稀疏性的度量，能够在去除噪声的同时很好地刻画图像的边缘细节信息，因此常用于模型（3-1）的求解。

联立最大后验概率条件下的保真项，得到如（3-2）所示的图像复原优化模型[11]：
12
2x *0.5m in τ+-Ax y x （3-2） 其中2
2Ax y -为数据保真项，是观测数据与真实数据误差的度量，1x 为正则化约束项，衡量细节保持的能力，τ为正则化参数，用来平衡保真项与正则化项所起的作用。

2007年，Figueiredo 等[7]将梯度投影法成功地应用于问题（3-2）的求解。

首先利用凸松弛变换，将21l l -模型转化为凸约束二次规划问题，并分别采用标准线搜索和Barzilai Borwein 搜索[2-6]提出了GPSR Basic 算法和GPSR BB 算法。

梯度投影算法的下降方向为负梯度方向，因此只需要计算目标函数的一阶导数，简单快速，同时，算法可在保证重建质量的同时提高收敛速度。

3.1 梯度投影（GPSR ）算法
图像复原模型（3-2）中一范数1x 不可微，直接求导困难，因此，文献[7]利用松弛变换将图像复原问题（3-2）转为有约束的凸二次规划问题（BCQP ），最后利用已经有的凸优化方法求解[9]。

首先，将变量x 化分为正、负两部分，在公式上表达式为：
0,0≥≥-=i i i
i i v u v u x （3-3） 这里的},0max {u i i x =， },0m ax {v i i x -=，其中一范数项可以表示为 v u x
T n T n 111+= 其中的]1.......11[1n =T。

将式（3-1）代入到图像复原模型（3-2）中。

然后，获得转化后的BCQP 问题
0..2
1)z (min ≥+=z t s Qz Pz z f T z （3-4） 式中的()R R Q W W W W P A y R A A W v u T
T T -+=⎥⎦⎤⎢⎣⎡--===⎥⎦⎤⎢⎣⎡=n 21,,,,x τ。

BCQP 问题（3-2）中的目标函数)(z f 为二次函数，并且在可行域范围内连续可微，所以BCQP 问题（3-2）是可行范围Ω上的有约束优化问题（3-5），梯度下降法、共轭梯度法、梯度投影法等凸优化理论的求解方法均可方便地对有约束优化问题（3-5）进行求解。

表达式为：
{}Ω∈z z f in z
:)(m （3-5） 所以，Figueiredo 等人提出解决此问题的Gradient Projection for Sparse Reconstruction （GPSR ）算法，即梯度投影算法，
3.2 GPSR Basic 算法
首先，在负梯度方向上获得一个正投影，其表达式为：
+-∇∂-=))(z (1k 1-k 1-k z f s k （3-6） 其中，k ∂ 为下降系数，其取值范围为k ∂>0。

其次，采用梯度投影法获得迭代方向，其表达式为：
M=11---k k z s
从而得到算法迭代表达式为：
M z k k k λ+=-1z （3-7） 其中，k λ为迭代步长。

再者，采用标准的线搜索，即采用反步搜索方法，其表达式为：
)))(z (()()()))(z ((1k 1-k 1k 1-k 111+---+--∇∂--∇-≤∇∂-k k z f z z f z f z f f k T k k μ （3-8）
进行回溯，确保迭代目标函数的单调性。

如果满足表达式（3-8），就是所求的解；若不满足，则令1k +→k ,返回至公式（3-6），再次进行计算。

综上，我们得到GPSR Basic 图像复原算法的具体流程为：
GPSR Basic 图像复原算法的具体步骤
Step 0 初始化0z ， 设置k 参数为0，设置参数∈β（0,1），∈μ（0，0.5），()k g Step 1 计算()k g 的值：
当0>k i z 或者()0)(<∇i k z f 时，其计算公式为：
()i k z f )(g k i ∇=
否则，0g k i =。

设定下降系数k ∂的最大值m ax ∂与最小值min ∂，并使),,max 0min 0∂∂∂=
∂（。

Step 2 反步搜索公式：
)))(z (()()()))(z ((1k 1-k 1k 1-k 111+---+--∇∂--∇-≤∇∂-k k z f z z f z f z f f k T k k μ Step 3 截止条件：如果满足表达式（3-8），就是所求的解；若不满足，
则令1k +→k ,返回至公式（3-6），再次进行计算。

GPSR Basic 图像复原算法的主程序：
while keep_going
x_previous = x;
% 计算梯度
temp = AT(resid_base);
term = temp - Aty;
gradu = term + tau;
Gradv = -term + tau;
%计算下降方向
%du = -gradu; dv = -gradv; dx = du-dv; dz=(du,dv)
dx = gradv-gradu;
old_u = u;
old_v = v;
du = max(u-lambda*gradu,0.0) - u;
u_new = u + du;
dv = max(v-lambda*gradv,0.0) - v;
v_new = v + dv;
dx = du-dv;
x_new = x + dx;
%反步搜索
if f_new <= f + mu * (gradu'*du + gradv'*dv)
%disp('OK')
break
end
lambda = lambda * lambda_backtrack;
fprintf(1,' reducing lambda to %6.2e\n', lambda)
end
%更新迭代
u = u_new;
v = v_new;
prev_f = f;
f = f_new;
uvmin = min(u,v);
u = u - uvmin;
v = v - uvmin;
x = u-v;
iter=iter+1
3.3 GPSR BB 算法
GPSR BB 图像复原算法是在GPSR Basic 图像复原算法的基础上进行改进的，采用了Barzila Borwein 的近似方式来获得参数的自适应设置[10]，其详细的步骤如下：
首先，计算下降方向，计算表达式为：
()()k k k k k f z z z z -∇∂-=+ （3-9）
其次为线搜索，在这过程中获得参数k λ，使得计算公式（3-10）的值最小，计算公式（3-10）的表达式为：
)z z (k k k f λ+ （3-10） 其中∈k λ[0,1]，所以，就得到新的迭代公式，表达式为：
k k k k r z z λ+= （3-11）
在者，更新参数∂，首先，计算表达式：
k T k k P δδγ)(=
如果k γ=0，则令1+∂k =m ax ∂，不然计算表达式就为:
⎪⎭
⎪⎬⎫⎪⎩⎪⎨⎧∂∂=∂+max 22m 1k ,z m k k in id γ， 最后，判断截止的条件，如果满足公式（3-12）
0,)f(z )
f(z -)f(z k k 1k >≤+δδ （3-12）
则 1z +k 为问题的解，若不满足，则令1k +→k ,返回至公式（3-9），再次进行计算。

综上，我们得到GPSR BB 图像复原算法的简单流程为：
GPSR BB 图像复原算法的简单步骤
Step 0 初始化0z ， 设置k 参数为0， 设定下降系数k ∂的最大值m ax ∂与
最小值min ∂，并使],[max min 0∂∂∈∂
Step 1 迭代公式：k k k k r z z λ+=
Step 2 更新参数∂。

如果k γ=0，则令1+∂k =m ax ∂，否则
⎪⎭
⎪⎬⎫⎪⎩⎪⎨⎧∂∂=∂+max 22m 1k ,z m k k in id γ，
Step 3 截止条件：如果满足表达式（3-12），就是所求的解；若不满足， 则令1k +→k ,返回至公式（3-9），再次进行计算。

GPSR-BB 算法的主流程与GPSR-Basic 算法类似，这里仅列出自适应参数设计的代码如下：
% 计算自适应参数
auv = A(dx);
dGd = auv(:)'*auv(:);
dd = du(:)'*du(:) + dv(:)'*dv(:);
if dGd <= 0
fprintf(1,' dGd=%12.4e, nonpositive curvature detected\n', dGd);
alpha = alphamax;
else
alpha = min(alphamax,max(alphamin,dd/dGd));
end
从本质上来说， GPSR 图像复原算法的下降方向是负梯度方向，因此导致在有复原信噪比比较低的源高图像时，对收敛性较差、局部最优值太过于敏感等缺点。

所以在下文中，介绍一种改进的梯度下降方向——动量梯度下降方向。