1. 目标
本文目标是描述了如何使用统一向量指令（universal intrinsics）功能对 C++ 代码进行向量化以获得更快运行时。将简要介绍SIMD 内置函数以及如何使用宽寄存器，然后是有关使用宽寄存器的基本操作的教程。
2. 理论
在本节中，将简要介绍一些概念，以更好地帮助理解该功能。
2.1 内置函数
内置函数是由编译器单独处理的函数。这些函数通常经过优化以尽可能以最有效的方式执行，通常比正常实现运行得更快。但是，由于这些函数依赖于编译器，因此很难编写可移植的应用程序。
2.2 SIMD
SIMD 代表单指令多数据。SIMD 指令 允许处理器向量化计算。存储数据的寄存器可以是128 位、256 位或512位宽。每个寄存器存储多个相同数据类型的值。寄存器的大小和每个值的大小决定了总共存储的值的数量。
根据CPU 支持的指令集，可以使用不同的寄存器。要了解更多信息，请看这里
3. universal intrinsics 统一向量指令
OpenCV 统一向量指令 提供了对 SIMD 向量化方法的抽象，并允许用户使用统一的内置函数而无需编写系统特定代码。
OpenCV Universal Intrinsics 支持以下指令集：
-
支持各种类型的128 位寄存器，适用于各种架构，包括
- x86(SSE/SSE2/SSE4.2),
- ARM（NEON），
- PowerPC(VSX),
- MIPS（MSA）。
-
支持x86(AVX2) 256 位寄存器
-
支持x86(AVX512)512位寄存器
现在将介绍可用的数据结构和函数：
- 寄存器结构
- 加载和存储
- 数学运算
- Reduce and Mask
3.1 寄存器结构
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
可用的数据类型和大小：
| 类型 | 字节数 |
|---|---|
| uint | 8, 16, 32, 64 |
| int | 8, 16, 32, 64 |
| float | 32, 64 |
- 固定大小的寄存器：这些结构具有固定的位大小并保存固定数量的值。需要知道系统支持哪些 SIMD 指令集并选择兼容的寄存器。仅当需要精确的位长度时才使用这些。
每个结构都遵循约定：
v_[type of value][size of each value in bits]x[number of values]
假设我们要
- 在128 位寄存器保存32位有符号整数。由于寄存器大小已知，可以计算出寄存器可以容纳的数量（128/32 = 4）
v_int32x8 reg1 // holds 8 32-bit signed integers.
- 在512 位寄存器中保存64 位浮点数：
v_float64x8 reg2 // reg2.nlanes = 8
3.2 加载和存储操作
现在知道了寄存器是如何工作的，下边看看操作寄存器的函数。
3.2.1 Load
Load函数允许将值加载到寄存器中。
- 构造函数- 当声明一个寄存器结构时，可以提供一个内存地址，寄存器将从那里获取连续值，或者显式提供多个参数的值（显式多个参数仅适用于常量大小的寄存器）：
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()功能。
3.2.2 存储
存储函数允许将寄存器中的值存储到特定的内存位置。
- 要将寄存器中的值存储到内存位置，您可以使用v_store() 函数：
float ptr[4];
v_store(ptr, reg); // store the first 128 bits(interpreted as 4x32-bit floats) of reg into ptr.
-
笔记
确保ptr与 register 具有相同的类型。还可以在执行操作之前将寄存器转换为正确的类型。简单地将指针类型转换为特定类型将导致对数据的错误解释。
3.3 二元和一元运算符
universal intrinsics 统一向量指令集函数提供了 对应元素(element wise)的二元和一元运算。
- 算术：可以按元素对两个寄存器进行加、减、乘和除。寄存器必须具有相同的宽度并保持相同的类型。将两个寄存器相乘，例如：
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（0xff 表示 8 位，0xffff 表示 16 位等），而假值则返回转换为零的位。
// 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_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 位寄存器。
3.4 Reduce and Mask
- Reduce操作：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
-
Mask操作：Mask操作允许在宽寄存器中复制条件。这些包括：
-
v_check_all() - 返回一个布尔值，如果寄存器中的所有值都小于零，则为真。
-
v_check_any() - 返回一个布尔值，如果寄存器中的任何值小于零，则为真。
-
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.
*/
4.范例代码
在下一节中，将对单通道的简单卷积函数进行矢量化，并将结果与标量实现进行比较。
- 笔记
并非所有算法都通过手动矢量化得到改进。事实上，在某些情况下，编译器可以自动向量化代码，从而为标量实现产生更快的结果。
可以从上一个教程中了解有关卷积的更多信息。使用与上一教程相同的简单实现，并将其与矢量化版本进行比较。
完整的教程代码在这里。
4.1 向量化卷积
我们将首先实现一维卷积，然后对其进行矢量化。2-D 矢量化卷积将跨行执行 1-D 卷积以产生正确的结果。
4.1.1 一维卷积：标量
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);
}
}
- 首先创建变量并填充 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);
- 对于主循环，选择索引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);
}
4.1.2 一维卷积：向量
现在来看看一维卷积的向量化版本。
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 = v_float32().nlanes;
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 = vx_load(ans + i) + kernel_wide * window;
v_store(ans + i, sum);
}
for (; i < len; i++)
{
*(ans + i) += sptr[i + k]*kptr[k];
}
}
}
- 在上边例子中，卷积核是一个浮点数。由于卷积核的数据类型最大，需要将 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);
-
现在，对于卷积核中的每一列，计算该值与所有长度为的窗口
step向量的标量积。并将这些值添加到 ans 中已经存储的值中int step = v_float32().nlanes; 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 = vx_load(ans + i) + kernel_wide * window; v_store(ans + i, sum); } for (; i < len; i++) { *(ans + i) += sptr[i + k]*kptr[k]; } }
- 声明一个指向 src_32 和卷积核的指针，并为每个卷积核元素运行一个循环
```
int step = v_float32().nlanes;
float *sptr = src_32.ptr<float>(row), *kptr = kernel.ptr<float>(rowk);
for (int k = 0; k < ksize; k++)
{
```
- 用当前卷积核元素加载一个寄存器。一个窗口从0移动到len-step并且它与 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 = vx_load(ans + i) + kernel_wide * window;
v_store(ans + i, sum);
}
```
- 由于长度可能不能被步长整除，我们直接处理剩余的值。尾值的数量将始终小于步长，不会显着影响性能。将所有值存储到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-D 卷积的中间步骤时使用
4.1.3 二维卷积
假设卷积核有ksize行。为了计算特定行的值，计算前一个 ksize/2 行和下一个ksize/2行与相应的内核行的一维卷积。 最终值只是单个 1-D 卷积的总和
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 = v_float32().nlanes;
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 = 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);
}
- 首先初始化变量，并在src矩阵的上方和下方填充边框。左侧和右侧由一维卷积函数处理。
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 = v_float32().nlanes;
- 对于每一行，计算它上面和下面的行的一维卷积。然后将值添加到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 = 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];
}
}
- 最终将dst矩阵转换为8 位
unsigned char矩阵
const int alpha = 1;
dst.convertTo(dst, CV_8UC1, alpha);
5. 结果
在本文中，使用了水平梯度卷积核。而且为两种方法获得了相同的输出图像。
运行速度的改进因 CPU 中可用的 SIMD 功能而异。
