4

我正在使用 AVX2 优化视差估计算法的“赢家通吃”部分。我的标量例程是准确的，但在 QVGA 分辨率和 48 个视差下，运行时间在我的笔记本电脑上大约 14 毫秒，令人失望。我创建了 LR 和 RL 视差图像，但为了简单起见，我将只包含用于 RL 搜索的代码。

我的标量例程：

int MAXCOST = 32000;
for (int i = maskRadius; i < rstep-maskRadius; i++) {

    // WTA "RL" Search:
    for (int j = maskRadius; j+maskRadius < cstep; j++) {
        int minCost = MAXCOST;
        int minDisp = 0;
        for (int d = 0; d < numDisp && j+d < cstep; d++) {
            if (asPtr[(i*numDisp*cstep)+(d*cstep)+j] < minCost) {
                minCost = asPtr[(i*numDisp*cstep)+(d*cstep)+j];
                minDisp = d;
            }
        }
        dRPtr[(i*cstep)+j] = minDisp;
    }
}

我尝试使用 AVX2：

int MAXCOST = 32000;
int* dispVals = (int*) _mm_malloc( sizeof(int32_t)*16, 32 );

for (int i = maskRadius; i < rstep-maskRadius; i++) {

    // WTA "RL" Search AVX2:
    for( int j = 0; j < cstep-16; j+=16) {

        __m256i minCosts = _mm256_set1_epi16( MAXCOST );
        __m128i loMask   = _mm_setzero_si128();
        __m128i hiMask   = _mm_setzero_si128();

        for (int d = 0; d < numDisp && j+d < cstep; d++) {
            // Grab 16 costs to compare
            __m256i costs = _mm256_loadu_si256((__m256i*) (asPtr[(i*numDisp*cstep)+(d*cstep)+j]));

            // Get the new minimums
            __m256i newMinCosts = _mm256_min_epu16( minCosts, costs );

            // Compare new mins to old to build mask to store minDisps
            __m256i mask   = _mm256_cmpgt_epi16( minCosts, newMinCosts );
            __m128i loMask = _mm256_extracti128_si256( mask, 0 );
            __m128i hiMask = _mm256_extracti128_si256( mask, 1 );
            // Sign extend to 32bits
            __m256i loMask32 = _mm256_cvtepi16_epi32( loMask );
            __m256i hiMask32 = _mm256_cvtepi16_epi32( hiMask );

            __m256i currentDisp = _mm256_set1_epi32( d );
            // store min disps with mask
            _mm256_maskstore_epi32( dispVals, loMask32, currentDisp );    // RT error, why?
            _mm256_maskstore_epi32( dispVals+8, hiMask32, currentDisp );  // RT error, why?

            // Set minCosts to newMinCosts
            minCosts = newMinCosts;
        }

        // Write the WTA minimums one-by-one to the RL disparity image
        int index = (i*cstep)+j;
        for( int k = 0; k < 16; k++ ) {
            dRPtr[index+k] = dispVals[k];
        }
    }
}
_mm_free( dispVals );

视差空间图像 (DSI) 的大小为 HxWxD (320x240x48)，我将其水平布局以便更好地访问内存，这样每一行的大小都是 WxD。

视差空间图像具有每像素匹配成本。这与一个简单的盒子过滤器聚合在一起，以制作另一张完全相同大小的图像，但成本总和超过了一个 3x3 或 5x5 的窗口。这种平滑使结果更加“稳健”。当我使用 asPtr 访问时，我正在索引这个聚合成本图像。

此外，为了节省不必要的计算，我一直在以掩码半径偏移的行开始和结束。这个掩码半径是我的人口普查掩码的半径。我可以做一些花哨的边界反射，但它更简单更快，只是为了不打扰这个边界的差异。这当然也适用于开始和结束列，但是当我强制我的整个算法只在列是 16 的倍数（例如 QVGA：320x240）的图像上运行时，在这里搞乱索引并不好，这样我就可以简单地索引并使用 SIMD 命中所有内容（无残留标量处理）。

此外，如果您认为我的代码一团糟，我鼓励您查看高度优化的 OpenCV 立体算法。我发现它们是不可能的，并且几乎没有使用它们。

我的代码编译但在运行时失败。我正在使用 VS 2012 Express Update 4。当我使用调试器运行时，我无法获得任何见解。我对使用内在函数相对较新，所以我不确定调试时应该看到哪些信息、寄存器数量、__m256i 变量是否应该可见等。

听从下面的评论建议，我通过使用更智能的索引将标量时间从 ~14 提高到 ~8。我的 CPU 是 i7-4980HQ，我在同一个文件的其他地方成功使用了 AVX2 内部函数。

信息图片

4

2 回答 2

2

我仍然没有发现问题，但我确实看到了一些您可能想要更改的内容。但是，您没有检查 的返回值_mm_malloc。如果它失败了，那就可以解释了。（也许它不喜欢分配 32 字节对齐的内存？）

如果您在内存检查器或其他东西下运行代码，那么它可能不喜欢从未初始化的内存中读取dispVals. （_mm256_maskstore_epi32即使掩码是全一，也可能算作读取-修改-写入。）

在调试器下运行您的代码并找出问题所在。“运行时错误”不是很有意义。

_mm_set1*功能很慢。 VPBROADCASTD需要它的源在内存或向量 reg，而不是 GP reg，因此编译器可以movd从 GP reg 到向量 reg 然后广播，或者存储到内存然后广播。无论如何，这样做会更快

const __m256i add1 = _mm256_set1_epi32( 1 );
__m256i dvec = _mm256_setzero_si256();
for (d;d...;d++) {
    dvec = _mm256_add_epi32(dvec, add1);
}

其他东西：
如果您不将内部循环的每次迭代都存储到内存中，这可能会运行得更快。使用混合指令 ( _mm256_blendv_epi8) 或类似指令来更新与最小成本相关的位移向量。混合 = 带有注册目标的蒙版移动。

此外，您的位移值应该适合 16b 整数，因此在找到它们之前不要将它们符号扩展为 32b。英特尔 CPU 可以动态将 16b 内存位置符号扩展到 gp 寄存器中，而不会降低速度（movsz与 一样快mov），所以很有可能。只需将您的dRPtr数组声明为uint16_t. 然后，您根本不需要向量代码中的符号扩展内容（更不用说在您的内部循环中了！）。希望_mm256_extracti128_si256( mask, 0 )编译为空，因为您想要的 128 已经是 low128，所以只需使用 reg 作为 src 的 src vmovsx，但仍然如此。

您还可以通过不先加载来保存指令（和融合域 uop）。（除非编译器足够聪明，不会忽略vmovdqu和使用vpminuw内存操作数，即使您使用了 load 内在函数）。

所以我在想这样的事情：

// totally untested, didn't even check that this compiles.
for(i) { for(j) {
// inner loop, compiler can hoist these constants.
const __m256i add1 = _mm256_set1_epi16( 1 );
__m256i dvec = _mm256_setzero_si256();
__m256i minCosts = _mm256_set1_epi16( MAXCOST );
__m256i minDisps = _mm256_setzero_si256();

for (int d=0 ; d < numDisp && j+d < cstep ;
     d++, dvec = _mm256_add_epi16(dvec, add1))
{
    __m256i newMinCosts = _mm256_min_epu16( minCosts, asPtr[(i*numDisp*cstep)+(d*cstep)+j]) );
    __m256i mask   = _mm256_cmpgt_epi16( minCosts, newMinCosts );
    minDisps = _mm256_blendv_epi8(minDisps, dvec, mask); // 2 uops, latency=2
    minCosts = newMinCosts;
}

// put sign extension here if making dRPtr uint16_t isn't an option.
int index = (i*cstep)+j;
_mm256_storeu_si256 (dRPtr + index, __m256i minDisps);
}}

minCosts0拥有两个并行依赖链： /minDisps0和minCosts1/ minDisps1，然后在最后将它们组合起来， 您可能会获得更好的性能。minDisps是一个循环携带的依赖项，但循环只有 5 条指令（包括vpadd，它看起来像循环开销，但不能通过展开来减少）。它们解码为 6 微指令（blendv 为 2），加上循环开销。它应该在 haswell 上以 1.5 个周期/迭代运行（不计算循环开销），但 dep 链会将其限制为每 2 个周期进行一次迭代。（假设展开以摆脱循环开销）。并行执行两个 dep 链可以解决此问题，并且与展开循环具有相同的效果：减少循环开销。

嗯，实际上在 Haswell 上，

  • pminuw可以在 p1/p5 上运行。（以及 p2/p3 上的负载部分）
  • pcmpgtw可以在 p1/p5 上运行
  • vpblendvbp5 是 2 微秒。
  • padduw可以在 p1/p5 上运行
  • movdqa reg,reg可以在 p0/p1/p5 上运行（并且可能根本不需要执行单元）。展开应该消除 的任何开销minCosts = newMinCosts，因为编译器可以newMinCosts从右寄存器中最后展开的循环体结束，用于下一次迭代的第一个循环体。
  • fused sub/ jge（循环计数器）可以在 p6 上运行。（在 dvec 上使用PTEST+jcc会更慢）。add/sub未与 . 融合时可以在 p0/p1/p5/p6 上运行jcc。

好的，所以实际上循环每次迭代需要 2.5 个周期，受只能在 p1/p5 上运行的指令的限制。展开 2 或 4 将减少循环/movdqa开销。由于 Haswell 可以每个时钟发出 4 个微指令，因此它可以更有效地将微指令排队以进行乱序执行，因为循环不会有超多的迭代次数。（48 是你的例子。）有很多 uops 排队会让 CPU 在离开循环后有事可做，并隐藏缓存未命中的任何延迟等。

_mm256_min_epu16( PMINUW) 是另一个循环携带的依赖链。将它与内存操作数一起使用会导致 3 或 4 个周期的延迟。但是，一旦知道地址，指令的加载部分就可以开始，因此将加载折叠到修改操作中以利用微融合不会使 dep 链比使用单独的加载更长或更短。

有时您需要对未对齐的数据使用单独的加载（AVX 删除了内存操作数的对齐要求）。我们受执行单元的限制比 4 uop / 时钟问题限制更多，因此使用专用加载指令可能没问题。

insn 端口/延迟的来源。

于 2015-06-08T05:56:06.620 回答
2

在进行特定于平台的优化之前，可以执行许多可移植的优化。提取循环不变量，将索引乘法转换为增量加法等...

这可能不准确，但可以大致了解：

int MAXCOST = 32000, numDispXcstep = numDisp*cstep;
for (int i = maskRadius; i < rstep - maskRadius; i+=numDispXcstep) {
    for (int j = maskRadius; j < cstep - maskRadius; j++) {
        int minCost = MAXCOST, minDisp = 0;
        for (int d = 0; d < numDispXcstep - j; d+=cstep) {
            if (asPtr[i+j+d] < minCost) {
                minCost = asPtr[i+j+d];
                minDisp = d;
            }
        }
        dRPtr[i/numDisp+j] = minDisp;
    }
}

一旦你这样做了，你就会明白实际发生了什么。看起来“i”是最大的步骤，其次是“d”，“j”实际上是对顺序数据进行操作的变量。...下一步将相应地重新排序循环，如果您仍需要进一步优化，请应用特定于平台的内在函数。

于 2015-06-08T08:36:03.157 回答