一种基于CUDA的快速图像类比合成方法

摘要

一种基于CUDA的快速图像类比合成方法，包括如下步骤：1)生成源图像A，滤波后的图像A’，源图像B的高斯金字塔；2)当前合成第L层，对金字塔从粗糙层到精细层执行以下过程：2.1)对B′

说明书

技术领域

本发明涉及图像处理领域，尤其是一种图像类比合成方法。

背景技术

纹理合成可以从小区域纹理样本自动生成视觉上相似和连续的大样本纹理。 基于样图的纹理合成方法是继纹理映射、过程纹理合成等方法后发展起来的一种 新的纹理拼贴技术，避免了以往方法中的接缝、走样、反复调整参数等问题。目 前，对基于样图的纹理合成的研究可大致分为二方面。一方面，希望在合成中加 入更多控制因素来约束纹理生成，产生更具多样性的纹理。另一方面，希望加快 合成速度，使合成趋于实时。

图像类比框架借用了类比思想，提供了一种非常自然的方法来定义图像间的 转换，使用户可以控制生成结果，产生多种图像特效。图像类比框架适用于各种 应用，如纹理合成、纹理转移、超分辨率、艺术风格学习，texture-by-number(TBN) 等。传统方法采用两种搜索策略来查找具有最相似邻域的候选点：最相似邻域查 找ANN(approximate-nearest-neighbor)获得全局最优的匹配点，经典的Ashikhmin 算法获得邻域相关的待选点。但为了获得一个匹配点，ANN需要上万次的迭代， 迭代次数随源图像分辨率的提高线性增加，Ashikhmin算法也根据邻域大小需要 进行多次计算才能获得一个待选点。合成一张图像要花费几十分钟甚至更长的时 间。

目前的研究中，研究者采用各种方法来加快合成速度，提高合成质量。如采 用粒子群优化算法来加快最相似邻域查找过程，提高合成速度。对TBN应用， 有研究者提出用Random Walk算法进行最相似邻域搜索，无需预处理，提高搜 索速度。对于艺术风格学习，研究者提出了适用于不同绘画风格的基于笔刷的算 法，可以生成油画、水粉画、钢笔画等风格效果，这些算法大多是模拟一种特定 的绘画材质，有的还需要人工交互。还有研究将纹理优化技术用于图像类比合成 来提高合成质量，但合成速度还是较慢。

计算机中多核/众核体系架构的引入，使得计算能力成倍增加，GPU的并行 处理能力也越来越强大。Nvidia公司推出的CUDA开发工具使程序员能在C语 言环境下完成对GPU并行处理架构的调用，能更有效的把计算问题映射到硬件 架构上。研究者提出利用GPU来进行并行可控纹理合成，但其主要针对同质纹 理合成，不适用于复合纹理的合成。

现有方法大多针对该框架中的某一种具体应用来加快合成速度或提高合成 质量，还没有一种快速的方法能较高质量的适用于整个图像类比框架的合成。

发明内容

为了克服已有图像类比合成方法的合成速度较慢、合成质量较差的不足，本 发明提供一种提升合成速度、提高合成质量的基于CUDA的快速图像类比合成方 法。

本发明解决其技术问题所采用的技术方案是：

一种基于CUDA的快速图像类比合成方法，所述合成方法包括如下步骤：

1)生成源图像A，滤波后的图像A’，源图像B的高斯金字塔；

2)当前合成第L层，对金字塔从粗糙层到精细层执行以下过程：

2.1)对B′_L+1超采样得到B′_L；

2.2)对B′_L中的每个像素点q，采用基于GPU的全局搜索策略找到全局搜 索最优解q_global；

2.3)对B′_L中的每个像素点q，采用基于GPU的k-coherence搜索策略找 到一致性搜索最优解q_kcoh；

2.4)如果全局搜索得到的结果更好，则将q_Best赋值为q_global，否则将q_Best赋 值为q_kcoh，把q_Best在A′_L中的对应坐标保存在B′_L(q)中；

2.5)循环执行步骤2.2)～2.4)根据设定的次数，得到滤波后的图像B’；

3)对滤波后的图像B’做合成放大，得到高分辨率的EHB’。

进一步，所述步骤2.1)中，将B′_L+1与纹理存储器绑定，每个线程执行的操作 为：得到B′_L(p)的坐标(u，v)，取得B′_L+1中(u/2，v/2)点颜色值color，B′_L(p) 的值为color*2+(u％2，v％2)。

再进一步，所述步骤2.2)中，输入为A，A’，B金字塔L层对应的A_L，A′_L，B_L图 像，输出为B′_L，B′_L中每个像素点的RG通道保存了由全局搜索得到的最相似点 的坐标q_global，同时保留了其对应的距离，对B′_L中的一点q，根据q点在B_L，B′_L的 邻域在A_L，A′_L中找到和其最相似的点。

所述步骤2.2)中，所述全局搜索策略的过程为：设R＝{r₁，r₂，...，r_m}为参考集 合，Q＝{q₁，q₂，...，q_n}为查询集合。对于纹理合成来说，m为样图中的点个数，n 为查询点个数。N()代表邻域向量，N(r)和N(q)为收集图像中邻域得到的向量。公 式(1)表示，对于查询集合Q中的每个点q，在R中根据L₂度量准则找到k个与其最相 似点的集合C_1..k(q)，也就是对于每个点q，在集合R中找到和其邻域距离最小的前 k个点，

$\left(\begin{array}{cc}{C}_{1..k}\left(q\right)=\arg {\min }_{r\in R}{\left|\right|N\left(r\right)-N\left(q\right)\left|\right|}_2& (r\in R,q\in Q)\end{array}\right)---\left(1\right)$

分解该运算，得到公式(2)：

$C_{1..k}\left(q\right)={\left|\right|N\left(q\right)\left|\right|}_2^2+\arg {\min }_{r\in R}({\left|\right|N\left(r\right)\left|\right|}_2^2-2{\left|\right|N\left(r\right)·N\left(q\right)\left|\right|}_2)---\left(2\right)$

从公式(2)中可以看出不依赖于R中点，也不依赖于查询点 q，因此对于q点与R中点比较，只要计算部分即可。收 集邻域向量得到参考集合对应的矩阵RMatrix＝|N(r₁)，N(r₂)，...，N(r_n)|，查询集合 对应的矩阵QMatrix＝|N(q₁)，N(q₂)，...，N(q_n)|。对于全局搜索来说，k取1，RMatrix 根据A，A’中对应的点来获得，QMatrix根据B，B’中对应的点来获得。

利用CUDA架构快速地进行此运算，在CUDA架构下程序被分为host端和 device端，host端执行CPU上的操作，device端执行显卡上的操作，内核函数Kernel 为CUDA运行在device端的并行计算函数，因此对Q进行分块，分批处理查询集合 中的点，动态确定当前可分配的最大Q分块数；对于一个Q分块，计算分块中的点 和参考集合中所有像素点间邻域的距离和坐标，采用插入排序法找到具有最小距 离的前k个点的坐标集合，

kernel函数计算步骤如下：

Step1：收集RMatrix和QMatrix，在device端分配存储空间，把RMatrix、 QMatrix拷贝给device端；在host端分配空间indexMatrix，存储计算得到的最相 似点坐标

Step2：查找显存中剩余空间，确定每次可同时搜索的Q集合大小 maxQueryNum，

Step3：根据maxQueryNum，对Q分块，对Q的每一个分块：

Step3.1：计算距离${\left|\right|N\left(r\right)\left|\right|}_2^2-2{\left|\right|N\left(r\right)·N\left(q\right)\left|\right|}_2$

Step3.2：找到距离最小的前k个最相似点集合C_1..k(q)

Step 3.3将C_1..k(q)中点的坐标拷贝到indexMatrix对应位置的内存中，

所述步骤2.3)中，输入同样为A，A’，B金字塔L层对应的A_L，A′_L，B_L图 像，输出为B′_L。B′_L中每个像素点的RG通道保存了由k-coherence搜索得到 的最相似点的坐标q_kcoh。在用CUDA实现时，线程网格的分割为与图像大小 相对应，每个线程块处理一个像素元素，它拥有的线程数为一个像素需要查找 的coherence候选集个数，通过一个预处理过程，对于A和A’中的每个点， 找到k个与其最相似的点，把这些点的坐标保存在纹理中。

kernel函数计算步骤如下：

Step1.在共享存储器中，声明类型为iht4的distance数组，xy分量存储坐标， z分量存储距离，声明类型为uchar2的kcoh数组，把预处理得到的纹理中像素对 应的坐标载入kcoh，

Step2.每个线程访问共享内存，取出kcoh中的坐标，根据CUDA内建函数 blockIdx.x，blockIdx.y获取待合成像素点位置，基于度量准则，计算两点邻域的 之间的dist

Step3：线程块中的每个线程将计算得到的dist及其对应的样图坐标保存到 distance数组中

Step4：当一个线程块中所有线程都计算完毕时，用树形规约算法找到具有最 小dist的样图坐标，即为所求的最相似点坐标，

采用用树形规约算法搜索最小值，每次比较时，存储在前半部分线程中的数 据和存储在后半部分线程中的数据相比较，取较小值更新前半部分的数据，每次 迭代后减半，直到最后剩下一个值，即为所求的最小值，数据单元之间的比较由 多线程并行执行，在CUDA架构下，使用cudaMallocHost函数直接在主机端内存 中取一块页面锁定内存。

所述步骤2.4)中，根据coherenceFac比较步骤2.2)，2.3)计算得到的结果， 得到q_Best的取值，比较准则如下：若d(q_kcoh，q)＞d(q_global，q)(1+coherencefac*2^L)， 则q_Best＝q_global，否则q_Best＝q_kcoh，

两个邻域的相似性度量准则定义如下，同时比较A，A’和B，B’邻域的特征向量：

d(p，q)＝||F_A(p)-F_B(q)||_L2+||F_A′(p)-F_B′(q)||_L2，

在RGB空间下，选取RGB颜色值作为特征向量；在YIQ空间下，选取亮度Y 作为特征向量，将合成得到的新Y值与原有的IQ值组合还原成RGB图像；对于艺术 风格学习应用来说，使用亮度分布的平均值和方差对图像的Y通道进行线性变换， 对原有的YIQ空间进行转化(remap)：

$Y\left(p\right)\Leftarrow \frac{{\sigma }_B}{{\sigma }_A}(Y\left(p\right)-{\mu }_A)+{\mu }_B$

其中Y(P)是图像A中P点的亮度(Y)值，μ_A和μ_B是图像A的亮度(Y)平均值，σ_A和σ_B是图像A的亮度(Y)方差，对A’做同样的映射以保持一致性。

所述步骤3)中，对于A′图，有一对对应的高分辨率EHA′图和降采样得到的 低分辨率A′图，对输出坐标图中的像素点，从EHA′中拷贝和其具有相同的块结 构的颜色值，放大细节并创建高分辨率合成图EHB′。

所述步骤3)中在RGB空间下，放大同时作用于RGB三个通道；在YIQ空间下， 放大只对Y通道进行。同时需要一张高分辨率的EHB图，然后把放大后得到的Y 通道值与EHB图中的IQ通道值组合得到高分辨率的颜色图EHB′。

本发明的技术构思为：基于CUDA的并行图像类比框架，同时采用适用于 GPU的最相似邻域查找策略来进行全局搜索和适用于GPU的k-coherence搜索进 行一致性搜索，根据一致性参数的不同，平衡两种搜索方法得到的结果，来决定 最终的待合成点。两个搜索过程都能并行地对每个像素点进行最相似邻域查找。 在合成结束后，采用合成放大的方法来提高合成图的分辨率。本方法在合成速度 和质量上都取得了令人满意的结果。

基于CUDA的快速图像类比合成算法利用图像金字塔分层合成，输入为三张 图：源图像A，滤波后的图像A’，源图像B，输出为滤波后的图像B’。算法能 学习得到两幅具有相关性的图像A和A’之间的关系，并将这种学习得到的关系作 用于另一幅图B上，生成B’，使得B到B’的变化在视觉上满足A到A’的变化过 程。假设p为A，A’中相同坐标位置的一点，q为B，B’中相同坐标位置的一点。 合成的每一层中依次执行超采样(upsample)和校正(correction)二个步骤。超 采样过程用来获得更高分辨率的图像。校正过程用来对待合成像素在样图中查找 最相似点，最相似点的坐标由全局搜索和k-coherence搜索共同决定。搜索时对于 待合成点q记录的都是在样图中相应点的坐标值。coherencefac为一致性参数。对 于校正过程，分为全局搜索和k-coherence搜索，全局搜索能查找到在L2准则下 的最匹配像素，但由于L2准则不能很好的度量感知相似性(perceptual similarity)， 因此加入k-coherence搜索进行调节。

本发明的有益效果主要表现在：用CUDA框架来进行并行图像类比合成，易 于实现。本方法同时采用全局搜索和k-coherence搜索在GPU上查找最相似邻域， 基于GPU全局最相似邻域查找策略与ANN查找方法相比，速度有较大提高。与 现有方法相比，本方法适用于图像类比框架中的各种应用，但大大加快了合成速 度，且合成质量也令人满意。

该方法的主要优点如下：

(1)思路新颖。利用CUDA架构进行图像类比算法设计，突破传统方法的时间 限制，具有较大创新意义。

(2)实现简单。在C语言环境下完成对GPU并行处理架构的调用，代码编写容 易。

(3)运行速度快，质量高。该方法充分利用GPU的并行性，合成图像速度快。 全局搜索和k-coherence搜索的结合大大提高了合成质量。

(4)方便易行。用户无需输入复杂的参数，只需提供输入图像，选择所需应用， 系统自动给出参数，为用户提供合成结果。

(5)适用范围广：本方法适用于图像类比框架中的各种应用，合成效果好。

附图说明

图1为不同应用所需的输入图像、采用的颜色空间、coherenceFac建议取值 示意图。

图2为基于CUDA的快速图像类比合成方法的流程图。

具体实施方式

下面结合附图对本发明作进一步描述。

参照图1和图2，一种基于CUDA的快速图像类比合成方法，包括如下步骤：

1)生成源图像A，滤波后的图像A’，源图像B的高斯金字塔

2)当前合成第L层，对金字塔从粗糙层到精细层循环做：

2.1)对B′_L+1超采样得到B′_L

2.2)对B′_L中的每个像素点q，采用基于GPU的全局搜索策略找到全局搜 索最优解q_global

2.3)对B′_L中的每个像素点q，采用基于GPU的k-coherence搜索策略找 到一致性搜索最优解q_kcoh

2.4)如果全局搜索得到的结果更好，则将q_Best赋值为q_global，否则将q_Best赋 值为q_kcoh。把q_Best在A′_L中的对应坐标保存在B′_L(q)中。步骤2.2)，2.3)，2.4)可循 环进行多次。

3)对B’做合成放大，得到高分辨率的EHB’。

本算法在实现时利用CUDA编程模型中线程、线程块、线程块网格的层次结 构，合理划分线程块、线程块网格，合理分配存储空间减少内存访问时间来加快 合成速度。

所述步骤2.1)中，要对已有的B′_L+1超采样得到B′_L。根据图片大小分配线程 网格和线程块。将B′_L+1与纹理存储器绑定。每个线程执行的操作为：得到B′_L(p) 的坐标(u，v)，取得B′_L+1中(u/2，v/2)点颜色值color，B′_L(p)的值为 color*2+(u％2，v％2)。

所述步骤2.2)中，输入为A，A’，B金字塔L层对应的A_L，A′_L，B_L图像，输 出为B′_L，B′_L中每个像素点的RG通道保存了由全局搜索得到的最相似点的坐标 q_global，同时保留了其对应的距离，以便和2.3)得到的结果进行比较。对B′_L中的 一点q，根据q点在B_L，B′_L的邻域在A_L，A′_L中找到和其最相似的点。该问题可以 归纳为高维度的最近k邻域查找问题。

设R＝{r₁，r₂，...，r_m}为参考集合，Q＝{q₁，q₂，...，q_n}为查询集合。对于纹理合成 来说，m为样图中的点个数，n为查询点个数。N()代表邻域向量，N(r)和N(q)为收 集图像中邻域得到的向量。公式(1)表示，对于查询集合Q中的每个点q，在R中根 据L₂度量准则找到k个与其最相似点的集合C_1..k(q)，也就是对于每个点q，在集合 R中找到和其邻域距离最小的前k个点，

$\left(\begin{array}{cc}{C}_{1..k}\left(q\right)=\arg {\min }_{r\in R}{\left|\right|N\left(r\right)-N\left(q\right)\left|\right|}_2& (r\in R,q\in Q)\end{array}\right)---\left(1\right)$

分解该运算，得到公式(2)：

$C_{1..k}\left(q\right)={\left|\right|N\left(q\right)\left|\right|}_2^2+\arg {\min }_{r\in R}({\left|\right|N\left(r\right)\left|\right|}_2^2-2{\left|\right|N\left(r\right)·N\left(q\right)\left|\right|}_2)---\left(2\right)$

从公式(2)中可以看出不依赖于R中点，也不依赖于查询点 q，因此对于q点与R中点比较，只要计算部分即可。收 集邻域向量得到参考集合对应的矩阵RMatrix＝|N(r₁)，N(r₂)，...，N(r_n)|，查询集合 对应的矩阵QMatrix＝|N(q₁)，N(q₂)，...，N(q_n)|。对于全局搜索来说，k取1，RMatrix 根据A，A’中对应的点来获得，QMatrix根据B，B’中对应的点来获得。

利用CUDA架构快速地进行此运算，在CUDA架构下程序被分为host端和 device端，host端执行CPU上的操作，device端执行显卡上的操作，内核函数Kernel 为CUDA运行在device端的并行计算函数，因此对Q进行分块，分批处理查询集合 中的点，动态确定当前可分配的最大Q分块数；对于一个Q分块，计算分块中的点 和参考集合中所有像素点间邻域的距离和坐标，采用插入排序法找到具有最小距 离的前k个点的坐标集合，

kernel函数计算步骤如下：

Step1：收集RMatrix和QMatrix，在device端分配存储空间，把RMatrix、 QMatrix拷贝给device端；在host端分配空间indexMatrix，存储计算得到的最相 似点坐标

Step2：查找显存中剩余空间，确定每次可同时搜索的Q集合大小 maxQueryNum，

Step3：根据maxQueryNum，对Q分块，对Q的每一个分块：

Step3.1：计算距离${\left|\right|N\left(r\right)\left|\right|}_2^2-2{\left|\right|N\left(r\right)·N\left(q\right)\left|\right|}_2$

Step3.2：找到距离最小的前k个最相似点集合C_1..k(q)

Step 3.3将C_1..k(q)中点的坐标拷贝到indexMatrix对应位置的内存中，

所述步骤2.3)中：输入同样为A，A’，B金字塔L层对应的A_L，A′_L，B_L图像， 输出为B′_L，B′_L中每个像素点的RG通道保存了由k-coherence搜索得到的最相 似点的坐标q_kcoh。在用CUDA实现时，线程网格的分割为与图像大小相对应，每 个线程块处理一个像素元素，它拥有的线程数为一个像素需要查找的coherence 候选集个数。通过一个预处理过程，对于A和A’中的每个点，找到k个与其最相 似的点，把这些点的坐标保存在纹理中。

计算步骤如下(kernel函数)：

Step1.在共享存储器中，声明类型为iht4的distance数组，xy分量存储坐标， z分量存储距离，声明类型为uchar2的kcoh数组，把预处理得到的纹理中像素对 应的坐标载入kcoh，

Step2.每个线程访问共享内存，取出kcoh中的坐标，根据CUDA内建函数 blockIdx.x，blockIdx.y获取待合成像素点位置，基于度量准则，计算两点邻域的 之间的dist

Step3：线程块中的每个线程将计算得到的dist及其对应的样图坐标保存到 distance数组中

Step4：当一个线程块中所有线程都计算完毕时，用树形规约算法找到具有最 小dist的样图坐标，即为所求的最相似点坐标，

树形规约算法可有效解决存取共享内存带宽冲突的问题。用树形规约算法搜 索最小值，每次比较时，存储在前半部分线程中的数据和存储在后半部分线程中 的数据相比较，取较小值更新前半部分的数据。每次迭代后减半，直到最后剩下 一个值，即为所求的最小值。数据单元之间的比较由多线程并行执行，速度较快。 在CUDA架构下，当从设备端复制数据到主机端内存时，数据一般先由直接内存 访问(DMA)模式复制到一块页面锁定内存，再由页面锁定内存复制到主机端内 存。为避免重复的内存赋值，优化主机端和设备端之间的内存传输，可使用 cudaMallocHost函数直接在主机端内存中取一块页面锁定内存，该内存可以被设 备端直接访问。页面锁定内存强制让系统在物理内存中完成内存申请和释放，不 参与页交换，比可分页内存可在更高带宽下读写，能有效地减少数据写回时间， 提高执行效率。

所述步骤2.4)中，根据coherenceFac比较步骤2.2)，2.3)计算得到的结果，得 到q_Best的取值。比较准则如下：若d(q_kcoh，q)＞d(q_global，q)(1+coherencefac*2^L)，则 q_Best＝q_global，否则q_Best＝q_kcoh。

两个邻域的相似性度量准则定义如下，同时比较A，A’和B，B’邻域的特征向量：

d(p，q)＝||F_A(p)-F_B(q)||_L2+||F_A′(p)-F_B′(q)||_L2。

对于不同应用，可选取不同的特征F(·)和coherenceFac。对于特征选取，在不 同颜色空间中计算邻域距离时F(·)的选取也不同。在RGB空间下，选取RGB颜色 值作为特征向量。在YIQ空间下，选取亮度Y作为特征向量，将合成得到的新Y值 与原有的IQ值组合还原成RGB图像。图1给出了对于不同应用输入样图、颜色空 间及coherenceFac的选取建议。

对于艺术风格学习应用来说，只把图像转到YIQ空间仍会产生图像邻域的直 方图重叠，使用亮度分布的平均值和方差对图像的Y通道进行线性变换，对原有 的YIQ空间进行转化(remap)：

$Y\left(p\right)\Leftarrow \frac{{\sigma }_B}{{\sigma }_A}(Y\left(p\right)-{\mu }_A)+{\mu }_B$

其中Y(P)是图像A中P点的亮度(Y)值，μ_A和μ_B是图像A的亮度(Y)平均值。σ_A和σ_B是图像A的亮度(Y)方差。对A’做同样的映射以保持一致性。

所述步骤3)中：采用合成放大(synthesis magnification)来快速地根据已有 的高分辨率样图放大得到高分辨率的合成图。对于A′图，有一对对应的高分辨率 EHA′图和降采样得到的低分辨率A′图。对输出坐标图中的像素点，从EHA′中 拷贝和其具有相同的块结构的颜色值，放大细节，创建高分辨率合成图EHB′。 合成放大对RGB空间和YIQ空间都适用。在RGB空间下，放大同时作用于RGB三 个通道。在YIQ空间下，放大只对Y通道进行。同时需要一张高分辨率的EHB图， 然后把放大后得到的Y通道值与EHB图中的IQ通道值组合得到高分辨率的颜色 图EHB′。