使用SIMD进行HOG优化

Question

有几种使用SIMD指令优化HOG描述符计算的尝试：OpenCV、Dlib和SIMD。它们都使用标量代码将结果幅值添加到HOG直方图中：

float histogram[height/8][width/8][18];
float ky[height], kx[width];
int idx[size];
float val[size]; 

for(size_t i = 0; i < size; ++i)
{
    histogram[y/8][x/8][idx[i]] += val[i]*ky[y]*kx[x];
    histogram[y/8][x/8 + 1][idx[i]] += val[i]*ky[y]*kx[x + 1];
    histogram[y/8 + 1][x/8][idx[i]] += val[i]*ky[y + 1]*kx[x];
    histogram[y/8 + 1][x/8 + 1][idx[i]] += val[i]*ky[y + 1]*kx[x + 1];
}

在那里，大小的值取决于实现，但一般意义相同。

我知道使用SIMD计算直方图的问题并没有简单有效的解决方案。但在这种情况下，我们有小尺寸（18）的直方图。它能帮助SIMD优化吗？

Answer 1

通过使用SIMD计算所有（展平的）直方图索引和bin增量，可以进行部分优化。然后在标量循环中处理这些。您可能还想将其剥离，以便一次处理一行，以便在缓存中保留临时bin索引和增量。由于使用了临时中间缓冲区，这似乎效率低下，但在实践中，我发现在类似的场景中总体上有了有益的提高。

uint32_t i = 0;

for (y = 0; y < height; ++y)   // for each row
{
    uint32_t inds[width * 4];  // flattened histogram indices for this row
    float vals[width * 4];     // histogram bin increments for this row

    // SIMD loop for this row - calculate flattened histogram indices and bin
    // increments (scalar code shown for reference - converting this loop to
    // SIMD is left as an exercise for the reader...)

    for (x = 0; x < width; ++x, ++i)
    {
        indices[4*x]   = (y/8)*(width/8)*18+(x/8)*18+idx[i];
        indices[4*x+1] = (y/8)*(width/8)*18+(x/8 + 1)*18+idx[i];
        indices[4*x+2] = (y/8+1)*(width/8)*18+(x/8)*18+idx[i];
        indices[4*x+3] = (y/8+1)*(width/8)*18+(x/8 + 1)*18+idx[i];

        vals[4*x]   = val[i]*ky[y]*kx[x];
        vals[4*x+1] = val[i]*ky[y]*kx[x+1];
        vals[4*x+2] = val[i]*ky[y+1]*kx[x];
        vals[4*x+3] = val[i]*ky[y+1]*kx[x+1];
    }

    // scalar loop for this row

    float * const histogram_base = &histogram[0][0][0]; // pointer to flattened histogram

    for (x = 0; x < width * 4; ++x) // for each set of 4 indices/increments in this row
    {
        histogram_base[indices[x]] += vals[x];  // update the (flattened) histogram
    }

}

Answer 2

我找到了解决办法。它是一个临时缓冲区。首先，我们将直方图和临时缓冲区相加（这个操作可以矢量化）。然后，我们将缓冲区和添加到输出直方图中（此操作也可以矢量化）：

float histogram[height/8][width/8][18];
float ky[height], kx[width];
int idx[size];
float val[size]; 
float buf[18][4];

for(size_t i = 0; i < size; ++i)
{
    buf[idx[i]][0] += val[i]*ky[y]*kx[x];
    buf[idx[i]][1] += val[i]*ky[y]*kx[x + 1];
    buf[idx[i]][2] += val[i]*ky[y + 1]*kx[x];
    buf[idx[i]][3] += val[i]*ky[y + 1]*kx[x + 1];
}

for(size_t i = 0; i < 18; ++i)
{
    histogram[y/8][x/8][i] += buf[i][0];
    histogram[y/8][x/8 + 1][i] += buf[i][1];
    histogram[y/8 + 1][x/8][i] += buf[i][2];
    histogram[y/8 + 1][x/8 + 1][i] += buf[i][3];
}