CN103839233A - 一种相机抖动造成的模糊图像复原方法 - Google Patents

Abstract

本发明涉及一种相机抖动模糊图像复原方法。相机抖动模糊图像复原过程中，模糊模型设计不当会导致复原效果不佳，计算效率低和内存占用量骤增，因此发展更好的相机抖动模糊图像复原方法已成为该领域亟待解决的问题。本发明设计了广义加性卷积模型，并基于该模型对相机抖动模糊图像进行复原：首先，估计出相机抖动动作路径；其次，计算出所有的切片状路径和纤维状路径，用贪心算法安排两种路径各自占有的比例；然后，用基于混合GACM的APG算法进行非盲复原。本发明复原视觉效果较好，兼顾高效快速，内存占用较少的特点，本发明适用于对各类相机抖动模糊图像复原。

Description

一种相机抖动造成的模糊图像复原方法

技术领域

本发明涉及的是一种图像复原方法，特别是针对相机抖动造成的模糊图像的复原方法，属于图像处理领域。

背景技术

图像复原是计算机视觉和图像处理领域一个非常基础的问题。由于散焦，相机抖动，图像中目标的动作等各种因素造成的图像模糊是难以避免的。在这样一个数字成像设备已经被广泛使用的时代，毫无疑问，对模糊图像的复原技术是一个非常热门的话题。目前为止，开发一个既精确，又拥有较强鲁棒性的图像复原技术依然是一个挑战性十足的问题。早期的科学研究都假设模糊图像中各个像素点上的模糊程度是相同的，基于这样的假设，很多方法也显示出较强的对模糊图像的复原能力。然而，近期的一系列研究分析结果证明：相机抖动造成的图像模糊是一种非均匀模糊，即，在模糊图像中，每个像素点上的模糊程度是空间变化的，这使得恢复相机抖动造成的模糊图像成为一个比均匀模糊图像去卷积更有难度的问题。虽然近年来有一些方法致力于恢复相机抖动模糊图像，但是其中主流的一些方法，例如，将整幅图像分割成很多重叠的区域并利用均匀复原的方法进行局部模糊核的估计，这种策略忽略了图像全局信息对局部模糊核的约束，显然没有合理利用图像信息，也无法得到一个令人满意的复原结果；另外一些基于映射动作路径模型的方法，通过将相机抖动造成的模糊定义为清晰图像经过一系列单应变换后结果的加权和，虽然能得到较好的复原效果，但是，基于这种模型的复原方法，计算效率低成了致命伤，虽然有很多方法改进了算法，并且提高了计算效率，但是这些方法又引进了另一种限制这种方法被广泛应用的弊端：计算过程中需要存储高维的稀疏矩阵，从而使整个复原方法的内存占用量骤增。总而言之，超负荷的计算量和严重的内存占用量都源自于复原过程中所用到的非均匀模糊模型。

发明内容

本发明的目的在于，研究开发一个合理的相机抖动模糊图像复原方法，其核心是设计一个模拟相机抖动模糊的模型，该模型物理概念清楚，并且能够高效计算，获得复原结果的几何模型，为后续采用复原后的清晰图像进行高级图像处理的工作打下基础。

基于上述目的，本发明技术方案的创新点在于将相机抖动造成的模糊模拟为是减少了数量的卷积变换结果的加权和，这种数量的减少是本发明特殊的技术实现的，从根本上降低了复原所需的时间，提高了复原效率。其中涉及到的卷积运算是数字图像处理领域常用的算子，在频域实现卷积算子，具有计算复杂度低(运算速度快)的优势，从而整体上降低复原方法的运行时间；另外，利用极变换将图像变换到极坐标下再进行对图像的旋转操作，结束后再将其用极坐标变换的共轭算子变换到普通空间，从而避免了复原过程中对高维稀疏矩阵的存储，降低了运算过程中的内存占用量。在本发明中建立了混合GACM-Generalized additive convolution model来模拟相机抖动模糊，在复原方法中，模糊模型的计算复杂度和内存占用量以及模拟的准确程度从根本上决定了复原方法的运算时间，内存占用量和复原效果。目前被广泛用来进行相机抖动模糊图像复原的算法都采用一种基于映射动作路径的模型，通过将相机抖动造成的模糊定义为清晰图像经过一系列单应变换后结果的加权和。相对于这种几何模型，混合GACM通过利用FFT变换及其逆变换能快速实现空间卷积运算的优势和利用极坐标变换避免了存储高维稀疏矩阵从而大幅降低计算过程中的内存占用量的优势提升了整体复原方法的效率。混合GACM在算法上相对于映射动作路径模型算法有以下的改进：

1)通过证明每个相机动作都可以分解为相机先进行平移操作，再进行旋转操作这一结论，GACM将3D空间(沿着x轴，y轴平移，绕z轴旋转，其中x，y轴为与相机镜头平面平行的坐标系所在的轴，z轴为垂直于相机镜头平面的轴)中相机抖动动作路径分解为两种不同的典型运动路径——切片和纤维，这两种动作路径的合理组合能够重现原始动作路径；

2)GACM分解出来的切片路径代表的是以给定的角度绕z轴旋转的相机抖动动作在平面内的平移动作构成的切片状路径。因此，切片路径的模糊模拟过程中只包含一次相机旋转动作，GACM通过极坐标变换实现；而大量的平移动作可以通过计算其平移程度的加权和合成为一次平移操作，通过一组FFT变换来快速实现，至此，整个切片路径中的模糊模拟得以精确实现；

3)GACM分解出来的纤维路径代表的是以给定的平移幅度沿x，y轴平移的相机抖动动作在平面内的旋转动作构成的纤维状路径。因此，纤维路径的模糊模拟过程中只包含一次相机平移动作，通过一组FFT变换来快速实现；而大量的旋转动作可以通过极坐标变换后在极坐标下快速实现，同时避免了映射动作路径模型中需要存储高维稀疏矩阵的缺点，降低了计算过程中的内存占用量；

4)然而单独使用基于切片的GACM或者是单独使用基于纤维的GACM都无法保证最大程度上提升方法的效率。本发明定义一种贪心算法来确定混合GACM中纤维和切片所占有的比例；

混合GACM的提出有严格的理论论证作为基础，拥有完美的物理角度的解释，因此，它可以精确的模拟相机抖动动作造成的图像模糊，提升了复原效果；与以上介绍的切片路径实现过程中的时间优势和纤维路径实现中的内存消耗优势使得本发明能够实现高效的相机抖动模糊图像复原。

本发明的具体内容如下：

建立一种新的相机抖动模糊模型——混合GACM，并利用该模型进行相机抖动模糊图像盲复原。复原工作是针对所有相机抖动模糊的自然图像进行的，在的验证实验中用到的图像都是现有的方法中常用到的相机抖动过程中拍摄到的模糊自然图像。首先，定义GACM为通过g_θ变换的向量x与模糊核的卷积结果再经过f_θ变换结果的累加。其中这两种变换的个数为C，C的数值越大，恢复图像所用的时间就越多，反之亦然。通过证明每个相机抖动动作都可以无误差的分解为两种子运动的合并，即，先对图像进行对应的平移操作，再进行旋转操作，GACM将整个相机抖动动作路径分解成了两种典型的路径——切片和纤维路径，其中切片路径中，所有的动作都拥有同样的绕z轴旋转的角度，在基于切片的GACM中，可以定义g_θ为恒等映射，定义f_θ为旋转变换，此时，在一个切片中的所有平移变换核函数的加权和构成一个核函数k_i；纤维路径中，所有的动作都拥有同样的相机平面内的平移幅度，在基于纤维的GACM中，可以定义g_θ为先平移输入图像，再对其进行极坐标变换映射到极坐标空间，定义f_θ为极坐标变换的逆变换，此时，在一个纤维中的所有旋转变换在极坐标空间内只是简单的加权运算，对极坐标变换和其逆变换，预先存储一组与图像大小相同的查找表来记录它们相应的坐标；而且其中的逆变换采用严格的极坐标变换的伴随算子，目的是避免由于插值和离散化造成的不连续性。本发明的目标在于快速高效的完成整个相机抖动模糊图像复原工作，而提高速度的关键在于降低g_θ变换或f_θ变换个数C。给定一个相机动作路径，纯粹的把整个路径都分解成是多个切片或多个纤维都无法保证得到足够小的C值，因此，本发明设计出一种贪心算法来将相机动作路径分解成切片和纤维的混合集合，以此来最大程度上降低C值。GACM可以和目前很多种图像先验条件和优化算法合作实现模糊图像复原。本发明中，采用APG-Accelerated proximal gradient算法来获得最终的复原的图像。

为实现以上内容，本发明的具体步骤如下：

1.估计出相机抖动动作路径以及动作路径中每个动作的权值。用共轭梯度法求解下面的优化问题解出W，即动作路径中每个动作的权值：

$\underset{W}{\min }{\left|\right|\underset{\mathrm{\θ}\∈P}{\mathrm{\Σ}}{w}_{\mathrm{\θ}}{K}_{\mathrm{\θ}}L-B\left|\right|}_2^2+{\mathrm{\Φ}}_1\left(W\right)---\left(1\right)$

其中，P是动作路径中所有动作构成的集合，B∈R^n×n是模糊图像，L∈R^n×n是初始的复原的图像，在这里预定为是已知量，Φ₁(W)是动作权值向量W的正则化项，是为了使得到的W更加接近真实值，本发明选择Φ₁(W)为W的l₂-范数；

首先定义加性卷积模型把相机抖动造成的非均匀模糊模拟为数目较小的卷积变换结果之和如下：

$B=\underset{i=1}{\overset{C}{\mathrm{\Σ}}}{f}_{\mathrm{\θ}}(g_{\mathrm{\θ}}\left(L\right)\⊗k_{\mathrm{\θ}})---\left(2\right)$

其中L和B是两个n×n大小的图像，

定义的是卷积算子，f_θ和g_θ是为了实现动作θ＝(θ_z，t_x，t_y)的关于像素的算子，k_θ是动作θ＝(θ_z，t_x，t_y)相关的滤波器，C是相机抖动路径中动作的个数，但通过使用切片状路径和纤维状路径集合，C的值会比相机抖动动作路径中的动作量少；加性卷积模型是广义加性卷积模型的一个特例，只要将上式中的f_θ定义为f_θ(x)＝α_θοx，把g_θ定义为g_θ(x)＝x；加性卷积模型可以用来做非盲非均匀模糊图像复原工作，其中非均匀模糊图像的范畴很广，例如，散焦图像，相机抖动造成的模糊图像等；广义加性卷积模型中，取适当的f_i和g_i，可以用来对相机抖动造成的非均匀模糊进行复原。

相机抖动动作路径表示为集合P＝{θ＝(θ_z，t_x，t_y)}，相机动作路径上的每个动作都是三维的，并且可以用以下形式的参数表示θ＝(θ_z，t_x，t_y)，其中θ_z是绕z轴进行旋转的角度，t_x和t_y分别是沿着x轴和y轴平移的量，然后根据证明相机抖动路径中每一个动作都等同于两个子动作的合体，即先对图像进行平移操作，然后对结果图像进行旋转：

假设路径P＝{θ＝(θ_z，t_x，t_y)}中动作个数为n个，每个动作我们记为{θ_j|j＝1，…n}，则上面这个公式中

和

是相机的3D动作路径中的动作{θ_j|j＝1，…n}，θ_j分解后的平移动作θ_j，t，θ_j分解后的旋转动作θ_j，r对应的运动变换矩阵，这些结果为GACM的建立提供基础；

先对图像进行平移操作，然后对结果图像进行旋转具体实现为

其中，

为平移动作θ_j，t对应的模糊核，是

的矩阵表示，

与L卷积就是对L进行平移操作，

为旋转动作θ_j，r对应的变换算子，

是算子的矩阵表示，

对做变换是对

做旋转操作。

变换矩阵K_θ的形成需要单应性矩阵H_θ来实现。对3D相机运动子空间中的每个动作θ_j＝(θ_z，j，t_x，j，t_y，j)，单应性可以定义为：

$\begin{array}{c}{H}_{{\mathrm{\θ}}_j}={\mathrm{CM}}_{{\mathrm{\θ}}_j}{C}^{-1}=C\left(\begin{array}{c}{R}_{{\mathrm{\θ}}_{z,j}}\end{array}+t_j\left[\begin{array}{ccc}0& 0& 1\end{array}\right]\right)C^{-1}\\ =C\left[\begin{array}{ccc}\cos \left({\mathrm{\θ}}_{z,j}\right)& -\sin \left({\mathrm{\θ}}_{z,j}\right)& t_{x,j}\\ \sin \left({\mathrm{\θ}}_{z,j}\right)& \cos \left({\mathrm{\θ}}_{z,j}\right)& t_{y,j}\\ 0& 0& 1\end{array}\right]C^{-1}\end{array}---\left(3\right)$

其中，t_j＝[t_x，j，t_y，j，0]^T是平移向量，且有：

$t_j\left[\begin{array}{ccc}0& 0& 1\end{array}\right]=\left[\begin{array}{c}{t}_{x,j}\\ t_{x,j}\\ 0\end{array}\right]\left[\begin{array}{ccc}0& 0& 1\end{array}\right]=\left[\begin{array}{ccc}0& 0& t_{x,j}\\ 0& 0& t_{x,j}\\ 0& 0& 0\end{array}\right],R_{{\mathrm{\θ}}_{z,j}}=\left[\begin{array}{ccc}\cos \left({\mathrm{\θ}}_{z,j}\right)& -\sin \left({\mathrm{\θ}}_{z,j}\right)& 0\\ \sin \left({\mathrm{\θ}}_{z,j}\right)& \cos \left({\mathrm{\θ}}_{z,j}\right)& 0\\ 0& 0& 1\end{array}\right]$

同时，假设相机内部参数是已知的，相机校准矩阵定义如下，

$C=\left|\begin{array}{ccc}f& 0& x_0\\ 0& f& y_0\\ 0& 0& 1\end{array}\right|---\left(4\right)$

其中，f是相机传感器分辨率与相机焦距的乘积除以相机景深得到的参数f＝l·s_max/ccd_max，单位为像素，其中，l是相机焦距长度，单位为毫米，s_max为相机拍摄到的图像的最大尺寸，ccd_max为相机传感器的尺寸最大值，(x₀，y₀)是图像中心像素点的坐标。

给定对应于动作θ_j＝(θ_z，j，t_x，j，t_y，j)的单应性矩阵

可以构造相应的将图像变换成模糊图像的变换矩阵

构造矩阵

的具体细节如下：

可以理解为是像素点到像素点的映射算子，即，对任意在位置[`，1]T的像素点，定义的是将图像x中处于

像素位置的像素值赋值到模糊图像上[x，y，1]位置的像素点。如果令θ_j＝(θ_z，j，t_x，j，t_y，j)，θ_j，r＝(θ_z，j，0，0)，和θ_j，t＝(0，t_x，j，t_y，j)，根据上述描述，定义

和

2.计算出所有的切片状路径。建立广义加性卷积模型，将相机抖动动作路径分解为两种特殊的相机动作路径集合，这一部分首先介绍如何分解出第一种路径集合：切片状路径集合。

1)在实际的相机抖动路径中，切片状路径集合

定义为S_θ＝{θ_j＝{θ_z，j，t_x，j，t_y，j}：θ_j∈P且θ_z，j＝θ}，即所有θ_z分量相同的相机抖动动作θ＝(θ_z，t_x，t_y)组成的集合，因为这些动作的旋转分量都相同，所以用一个旋转变换

就可以实现所有的旋转动作，如附图3所示；

2)定义f_θ(x)＝R_θ(x)，

则基于切片的广义加性卷积模型为：

$\mathrm{Kx}={\mathrm{\Σ}}_{\mathrm{\θ}}{R}_{\mathrm{\θ}}(k_{S_{\mathrm{\θ}}}\⊗x)---\left(5\right)$

其中，

的含义为切片S_θ对应的所有平移动作组成的平移模糊核，R_θ是相机抖动动作θ＝(θ_z，t_x，t_y)对应的旋转操作；

3)则给定一个切片S_θ，由切片S_θ中所有动作造成的相机抖动模糊可以用以下公式来模拟：

$\begin{array}{c}{K}_{\mathrm{\θ}}x=\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{K}_{\mathrm{\θ},r}{K}_{{\mathrm{\θ}}_{j,t}}x=\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{R}_{\mathrm{\θ}}(k_{t_j}\⊗x)\\ =R_{\mathrm{\θ}}(\left(\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{k}_{t_j}\right)\⊗x)=R_{\mathrm{\θ}}(k_{S_{\mathrm{\θ}}}\⊗x)\end{array}---\left(6\right)$

其中，

定义了关于z-轴旋转角度定义的平移模糊核，即切片上所有的动作拥有同一个旋转角度，但是每个动作的平移模糊核都不同，根据上式，这些不同的模糊核可以累加作为一个模糊核来运算，上式中第3个等式成立的原因是卷积算子和变换R_θ都是线性算子，简洁起见，这里不做证明。

4)不失一般性，将给定的相机动作路径集合P中的所有动作进行分类，P＝∪_θ{S_θ}，即，将整个路径都分解为无交集的切片集合的并集。

那么，公式(2)可以写成等同的形式：

$\mathrm{Kx}={\mathrm{\Σ}}_{\mathrm{\θ}}{R}_{\mathrm{\θ}}(k_{S_{\mathrm{\θ}}}\⊗x)---\left(7\right)$

同样的，我们可以定义Kx的伴随算子K^Ty如下：

$K^{T}y={\mathrm{\Σ}}_{\mathrm{\θ}}{\tilde{k}}_{S_{\mathrm{\θ}}}\⊗R_{\mathrm{\θ}}^T\left(y\right)---\left(8\right)$

其中

是

的伴随算子，其构成方法就是将

算子进行上下左右翻转，

定义为R_θ算子离散形式的严格定义的伴随算子。

5)为了实现R_θ的运算，可以直接利用Matlab语言中的imrotate函数来实现该操作。但是为了提高效率，采用了一种查表法(LUT：look up table)，即，在对离散角度进行旋转的前后，我们都存储一个坐标向量，这种方法的缺陷是我们需要提前计算和存储这n_z个查找表。在连续形式下，R_-θ是R_θ算子的伴随算子，但是在离散情况下，将整个路径离散化和计算过程中的差值运算造成的误差不可忽视。所以，相比于R_-θ算子，在整个运算过程中本发明用严格定义的伴随算子来避免这个问题，具体实现方法见说明书。

3.计算出所有的纤维状路径。这一部分介绍如何从相机抖动动作路径分解出纤维状动作路径集合。

1)在实际的相机抖动路径中，所有以给定的平移幅度沿x，y轴平移的相机抖动动作在平面内的旋转动作构成纤维状路径，定义为纤维；

2)定义公式(2)中f_θ(x)＝IPT(x)，则给定一个纤维F_t，由纤维中所有动作造成的相机抖动模糊可以用以下公式来模拟：

$K_{t}x=\underset{{\mathrm{\θ}}_j\∈F_t}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{K}_{\mathrm{\θ},r}{K}_{{\mathrm{\θ}}_{j,t}}x=\underset{{\mathrm{\θ}}_j\∈F_t}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{R}_{{\mathrm{\θ}}_{z,j}}\left(K_{t}x\right)=\mathrm{IPT}(w_t\⊗\mathrm{PT}\left(K_{t}x\right))---\left(9\right)$

其中，PT(·)和IPT(·)算子是极坐标变换及其逆变换，w_t是纤维f对应的动作集合中所有动作同样拥有的单一平移模糊核，在极坐标变换中，用分解旋转角度时一样的间隔来分段旋转路径，则上式中的基滤波器w_t可以定义如下：

$w_t=[w_{t,{\mathrm{\θ}}_1},w_{t,{\mathrm{\θ}}_2},...,w_{t,{\mathrm{\θ}}_{\mathrm{nf}}}]---\left(10\right)$

其中，θ₁是最小的旋转角度，θ_nf是最大的旋转角度。该纤维f中共有nf个动作，向量中每个元素

为纤维F_t中的动作θ_t在曝光过程中所占用的时间；

变换

的定义是：预先存储对应于一系列旋转操作w_t(称为一系列是指w_t是个向量)的变换前的图像(·)的坐标位置组成的向量和旋转变换后的图像的坐标位置组成的向量，即，如果该纤维中总共有m个动作，那么预先存储2×m个坐标位置向量，旋转操作是通过对原始图像用这些坐标位置向量进行读取操作实现的。

对一般的相机运动路径上的动作集合P中所有动作进行分类P＝∪_t{F_t}，即，将整个路径都分解为无交集的纤维集合的并集。那么，公式(2)可以写成等同的形式 $\mathrm{Kx}={\mathrm{\Σ}}_t\mathrm{IPT}(w_t\⊗\mathrm{PT}\left(K_{t}x\right))$

那可以定义Kx的伴随算子K^Ty如下： $K^{T}y={\mathrm{\Σ}}_t{K}_{-t}\left({\mathrm{PT}}^T({\tilde{w}}_t\⊗{\mathrm{IPT}}^T\left(x\right))\right)$

其中，

是算子w_t的伴随算子，PT^T和IPT^T分别是算子PT和IPT的伴随算子。

4.用贪心算法安排切片状路径和纤维状路径各自占有的比例。构建混合GACM：单独使用基于切片的广义加性卷积模型或者是单独使用基于纤维的广义加性卷积模型都无法保证最小化C，在广义加性卷积模型中，C越小，整个复原过程占用的时间就越少：

$\mathrm{KL}={\mathrm{\Σ}}_{s=1}^{n_s}{R}_{{\mathrm{\θ}}_s}(k_{s,t}\⊗L)+{\mathrm{\Σ}}_{f=1}^{n_f}\mathrm{IPT}(w_f\⊗\mathrm{PT}(k_{t_f}\⊗L))$

其中，n_s指的是相机抖动路径中利用下面的贪心算法计算出来的切片状路径的个数，相似的，n_f指的是利用下面的贪心算法计算出来的纤维状路径的个数。

1)针对混合模型中纤维和切片所占有的比例，定义一种贪心算法，首先输入：动作路径P和每个动作的权值(即每个动作在整个曝光过程中所占的时间）w；

2)初始化所有切片的权值

初始的误差设置为 $e=e_0={\mathrm{\Σ}}_{{\mathrm{\θ}}_z}\mathrm{WS}\left({\mathrm{\θ}}_z\right),$ 所有纤维的权值 $\mathrm{WF}(t_x,t_y)={\mathrm{\Σ}}_{{\mathrm{\θ}}_z}w({\mathrm{\θ}}_z,t_x,t_y),$ n_s＝0，n_f＝0，ε，0＜α≤0.9，初始的纤维集合和切片集合定义为空集；

3)当e/e₀＞ε时执行步骤4)-7)，如果不满足条件时，执行步骤(8)；

4)找出纤维中权值最大的那个纤维对应的平移动作以及相应的权值，

和找出切片中权值最大的那个切片对应的旋转动作以及相应的权值， ${\widehat{\mathrm{\θ}}}_z=\arg \underset{{\mathrm{\θ}}_z}{\max }\mathrm{WS}\left({\mathrm{\θ}}_z\right),w_s=\mathrm{WS}\left({\widehat{\mathrm{\θ}}}_z\right);$

5)更新切片集合：如果w_s＞αw_f，将第i_s个切片添加到集合S_s中，使

赋值 $W_{t_x,t_y}(j,k)=W_{t_x,t_y}(j,k)-w_{t_s,j,k},$ 且n_s＝n_s+1；

6)更新纤维集合：如果αw_f≥w_s，将第(j_s，k_s)个纤维添加到集合S_f中，并且使 $W_{t_x,t_y}(j_f,k_f)=0,$ 修改 $W_{{\mathrm{\θ}}_z}\left(i\right)=W_{{\mathrm{\θ}}_z}\left(i\right)-w_{i,j_f,k_f},$ n_f＝n_f+1；

7)更新阈值： $e/e_0={\mathrm{\Σ}}_t{W}_{{\mathrm{\θ}}_z}\left(i\right)/{\mathrm{\Σ}}_i{W}_{{\mathrm{\θ}}_z}\left(i\right);$

8)如果e/e₀≤ε，则结束，输出：纤维集合S_f＝{(t_i，F_i，w_i)：i＝1，...，n_f}和切片集合S_s＝{(q_j，S_j，k_j)：j＝1，...，n_s}，其中，n_f是集合S_j中纤维的个数，n_s是集合S_s中切片的个数，它们决定了切片和纤维路径在分解后的路径中的比例，t_i指的是第i个纤维F_i中的平移幅度，w_i指的是该纤维中所有旋转动作的参数，包括旋转角度和该动作对应的权值；q_j指的是第j个切片S_j中的旋转幅度，k_j指的是该切片中所有平移动作模糊核的加权和。

5.用基于混合GACM的APG算法进行非盲去模糊。根据步骤1计算出的动作路径中每个动作的权值W，解决下述的优化问题来对相机抖动图像进行非盲复原：

$\underset{L}{\min }{\left|\right|{\mathrm{\Σ}}_{s=1}^{n_s}{R}_{{\mathrm{\θ}}_s}(k_{s,t}\⊗L)+{\mathrm{\Σ}}_{f=1}^{n_f}\mathrm{IPT}(w_f\⊗\mathrm{PT}(k_{t_f}\⊗L))-B\left|\right|}_2^2+{\mathrm{\Φ}}_2\left(L\right)---\left(12\right)$

其中Φ₂(L)是复原后图像L的正则化项，该正则化项用的是梯度的信息，是为了使得到的L保留和恢复更清晰的边缘信息；

6.最后，判断||L_opt-L||≤1e-4是否成立，其中L_opt为当前得到的优化结果，L为上一次迭代得到的优化结果，如果不成立，则判断到此为止的迭代次数是否已经超过最大迭代次数20，如果没超过，就继续循环执行1-5步骤，如果||L_opt-L||≤1e-4成立或者目前的迭代次数已经大于20，则结束迭代，当前解出的L_opt就是该复原过程的最终解；

本发明的优点在于，GACM可以精确的模拟相机抖动造成的模糊图像，通过利用卷积变换和定义的f_i和g_i计算复杂度较低这个特性，降低复原过程的时间复杂度，整个算法运行过程在20次迭代内必定收敛，故设定最大迭代次数为20，通过与现行的其他方法做比较，发现本发明方法所用的时间是现有的最好的方法所用时间的1/3；GACM避免了目前的方法中将相机抖动每个动作代表的变换矩阵用高维的稀疏矩阵存储的缺陷，很大幅度上降低复原过程中内存的占用。本发明利用收敛速度比较快的APG优化算法，可以保证对相机抖动造成的模糊图像复原有非常好的结果，在估计出的清晰图像和上一次估计出的清晰图像之误差小于10^-4时，复原视觉效果已经是目前最好的结果，比较令人满意。本发明同时兼顾高效快速的特点和降低所使用的计算机内存，运行过程中占用的内存的最高值是现有的最好的方法所用内存的1/6，也就是说，本发明的效率很高。本研究适用于对任何相机抖动模糊图像复原，更重要的是GACM可以很好地与很多优化算法和先验项进行融合并高效快速的得到复原结果，克服了其他相机抖动模糊模型只能跟个别的优化算法结合的缺点。

附图说明

图1是该发明的具体实施流程；

图2是模拟的相机抖动动作路径；

图3所有以给定的角度绕z轴旋转的相机抖动动作在平面内的平移动作构成切片状路径——切片；

图4所有以给定的平移幅度沿x，y轴平移的相机抖动动作在平面内的旋转动作构成纤维状路径——纤维；

图5利用设计的贪心算法来对切片和纤维路径在分解后的相机抖动动作路径中的比例；

具体实施方式

本发明的具体实施流程见图1，下面结合附图，对本发明的具体实施方式作进一步描述：

1.首先，由给定的L∈R^n×n，估计出相机抖动动作路径，即，集合P＝{θ＝(θ_z，t_x，t_y)}，如附图2所示是一个一维的曲线，相机动作路径上的每个动作都是三维的，并且可以用以下形式的参数表示θ＝(θ_z，t_x，t_y)，其中θ_z是绕z轴进行旋转的角度，t_x和t_y分别是沿着x轴和y轴平移的量，用共轭梯度法求解下面的优化问题解出W，即动作路径中每个动作的权值：

$\underset{W}{\min }{\left|\right|\underset{\mathrm{\θ}\∈P}{\mathrm{\Σ}}{w}_{\mathrm{\θ}}{K}_{\mathrm{\θ}}L-B\left|\right|}_2^2+{\mathrm{\Φ}}_1\left(W\right)---\left(13\right)$

其中

是相应的相机抖动动作对应的模糊变换矩阵，W是所有w_θ组成的向量，其维数与集合P中元素的个数相同。w_θ为动作θ＝(θ_z，t_x，t_y)在曝光过程中占用的时间，Φ₁(W)是动作权值向量W的正则化项，是为了使得到的W更加接近真实值，选择Φ₁(W)为W的l₂-范数。这里必须要强调的是，实际的相机动作路径可以构成一条1维的动作路径，大部分动作的权值为零。所以，我们只需要考虑其中一些动作权值大于零的动作构成的子集P＝{θ：w_θ＞0}。

满足图1中的流程图中判断步骤的条件，即小于阈值或大于最大迭代次数之前，步骤1-6是循环迭代进行的；

假设附图2中的相机抖动路径就是通过第一步中的优化算法计算出来的，如表1，附图5中按照从左到右的顺序来区分每个动作的权值

在图中出现的顺序，下面举例说明。

表1相机抖动动作路径参数

θ∈P	θ15	θ16	θ17	θ18	θ19	θ20	θ21	θ22	θ23	θ24	θ25	θ26	θ27
														θ2/弧度	0.005	0.005	0.005	0.005	0.004	0.003	0.002	0.001	0.001	0.001	0.001	0.001	0.001
tx/毫米	0.002	0.003	0.004	0.005	0.005	0.005	0.005	0.005	0.004	0.004	0.004	0.004	0.005
														ty/毫米	0.002	0.002	0.002	0.002	0.002	0.002	0.002	0.002	0.002	0.003	0.004	0.005	0.005
w/秒	0.0589	0.0437	0.0332	0.0289	0.0560	0.0589	0.0616	0.0603	0.0458	0.0062	0.0362	0.0115	0.0536

如表1所示，假设相机曝光时间设置为1秒，相机焦距为20毫米，传感器尺寸最大为16毫米拍得的图像大小为256*256大小，那么这27个相机抖动动作构成一个相机抖动路径，其在相机3维空间中的表示见附图2，这27个动作在曝光过程中所停留的时间即为表格最下方的w值。从表1也可以清楚的看到，每个动作θ∈P都由1个旋转动作θ_z和两个平移动作t_x和t_y组成。

变换矩阵K_θ的形成需要单应性矩阵H_θ来实现。对3维相机运动子空间中的每个动作θ_j＝(θ_z，j，t_x，j，t_y，j)，单应性可以定义为：

$\begin{array}{c}{H}_{{\mathrm{\θ}}_j}={\mathrm{CM}}_{{\mathrm{\θ}}_j}{C}^{-1}=C\left(\begin{array}{c}{R}_{{\mathrm{\θ}}_{z,j}}\end{array}+t_j\left[\begin{array}{ccc}0& 0& 1\end{array}\right]\right)C^{-1}\\ =C\left[\begin{array}{ccc}\cos \left({\mathrm{\θ}}_{z,j}\right)& -\sin \left({\mathrm{\θ}}_{z,j}\right)& t_{x,j}\\ \sin \left({\mathrm{\θ}}_{z,j}\right)& \cos \left({\mathrm{\θ}}_{z,j}\right)& t_{y,j}\\ 0& 0& 1\end{array}\right]C^{-1}\end{array}---\left(14\right)$

其中，t_j＝[t_x，j，t_y，j，0]^T是平移向量，且有：

$t_j\left[\begin{array}{ccc}0& 0& 1\end{array}\right]=\left[\begin{array}{c}{t}_{x,j}\\ t_{x,j}\\ 0\end{array}\right]\left[\begin{array}{ccc}0& 0& 1\end{array}\right]=\left[\begin{array}{ccc}0& 0& t_{x,j}\\ 0& 0& t_{x,j}\\ 0& 0& 0\end{array}\right],R_{{\mathrm{\θ}}_{z,j}}=\left[\begin{array}{ccc}\cos \left({\mathrm{\θ}}_{z,j}\right)& -\sin \left({\mathrm{\θ}}_{z,j}\right)& 0\\ \sin \left({\mathrm{\θ}}_{z,j}\right)& \cos \left({\mathrm{\θ}}_{z,j}\right)& 0\\ 0& 0& 1\end{array}\right]$

同时，我们假设相机内部参数是已知的，相机校准矩阵定义如下，

$C=\left|\begin{array}{ccc}f& 0& x_0\\ 0& f& y_0\\ 0& 0& 1\end{array}\right|$

其中，f是相机传感器分辨率与相机焦距的乘积除以相机景深得到的参数f＝l·s_max/ccd_max，单位为像素，其中，l是相机焦距长度，单位为毫米，s_max为相机拍摄到的图像的最大尺寸，ccd_max为相机传感器的尺寸最大值，(x₀，y₀)是图像中心像素点的坐标。

那么，举个例子，对于动作θ¹，单应性矩阵可以定义为：

$\begin{array}{c}{H}_{{\mathrm{\θ}}^1}=\left|\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right|\left[\begin{array}{ccc}\cos \left(0.02\right)& -\sin \left(0.002\right)& 0.001\\ \sin \left(0.002\right)& \cos \left(0.002\right)& 0.008\\ 0& 0& 1\end{array}\right]{\left|\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right|}^{-1}\\ =\left[\begin{array}{ccc}1& -0.0020& 0.5763\\ 0.0020& 1& 2.3043\\ 0& 0& 1\end{array}\right]\end{array}$

这里需要强调的是，

可以理解为是像素点到像素点的映射算子，即，对任意在位置[x，y，1]^T的像素点，

定义的是将图像x中处于

像素位置的像素值赋值到模糊图像

上[x，y，1]位置的像素点。

$H_{{\mathrm{\θ}}_j}^{-1}{\left[\mathrm{1,1,1}\right]}^T=\left[\begin{array}{ccc}1& 0.0020& -0.5809\\ -0.0020& 1& -2.3031\\ 0& 0& 1\end{array}\right]*\left[\begin{array}{c}1\\ 1\\ 1\end{array}\right]=\left[\begin{array}{c}0.4211\\ -1.3051\\ 1\end{array}\right]$

那么，举个例子，对于模糊图像上位于[1，1，1]的像素点，它的像素值为原始图像x中位于[0.4211  -13051  1]位置的像素值。显然，原始图像中不存在这样的点，对这种情况，本发明用双线性插值来计算这个像素点的像素值，超出图像的像素点的像素值通过预先假设的边界值进行运算，如果超出图像边界过多，则用0值代替。

2.其次，计算出所有的切片状路径。可以证明相机抖动动作路径中每一个动作对图像的变换都等同于两个子动作轮流对图像变换的结果，它们分别是先对图像进行平移操作，然后对结果图像进行旋转，用公式表示是

其中和

是相机的3D动作路径中的动作θ_j，平移动作θ_j，t，旋转动作θ_j，r对应的运动变换矩阵，该过程具体实现时

其中，

为平移动作θ_j，t对应的模糊核，是

的矩阵表示，

与L卷积就是对L进行平移操作，

为旋转动作

对应的变换算子，是

算子的矩阵表示，

对

做变换是对

做旋转操作。根据定义的广义加性卷积模型把相机抖动造成的非均匀模糊模拟为数目较小的卷积变换结果之和如下：

$y=\underset{i=1}{\overset{C}{\mathrm{\Σ}}}{f}_i(g_l\left(x\right)\⊗k_i)---\left(15\right)$

通过证明相机抖动路径中每一个动作都等同于两个子动作的合体，即先对图像进行平移操作，然后对结果图像进行旋转，把相机抖动动作路径分解为两种特殊的相机动作路径：

(1)计算出所有的切片状路径集合。在实际的相机抖动路径中，所有以给定的角度绕z轴旋转的相机抖动动作在平面内的平移动作构成切片状路径集合，定义为切片；

公式(2)中定义f_θ(x)＝R_θ(x)，g_θ(x)＝x，则给定一个切片S_θ，由切片S_θ中所有动作造成的相机抖动模糊可以用以下公式来模拟：

$\begin{array}{c}{K}_{\mathrm{\θ}}x=\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{K}_{\mathrm{\θ},r}{K}_{{\mathrm{\θ}}_{j,t}}x=\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{R}_{\mathrm{\θ}}(k_{t_j}\⊗x)\\ =R_{\mathrm{\θ}}(\left(\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{k}_{t_j}\right)\⊗x)=R_{\mathrm{\θ}}(k_{S_{\mathrm{\θ}}}\⊗x)\end{array}---\left(16\right)$

其中，定义了关于z-轴旋转角度定义的平移模糊核，即切片上所有的动作拥有同一个旋转角度，但是每个动作的平移模糊核都不同，根据上式，这些不同的模糊核可以累加作为一个模糊核来运算，上式中第3个等式成立的原因是卷积算子

和变换R_θ都是线性算子，简洁起见，这里不做证明。

不失一般性，将给定的相机动作路径集合P中的所有动作进行分类，P＝∪_θ{S_θ}，即，将整个路径都分解为无交集的切片集合的并集。

那么，公式(2)可以写成等同的形式：

同样的，我们可以定义Kx的伴随算子K^Ty如下：

其中

是

的伴随算子，其构成方法就是将

算子进行上下左右翻转，

定义为R_θ算子离散形式的严格定义的伴随算子。

为了实现R_θ的运算，可以直接利用Matlab语言中的imrotate函数来实现该操作。但是为了提高效率，采用了一种查表法(LUT：look up table)，即，在对离散角度进行旋转的前后，我们都存储一个坐标向量，这种方法的缺陷是我们需要提前计算和存储这n_z个查找表。在连续形式下，R_-θ是R_θ算子的伴随算子，但是在离散情况下，将整个路径离散化和计算过程中的差值运算造成的误差不可忽视。所以，相比于R_-θ算子，在整个运算过程中本发明用严格定义的伴随算子

来避免这个问题。

(2)计算出所有的纤维状路径集合。在实际的相机抖动路径中，如图4所示，所有以给定的平移幅度沿x，y轴平移的相机抖动动作在平面内的旋转动作构成纤维状路径，定义为纤维，公式(2)中定义f_θ(x)＝IPT(x)，g_θ(x)＝PT(K_tx)，则给定一个纤维F_t，由纤维F_t中所有动作造成的相机抖动模糊GACM可以用以下公式来模拟：

$K_{t}x=\underset{{\mathrm{\θ}}_j\∈F_t}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{K}_{\mathrm{\θ},r}{K}_{{\mathrm{\θ}}_{j,t}}x=\underset{{\mathrm{\θ}}_j\∈F_t}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{R}_{{\mathrm{\θ}}_{z,j}}\left(K_{t}x\right)=\mathrm{IPT}(w_t\⊗\mathrm{PT}\left(K_{t}x\right))---\left(17\right)$

其中，PT(·)和IPT(·)算子是极坐标变换及其逆变换.在极坐标变换中，我们用分解旋转角度时一样的间隔来分段旋转路径。则上式中的基滤波器w_t可以定义为 $w_t=[w_{t,{\mathrm{\θ}}_1},w_{t,{\mathrm{\θ}}_2},...,w_{t,{\mathrm{\θ}}_{\mathrm{nz}}}]$

其中，θ₁是最小的旋转角度，θ_nz是最大的旋转角度。该纤维f中共有nz个动作，向量中每个元素

为纤维F_t中的动作θ_i在曝光过程中所占用的时间；

变换

的定义是：预先存储对应于一系列旋转操作w_t(称为一系列是指w_t是个向量)的变换前的图像(·)的坐标位置组成的向量和旋转变换后的图像的坐标位置组成的向量，即，如果该纤维中总共有m个动作，那么预先存储2×m个坐标位置向量，旋转操作是通过对原始图像用这些坐标位置向量进行读取操作实现的。

对一般的相机运动路径上的动作集合P中所有动作进行分类P＝∪_t{F_t}，即，将整个路径都分解为无交集的纤维集合的并集。那么，公式(2)可以写成等同的形式： $\mathrm{Kx}={\mathrm{\Σ}}_t\mathrm{IPT}(w_t\⊗\mathrm{PT}\left(K_{t}x\right))$

那可以定义Kx的伴随算子K^Ty如下：

$K^{T}y={\mathrm{\Σ}}_t{K}_{-t}\left({\mathrm{PT}}^T({\tilde{w}}_t\⊗{\mathrm{IPT}}^T\left(x\right))\right)---\left(18\right)$

其中，

是算子w_t的伴随算子，PT^T和IPT^T分别是算子PT和IPT的伴随算子。

3.利用设计的贪心算法来计算出切片和纤维路径在分解后的相机抖动动作路径中的比例，即算法中的输出纤维集合S_f＝{(t_i，F_i，w_i)：i＝1，...，n_f}和切片集合S_s＝{(q_j，S_j，k_j)：j＝1，...，n_s}，其中，n_f是集合S_j中纤维的个数，n_s是集合S_s中切片的个数，它们决定了切片和纤维路径在分解后的路径中的比例，t_i指的是第i个纤维F_i中的平移幅度，w_i指的是该纤维中所有旋转动作的参数，包括旋转角度和该动作对应的权值；q_j指的是第j个切片S_j中的旋转幅度，k_j指的是该切片中所有平移动作模糊核的加权和：

1)针对混合模型中纤维和切片所占有的比例，定义贪心算法，首先输入：动作路径P和每个动作的权值w，即每个动作在整个曝光过程中所占的时间；

2)初始化所有切片的权值初始的误差设置为 $e=e_0={\mathrm{\Σ}}_{{\mathrm{\θ}}_z}\mathrm{WS}\left({\mathrm{\θ}}_z\right),$ 初始化所有纤维的权值 $\mathrm{WF}(t_x,t_y)={\mathrm{\Σ}}_{{\mathrm{\θ}}_z}w({\mathrm{\θ}}_z,t_x,t_y),$ n_s＝0，n_f＝0，即初始的纤维集合和切片集合定义为空集，ε为判断迭代是否结束的阈值，0＜α≤1为贪心算法中决定切片和纤维个数比例的参数。

以表1为例：首先将27个动作全部分解为无交集的切片集合的并集P＝∪_θ{S_θ}，可分为5个，即如下表2，将所有拥有相同旋转角度的动作分配到一个切片集合中其中，

包括的动作总共6个{θⁱ|i＝22，23，24，25，26，27}，包括的动作总共15个{θⁱ|i＝4，5，…，18}。

表2相机抖动动作路径切片纤维集合

以表1为例：其次将27个动作全部分解为无交集的纤维集合的并集P＝∪_θ{F_θ}，可分为19个，即如下表3(表3-1和表3-2分成了两个表格，统称为表3)将所有拥有相同旋转角度的动作分配到一个切片集合中：

表3-1相机抖动动作路径纤维集合

表3-2相机抖动动作路径纤维集合

3)当e/e₀＞ε，在本例中，

4)找出纤维中权值最大的那个纤维对应的平移动作以及相应的权值，

和

在上述的例子中，我们可以看到，权值最大的纤维

中的平移为

权值为w_f＝0.2657，找出切片中权值最大的那个切片对应的旋转动作以及相应的权值，

在上述的例子中，我们可以看到，权值最大的切片

中的旋转角度为权值为w_s＝0.4499；

5)更新切片集合：如果w_s＞αw_f，将第i_s个切片添加到集合S_s中，使

赋值

且n_s＝n_s+1，上述例子中，如果取α＝1，则有，w_s＞αw_f，所以我们将第5个切片添加到集合S_s中，见附图5中第一个分图，即Slice1，使

表2即如下表4所示，由于15个动作{θⁱ|i＝4，5，…，18}已经被分到第一个切片集合，所以这些动作所在的所有纤维权值中都要减去这些权值的权重(因为它们都已经完成分类)，赋值 $W_{t_x,t_y}\left(8\right)=W_{t_x,t_y}\left(8\right)-w_4,W_{t_x,t_y}\left(6\right)=W_{t_x,t_y}\left(6\right)-w_5,$ $W_{t_x,t_y}\left(4\right)=W_{t_x,t_y}\left(4\right)-w_6,......,W_{t_x,t_y}\left(18\right)=W_{t_x,t_y}\left(18\right)-w_{18}$ 等(表5-1和表5-2由于空间不足分成了两个表格，统称为表5，表5中由该步骤引起变化的权值都已加下划线，动作集合中被选中的动作已设置为空)，定义n_s＝1，即现在已经选出1个切片集合；故这次迭代中，跨过第6)步，直接进行第7)步。

表4相机抖动动作路径切片纤维集合

表5-1相机抖动动作路径纤维集合

表5-2相机抖动动作路径纤维集合

6)更新纤维集合：如果αw_f≥w_s，将第(j_s，k_s)个纤维添加到集合S_f中，并且使 $W_{t_x,t_y}(j_f,k_f)=0,$ 修改 $W_{{\mathrm{\θ}}_z}\left(i\right)=W_{{\mathrm{\θ}}_z}\left(i\right)-w_{i,j_f,k_f},$ n_f＝n_f+1；

7)更新阈值：

本例中，e/e₀＝0.5501；

8)如果e/e₀≤ε，则结束，输出：纤维集合S_f＝{(t_i，F_i，w_i)：i＝1，...，n_f}和切片集合S_s＝{(q_j，S_j，k_j)：j＝1，...，n_s}，如果不满足终止条件，则继续迭代。

附图5是通过步骤3中计算好的切片和纤维来构造混合GACM，其中每个动作的编号顺序是从左到右，是进行分类的根据。本例中，如附图5所示，总共迭代4次结束，选出2个切片集合，两个纤维集合。最后迭代的结果用表格表示是：

表6相机抖动动作路径切片纤维集合

表7相机抖动动作路径纤维集合

$\mathrm{KL}={\mathrm{\Σ}}_{i=1}^{n_f}\mathrm{IPT}(w_i\⊗\mathrm{PT}(k_{t_i}\⊗x))+{\mathrm{\Σ}}_{j=1}^{n_s}{R}_{{\mathrm{\θ}}_j}(k_j\⊗x)---\left(19\right)$

按照以上例子中的情形，上式应该写为：

$\mathrm{KL}={\mathrm{\Σ}}_{i=1}^2\mathrm{IPT}(w_i\⊗\mathrm{PT}(k_{t_i}\⊗x))+{\mathrm{\Σ}}_{j=1}^2{R}_{{\mathrm{\θ}}_j}(k_j\⊗x)---\left(20\right)$

其中，

为切片集合

中3个动作θ¹，θ²，θ³共有的平移模糊核，平移参数为(0.001，0.008)，根据公式(14)，则有：

$K_{t_1}=\left|\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right|\left[\begin{array}{ccc}1& 0& 0.001\\ 0& 1& 0.008\\ 0& 0& 1\end{array}\right]{\left[\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right]}^{-1}=\left[\begin{array}{ccc}1& 0& 0.3200\\ 0& 1& 2.5600\\ 0& 0& 1\end{array}\right]$

这个模糊核可以不用在图像每个像素点进行逐一寻找平移结果，而是通过公式中的卷积就能快速实现整个图像的平移变换，而对于

过程的实现，是不能利用卷积算子的，因为旋转变换对每个像素点的作用都不同，而是与被变换的像素点位置有关，故本发明通过存储查找表的方法来准确实现该过程：对于

中动作θ¹，θ²，θ³，所有的旋转操作θ_z有3种：

那么对于像素点位置[x，y，1]＝[25，48，1]的像素，本发明的方法采用这样的策略，对于3种旋转，记录3种映射后的坐标值：

$K_{{\mathrm{\θ}}_z^1}=\left|\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right|\left[\begin{array}{ccc}\cos \left(0.0522\right)& -\sin \left(0.0522\right)& 0\\ \sin \left(0.0522\right)& \cos \left(0.0522\right)& 0\\ 0& 0& 1\end{array}\right]{\left[\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right]}^{-1}=\left[\begin{array}{ccc}0.9986& -0.0522& 6.8529\\ 0.0522& 0.9986& -6.5042\\ 0& 0& 1\end{array}\right]$

$K_{{\mathrm{\θ}}_z^2}=\left|\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right|\left[\begin{array}{ccc}\cos \left(0.0477\right)& -\sin \left(0.0477\right)& 0\\ \sin \left(0.0477\right)& \cos \left(0.0477\right)& 0\\ 0& 0& 1\end{array}\right]{\left[\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right]}^{-1}=\left[\begin{array}{ccc}0.9989& -0.0477& 6.2489\\ 0.0477& 0.9989& -5.9577\\ 0& 0& 1\end{array}\right]$

$K_{{\mathrm{\θ}}_z^3}=\left|\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right|\left[\begin{array}{ccc}\cos \left(0.0603\right)& -\sin \left(0.0603\right)& 0\\ \sin \left(0.0603\right)& \cos \left(0.0603\right)& 0\\ 0& 0& 1\end{array}\right]{\left[\begin{array}{ccc}320& 0& 128\\ 0& 320& 128\\ 0& 0& 1\end{array}\right]}^{-1}=\left[\begin{array}{ccc}0.9982& -0.0603& 7.9464\\ 0.0603& 0.9982& -7.4811\\ 0& 0& 1\end{array}\right]$

$K_{{\text{\θ}}_z^1}\left(\begin{array}{c}25\\ 48\\ 1\end{array}\right)=\begin{array}{c}\left(\begin{array}{c}29.3144\\ 42.7348\\ 1\end{array}\right),K_{{\mathrm{\θ}}_z^2}\left(\begin{array}{c}25\\ 48\\ 1\end{array}\right)\end{array}=\left(\begin{array}{c}28.9317\\ 43.1798\\ 1\end{array}\right),K_{{\mathrm{\θ}}_z^3}\left(\begin{array}{c}25\\ 48\\ 1\end{array}\right)=\left(\begin{array}{c}30.0083\\ 41.9383\\ 1\end{array}\right)$

本发明的方法是记录三个像素位置[29.3144，42.7348]，[28.9317，43.1798]，[30.0083，41.9383]，那么具体的实现过程如下：这3种变换形成的新的图像像素位置[x，y，1]＝[25，48，1]的像素值为原始图像像素位置[29.3144，42.7348]，[28.9317，43.1798]，[30.0083，41.9383]的像素值的加权和，权值为之前定义的W，这些像素位置大部分并非整数，这就需要对这些像素值进行双线性插值获得。该技术为数字图像处理领域常用的技术，此处不做详细解释。需要说明的是，上述过程中存储的变换前像素位置和变换后像素位置即为查找表。

4.然后用APG算法进行复原的图像结果：

$L_{\mathrm{opt}}=\underset{L}{\min }{\left|\right|{\mathrm{\Σ}}_{s=1}^{n_s}{R}_{{\mathrm{\θ}}_s}(k_{s,t}\⊗L)+{\mathrm{\Σ}}_{f=1}^{n_f}\mathrm{IPT}(w_f\⊗\mathrm{PT}(k_{t_f}\⊗L))-B\left|\right|}_2^2+{\mathrm{\Φ}}_2\left(L\right)---\left(21\right)$

其中，Φ_２(L)是复原后图像L的正则化项，该正则化项用的是梯度的信息，是为了使得到的L保留和恢复更清晰的边缘信息。复原过程中会涉及到Kx的伴随算子

如公式(18)，伴随算子在具体实施过程中用到了之前提到的预先存储的查找表，伴随算子的具体实施办法就是将变换后的坐标位置进行关于原始位置的中心对称，下面举例说明：

步骤3的例子中，3个变换将原始像素点[x，y，1]＝[25，48，1]映射到下面3个新位置：[29.3144，42.7348]，[28.9317，43.1798]，[30.0083，41.9383]，取整后即为[29，43]，[29，43]，[30，42]，那么其伴随算子的作用就是将[x，y，1]＝[25，48，1]点映射到这3个位置关于位置[25，48]的中心对称点：[21，51]，[[21，51]]，[20，54]。这里只提供该例，实际操作时每个原始点对应的是4个线性插值的新像素点，这里不做具体解释。

5.最后，判断||L_opt-L||≤1e-4是否成立，如果不成立，则判断到此为止的迭代次数是否已经超过最大迭代次数20，如果没超过，就继续循环执行1-5步骤，如果||L_opt-L||≤1e-4成立或者目前的迭代次数已经大于20，则结束迭代，当前解出的L_opt就是该复原过程的最终解，如流程图图1所示；

表8对于不同大小的图像，基于广义加性卷积(GACM)模型的相机抖动模糊图像复原方法所占用的时间和计算机内存；观察该统计结果，可以发现本发明在给相机抖动造成的模糊图像进行复原时，所用的时间非常少(为现有的最好的方法，如基于几何模型的方法和一些基于EFF的非几何模型的方法，所用时间的1/3)，并且运行过程中占用的内存也保持非常稳定的低值(现有的最好的方法所用内存的1/6)。

表8图像复原过程中时间和内存占用情况

图像编号	图像大小	时间(s)	内存(Gb)
				#1	512*768	945.67	2.06

#2	512*768	973.15	2.06
				#3	512*768	901.61	2.05
#4	512*768	530.69	2.06
				#5	512*768	885.72	2.05
#6	512*768	927.99	0.98
				#7	406*679	531.34	1.43
#8	512*768	904.92	1.96

最终结合流程图1可将本实时评价系统总结为如下步骤：

a)估计出相机抖动动作路径及每个动作的权值；

b)计算出所有的切片状路径；

c)计算出所有的纤维状路径；

d)用贪心算法安排切片状路径和纤维状路径各自占有的比例；

e)用基于GACM计算出的混合GACM的APG算法进行复原；

f)最终采用阈值和最大迭代次数判断该算法是否达到收敛条件。

本发明说明书中未作详细描述的内容属于本领域专业技术人员公知的现有技术。

Claims (3)

1.一种相机抖动造成的模糊图像复原方法，其特征在于：首先，估计出相机抖动动作路径；其次，计算出所有的切片状路径和纤维状路径，用贪心算法来计算切片状路径和纤维状路径各自占有的比例；然后，用基于混合广义加性卷积模型GACM的APG算法进行非盲复原；具体包括以下步骤：

1)首先，估计出相机抖动动作路径及每个动作的权值，即集合P＝{θ＝(θ_z，t_x，t_y)}和向量W，相机动作路径上的每个动作都是三维的，并且可以用以下形式的参数表示θ＝(θ_z，t_x，t_y)，其中θ_z是图像绕z轴旋转的角度，t_x和t_y分别是沿着x轴和y轴平移的量，由初始的假设得到的清晰图像L∈R^n×n，用共轭梯度法求解下面的优化问题解出W，即动作路径中每个动作的权值：

$\underset{W}{\min }{\left|\right|\underset{\mathrm{\θ}\∈P}{\mathrm{\Σ}}{w}_{\mathrm{\θ}}{K}_{\mathrm{\θ}}L-B\left|\right|}_2^2+{\mathrm{\Φ}}_1\left(W\right)$

其中

是相机抖动动作对应的模糊变换矩阵，W是所有w_θ组成的向量，其维数与集合P中元素的个数相同；w_θ为动作θ＝(θ_z，t_x，t_y)在曝光过程中占用的时间，Φ₁(W)是动作权值向量W的正则化项，为W的l₂-范数；

2)其次，计算出所有的切片状路径和纤维状路径；

相机抖动造成的非均匀模糊定义为如下广义加性卷积模型GACM：

$B=\underset{i=1}{\overset{C}{\mathrm{\Σ}}}{f}_{\mathrm{\θ}}(g_{\mathrm{\θ}}\left(L\right)\⊗k_{\mathrm{\θ}})$

其中L和B是两个n×n大小的图像，

定义的是卷积算子，f_θ和g_θ是为了实现动作θ＝(θ_z，t_x，t_y)的关于像素的算子，k_θ是动作θ＝(θ_z，t_x，t_y)相关的滤波器，C是相机抖动路径中动作的个数；

相机抖动动作路径中每一个动作θ对图像的变换都等同于两个子动作轮流对图像变换的结果，它们分别是先对图像进行平移操作，然后对结果图像进行旋转，用公式表示是

假设路径P＝{θ＝(θ_z，t_x，t_y)}中动作个数为n个，每个动作记为{θ_j|j＝1，…n}，则之前所述

和

是相机的3D动作路径中的动作{θ_j|j＝1，…n}，θ_j分解后的平移动作θ_j，t和θ_j分解后的旋转动作θ_j，r对应的运动变换矩阵，该过程具体实现为

其中，

为平移动作θ_j，t对应的模糊核，

是

的矩阵表示，与L卷积就是对L进行平移操作，

为旋转动作θ_j，r对应的变换算子，

是

(·)算子的矩阵表示，对做变换是对

做旋转操作。

变换矩阵K_θ的形成需要单应性矩阵H_θ来实现，对3D相机运动子空间中的每个动作θ_j＝(θ_z，j，t_x，j，t_y，j)，单应性可以定义为：

$\begin{array}{c}{H}_{{\mathrm{\θ}}_j}={\mathrm{CM}}_{{\mathrm{\θ}}_j}{C}^{-1}=C\left(\begin{array}{c}{R}_{{\mathrm{\θ}}_{z,j}}\end{array}+t_j\left[\begin{array}{ccc}0& 0& 1\end{array}\right]\right)C^{-1}\\ =C\left[\begin{array}{ccc}\cos \left({\mathrm{\θ}}_{z,j}\right)& -\sin \left({\mathrm{\θ}}_{z,j}\right)& t_{x,j}\\ \sin \left({\mathrm{\θ}}_{z,j}\right)& \cos \left({\mathrm{\θ}}_{z,j}\right)& t_{y,j}\\ 0& 0& 1\end{array}\right]C^{-1}\end{array}$

其中，t_j＝[t_x，j，t_y，j，0]^T是平移向量，且有：

$t_j\left[\begin{array}{ccc}0& 0& 1\end{array}\right]=\left[\begin{array}{c}{t}_{x,j}\\ t_{x,j}\\ 0\end{array}\right]\left[\begin{array}{ccc}0& 0& 1\end{array}\right]=\left[\begin{array}{ccc}0& 0& t_{x,j}\\ 0& 0& t_{x,j}\\ 0& 0& 0\end{array}\right],R_{{\mathrm{\θ}}_{z,j}}=\left[\begin{array}{ccc}\cos \left({\mathrm{\θ}}_{z,j}\right)& -\sin \left({\mathrm{\θ}}_{z,j}\right)& 0\\ \sin \left({\mathrm{\θ}}_{z,j}\right)& \cos \left({\mathrm{\θ}}_{z,j}\right)& 0\\ 0& 0& 1\end{array}\right]$

假设相机内部参数是已知的，相机校准矩阵定义如下，

$C=\left|\begin{array}{ccc}f& 0& x_0\\ 0& f& y_0\\ 0& 0& 1\end{array}\right|$

其中，f是相机传感器分辨率与相机焦距的乘积除以相机景深得到的参数f＝l·s_max/ccd_max，单位为像素，其中，l是相机焦距长度，单位为毫米，s_max为相机拍摄到的图像的最大尺寸，ccd_max为相机传感器的尺寸最大值，(x₀，y₀)是图像中心像素点的坐标。

给定对应于动作θ_j＝(θ_z，j，t_x，j，t_y，j)的单应性矩阵

可以构造将图像变换成模糊图像的变换矩阵

构造矩阵的具体步骤如下：

是像素点到像素点的映射算子，即对任意在位置[x，y，1]^T的像素点，将图像x中处于

像素位置的像素值赋值到模糊图像

位置的像素点；如果令θ_j，r＝(θ_z，j，0，0)，和θ_j，t＝(0，t_x，j，t_y，j)，类似的，定义

为对应于动作θ_j，r＝(θ_z，j，0，0)的单应性矩阵，

$H_{{\mathrm{\θ}}_{j,r}}={\mathrm{CM}}_{{\mathrm{\θ}}_j}{C}^{-1}=C\left[\begin{array}{ccc}\cos \left({\mathrm{\θ}}_{z,j}\right)& -\sin \left({\mathrm{\θ}}_{z,j}\right)& 0\\ \sin \left({\mathrm{\θ}}_{z,j}\right)& \cos \left({\mathrm{\θ}}_{z,j}\right)& 0\\ 0& 0& 1\end{array}\right]C^{-1}$

为对应于动作θ_j，t＝(0，t_x，j，t_y，j)的单应性矩阵，

$H_{{\mathrm{\θ}}_{j,t}}={\mathrm{CM}}_{{\mathrm{\θ}}_j}{C}^{-1}=C\left[\begin{array}{ccc}1& 0& t_{x,j}\\ 0& 1& t_{y,j}\\ 0& 0& 1\end{array}\right]C^{-1}$

将图像x中处于

像素位置的像素值赋值到模糊图像上[x，y，1]位置的像素点；

将图像x中处于

像素位置的像素值赋值到模糊图像

上[x，y，1]位置的像素点；同时，相机抖动路径中每一个动作都等同于两个子动作的合体，即先对图像进行平移操作，然后对结果图像进行旋转：

定义f_θ(x)＝R_θ(x)，g_θ(x)＝x，其中，R_θ是相机抖动动作θ＝(θ_z，t_x，t_y)对应的旋转操作，切片状路径集合S_θ定义为S_θ＝{θ_j＝{θ_z，j，t_x，j，t_y，j}：θ_j∈P且θ_z，j＝θ}，即所有θ_z分量相同的相机抖动动作θ＝(θ_z，t_x，t_y)组成的集合，因为这些动作的旋转分量都相同，所以用一个旋转变换R_θ就可以实现所有的旋转动作；

定义f_t(x)＝IPT(x)，g_t(x)＝PT(K_tx)，其中，PT(·)和IPT(·)算子是极坐标变换及其逆变换，K_t是相应的纤维对应的平移模糊矩阵，纤维状路径集合F_t定义为F_t＝{θ_j＝{θ_z，j，t_x，j，t_y，j}：θ_j∈P and(t_x，j，t_y，j)＝t}，即每一个t_x和t_y分量都相同的相机抖动动作组成的集合，因为这些动作的平移分量都相同，所以用一个平移模糊矩阵K_t就可以实现所有的平移动作；

相机抖动动作路径分解为切片状路径和纤维状路径这两种相机动作路径，然后对结果图像进行旋转，建立了混合广义加性卷积模型：

$\mathrm{KL}={\mathrm{\Σ}}_{s=1}^{n_s}{R}_{{\mathrm{\θ}}_s}(k_{s,t}\⊗L)+{\mathrm{\Σ}}_{f=1}^{n_f}\mathrm{IPT}(w_f\⊗\mathrm{PT}(k_{t_f}\⊗L))$

其中，L为清晰图像，在极坐标变换中，用分解旋转角度时一样的间隔来分段旋转路径，则基滤波器w_f可以定义为

其中，θ_i是纤维f中的第i个动作，该纤维f中共有n_f个动作，向量中每个元素为对应的纤维中的动作在曝光过程中所占用的时间，n_s指的是利用贪心算法计算出来的相机抖动路径中切片状路径集合的个数，相似的，n_f指的是利用贪心算法计算出来的纤维状路径集合的个数；

指的是从第一个切片状路径到第n_s个切片状路径将所有路径模糊的图像累加，每个路径的模糊过程是先对清晰图像做平移操作，平移的模糊核k_s，t为相应切片对应的所有平移动作组成的模糊核，再对该模糊结果进行旋转操作，旋转尺度为第s个切片所对应的唯一的旋转动作分量θ_s；

指的是从第一个纤维状路径到第n_f个纤维状路径，将所有路径模糊的图像累加，每个路径的模糊过程是先对清晰图像做平移操作，平移的模糊核为相应的纤维动作集合对应的单一平移模糊核，再通过PT(·)变换将该模糊结果映射到极坐标空间中，在极坐标空间中进行旋转的权值为w_f，结束后再用极坐标变换的逆变换IPT(·)将结果图像映射到原来的空间；

3)设计贪心算法来计算切片状路径和纤维状路径各自占有的比例；为了保证广义加性卷积模型GACM中C最小化，使得整个复原过程占用的时间越少，建立贪心算法，建立切片状路径组成的集合和纤维状路径组成的集合，并统计相应集合中动作子集合的个数n_s和n_f；

4)用基于混合GACM的APG算法完成非盲复原：固定动作路径中每个动作的权值，解决下述的优化问题来对相机抖动图像进行复原

$\underset{L}{\min }{\left|\right|{\mathrm{\Σ}}_{s=1}^{n_s}{R}_{{\mathrm{\θ}}_s}(k_{s,t}\⊗L)+{\mathrm{\Σ}}_{f=1}^{n_f}\mathrm{IPT}(w_f\⊗\mathrm{PT}(k_{t_f}\⊗L))-B\left|\right|}_2^2+{\mathrm{\Φ}}_2\left(L\right)$

其中，B为需要被复原的模糊图像，Φ₂(L)是复原后图像L的正则化项，该正则化项用的是梯度的信息，是为了使得到的L保留和恢复更清晰的边缘信息。

5)最后，判断||L_opt-L||≤10^-4是否成立，其中L_opt为当前得到的优化结果，L为上一次迭代得到的优化结果，如果不成立，则判断当时的迭代次数是否已经超过最大迭代次数20，如果没超过，就继续循环执行1)-4)步骤，如果||L_opt-L||≤1e-4成立或者目前的迭代次数已经大于20，则结束迭代，当前解出的L_opt就是该复原过程的最终解。

2.根据权利要求1所述相机抖动造成的模糊图像复原方法，其特征在于：所述的相机抖动动作路径分解为切片状路径和纤维状路径这两种相机动作路径，其具体步骤为：

1)在实际的相机抖动路径中，所有以给定的角度绕z轴旋转的相机抖动动作在平面内的平移动作构成切片状路径集合，定义为切片；

2)定义f_θ(x)＝R_θ(x)，g_θ(x)＝x，则给定一个切片S_θ，由切片S_θ中所有动作造成的相机抖动模糊可以模拟为：

$K_{\mathrm{\θ}}x=\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{K}_{\mathrm{\θ},r}{K}_{{\mathrm{\θ}}_{j,t}}x=\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{R}_{\mathrm{\θ}}(k_{t_j}\⊗x)=R_{\mathrm{\θ}}(\left(\underset{{\mathrm{\θ}}_j\∈S_{\mathrm{\θ}}}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{k}_{t_j}\right)\⊗x)=R_{\mathrm{\θ}}(k_{S_{\mathrm{\θ}}}\⊗x)$

其中，

定义了关于z-轴旋转角度定义的平移模糊核，即切片上所有的动作拥有同一个旋转角度，但是每个动作的平移模糊核都不同，根据上式，这些不同的模糊核可以累加作为一个模糊核来运算；

将给定的相机动作路径集合P中的所有动作进行分类，将整个路径都分解为无交集的切片集合的并集，

P＝{θ＝(θ_z，t_x，t_y)}＝∪_θ{S_θ}

那么，清晰图像模糊变换的过程可以写成等同的形式：

$\mathrm{Kx}={\mathrm{\Σ}}_{\mathrm{\θ}}{R}_{\mathrm{\θ}}(k_{S_{\mathrm{\θ}}}\⊗x)$

同样的，可以定义Kx的伴随算子K^Ty如下：

$K^{T}y={\mathrm{\Σ}}_{\mathrm{\θ}}{\tilde{k}}_{S_{\mathrm{\θ}}}\⊗R_{\mathrm{\θ}}^T\left(y\right)$

其中是的伴随算子，其构成方法就是将

算子进行上下左右翻转，

定义为R_θ算子离散形式的严格定义的伴随算子；

3)在实际的相机抖动路径中，所有以给定的平移幅度沿x，y轴平移的相机抖动动作在平面内的旋转动作构成纤维状路径，定义为纤维；即在GACM中定义f_t(x)＝IPT(x)，g_t(x)＝PT(K_tx)，给定一个纤维F_t，由纤维F_t中所有动作造成的相机抖动模糊可以用以下公式来模拟：

$K_{t}x=\underset{{\mathrm{\θ}}_j\∈F_t}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{K}_{\mathrm{\θ},r}{K}_{{\mathrm{\θ}}_{j,t}}x=\underset{{\mathrm{\θ}}_j\∈F_t}{\mathrm{\Σ}}{w}_{{\mathrm{\θ}}_j}{R}_{{\mathrm{\θ}}_{z,j}}\left(K_{t}x\right)=\mathrm{IPT}(w_t\⊗\mathrm{PT}\left(K_{t}x\right))$

其中，PT(·)和IPT(·)算子是极坐标变换及其逆变换.在极坐标变换中，我们用分解旋转角度时一样的间隔来分段旋转路径；则上式中的基滤波器w_t可以定义为

其中，θ₁是最小的旋转角度，θ_nz是最大的旋转角度。该纤维f中共有nz个动作，向量中每个元素为纤维F_t中的动作θ_t在曝光过程中所占用的时间；

变换

的定义是：预先存储对应于一系列旋转操作w_t的变换前的图像(·)的坐标位置组成的向量和旋转变换后的图像的坐标位置组成的向量，即，如果该纤维中总共有m个动作，那么预先存储2×m个坐标位置向量，旋转操作是通过对原始图像用这些坐标位置向量进行读取操作实现的；

对一般的相机运动路径上的动作集合P中所有动作进行分类，将整个路径都分解为无交集的纤维集合的并集P＝{θ＝(θ_z，t_x，t_y)}＝∪_t{F_t}，那么，GACM可以写成等同的形式：

$\mathrm{Kx}={\mathrm{\Σ}}_t\mathrm{IPT}(w_t\⊗\mathrm{PT}\left(K_{t}x\right))$

那可以定义Kx的伴随算子K^Ty如下：

其中，是算子w_t的伴随算子，PT^T和IPT^T分别是算子PT和IPT的伴随算子。

3.根据权利要求1所述的贪心算法来计算切片状路径和纤维状路径各自占有的比例，其特征在于：所述的贪心算法的实现步骤为：

1)针对混合模型中纤维和切片所占有的比例，首先输入：动作路径P和每个动作的权值w，即每个动作在整个曝光过程中所占的时间；

2)初始化所有切片的权值

初始的误差设置为 $e=e_0={\mathrm{\Σ}}_{{\mathrm{\θ}}_z}\mathrm{WS}\left({\mathrm{\θ}}_z\right),$ 初始化所有纤维的权值 $\mathrm{WF}(t_x,t_y)={\mathrm{\Σ}}_{{\mathrm{\θ}}_z}w({\mathrm{\θ}}_z,t_x,t_y),$ n_s＝0，n_f＝0，即初始的纤维集合和切片集合定义为空集，ε为判断迭代是否结束的阈值，0＜α≤1为贪心算法中决定切片和纤维个数比例的参数；

3)当e/e₀＞ε，初始值为1；

4)找出纤维中权值最大的那个纤维对应的平移动作以及相应的权值，和

找出切片中权值最大的那个切片对应的旋转动作以及相应的权值， ${\widehat{\mathrm{\θ}}}_z=\arg \underset{{\mathrm{\θ}}_z}{\max }\mathrm{WS}\left({\mathrm{\θ}}_z\right),w_s=\mathrm{WS}\left({\widehat{\mathrm{\θ}}}_z\right);$

5)更新切片集合：如果w_s＞αw_f，将第i_s个切片添加到集合S_s中，使

赋值 $W_{t_x,t_y}(j,k)=W_{t_x,t_y}(j,k)-w_{t_s,j,k},$ 且n_s＝n_s+1；

6)更新纤维集合：如果αw_f≥w_s，将第(j_s，k_s)个纤维添加到集合S_f中，并且使 $W_{t_x,t_y}(j_f,k_f)=0,$ 修改 $W_{{\mathrm{\θ}}_z}\left(i\right)=w_{{\mathrm{\θ}}_z}\left(i\right)-w_{{i,j}_f,k_f},$ n_f＝n_f+1；

7)更新阈值： $e/e_0={\mathrm{\Σ}}_i{W}_{{\mathrm{\θ}}_z}\left(i\right)/{\mathrm{\Σ}}_i{W}_{{\mathrm{\θ}}_z}\left(i\right);$

8)如果e/e₀≤ε，则结束，输出：纤维集合S_f＝{(t₁，F_t，w_t)：i＝1，...，n_f}和切片集合S_s＝{(_qj，S_j，k_j)：j＝1，...，n_s}，其中，n_f是集合S_j中纤维的个数，n_s是集合S_s中切片的个数，它们决定了切片和纤维路径在分解后的路径中的比例，t_i指的是第i个纤维F_i中的平移幅度，w_i指的是该纤维中所有旋转动作的参数；q_j指的是第j个切片S_j中的旋转幅度，k_j指的是该切片中所有平移动作模糊核的加权和。

Images