OpenCV  4.10.0
开源计算机视觉
加载...
搜索...
没有匹配项
使用通用内在函数矢量化您的代码

前一个教程: 如何使用 OpenCV parallel_for_对代码进行并行处理

兼容性OpenCV >= 3.0

目标

本教程的目标是提供一个指南,利用通用内在函数对 C++ 代码进行矢量化处理,以实现更快的运行时间。我们首先简要介绍SIMD 内在函数以及如何使用宽寄存器,然后完成一个关于使用宽寄存器的基本操作的教程。

理论

在本部分,我们简要介绍几个概念,以便更好地理解该功能。

内在函数

内在函数由编译器单独处理。这些函数通常经过优化,以尽可能高效地执行,因此运行速度比一般实现要快。然而,由于这些函数依赖于编译器,因此编写可移植应用程序变得困难。

SIMD

SIMD 的全称是单指令,多数据。SIMD 内在函数允许处理器对计算进行矢量化。数据存储在称为寄存器的地方。寄存器可能为128 位256 位512 位宽。每个寄存器存储相同数据类型多个值。寄存器的宽度和每个值的大小决定了总共存储的值的数量。

根据 CPU 支持的指令集,您可能可以使用不同的寄存器。如需了解更多信息,请查看 此处

通用内在函数

OpenCV 的通用内在函数提供了一个对 SIMD 矢量化方法的抽象,允许用户使用内在函数,而不需要编写特定系统的代码。

OpenCV 通用内在函数支持以下指令集

  • 对各种体系结构(包括以下几个)实现了支持各种类型的128 位寄存器
    • x86(SSE/SSE2/SSE4.2),
    • ARM(NEON),
    • PowerPC(VSX),
    • MIPS(MSA)。
  • x86(AVX2) 上支持256 位寄存器
  • x86(AVX512) 支持512 位寄存器

我们现在将介绍可用结构和函数

  • 寄存器结构
  • 加载和存储
  • 数学运算
  • 归约和掩码

寄存器结构

Universal Intrinsics 集将每个寄存器实现为一个基于特定 SIMD 寄存器的结构。所有类型均包含枚举 nlanes,该枚举给出了类型可以容纳的确切值数量。这消除了在实现过程中对值数量进行硬编码的需要。

注意
每个寄存器结构都属于 cv 命名空间。

两种类型的寄存器

  • 可变大小的寄存器:这些结构没有固定大小,其确切位长在编译期间根据可用的 SIMD 功能推断出来。因此,nlanes 枚举的值是在编译时确定的。

    每个结构遵循以下约定

    v_[type of value][size of each value in bits]
    

    例如,v_uint8 容纳 8 位无符号整数v_float32 容纳 32 位浮点值。然后,我们声明一个寄存器,就像我们在 C++ 中声明任何对象一样

    根据可用的 SIMD 指令集,特定寄存器将容纳不同数量的值。例如:如果您的计算机最多支持 256 位寄存器,

    • v_uint8 将容纳 32 个 8 位无符号整数
    • v_float64 将容纳 4 个 64 位浮点数(双精度数)
        v_uint8 a;                            // a is a register supporting uint8(char) data
        int n = a.nlanes;                     // n holds 32
      
      可用的数据类型和大小
      类型位大小
      uint8, 16, 32, 64
      int8, 16, 32, 64
      float32, 64
  • 固定大小的寄存器:这些结构具有固定的位大小,并容纳恒定数量的值。我们需要了解系统支持的 SIMD 指令集,并选择兼容的寄存器。仅在需要确切位长时才使用它们。

    每个结构遵循惯例

    v_[type of value][size of each value in bits]x[number of values]
    

    假设我们要存储

    • 32 位(以位为单位的大小)有符号整数在128 位寄存器中。由于寄存器大小已经知道,我们可以找出寄存器中的数据点数128/32 = 4
        v_int32x8 reg1                       // holds 8 32-bit signed integers.
      
    • 512 位寄存器中的 64 位浮点数
        v_float64x8 reg2                     // reg2.nlanes = 8
      

加载和存储操作

既然我们知道了寄存器的工作原理,让我们来看看用于向这些寄存器填充值的函数。

  • 加载:加载函数允许您向寄存器中加载值。
    • 构造函数 - 在声明寄存器结构时,我们可以提供一个寄存器将从中选取连续值,或者提供多个参数作为显式值(显式多个参数仅对固定大小的寄存器可用)
        float ptr[32] = {1, 2, 3 ..., 32};   // ptr is a pointer to a contiguous memory block of 32 floats
      
        // Variable Sized Registers //
        int x = v_float32().nlanes;          // set x as the number of values the register can hold
      
        v_float32 reg1(ptr);                 // reg1 stores first x values according to the maximum register size available.
        v_float32 reg2(ptr + x);             // reg stores the next x values
      
        // Constant Sized Registers //
        v_float32x4 reg1(ptr);               // reg1 stores the first 4 floats (1, 2, 3, 4)
        v_float32x4 reg2(ptr + 4);           // reg2 stores the next 4 floats (5, 6, 7, 8)
      
        // Or we can explicitly write down the values.
        v_float32x4(1, 2, 3, 4);
      

    • 加载函数 - 我们可以使用加载方法并提供数据的内存地址
        float ptr[32] = {1, 2, 3, ..., 32};
        v_float32 reg_var;
        reg_var = vx_load(ptr);              // loads values from ptr[0] upto ptr[reg_var.nlanes - 1]
      
        v_float32x4 reg_128;
        reg_128 = v_load(ptr);               // loads values from ptr[0] upto ptr[3]
      
        v_float32x8 reg_256;
        reg_256 = v256_load(ptr);            // loads values from ptr[0] upto ptr[7]
      
        v_float32x16 reg_512;
        reg_512 = v512_load(ptr);            // loads values from ptr[0] upto ptr[15]
      
      注意
      加载函数假设数据未对齐。如果您的数据已对齐,则可以使用 vx_load_aligned() 函数。
  • 存储:存储函数允许您将寄存器中的值存储到特定内存位置。
    • 要将寄存器中的值存储到内存中,可以使用v_store()函数
        float ptr[4];
        v_store(ptr, reg); // store the first 128 bits(interpreted as 4x32-bit floats) of reg into ptr.
      

      注意
      确保ptr与寄存器类型一致。你也可以在执行操作之前将寄存器转换为适当的类型。只需将指针类型转换到一个具体类型,就会导致数据错误解释。

二元和一元运算符

通用内置命令集提供了元素级二元和一元运算。

  • 算术:例如,我们可以逐元素加减乘除两个寄存器。这些寄存器必须宽度相同且类型相同。要乘以两个寄存器:
      v_float32 a, b;                          // {a1, ..., an}, {b1, ..., bn}
      v_float32 c;
      c = a + b                                // {a1 + b1, ..., an + bn}
      c = a * b;                               // {a1 * b1, ..., an * bn}
    

  • 位运算和移位:我们可以位左移或位右移每个寄存器的每个元素的位元素。我们还可以在两个寄存器之间逐元素使用位运算 &、|、^ 和 ~
      v_int32 as;                              // {a1, ..., an}
      v_int32 al = as << 2;                    // {a1 << 2, ..., an << 2}
      v_int32 bl = as >> 2;                    // {a1 >> 2, ..., an >> 2}
    
      v_int32 a, b;
      v_int32 a_and_b = a & b;                 // {a1 & b1, ..., an & bn}
    

  • 比较运算符:我们可以使用 <、>、<=、>=、== 和 != 等运算符比较两个寄存器之间的值。由于每个寄存器包含多个值,因此我们无法对这些运算进行单一的布尔运算。而对于真值,则全部位转换到 1(对于 8 位是 0xff、对于 16 位是 0xffff 等),而假值则返回转换为 0 的位。
      // let us consider the following code is run in a 128-bit register
      v_uint8 a;                               // a = {0, 1, 2, ..., 15}
      v_uint8 b;                               // b = {15, 14, 13, ..., 0}
    
      v_uint8 c = a < b;
    
      /*
          let us look at the first 4 values in binary
    
          a = |00000000|00000001|00000010|00000011|
          b = |00001111|00001110|00001101|00001100|
          c = |11111111|11111111|11111111|11111111|
    
          If we store the values of c and print them as integers, we will get 255 for true values and 0 for false values.
      */
      ---
      // In a computer supporting 256-bit registers
      v_int32 a;                               // a = {1, 2, 3, 4, 5, 6, 7, 8}
      v_int32 b;                               // b = {8, 7, 6, 5, 4, 3, 2, 1}
    
      v_int32 c = (a < b);                     // c = {-1, -1, -1, -1, 0, 0, 0, 0}
    
      /*
          The true values are 0xffffffff, which in signed 32-bit integer representation is equal to -1.
      */
    

  • 最小值/最大值运算:我们可以使用v_min()v_max()函数来返回包含两个寄存器的最小值或最大值的寄存器
      v_int32 a;                               // {a1, ..., an}
      v_int32 b;                               // {b1, ..., bn}
    
      v_int32 mn = v_min(a, b);                // {min(a1, b1), ..., min(an, bn)}
      v_int32 mx = v_max(a, b);                // {max(a1, b1), ..., max(an, bn)}
    

注意
64 位整数不可用比较运算和最小值/最大值运算。位移运算、逻辑运算仅可用于整数。位移运算仅可用于 16、32 和 64 位寄存器。

约简和掩码

  • 约简运算:v_reduce_min()v_reduce_max()v_reduce_sum() 返回单个值,表示整个寄存器的最小值、最大值或和
      v_int32 a;                                //  a = {a1, ..., a4}
      int mn = v_reduce_min(a);                 // mn = min(a1, ..., an)
      int sum = v_reduce_sum(a);                // sum = a1 + ... + an
    

  • 掩码运算:掩码运算允许我们在宽寄存器中复制条件。这些运算包括:
    • v_check_all() - 返回一个布尔值,如果寄存器中的全部值小于 0,则为真。
    • v_check_any() - 返回一个布尔值,如果寄存器中的任何值小于 0,则为真。
    • v_select() - 返回一个寄存器,基于掩码混合两个寄存器。
        v_uint8 a;                           // {a1, .., an}
        v_uint8 b;                           // {b1, ..., bn}
      
        v_int32x4 mask:                      // {0xff, 0, 0, 0xff, ..., 0xff, 0}
      
        v_uint8 Res = v_select(mask, a, b)   // {a1, b2, b3, a4, ..., an-1, bn}
      
        /*
            "Res" will contain the value from "a" if mask is true (all bits set to 1),
            and value from "b" if mask is false (all bits set to 0)
      
            We can use comparison operators to generate mask and v_select to obtain results based on conditionals.
            It is common to set all values of b to 0. Thus, v_select will give values of "a" or 0 based on the mask.
        */
      

演示

在下面一节中,我们将对单通道的简单卷积函数进行矢量化处理,并将结果与标量实现进行比较。

注意
并不是所有的算法都能通过手动矢量化处理得到优化。事实上,某些情况下,编译器会对代码自动矢量化,因此产生更快的标量实现结果。

你可以在之前的教程中了解更多有关卷积的内容。我们使用之前的教程中的简单实现,并将其与矢量化版本进行比较。

完整教程代码 在此

对卷积进行矢量化处理

我们首先实现一维卷积,然后将其矢量化。2D 向量化卷积将执行行上的 1D 卷积以生成正确的结果。

1D 卷积:标量

void conv1d(Mat src, Mat &dst, Mat kernel)
{
int len = src.cols;
dst = Mat(1, len, CV_8UC1);
int sz = kernel.cols / 2;
copyMakeBorder(src, src, 0, 0, sz, sz, BORDER_REPLICATE);
for (int i = 0; i < len; i++)
{
double value = 0;
for (int k = -sz; k <= sz; k++)
value += src.ptr<uchar>(0)[i + k + sz] * kernel.ptr<float>(0)[k + sz];
dst.ptr<uchar>(0)[i] = saturate_cast<uchar>(value);
}
}
unsigned char uchar
Definition interface.h:51
#define CV_8UC1
Definition interface.h:88
  1. 我们首先设置变量,并在 src 矩阵的两边创建边框以处理边缘情况。
    int len = src.cols;
    dst = Mat(1, len, CV_8UC1);
    int sz = kernel.cols / 2;
    copyMakeBorder(src, src, 0, 0, sz, sz, BORDER_REPLICATE);
  2. 对于主循环,我们选择一个索引 i,并使用 k 变量将其与内核一起向两侧偏移。我们将值存储在 value 中并将其添加到 dst 矩阵中。
    for (int i = 0; i < len; i++)
    {
    double value = 0;
    for (int k = -sz; k <= sz; k++)
    value += src.ptr<uchar>(0)[i + k + sz] * kernel.ptr<float>(0)[k + sz];
    dst.ptr<uchar>(0)[i] = saturate_cast<uchar>(value);
    }

1D 卷积:矢量

现在我们将研究 1D 卷积的矢量化版本。

void conv1dsimd(Mat src, Mat kernel, float *ans, int row = 0, int rowk = 0, int len = -1)
{
if (len == -1)
len = src.cols;
Mat src_32, kernel_32;
const int alpha = 1;
src.convertTo(src_32, CV_32FC1, alpha);
int ksize = kernel.cols, sz = kernel.cols / 2;
copyMakeBorder(src_32, src_32, 0, 0, sz, sz, BORDER_REPLICATE);
int step = VTraits<v_float32x4>::vlanes();
float *sptr = src_32.ptr<float>(row), *kptr = kernel.ptr<float>(rowk);
for (int k = 0; k < ksize; k++)
{
v_float32 kernel_wide = vx_setall_f32(kptr[k]);
int i;
for (i = 0; i + step < len; i += step)
{
v_float32 window = vx_load(sptr + i + k);
v_float32 sum = v_add(vx_load(ans + i), v_mul(kernel_wide, window));
v_store(ans + i, sum);
}
for (; i < len; i++)
{
*(ans + i) += sptr[i + k]*kptr[k];
}
}
}
#define CV_32FC1
定义 interface.h:118
  1. 在本例中,内核是 float。由于内核的数据类型最大,因此我们将 src 转换为 float32,形成 *src_32*。我们还会像处理原始情况一样制作边框。
    Mat src_32, kernel_32;
    const int alpha = 1;
    src.convertTo(src_32, CV_32FC1, alpha);
    int ksize = kernel.cols, sz = kernel.cols / 2;
    copyMakeBorder(src_32, src_32, 0, 0, sz, sz, BORDER_REPLICATE);
  2. 现在,对于 *kernel* 中的每一列,我们计算该值与长度为 step 的所有 *window* 向量的标量积。我们将这些值添加到已存储在 ans 中的值
    int step = VTraits<v_float32x4>::vlanes();
    float *sptr = src_32.ptr<float>(row), *kptr = kernel.ptr<float>(rowk);
    for (int k = 0; k < ksize; k++)
    {
    v_float32 kernel_wide = vx_setall_f32(kptr[k]);
    int i;
    for (i = 0; i + step < len; i += step)
    {
    v_float32 window = vx_load(sptr + i + k);
    v_float32 sum = v_add(vx_load(ans + i), v_mul(kernel_wide, window));
    v_store(ans + i, sum);
    }
    for (; i < len; i++)
    {
    *(ans + i) += sptr[i + k]*kptr[k];
    }
    }
    • 我们声明一个指向 src_32 和 kernel 的指针,并为每个内核元素运行一个循环
      int step = VTraits<v_float32x4>::vlanes();
      float *sptr = src_32.ptr<float>(row), *kptr = kernel.ptr<float>(rowk);
      for (int k = 0; k < ksize; k++)
      {
    • 我们使用当前的内核元素加载寄存器。一个窗口从 *0* 移动到 *len - step*,并将它的 product 与 kernel_wide 阵列相乘并添加到存储在 *ans* 中的值。我们将值存储回 *ans*。
      v_float32 kernel_wide = vx_setall_f32(kptr[k]);
      int i;
      for (i = 0; i + step < len; i += step)
      {
      v_float32 window = vx_load(sptr + i + k);
      v_float32 sum = v_add(vx_load(ans + i), v_mul(kernel_wide, window));
      v_store(ans + i, sum);
      }
    • 由于长度不可能被 steps 整除,因此我们将直接处理剩余值。*tail* 值的数量将始终小于 *step*,并且不会对性能造成显着影响。我们将所有值存储到 float 指针 *ans*。我们还可以将它们直接存储在 `Mat` 对象中。
      for (; i < len; i++)
      {
      *(ans + i) += sptr[i + k]*kptr[k];
      }
    • 以下是一个迭代示例
        For example:
        kernel: {k1, k2, k3}
        src:           ...|a1|a2|a3|a4|...
      
      
        iter1:
        for each idx i in (0, len), 'step' idx at a time
            kernel_wide:          |k1|k1|k1|k1|
            window:               |a0|a1|a2|a3|
            ans:               ...| 0| 0| 0| 0|...
            sum =  ans + window * kernel_wide
                =  |a0 * k1|a1 * k1|a2 * k1|a3 * k1|
      
        iter2:
            kernel_wide:          |k2|k2|k2|k2|
            window:               |a1|a2|a3|a4|
            ans:               ...|a0 * k1|a1 * k1|a2 * k1|a3 * k1|...
            sum =  ans + window * kernel_wide
                =  |a0 * k1 + a1 * k2|a1 * k1 + a2 * k2|a2 * k1 + a3 * k2|a3 * k1 + a4 * k2|
      
        iter3:
            kernel_wide:          |k3|k3|k3|k3|
            window:               |a2|a3|a4|a5|
            ans:               ...|a0 * k1 + a1 * k2|a1 * k1 + a2 * k2|a2 * k1 + a3 * k2|a3 * k1 + a4 * k2|...
            sum =  sum + window * kernel_wide
                =  |a0*k1 + a1*k2 + a2*k3|a1*k1 + a2*k2 + a3*k3|a2*k1 + a3*k2 + a4*k3|a3*k1 + a4*k2 + a5*k3|
      
注意
函数参数还包括 *row*、*rowk* 和 *len*。当将函数用作 2 维卷积的中间步骤时,将使用这些值。

2 维卷积

假设我们的内核有 *ksize* 行。要计算特定行的值,我们计算前 *ksize/2* 行和后 *ksize/2* 行与相应内核行的 1 维卷积。最终值是各个 1 维卷积的和。

void convolute_simd(Mat src, Mat &dst, Mat kernel)
{
int rows = src.rows, cols = src.cols;
int ksize = kernel.rows, sz = ksize / 2;
dst = Mat(rows, cols, CV_32FC1);
copyMakeBorder(src, src, sz, sz, 0, 0, BORDER_REPLICATE);
int step = VTraits<v_float32x4>::vlanes();
for (int i = 0; i < rows; i++)
{
for (int k = 0; k < ksize; k++)
{
float ans[N] = {0};
conv1dsimd(src, kernel, ans, i + k, k, cols);
int j;
for (j = 0; j + step < cols; j += step)
{
v_float32 sum = v_add(vx_load(&dst.ptr<float>(i)[j]), vx_load(&ans[j]));
v_store(&dst.ptr<float>(i)[j], sum);
}
for (; j < cols; j++)
dst.ptr<float>(i)[j] += ans[j];
}
}
const int alpha = 1;
dst.convertTo(dst, CV_8UC1, alpha);
}
  1. 我们首先初始化变量,并在 *src* 矩阵之上和之下创建一个边框。左侧和右侧由 1 维卷积函数处理。
    int rows = src.rows, cols = src.cols;
    int ksize = kernel.rows, sz = ksize / 2;
    dst = Mat(rows, cols, CV_32FC1);
    copyMakeBorder(src, src, sz, sz, 0, 0, BORDER_REPLICATE);
    int step = VTraits<v_float32x4>::vlanes();
  2. 对于每一行,我们计算其上下行的 1 维卷积。然后,我们将这些值添加到 *dst* 矩阵中。
    for (int i = 0; i < rows; i++)
    {
    for (int k = 0; k < ksize; k++)
    {
    float ans[N] = {0};
    conv1dsimd(src, kernel, ans, i + k, k, cols);
    int j;
    for (j = 0; j + step < cols; j += step)
    {
    v_float32 sum = v_add(vx_load(&dst.ptr<float>(i)[j]), vx_load(&ans[j]));
    v_store(&dst.ptr<float>(i)[j], sum);
    }
    for (; j < cols; j++)
    dst.ptr<float>(i)[j] += ans[j];
    }
    }
  3. 最后,将dst矩阵转换为8 位unsigned char矩阵
    const int alpha = 1;
    dst.convertTo(dst, CV_8UC1, alpha);

结果

在本教程中,我们使用了水平梯度核。两者的输出图像相同。

运行时优化程度因 CPU 可用 SIMD 能力而异。