0 引言

图像修复是计算机视觉和图形图像处理中的一个研究热点，被广泛应用于古文物的保护、影视特技的制作、目标物体的移除、视频图像差错隐藏等视觉处理领域。图像修复技术主要是对图像中丢失或者损坏的部分，根据图像自身信息，按照一定的修复原则补全其损坏的图像信息，使得处理后的图像接近或达到人们视觉上的连通性^[1]。

图像修复问题主要分为两类:一类是用于修复小尺度缺损的数字图像修复，比如划痕、折痕或者污点，主要的解决思路为基于变分非线性偏微分方程(Partial Differential Equation, PDE)的图像修复模型, 即最早由Bertalmio等^[2]提出的基于偏微分方程的图像修复模型, 同时利用三阶PDE模型提出了BSCB(Bertalmio-Sapiro-Caselles-Bellester)模型。该模型是模拟专业人员的手工修复方式，将破损区域附近的图像信息，按照破损边缘等照度线的方向，向破损区域内部进行扩散。之后，延续BSCB的算法思想，Chan等^[3]提出了改进的基于曲率驱动扩散(Curvature Driven Diffusion, CDD)模型，该模型通过引入曲率的概念改善了全变分(Total Variation, TV)模型。但是以上基于PDE的修复模型对于非纹理图像有较好的修复效果，对于具有大破损区域的图像进行修复时常常出现结构线断裂和图像模糊的现象。另一类是用于填充补全图像中大破损区域的图像修复，比如目标移除，主要的解决思路为基于纹理合成的图像修复、基于稀疏表示的图像修复、基于图像分解的图像修复以及基于参照样例的图像修复。Bornard等^[4]提出了以像素点为合成单元的纹理合成修复算法，具有显著的局部性，在复杂纹理区域中容易导致无法匹配的问题。因此，Criminisi等^[5]提出了一种基于优先级的图像修复算法，通过优先级计算出破损区域边缘最先修复的待修复块，根据待修复块中未破损区域信息，在整幅图像的完好区域寻找匹配的样本块，利用最优样本块填充待修复块，完成修复。然而，在修复过程中因为优先级的计算顺序和样本块搜索的问题，导致纹理的错误延伸和结构线断裂。之后，Sun等^[6]提出通过人工手绘使一些重要的曲线或者线段从完好区域延伸至破损区域，然后利用完好区域结构线附近的图像信息完成填充。该算法对于修复结构边缘图像具有良好的效果。Elad等^[7]以信号的稀疏模型为基础，对图像的纹理和结构部分分别进行稀疏表示, 然后进行重构。Xu等^[8]提出一种基于结构稀疏的传播方法, 降低了块传播的贪婪风险。Bertalmio等^[9]是将图像分解为结构和纹理后各自进行修复。Tang等^[10]提出了一种基于参照图像的图像补全算法，该算法首先从图像库中选取最佳参照图像，再利用参照图像中的结构信息重建原始图像中的缺失结构，最后通过纹理合成完成图像修复。以上将纹理和结构部分分开进行修复的算法，虽然能避免结构线断裂的产生，但是容易导致样本块的错误累积和块效应现象。Liu等^[11]提出一种新的基于样本的图像修复算法，该算法通过引入结构张量理论构造局部测量功能, 获得最优匹配块填充破损区域从而修复图像。近三年，提出了一种基于统计的图像修复算法，He等^[12]采用相似块统计的修复方法，通过相似块的匹配并且获得它们的相对位置参数，即偏移量，根据偏移量的统计结果为修复破损区域提供可靠信息。该算法虽然能够提高修复效率，但是具有一定的局限性，适用于相似块多的图像。

针对图像修复中出现的纹理块的错误延伸、错误累积、结构线断裂以及修复效率等问题，本文提出一种通过划分特征子区域，局部寻找最佳样本块的修复算法。所提算法主要针对含有丰富纹理信息和复杂结构信息的大破损区域的图像进行修复，其主要的解决思路是以纹理合成中Criminisi算法为基础，Criminisi算法最主要是求破损区域边缘上的所有待修复块的优先级值和在完好区域中搜索填充待修复块的样本块。Criminisi算法中优先级值的计算方法简单，但是用该优先级的计算公式并不能很好地突出结构信息，会导致修复后的图像出现结构断裂现象; 而且采用单个样本块复制来填充待修复区域会造成已修复区域出现单个块重复复制，并容易产生错误匹配; 以及Criminisi算法的搜索范围是全局搜索，使得该算法的运行速度较慢。本文在Criminisi算法的基础上对优先级和选取样本块上都作出了改进，并且根据文献[13-15]所提到的有关图像特征提取的描述，可知图像的不同区域包含不同的图像特征并且可以进行特征的提取。本文就是根据把不同的图像特征进行提取，用提取后的特征值进行分类，划分出具有不同特征的子区域。

1 算法描述

本文所提的划分特征子区域算法主要思想是通过划分特征子区域，利用约束条件在各个对应的特征子区域中搜索样本块集，根据权值分配法最终合成一个最佳样本块。该算法主要分为两个部分:第一个部分是预处理部分，用于划分特征子区域; 第二个部分为修复部分，用于对破损区域进行修复。该算法具体步骤如下：

第一部分：

1) 将破损图像划分成小图像块Patch，并分别对每个图像块进行特征提取。

2) 对提取的特征值进行统计分类，根据不同特征集对图像进行子区域划分。

第二部分：

1) 对破损区域的边缘待修复块用带结构因子的优先级计算公式计算其优先级顺序，并筛选出最高优先级的待修复的目标块。

2) 在与该目标块相对应的子区域中的完好部分进行局部搜索构造候选样本块集。先根据目标块与样本块自身约束条件选取初始候选样本块集；再根据目标块的邻域窗口中确定的最佳相似块与样本块之间的约束条件确定最佳样本块集；最后，对最佳样本块集中的所有样本块进行权值分配从而合成一个最佳样本块。

3) 把最佳样本块复制至目标块，并沿着新形成的边缘更新各像素点的优先级，重复第二部分的各个步骤直至破损区域全部修复完成。

2 修复方法 2.1 预处理——划分特征子区域

本文进行预处理的目的是通过图像特征分类划分特征子区域局部搜索候选样本块集，不仅提升了图像修复的速度，而且提高了其匹配精度。对一幅破损图像进行预处理，该子区域划分示意图如图 1所示。

首先，把输入的破损图像I划分成大小为n×n的小图像块Patch(如图 1(a)中虚线所分割出来的Patch)。

其次，本文选取Sobel算子提取边缘结构特征，选取局部二值模式(Local Binary Pattern, LBP)算子提取纹理特征，将Sobel算子与LBP算子作为卷积核与原始图像中的各个图像块Patch进行卷积运算从而得到特征矩阵，特征公式如下：

$ I\left( {x, y} \right)*\omega \left( {x, y} \right) \equiv \sum\limits_{s = 0}^a {\sum\limits_{t = 0}^b {\omega \left( {s, t} \right)} } I\left( {x - s, y - t} \right) $

(1)

其中：I(x, y)为图像上(x, y)处的像素值；ω(x, y)为卷积核；a、b分别为图像块Patch的宽度和长度。然后，采用Max-pooling技术对特征矩阵进行处理从而获得显著特征。

最后，将每个图像块的显著特征矩阵进行归一化计算出其特征值，并根据每个图像块的特征值与位置，对特征值进行分类统计。根据图像的局部相似性原理^[16]可知，相同或者相邻区域的特征值差距较小，而不同特征区域的特征值差别较大。那么，把特征值差距较小且位置相邻的特征值归为一类，组成一个特征集，从而依照不同的特征集把原图像划分为不同的子区域(如图 1(b)所示)，并区分为结构类子区域(S₁, S₂, S₃)和纹理类子区域(C₁, C₂, C₃, C₄, C₅)。

2.2 破损区域修复

本文从优先级计算和样本块选取两个方面对Criminisi算法^[17]进行了改进从而构造了一种新的基于样本块的修复方法。在优先级计算中引入了结构因子，增强结构信息的作用，避免修复后的图像出现结构断裂现象。在样本块选取中，本文是先对破损图像进行子区域的划分，把原算法的全局搜索改进为仅在对应的子区域的局部搜索。在选取候选样本时，先在目标块所对应的子区域中选取候选样本块组成初始候选样本块集；再在目标块的邻域中确定其最佳邻域相似块，并根据该最佳邻域相似块对初始候选样本块集进行筛选从而进一步缩小样本块集的大小并得到最佳样本块集。最后，为最佳样本块集里面的所有样本块分配不同的权值, 并合成这些样本块得到一个最佳样本块来填充目标块。

2.2.1 改进的优先级计算

优先级的计算主要是计算破损区域(待修复区域)中各个目标块(待修复块)的优先级，该目标块是一个通过选择破损区域边缘上的一个中心像素点作为中心，并且包含破损区域和完好区域信息的矩阵块。破损区域边缘上的某个目标块示意图如图 2所示。图 2中，I为一幅破损图像，Ω为破损区域(待修复区域)，∂Ω为破损区域的边缘，Ω为完好区域。假设Φ_P为以P点为中心、位于破损区域边缘∂Ω上的目标块(待修复块)，W_P是目标块的邻域窗口，其大小为图像块Patch大小的3倍，其已知像素区域为：

$ {\mathit{\boldsymbol{W}}_P}\left( q \right) = \left\{ {{q_i}:{q_i} \in {\mathit{\boldsymbol{W}}_p} \cap \mathit{\overline \Omega } } \right\} $

(2)

改进后的优先级计算公式为：

$ P' = {\alpha _W} \cdot {P'_P} $

(3)

其中, α_W为结构因子，是一个比值，其分子分母分别表示邻域窗口中结构类子区域和纹理类子区域中的已知像素，它的计算式为：

$ {\alpha _W} \equiv \frac{{\sum {\left[{{\mathit{\boldsymbol{W}}_P}\left( q \right) \cap {S_n}} \right]} }}{{\sum {\left[{{\mathit{\boldsymbol{W}}_P}\left( q \right) \cap {C_n}} \right]} }} + 1 $

(4)

其中：W_P(q)表示邻域窗口的已知像素点；S_n表示结构类子区域；C_n表示纹理类子区域。结构因子的分母表示邻域窗口内的已知像素点属于纹理类子区域的像素点个数之和，分子表示邻域窗口内的已知像素点属于结构类子区域的像素点个数之和。当分子为零时，说明该目标块的邻域窗口中不包含结构信息，为了防止其最后的优先级为零，在结构因子式子后加1，确保对其自身目标块的优先级不影响。结构因子用来放大结构信息在优先级中的重要性，如果在目标块邻域中的结构信息越丰富，那么该目标块自身的结构信息也丰富，相应地它就具有较高的优先级。这样保证了含有丰富结构信息的目标块可以得到优先修复，避免结构特征的断裂。

目标块P的原始优先级计算公式为：

$ {P'_P} = {C_P} \cdot {D_P} = \frac{{\sum {q \in {\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_p} \cap \mathit{\overline \Omega } } }}{{\left| {{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_P}} \right|}} \cdot \frac{{\left| {\nabla {I^ \bot } \cdot \eta } \right|}}{\alpha } $

(5)

式中包含两个计算项：C_P为置信度项，是对P点周围已知信息进行统计，P点周围的已知信息越多其优先级越大；D_P为数据项，表示该修复应该尽可能兼顾等照度方向和边缘法线方向，等照度线方向与边缘法线方向两者夹角越小，优先级越大。

2.2.2 合成最佳样本块

样本块的选取最主要是为了找出最符合目标块信息的样本块，本文构造了一种新的选取方法从而确定其最佳样本块集，某个目标块确定其最佳样本块集的示意图如图 3所示。

首先，通过优先级计算确定出边缘目标块的修复顺序，假设Φ_P是最优先修复的目标块，在它所对应的子区域中进行局部搜索得到其初始候选样本块集A={ψ_q}_q=1^N (如图 3(a)所示)。该集合A是由目标块自身与样本块的相似性确定，其计算公式为离差平方和(Sum of Squares of Deviations, SSD) ^[17]：

$ \begin{array}{l} d\left( {{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_p}, {\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_q}} \right) = \sum\limits_{a \in {\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_p}, b \in {\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_q}} {\left[{{{\left( {R\left( a \right)-R\left( b \right)} \right)}^2} + } \right.} \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\left. {{{\left( {G\left( a \right)-G\left( b \right)} \right)}^2} + {{\left( {B\left( a \right)-B\left( b \right)} \right)}^2}} \right] \end{array} $

(6)

其中，R(a)、G(a)、B(a)为三原色的三个通道。当两块之间的相似性小于某个阈值ε时(如式(7))，那么，就将该样本块作为候选样本块并将它放入初始候选样本块集。

$ \left| {d\left( {{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_p}, {\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_q}} \right)} \right| < \varepsilon $

(7)

其次，在该目标块的邻域窗口中确定一个最佳邻域相似块Φ′_P (如图 3(b)所示)，通过计算该目标块与邻域中的各相似块之间的离差平方和(如式(6)), 并选出离差平方和最小的相似块作为最佳邻域相似块，计算式如下：

$ {\mathit{\boldsymbol{ \boldsymbol{\varPhi} '}}_p} = \mathop {{\rm{arg}}\;{\rm{min}}}\limits_{{{\mathit{\boldsymbol{ \boldsymbol{\varPsi '} }}}_q} \in {\mathit{\boldsymbol{W}}_P}} \left( {{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_P}, {{\mathit{\boldsymbol{ \boldsymbol{\varPsi '} }}}_q}} \right) $

(8)

第三，根据最佳邻域相似块对初始候选样本块集A进行再次筛选以确定最佳样本块集B={ψ′_qM}_q=1^M (如图 3(c)所示)。根据图像的局部相似性原理可知，初始候选样本块集A中的各个样本块与最佳邻域相似块之间的差别不应该太大，比如小于阈值ε，如果两者之间的差距大于ε就将对应的候选样本块从初始候选样本块集中删除。经过这样的筛选步骤之后，就可以将初始候选样本块集压缩形成最终的最佳样本块集。

然后，根据权值分配法合成一个最佳样本块复制至目标块进行修复。根据权值分配公式合成一个最佳样本块$ {\mathit{\boldsymbol{ \boldsymbol{\hat \varPhi} }}_P}$，计算公式为：

$ {\mathit{\boldsymbol{ \hat \varPhi} }_P} = \frac{{{\beta _1}{{\mathit{\boldsymbol{ \boldsymbol{\varPsi '} }}}_{{q_1}}} + {\beta _2}{{\mathit{\boldsymbol{ \boldsymbol{\varPsi '} }}}_{{q_2}}} + ... + {\beta _M}{{\mathit{\boldsymbol{ \boldsymbol{\varPsi '} }}}_{{q_M}}}}}{M} $

(9)

其中：ψ′_qM为最佳样本块集B中的系列样本块；β为权值，它由样本块ψ′_qM与目标块Φ_P的相似度来确定，相似度采用离差平方和(SSD)方法进行求得，即当某一样本块ψ′_qi与目标块的SSD值越小，则该样本块分配的权值就越大，相反分配权值越小。权值的计算公式为：

$ {\beta _j} = \frac{{\sum\limits_{i = 1}^M {d\left( {{\mathit{\boldsymbol{ \varPhi }}_p}, {{\mathit{\boldsymbol{ \boldsymbol{\varPsi '} }}}_{{q_j}}}} \right)} }}{{d\left( {{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_p}, {{\mathit{\boldsymbol{ \boldsymbol{\varPsi '} }}}_{{q_j}}}} \right)}}{\rm{;}}\;\;\;\;\mathit{j} = 1, {\rm{ }}2, \ldots {\rm{ }}, \mathit{M} $

(10)

式中：分子是一个定值，代表了目标块与最佳样本块集中的各个样本块之间的相似度之和；分母表示目标块与某一样本块之间的相似度。

最后，把得到的最佳样本块复制至目标块，并沿着新形成的边缘更新各像素点的优先级，直至破损区域全部被修复。

3 实验结果与分析

为了验证所提算法的性能，进行以下实验仿真。实验环境说明：软件环境是操作系统为Windows 7 SP1，64位，编译环境为Visual Studio 2015;硬件环境是CPU为Intel赛扬处理器，主频为1.8 GHz，内存为4 GB。通过四组实验，将Criminisi算法、文献[8]算法和本文算法进行了主观和客观的对比。

3.1 主观对比

主观对比的结果常常因为评价者和评价内容的不同而产生较大差异，本文根据实验得到的图像进行主观对比，即针对同一幅破损图像使用不同的算法分别完成修复，对比得到的修复效果，不同图像的修复结果如图 4所示。预先设置算法的各项参数：小图像块Patch值、目标块矩阵值、阈值ε。

图 4(a)是经典蹦极图像，该修复主要是对蹦极人进行了移除。对比图 4(a)的结果可知：Criminisi算法修复的图像结果可以看到房顶处和水中都出现了错误的纹理延伸；文献[8]算法修复的图像房顶处出现的错误的纹理延伸较少；而本文的算法进行修复的图像克服了匹配上的错误，使修复的图像效果更好。

图 4(b)是对沙滩上的小狗进行移除。对比图 4(b)的结果可知：Criminisi算法修复的图像结果在对该小狗进行移除时把栈道上的纹理错误地匹配为了天空处的纹理块；而文献[8]算法修复的图像把栈道上鞋子的纹理块匹配为了天空上的纹理块，使得鞋子的图像信息被覆盖；但是本文算法通过子区域的划分，把匹配块的选取精确为同一个区域，因此克服了该错误，使视觉上达到良好的效果。

图 4(c)是把大象作为一个大破损区域进行修复。对比图 4(c)的结果可知：Criminisi算法修复的图像把草地的纹理错误的延伸至了水中；而文献[8]算法修复的图像把天空的纹理也错误延伸至了水中。两个算法所修复的图像纹理块都出现了错误传播，不符合实际情况。而本文算法通过对选取的匹配块进行双重约束，使匹配的纹理块更加准确，避免了该错误，使得修复效果更好，更加符合实际情况。

图 4(d)是对小雏菊上破损区域进行修复。对比图 4(d)的结果可知：Criminisi算法修复的图像小雏菊的花瓣上出现了严重的结构断裂现象；文献[8]算法修复的图像结构断裂的部分修复相对完整些，但是出现了冗余的结构块使得结构断裂部分重叠；而本文算法增加了对结构信息作用的放大，减少了图像上的结构断裂现象，并且解决了结构断裂部分重叠的问题，使视觉上达到良好的效果。

3.2 客观对比。

客观评价方法是指采用合理的客观评价算子对修复质量进行评价。修复算法的执行时间代表一个算法的效率高低，是一个很容易得到和量化的数据。本文选择峰值信噪比(Peak Signal-to-Noise Ratio, PSNR)作为客观评价算子对修复质量进行评价, 其计算公式为：

$ PSNR = 10\log \frac{{{{255}^2}}}{{\frac{1}{{M \times N}}\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^N {{{\left[{g\left( {i, j} \right)-u\left( {i, j} \right)} \right]}^2}} } }} $

(11)

其中：u(i, j)表示原始图像；图像的尺寸大小为M×N；g(i, j)表示修复后的图像。

本文用峰值信噪比(PSNR)(如表 1所示)和修复时间(如表 2所示)作为客观评价指标。从表 1中可以看出，本文算法的PSNR值大于Criminisi算法，说明本文算法所修复的图像失真更小。本文算法的PSNR值同样大于文献[8]算法，除了图 4(d)的修复结果，说明用本文算法所修复的图像虽然对具有过多的结构断裂的图像修复效果稍次于文献[8]算法，但是不会出现结构断裂部分重叠的问题，更加符合人们的视觉要求。从表 2中可以看出，本文算法的修复时间是最短的，明显少于Criminisi算法和文献[8]中的算法，说明本文算法的修复效率高。

综合主观和客观对比所得结果，本文所提的划分特征子区域的图像修复算法，无论是针对目标物体的移除还是大破损区域的修复，不仅能够基本不丧失修复精度，而且能够提高算法效率，普适性强。

4 结语

本文通过预处理进行了特征子区域的划分，对Criminisi算法中的优先级和样本块的选取进行了改进，提出了一种划分特征子区域的图像修复算法。所提算法有效地对结构和纹理信息进行区分，解决了Criminisi算法容易造成错误的修复顺序的问题，提高了匹配的精确度，解决了错误匹配的问题。通过上述实验主观和客观上的对比可以看出，本文算法普遍适用于修复具有大破损区域的图像，通过与Criminisi算法和文献[8]中的算法进行比较，本文算法的修复效果较为优秀，修复效率提升明显，并且修复后的图像也更加符合人们的视觉上的连通性。但是本文算法对于结构信息特别复杂的图像(如图 4(d)所示)，虽然修复效果比Criminisi算法好，但是和文献[8]中的算法比较, 本文算法效果较差。因此，对于这种结构线过多的断裂，本文算法并不能完全修复，需要近一步研究如何修复具有过多结构线断裂图像的方法。