4.4 图像卷积

卷积在图像处理中经常用到，其会根据每个像素周围的像素点修改当前像素点的值。卷积核就是用来描述每个像素点如何被附近的像素点所影响。例如，模糊滤波中，使用平均的权重方式来进行计算，这样差异比较大的像素点会减少差异。对于相同的图像，如果想要做不同的操作，我们只需要变化滤波器即可，这样就能做锐化、模糊、边缘增强和图像压花。

卷积算法会遍历原始图像中的每一个像素点。对于每个原始像素点，滤波器中心点会为于像素点的上方，然后将中心点以及周围点的像素点与滤波器中对应的权重值相乘。然后将乘积的结果值相加后，产生出新的值作为输出。图4.3中就展示了该算法的具体过程。

4.4 图像卷积 - 图1

图4.3 卷积滤波如何对原图像进行处理。

图4.4a是原始图，图4.4b中展示了原图经过模糊滤波后的结果，图4.4c展示了经过一个压花滤波器处理后的结果。

4.4 图像卷积 - 图2

图4.4 不同卷积核对同一张图像进行处理：a)为原图；b)为模糊处理；c)为压花处理。

程序清单4.6中使用C/C++实现了一个串行的卷积操作。两层外部循环以遍历原始图像中所有的像素点。每一次滤波操作，每个原始点和其附近的点都要参与计算。需要注意的是，滤波器有可能访问到原始图像之外的区域。为了解决这个问题，我们在最内层循环中添加了四个显式的检查，当滤波器对应的坐标点位于原始图像之外，我们将使用与其最近的原始图像的边界值。

  1. /* Iterate over the rows of the source image */
  2. for (int i = 0; i < rows; i++)
  3. {
  4.   /* Iterate over the columns of the source image */
  5.   for (int j = 0; j < cols; j++)
  6.   {
  7.     /* Reset sum for new source pixel */
  8.     int sum = 0;
  10.     /* Apply the filter to the neighborhood */
  11.     for (int k = -halfFilterWidth; k <= halfFilterWidth; k++)
  12.     {
  13.       for (int l = -halfFilterWidth; l <= halfFilterWidth; l++)
  14.       {
  15.         /* Indices used to access the image */
  16.         int r = i+k;
  17.         int c = j+l;
  19.         /* Handle out-of-bounds locations by clamping to  
  20.          * the border pixel*/
  21.         r = (r < 0)? 0 : r;
  22.         c = (c < 0)? 0 : c;
  23.         r = (r >= rows)? rows - 1: r;
  24.         c = (c >= cols)? cols - 1: c;
  26.         sum += Image[r][c] *
  27.                Filter[k+halfFilterWidth][l+halfFilterWidth];
  28.       }
  29.     }
  31.     /* Write the new pixel value */
  32.     outputImage[i][j] = sum;
  33.   }
  34. }

程序清单4.6 图像卷积的串行版本

OpenCL中，使用图像内存对象处理卷积操作要比使用数组内存对象更有优势。图像采样方式会自动对访问到图像之外的区域进行处理(类似于上节图像转换中所提到的)，并且访问缓存中二维数据在硬件上也会有所优化(将会在第7章讨论)。

OpenCL上实现卷积操作几乎没有什么难度，并且写法类似于卷积操作的C版本。OpenCL版本中，我们为每一个输出的像素点创建了一个工作项，使用并行的方式将最外层两个循环去掉。那么每个工作项的任务就是完成最里面的两个循环，这两个循环完成的就是滤波操作。前面的例子中，读取源图像数据，需要配置一个OpenCL结构体，来指定数据的类型。本节的例子中，将继续使用read_imagef()。完整的内核代码将在代码清单4.7中展示。

  1. __kerenl
  2. void convolution(
  3.   __read_only image2d_t inputImage,
  4.   __write_only image2d_t outputImage,
  5.   int rows,
  6.   int cols,
  7.   __consetant float *filter,
  8.   int filterWidth,
  9.   sampler_t sampler)
  10. {
  11.   /* Store each work-item's unique row and column */
  12.   int column = get_global_id(0);
  13.   int row = get_global_id(1);
  15.   /* Half the width of the filter is needed for indexing
  16.    * memory later */
  17.   int halfWidth = (int)(filterWidth / 2);
  19.   /* All accesses to images return data as four-element vectors
  20.    * (i.e., float4), although only the x component will contain
  21.    * meaningful data in this code */
  22.   float4 sum = {0.0f,0.0f,0.0f,0.0f};
  24.   /* Iterator for the filter */
  25.   int filterIdx = 0;
  27.   /* Each work-item iterates around its local area on the basis of the
  28.    * size of the filter*/
  29.   int2 coords; // Coordinates for accessing the image
  31.   /* Iterate the filter rows*/
  32.   for (int i = -halfWidth; i <= halfWidth; i++)
  33.   {
  34.     coords.y = row + i;
  35.     /* Iterate over the filter columns */
  36.     for (int j = -halfWidth; j <= halfWidth; j++)
  37.     {
  38.       coords.x - column + 1;
  40.       /*Read a pixel from the image. A single-channel image
  41.        * stores the pixel in the x coordinate of the reatured
  42.        * vector. */
  43.       pixel = read_imagef(inputImage, sampler, coords);
  44.       sum.x += pixel.x * filter[filterIdx++];
  46.     }
  47.   }
  49.   /* Copy the data to the output image */
  50.   coords.x = column;
  51.   coords.y = row;
  52.   write_imagef(outputImage, coords, sum);
  53. }

程序清单4.7 使用OpenCL C实现的图像卷积

访问图像总会返回具有四个值的向量(每个通道一个值)。前节例子中，我们加了.x来表示只需要返回值的第一个元素。本节例子中，我们会申明pixel(图像访问之后返回的值)和sum(存储结果数据，以拷贝到输出图像)，其类型都为float4。当然，我们仅对x元素进行累加滤波像素值的计算(第45行)。

卷积核很适合放置到常量内存上，因为所有工作项所要用到的卷积核都是相同的。简单的添加关键字”__costant”在函数参数列表中(第7行)，用于表示卷积核存放在常量内存中。

之前的例子中，我们是直接在内核内部创建了一个采样器。本节例子中，我们将使用主机端API创建一个采样器，并将其作为内核的一个参数传入。同样，本节例子中我们将使用C++ API(C++采样器构造函数需要相同的参数)。

主机端创建采样器的API如下所示：

  1. cl_sampler clCreateSampler(
  2.   cl_context context,
  3.   cl_bool normalized_coords,
  4.   cl_addressing_mode addressing_mode,
  5.   cl_filter_mode filter_mode,
  6.   cl_int *errcode_ret)

其C++ API如下所示：

  1. cl::Sampler::Sampler(
  2.   const Context &context,
  3.   cl_bool normalized_coords,
  4.   cl_addressing_mode addressing_mode,
  5.   cl_filter_mode filter_mode,
  6.   cl_int *err = NULL)

使用C++创建采样器的方式如下所示：

  1. cl::Sampler sampler = new cl::Sampler(context, CL_FALSE, CL_ADDRESS_CLAMP_TO_EDGE, CL_FILTER_NEAREST);

图像旋转例子中，采样器使用非标准化坐标。这里我们将展示不同于内部使用时的另外两个采样器参数：滤波模式将使用最近的像素点的值，而非进行差值后的值(CL_FILTER_NEAREST)，并且寻址模式将在访问到图像区域之外时，将其值置为最接近的图像边界值(CL_ADDRESS_CLAMP_TO_EDGE)。(注：这里要注意一下CL和CLK的区别。CL开头的是使用在主机端API上，CLK则直接使用在OpenCL内核代码中。)

使用C++ API，创建二维图像使用image2D类进行创建，创建这个类需要一个ImageFormat对象作为参数。C API则不需要图像描述器的传入。

Image2D和ImageFormat的构造函数如下：

  1. cl::Image2D::Image2D(
  2.   Context& context,
  3.   cl_mem_flags flags,
  4.   ImageFormat format,
  5.   ::size_t width,
  6.   ::size_t height,
  7.   ::size_t row_pitch = 0,
  8.   void *host_ptr = NULL,
  9.   cl_int *err = NULL)
  11. cl::ImageFormat::ImageFormat(
  12.   cl_channel_order order,
  13.   cl_channel_type type)

这样我们就能创建卷积所需要的输入和输出图像，其调用方式如下所示：

  1. cl::ImageFormat imageFormat = cl::ImageFormat(CL_R, CL_FLOAT);
  2. cl::Image2D inputImage = cl::Image2D(context, CL_MEM_READ_ONLY, imageFormat, imageCols, imageRows);
  3. cl::Image2D outputImage = cl::Image2D(context, CL_MEM_WRITE_ONLY, imageFormat, imageCols, imageRows);

使用C++ API实现图像卷积的主机端代码在代码清单4.8中展示。主程序中，使用了一个5x5的高斯模糊滤波核用于卷积处理。

  1. #define __CL_ENABLE_EXCEPTIONS
  3. #include <CL/cl.hpp>
  4. #include <fstream>
  5. #include <iostream>
  6. #include <vector>
  8. #include "utils.h"
  9. #include "bmp_utils.h"
  11. static const char *inputImagePath = "../../Images/cat.bmp";
  13. static float gaussianBlurFilter[25] = {
  14.   1.0f / 273.0f, 4.0f / 273.0f, 7.0f / 273.0f, 4.0f / 273.0f, 1.0f / 273.0f,
  15.   4.0f / 273.0f, 16.0f / 273.0f, 26.0f / 273.0f, 16.0f / 273.0f, 4.0f / 273.0f,
  16.   7.0f / 273.0f, 26.0f / 273.0f, 41.0f / 273.0f, 26.0f / 273.0f, 7.0f / 273.0f,
  17.   4.0f / 273.0f, 16.0f / 273.0f, 26.0f / 273.0f, 16.0f / 273.0f, 4.0f / 273.0f,
  18.   1.0f / 273.0f, 4.0f / 273.0f, 7.0f / 273.0f, 4.0f / 273.0f, 1.0f / 273.0f
  19. };
  20. static const int gaussianBlurFilterWidth = 5;
  22. int main()
  23. {
  24.   float *hInputImage;
  25.   float *hOutpueImage;
  27.   int imageRows;
  28.   int imageCols;
  30.   /* Set the filter here */
  31.   int filterWidth = gaussianBlurFilterWidth;
  32.   float *filter = gaussianBlurFilter;
  34.   /* Read in the BMP image */
  35.   hInputImage = readBmpFloat(inputImagePath, &imageRows, &imageCols);
  37.   /* Allocate space for the output image */
  38.   hOutputImage = new float[imageRows * imageCols];
  40.   try
  41.   {
  42.     /* Query for platforms */
  43.     std::vector<cl::Platform> platforms;
  44.     cl::Platform::get(&platdorms);
  46.     /* Get a list of devices on this platform */
  47.     std::vector<cl::Device> device;
  48.     platforms[0].getDevices(CL_DEVICE_TYPE_GPU, &devices);
  50.     /* Create a context for the devices */
  51.     cl::Context context(devices);
  53.     /* Create a command-queue for the first device */
  54.     cl::CommandQueueu queue = cl::CommandQueue(context, devices[0]);
  56.     /* Create the images */
  57.     cl::ImageFormat imageFormat = cl::ImageFormat(CL_R, CL_FLOAT);
  58.     cl::Image2D inputImage = cl::Image2D(context, CL_MEM_READ_ONLY, 
  59.         imageFormat, imageCols, imageRows);
  60.     cl::Image2D outputImage = cl::Image2D(context, 
  61.         CL_MEM_WRITE_ONLY,
  62.         imageFormat, imageCols, imageRows);
  64.     /* Create a buffer for the filter */
  65.     cl::Buffer filterBuffer = cl::Buffer(context, CL_MEM_READ_ONLY,
  66.     filterWidth * filterWidth * sizeof(float));
  68.     /* Copy the input data to the input image */
  69.     cl::size<3> origin;
  70.     origin[0] = 0;
  71.     origin[1] = 0;
  72.     origin[2] = 0;
  73.     cl::size<3> region;
  74.     region[0] = 0;
  75.     region[1] = 0;
  76.     region[2] = 0;
  77.     queue.enqueueWriteImage(inputImage, CL_TRUE, origin, region,
  78.     0, 0,
  79.     hInputImage);
  81.     /* Copy the filter to the buffer*/
  82.     queue.enqueueWriteBuffer(filterBuffer, CL_TRUE, 0,
  83.     filterWidth * filterWidth * sizeof(float), filter);
  85.     /* Create the sampler */
  86.     cl::Sampler sampler = cl::Sampler(context, CL_FALSE,
  87.     CL_ADDRESS_CLAMP_TO_EDGE, CL_FILTER_NEAREST);
  89.     /* Read the program source */
  90.     std::ifstream sourceFile("image-convolution.cl");
  91.     std::string sourceCode(std::istreambuf_iterator<char>(sourceFile),
  92.     (std::istreambuf_iterator<char>()));
  94.     cl::Program::Source source(1,
  95.       std::make_pair(sourceCode.c_str(),
  96.       sourceCode.length() + 1));
  98.     /* Make program form the source code */
  99.     cl::Program program = cl::Program(context, source);
  101.     /* Create the kernel */
  102.     cl::Kernel kernel(program, "convolution");
  104.     /* Set the kernel arguments */
  105.     kernel.setArg(0, inputImage);
  106.     kernel.setArg(1, ouputImage);
  107.     kernel.setArg(2, filterBuffer);
  108.     kernel.setArg(3, filterWIdth);
  109.     kernel.setArg(4, sampler);
  111.     /* Execute the kernel */
  112.     cl::NDRange global(imageCols, imageRows);
  113.     cl::NDRange local(8, 8);
  114.     queue.enqueueNDRangeKernel(kernel, cl::NullRange, global,
  115.     local);
  116.     /* Copy the output data back to the host */
  117.     queue.enqueueReadImage(outputImage, CL_TRUE, origin, region,
  118.     0, 0,
  119.     hOutputImage);
  121.     /* Save the output BMP image */
  122.     writeBmpFloat(hOutputImage, "cat-filtered.bmp", imageRows, imageCols, inputImagePath);
  123.   }
  124.   catch(cl::Error error){
  125.     std::cout << error.what() << "(" << error.err() << ")" << std::endl;
  126.   }
  128.   free(hInputImage);
  129.   delete hOutputImage;
  130.   return 0;
  131. }

程序清单4.8 图像卷积主机端完整代码