开放科学（资源服务）标志码（OSID）：

0 概述

图像运动模糊是相机在曝光时被摄物体和相机间发生相对位移导致的一种图像退化现象^[1]，成像时相机和被摄物体间发生的相对位移会使得不同位置的像素发生叠加，从而导致成像质量下降。运动模糊作为最常见的图像退化现象之一，广泛存在于众多视觉处理任务中，而高质量的清晰图像是高效完成此类任务的关键，因此，运动模糊图像复原问题一直是国内外学者们关注的焦点。

根据退化图像的点扩散函数（Point Spread Function，PSF）是否已知，可以将图像复原算法分为图像非盲复原算法和图像盲复原算法两类。图像盲复原是在PSF完全未知的情况下进行图像复原，如基于正则化的方法和基于深度学习的方法。基于正则化的图像复原方法通过设立正则项约束，以引导复原模型在迭代的过程中向符合图像先验和真实清晰图像的方向进行优化。TU等^[2]对图像复原算法的边缘保留正则化项进行改进，减少了由运动模糊引起的边界误差，并通过模糊映射法和下采样插值技术有效提高了算法效率。为了解决多帧超分辨率（Multi-Frame Super-Resolution，MFSR）中的运动模糊问题，MA等^[3]设计一个EM框架来引导PSF估计和高分辨率图像重建，并加入L2正则化项来抑制噪声。基于深度学习的方法利用预先收集的大量训练数据隐式地学习一个从降质图像到清晰图像的映射。王伟鹏等^[4]提出一种结合卷积变换与细节优化的图像去雾算法，其采用卷积变换的思想和分段映射后处理操作，并结合大气散射模型实现雾化图像的细节提升。缪斯等^[5]提出一种可微分神经网络架构搜索方法，该方法可用于图像盲去模糊时神经网络的参数调节。此外，文献[6-8]首先从退化图像中提取出模糊参数，然后用估计出的PSF复原出清晰图像。

本文通过分析运动模糊图像的频谱特征，提出一种改进的PSF参数估计方法，以提高估计结果的精度和有效估计范围，并通过PSF参数估计和运动模糊图像复原实验来验证所提方法的有效性。

1 相关工作

PSF求解是图像复原的关键步骤。LEE等^[9]通过运动模糊图像的频域周期性估计PSF参数，并利用图像后处理技术改善重建图像的质量。ZHU等^[10]提出一种利用多幅图像间接估计PSF的方法，其通过结合多通道PSF正则项和L1图像稀疏正则化的最小化代价函数来粗略估计PSF，然后根据PSF的稀疏性和正确率将估计结果进行细化并输出解卷积。MOGHADDAM等^[11]引入模糊集概念，以对基于Radon变换的PSF参数估计算法进行改进。

由于实际中相机成像尺寸有限引起的边缘截断破坏了卷积的连续性，使得运动模糊图像的频谱图像中心位置会出现一道十字亮线，如图 1所示，这种十字亮线会对模糊参数的估算造成极大干扰，可能导致估算结果误差过大甚至不可信。

针对十字亮线干扰问题，研究人员提出了不同的解决方案^[12-13]。胡发焕等^[14]利用骨架化变换操作细化暗条纹，有效提高了条纹检测精度。乐翔等^[15]通过频谱分块切除亮线干扰，将暗亮条纹对比最明显的块图反色后再进行Radon变换以提取PSF参数。胡硕等^[16]对频谱图像进行自适应形态学滤波以消除亮线干扰，但是，该方法需要设定阈值以对中心亮条纹进行迭代腐蚀，导致其效率和估计精度有限。

2 运动模糊图像退化模型及其频谱特性分析

图像复原^[17]是一个病态的反问题，在数学定义上，图像退化是图像与PSF的卷积过程，图像复原则是这一过程的逆过程^[18-19]。若要恢复出清晰的真实图像，必须先建立图像退化模型。图 2所示为空不变系统中运动模糊图像的退化模型。

图像退化模型的时域表示为：

$ g(x, y)=f(x, y)\mathrm{*}h(x, y)+n(x, y) $

(1)

其中：g（x，y）表示退化图像；f（x，y）表示原始图像；h（x，y）表示引起图像模糊的PSF（退化函数）；n（x，y）表示噪声函数；“*”表示卷积运算。

当$ n(x, y)=0 $时，即图像除产生运动模糊外不受噪声干扰，根据卷积定理，式（1）的傅立叶变换（频域表达）可简化为：

$ G(u, v)=F(u, v)H(u, v) $

(2)

在实际运算时，图像分析是二维离散域的处理过程。设图像尺寸为M×N，则原始图像的傅立叶变换为：

$ F(u, v)=\frac{1}{MN}\sum\limits_{x=0}^{M-1}\sum\limits_{y=0}^{N-1}f(x, y){\mathrm{e}}^{-\mathrm{j}2\mathrm{\pi }\left(\frac{ux}{M}+\frac{vy}{N}\right)} $

(3)

运动模糊图像的PSF主要由模糊角度和模糊尺度这2个重要参数定义。假设在曝光时间T内，目标物体相对于摄像机在θ角（运动方向与水平轴正向的夹角）方向上匀速运动，产生的位移为L，则运动模糊图像的PSF可描述为：

$ h(x, y)=\left\{\begin{array}{ll}\frac{1}{L}, \sqrt{{x}^{2}+{y}^{2}}\le \frac{L}{2}, y=x\mathrm{t}\mathrm{a}\mathrm{n}\theta & \\ 0, \mathrm{其}\mathrm{他}& \end{array}\right. $

(4)

设a和b分别为曝光时间T内像素点在水平和垂直方向上的位移，则h（x，y）的傅立叶变换为：

$ \begin{array}{l}H(u, v)=\iint h(x, y){\mathrm{e}}^{-\mathrm{j}2\mathrm{\pi }(ux+vy)}\mathrm{d}x\mathrm{d}y=\\ \iint \frac{1}{L}{\mathrm{e}}^{-\mathrm{j}2\mathrm{\pi }(ux+vy)}\mathrm{d}x\mathrm{d}y={\int }_{0}^{T}{\mathrm{e}}^{-\mathrm{j}2\mathrm{\pi }\left[ux\right(t)+vy(t\left)\right]}\mathrm{d}t=\\ T\frac{\mathrm{s}\mathrm{i}\mathrm{n}\left(\mathrm{\pi }\right(ua+vb\left)\right)}{\mathrm{\pi }(ua+vb)}{\mathrm{e}}^{-\mathrm{j}\mathrm{\pi }(ua+vb)}\end{array} $

(5)

PSF的离散傅立叶变换为：

$ H(u, v)=\frac{T\mathrm{s}\mathrm{i}\mathrm{n}\left(\mathrm{\pi }\right(ua/N+vb/M\left)\right)}{\mathrm{\pi }(ua/N+vb/M)}{\mathrm{e}}^{-\mathrm{j}\mathrm{\pi }(ua/N+vb/M)} $

(6)

其中：u=0，1，…，M-1；v=0，1，…，N-1。

$ \left|H(u, v)\right|=T\frac{\mathrm{s}\mathrm{i}\mathrm{n}\left(\mathrm{\pi }\right(ua/N+vb/M\left)\right)}{\mathrm{\pi }(ua/N+vb/M)} $

(7)

由式（7）可知，$ \left|H(u, v)\right| $的表达式中包含sin函数，由sin函数性质可知，$ \left|H(u, v)\right| $在（ua/N+vb/M）=0时有最大值，$ \left|H(u, v)\right| $在（ua/N+vb/M）取非零值时为零。因此，模糊图像的频域上会呈现以下特征：模糊图像频谱图上呈现出亮暗相间的平行条纹，且条纹方向随着模糊角度的变化而相应改变。

设l_u和l_v分别为运动模糊长度L在垂直和水平2个方向上的分量，则：

$ L=\sqrt{{l}_{u}^{2}+{l}_{v}^{2}} \text{，} \mathrm{t}\mathrm{a}\mathrm{n}\theta =\frac{{l}_{u}}{{l}_{v}} $

(8)

设亮暗条纹相对水平轴正向的倾斜角度为φ，则：

$ \mathrm{t}\mathrm{a}\mathrm{n}\varphi =-\frac{{l}_{u}N}{{l}_{v}M} $

(9)

可得条纹的倾斜角度φ与运动模糊角度θ间的关系为：

$ \mathrm{t}\mathrm{a}\mathrm{n}\theta =-\mathrm{t}\mathrm{a}\mathrm{n}\varphi \cdot \frac{M}{N} $

(10)

由式（10）可知，通过检测模糊图像频谱图像中的亮暗条纹倾角φ可以得出运动模糊角度θ。设频谱图中相邻条纹的间距为d，则频谱图像中心第一个条纹的方程式为：

$ u\frac{{l}_{u}}{M}+v\frac{{l}_{v}}{N}=1 $

(11)

$ d=\frac{1}{\sqrt{\left({l}_{u}/M\right)+\left({l}_{v}/N\right)}} $

(12)

由式（8）、式（11）、式（12）推导可得：

$ L=\frac{M}{d\sqrt{\mathrm{c}\mathrm{o}{\mathrm{s}}^{2}\theta +（M/N{）}^{2}\mathrm{s}\mathrm{i}{\mathrm{n}}^{2}\theta }} $

(13)

由式（13）可知，当图像尺寸一定时，条纹间距d与模糊长度L呈反比，随着模糊长度的增加，条纹间距减小，条纹数量变多；反之，条纹间距增大，条纹数量变少。

3 运动模糊图像PSF参数估计方法

根据上文分析可知，通过运动模糊图像的频域周期性特征与模糊参数间的关系，可以估算图像的PSF参数。但是，由于存在外界限制和环境干扰，对运动模糊图像直接进行傅立叶变换获得的频谱图像往往会存在大量毛刺、噪点以及十字亮线干扰。为了解决这一问题，现有方法一般通过滤波或频谱分块的方式来消除十字亮线和冗余信息干扰，但是，这种方式存在重要信息丢失和残留干扰的问题。现有PSF参数估计方法主要存在以下2个方面的局限性：

1）根据Radon变换原理，在检测中心亮条纹获得条纹倾角时会因亮条纹本身的宽度而导致角度误差，如图 3（a）所示，从而使得参数估计结果不够精确。

2）由式（13）可知，当运动模糊图像的模糊长度变小时，相应地，频谱图中的条纹间距变大，条纹数量变少，此时对频谱图像进行分块得出的块图可能没有足够可用的条纹信息可供Radon变换检测，如图 3（b）所示。

综上所述，现有PSF参数估计方法在消除十字亮线和噪点、毛刺干扰时，出现了重要信息丢失和残留干扰的问题，因此，本文提出一种改进的PSF参数估计方法，该方法的具体步骤如下：

1）获取运动模糊图像的频谱图。图 4（a）所示为模糊长度L=30、模糊角度θ=120°的运动模糊图像，图像分辨率大小为720×720像素。对图 4（a）作傅立叶变换后得到其频谱图，将频谱图像动态范围压缩且居中，可以直观地观测到模糊图像的频谱条纹特征，如图 4（b）所示。

2）频谱图像增强处理。为了突出条纹特征，对频谱图像进行直方均衡化（如图 4（c））后再将频谱图像二值化。通过直方图均衡化拉伸像素值，增强亮暗条纹的对比度，突出条纹特征，便于后续条纹检测。然后，根据图像整体灰度值和最大类间方差法设定阈值^[20]，对频谱图像进行灰度阈值转换，将灰度图像转换为二值图像，如图 4（d）所示。二值图像的冗余信息更少，同时保留了必要的条纹结构，亮暗条纹的特征也更明显^[21]。

3）形态学变换。为了消除二值图像中十字亮线以及噪点、边缘毛刺等干扰，保留形态合适的条纹图像以供参数估计检测，对二值图像进行形态学开运算，如图 4（e）所示。形态学开运算可以去除比结构元素更小的明亮部分，为了控制图像开运算操作的精度以保留形态合适的中心亮条纹，需要选择合适的结构元素作形态学开运算^[22]。频谱图中的亮条纹形状大致呈长方形，因此，将结构元素的形状类型定义为矩形。同时，为了保证算法的鲁棒性，根据二值化图像自适应选择合适的结构元素大小，设置程序以二值化图像中心为原点，从上、下、左、右、右上、右下、左上、左下这8个不同方向寻找相邻的暗条纹，将找到的图像中心点与暗条纹间的最短距离作为定义结构元素大小的基准尺度。搜索相邻条纹的方式如图 5所示。

4）条纹图像边缘测定。为了避免亮条纹自身宽度导致的角度检测误差（图 3（b）），对开运算后的图像进行边缘检测。具体操作是对测定图像分别进行一次膨胀和腐蚀操作，然后将膨胀后的图像与腐蚀后的图像相减，获取到的边缘图像如图 4（f）所示。相比使用Canny、Soble等微分算子对图像进行边缘检测的方法，形态学方法对边缘方向性的依赖更小，能检测出更加完整的条纹边缘。

5）运动模糊角度估算。Radon变换是一种经典的特征提取技术，对图像作Radon变换是在特定方向上对图像沿直线作线性积分运算，通过这种方式将包含线条的图像变换为相应的含有峰值的参数场。对条纹边缘线条图像（图 4（f））进行1°~180°的Radon变换，提取出条纹倾角φ后由式（10）计算出模糊角度θ。图 4（g）和图 4（h）分别为直接检测频谱和本文方法的Radon变换结果。图 4（g）中亮点面积较大，Radon变换的角度检测结果存在较大的波动（大概在118°~122°），与直接检测的方法相比，本文方法的角度估算结果更精确。

6）运动模糊长度估算。为了降低模糊长度估算的误差，将图 4（b）所示的频谱图像按估算出的角度θ顺时针旋转后作垂直投影，可得到图 6所示的频谱图像垂直投影图，此时得出的垂直投影曲线在各个峰上的峰值点分别对应频谱图像中各条纹所在位置，投影图中相邻负峰值点间的像素距离即为暗条纹间距。为了保证测算精度，根据频谱投影结果，选取投影图中与中心主峰左右相邻的2个负峰值点间的总距离作为间距d值的估计度量（中心主峰值点距离为2d，取到的总距离为4d），求出条纹间距d后代入式（13）可得模糊长度L。

利用上述方法估算出模糊参数后，根据式（4）构造PSF，然后通过维纳滤波去卷积算法复原模糊图像。

4 实验结果及分析 4.1 PSF参数估计实验

为了验证本文PSF参数估计方法的有效性，将其与MOGHADDAM等提出的经典方法^[11]、胡发焕等^[14]提出的改进方法进行对比实验，所有算法均使用MATLAB R2018b编程，计算机操作系统为Windows10，CPU为Intel Core i5-5200U 2.20 GHz，内存为8 GB。

以图 4（a）为标准测试图像，该图像来源于公开图像数据集GOPRO^[23]，该数据集包括真实拍摄的1 111张模糊图像。PSF参数估计实验的模糊参数设置参照文献[14]，将模糊尺度L设置为10 px、20 px和30 px这3组，模糊角度θ在0°~180°范围内每间隔15°取一个测试角度。表 1所示为3种对比方法对同一运动模糊图像的PSF参数估计结果。

从表 1中的参数估计精确度对比结果可以看出，在不同的模糊条件下，本文方法的最大误差和平均误差均小于其他2种估计方法。本文利用二值频谱图像条纹特征自适应地控制形态学运算精度，尽可能地消除十字亮线和冗余噪声干扰，并且通过边缘测定来避免因中心条纹自身宽度导致的角度估计误差，此外，由于本文方法未对频谱进行分块处理，在模糊长度较小时处理后的待估频谱图像中仍保留足够的条纹线条，在进行模糊长度估算时对频谱图像进行直方图均衡处理，提高了亮暗条纹的对比度，由此检测的相邻条纹间距更加准确，从而保证了模糊长度估算精度。因此，本文方法在不同大小的运动模糊图像中都能得到更精确的估计结果。

在实验中发现，在模糊长度分别为10 px、20 px、30 px的3组参数估计实验中，当L=10 px时，参数估计结果的误差波动明显，且误差主要集中在模糊角度为75°、90°和105°附近。分析原因有以下2点：1）中心十字亮线的2条线分别呈180°和90°，当十字亮线效应严重时，残留部分会干扰相邻角度的检测结果；2）当模糊长度较小时，一些暗条纹会被中央亮条纹掩盖，从而导致对2个相邻暗条纹间的距离估计误差变大。

图 7所示为不同PSF参数估计方法在模糊长度为10 px时的模糊角度估计误差。从中可以看出：文献[11]方法的参数估计结果开始出现极大偏离值，参数估计结果绝度误差超过1，估计结果基本超出可信的有效范围；文献[14]方法的估计结果准确度明显优于文献[11]方法，误差波动也较小，但是在模糊角度为75°、90°和105°附近时，估计结果绝度误差仍大于0.5，影响了估计结果的可信度；本文方法在不同大小的运动模糊图像中都能得到更精确的估计结果，参数估计有效范围可达0~180°且误差更小，在参数估计误差波动明显的情况下仍能保证估计结果的可信度。总体而言，本文PSF参数估计方法在精度和有效估计范围方面均能取得良好的性能表现。

4.2 图像复原实验

维纳滤波法是一种经典且简单高效的线性复原算法，其以原始图像和复原图像之间的均方误差最小为约束条件对原始图像进行求解。根据PSF参数估计结果可以构造出运动模糊图像的点扩散函数，当确定了模糊图像的点扩散函数后，采用维纳滤波器可以快速复原模糊图像。维纳滤波器在对原始图像进行函数求解时表现为如下形式：

$ {e}^{2}=\mathrm{m}\mathrm{i}\mathrm{n}E\left\{\right[f(x, y)-g{(x, y)]}^{2}\} $

(14)

其中：e²为统计误差；E表示数学期望。

根据式（14）及图像退化模型，维纳滤波器可表示为：

$ F(u, v)=\frac{1}{H(u, v)}\frac{\left|H\right(u, v){|}^{2}}{\left|H\right(u, v){|}^{2}+S}G(u, v) $

(15)

其中：S为噪声功率谱与图像功率谱的比值。

图像复原结果可以直观地体现PSF参数估计方法的准确性。为了进一步验证本文方法的可信度以及泛化性能，在测试图像中加入GOPRO数据集中的更多图像，图像分辨率为720×720像素，运动模糊核为随机大小。图 8所示为不同方法所构建的点扩散函数结合维纳滤波的运动模糊图像复原实验结果，为了便于描述，将图 8中的图像由上至下分别命名为“汽车”“花纹”“人群”和“行人”。

对比图 8（c）、图 8（d）与图 8（e）可以看出，本文方法得出了更接近原始图像的清晰图像，图 8（c）、图 8（d）的一些复原图像中不仅没有完全消除模糊效果，且有严重的振铃效应和伪影残留，相比之下，本文方法的复原图像伪影残留更少，振铃效应更小，保留了更多的图像细节。

4.3 图像质量评价

为了客观评价复原图像的质量，本文引入峰值信噪比（Peak Signal-to-Noise Ratio，PSNR）和结构相似性（Structural Similarity，SSIM）2种客观评价指标。PSNR和SSIM的值越大，表明图像的保真度越高。PSNR根据图像对应像素点间的误差来进行图像质量评价，定义式为：

$ {P}_{\mathrm{P}\mathrm{S}\mathrm{N}\mathrm{R}}=10\mathrm{l}\mathrm{g}\left(\frac{{255}^{2}}{{M}_{\mathrm{M}\mathrm{S}\mathrm{E}}}\right) $

(16)

$ {M}_{\mathrm{M}\mathrm{S}\mathrm{E}}=\frac{1}{MN}\sum\limits_{y=1}^{M}\sum\limits_{y=1}^{N}{\left(I(x, y)-\widehat{I}(x, y)\right)}^{2} $

(17)

其中：M_MSE表示2幅图像间的均方误差；$ I(x, y) $和$ \widehat{I}(x, y) $分别是图像复原前后在像素点（x，y）处的灰度值；M、N分别为图像长度与宽度。

SSIM从亮度（Luminance）、对比度（Contrast）和结构（Structure）三方面度量图像间的相似性，相比PSNR等传统图像质量评价指标，SSIM更加符合人眼对图像的视觉感知。SSIM的取值范围为-1~1，定义式如下：

$ \mathrm{S}\mathrm{S}\mathrm{I}\mathrm{M}(I, J)=\frac{(2{\mu }_{I}{\mu }_{J}+{C}_{1})(2{\sigma }_{IJ}+{C}_{2})}{({\mu }_{I}^{2}+{\mu }_{J}^{2}+{C}_{1})({\sigma }_{I}^{2}+{\sigma }_{J}^{2}+{C}_{2})} $

(18)

其中：I、J为输入图像；$ \mu $为图像的均值；$ \sigma $为图像的标准差；$ {\sigma }_{IJ} $为2幅图像的协方差；C₁=（K₁L）²、C₂=（K₂L）²为常数；L为像素值的范围。本文设置K₁ =0.01，K₂ =0.03。

表 2所示为图 9中4个不同的运动模糊图像在不同方法下的复原图像质量评价结果，从中可以看出，本文方法取得了最大的PSNR和SSIM值，表明在PSNR和SSIM客观评价指标上，本文方法结合维纳滤波复原出的图像质量最优，证明本文方法精度更高，且对不同特征、不同模糊程度的图像的PSF参数估计结果稳定保持在可信范围内。

5 结束语

成像时相机抖动或目标物移动等会使图像产生运动模糊，从而导致图像清晰度大幅下降。为了提高运动模糊图像PSF参数估计结果的精度，本文提出一种改进的运动模糊参数估计方法，并利用该方法估算出的参数来复原运动模糊图像。PSF参数估计和图像复原实验结果表明，相比文献[11, 14]所提方法，该方法具有有效性且估计结果精度更高，有效估计范围更广，在相同条件下能复原出具有更少振铃和伪影的图像。在实际应用场景中，除了相机和场景之间产生线性运动造成的运动模糊外，往往还存在其他类型的模糊和噪声干扰因素，因此，下一步将寻找更合适的数学模型，以解决含有多个非线性模糊或大量噪声干扰的复杂图像复原问题。