如何在块复制期间矢量化范围检查？

Question

我有以下功能：

void CopyImageBitsWithAlphaRGBA(unsigned char *dest, const unsigned char *src, int w, int stride, int h,
    unsigned char minredmask, unsigned char mingreenmask, unsigned char minbluemask, unsigned char maxredmask, unsigned char maxgreenmask, unsigned char maxbluemask)
{
    auto pend = src + w * h * 4;
    for (auto p = src; p < pend; p += 4, dest += 4)
    {
        dest[0] = p[0]; dest[1] = p[1]; dest[2] = p[2];
        if ((p[0] >= minredmask && p[0] <= maxredmask) || (p[1] >= mingreenmask && p[1] <= maxgreenmask) || (p[2] >= minbluemask && p[2] <= maxbluemask))
            dest[3] = 255;
        else
            dest[3] = 0;
    }
}

它的作用是将一个 32 位位图从一个内存块复制到另一个内存块，当像素颜色落在某个颜色范围内时，将 Alpha 通道设置为完全透明。

如何在 VC++ 2017 中使用 SSE/AVX？现在它没有生成矢量化代码。如果无法自动执行此操作，我可以使用哪些功能自己执行此操作？

因为真的，我想测试字节是否在一个范围内可能是最明显有用的操作之一，但我看不到任何内置函数来处理它。

Answer 1

我不认为你会得到一个编译器来自动矢量化，就像你可以手动使用英特尔的内在函数那样。 （错误，以及 我 无论如何都可以手动完成：P）。

可能一旦我们手动对其进行向量化，我们就可以看到如何使用以这种方式工作的标量代码来手持编译器，但我们确实需要将字节元素打包比较为 0/0xFF，而且很难编写一些东西在 C 中，编译器会很好地自动向量化。默认整数提升意味着大多数 C 表达式实际上会产生 32 位结果，即使您使用 uint8_t，这通常会欺骗编译器将 8 位元素解包为 32 位元素，从而在自动4 倍的吞吐量损失（每个寄存器的元素更少），like in @harold's small tweak to your source。

---

SSE/AVX（在 AVX512 之前）对 SIMD 整数进行有符号比较，而不是无符号。但是您可以通过减去 128 将事物范围转移到有符号的 -128..127。XOR（加无进位）在某些 CPU 上的效率稍高一些，因此您实际上只需对 0x80 进行 XOR 即可翻转高位。但从数学上讲，你是从 0..255 无符号值中减去 128，得到 -128..127 有符号值。

甚至仍然可以实现(x-min) < (max-min) 的“无符号比较技巧”。 （例如，detecting alphabetic ASCII characters）。作为奖励，我们可以将范围偏移烘焙到该减法中。如果x<min，它会环绕并变成一个大于max-min 的大值。这显然适用于无符号，但它确实适用于 SSE/AVX2 有符号比较指令（使用范围移位 max-min）。 （此答案的先前版本声称此技巧仅在max-min < 128 时才有效，但事实并非如此。x-min 不能一直环绕并低于max-min，或者如果它开始则进入该范围以上max)。

此答案的早期版本具有使范围独占的代码，即不包括末端，因此即使 redmin=0 / redmax=255 也会排除 red=0 或 red= 的像素255.但我通过比较其他方式解决了这个问题（感谢@Nejc 和@chtz 的回答）。

@chtz 使用饱和加/减代替比较的想法非常酷。如果你安排的东西饱和意味着在范围内，它适用于一个包容的范围。 （并且您可以通过选择使所有 256 个可能的输入在范围内的最小值/最大值来将 Alpha 分量设置为已知值）。 这让我们避免将范围转移到有符号，因为无符号饱和是可用的

我们可以将 sub/cmp 范围检查与饱和技巧结合起来执行 sub（在越界低处回绕）/subs（只有在第一个 sub 没有回绕时才达到零)。然后我们不需要andnot 或or 来对每个组件进行两次单独的检查；我们已经在一个向量中有一个0 / 非零结果。

所以只需要两次操作就可以为我们可以检查的整个像素提供一个 32 位的值。如果所有 3 个 RGB 分量都在范围内，则该元素将具有特定值。 （因为我们已经安排 Alpha 组件也已经给出了一个已知值）。如果 3 个组件中的任何一个超出范围，它将具有其他值。

如果您以另一种方式执行此操作，则饱和意味着超出范围，那么您在该方向上有一个排他范围，因为您无法选择没有值达到 0 或达到 255 的限制。您可以始终使 alpha 分量饱和，以便在那里给自己一个已知值，而不管它对 RGB 分量意味着什么。独占范围会让您通过选择一个像素永远无法匹配的范围来滥用此功能，使其始终为假。 （或者如果有第三个条件，除了每个组件的最小/最大，那么也许你想要一个覆盖）。

---

显而易见的是使用 32 位元素大小的压缩比较指令 (_mm256_cmpeq_epi32 / vpcmpeqd) 生成0xFF 或0x00 （我们可以将其应用/混合到原始 RGB 像素值中）用于范围内/外。

// AVX2 core idea: wrapping-compare trick with saturation to achieve unsigned compare
__m256i tmp = _mm256_sub_epi8(src, min_values);       // wraps to high unsigned if below min
__m256i RGB_inrange = _mm256_subs_epu8(tmp, max_minus_min);  // unsigned saturation to 0 means in-range
__m256i new_alpha = _mm256_cmpeq_epi32(RGB_inrange, _mm256_setzero_si256());

// then blend the high byte of each element with RGB from the src vector
__m256i alpha_replaced = _mm256_blendv_epi8(new_alpha, src, _mm256_set1_epi32(0x00FFFFFF));  // alpha from new_alpha, RGB from src

请注意，SSE2 版本只需要一条 MOVDQA 指令即可复制src；同一个寄存器是每条指令的目的地。

还请注意，您可以使另一个方向饱和：add 然后adds（我认为是(256-max) 和(256-(min-max))）饱和到0xFF 以在范围内。这对于 AVX512BW 如果您使用零掩码 和固定掩码（例如，用于 alpha）或 可变掩码（对于某些其他条件） 来排除基于其他一些情况。 sub/subs 版本的 AVX512BW 零掩码会考虑组件在范围内，即使它们不在，这也很有用。

---

但是将其扩展到 AVX512 需要一种不同的方法：AVX512 比较产生一个位掩码（在掩码寄存器中），而不是向量，所以我们不能转而使用高字节分别对每个 32 位比较结果。

代替cmpeq_epi32，我们可以使用从左到右传播的减法进位/借位在每个像素的高字节中生成我们想要的值。

0x00000000 - 1 = 0xFFFFFFFF     # high byte = 0xFF = new alpha
0x00?????? - 1 = 0x00??????     # high byte = 0x00 = new alpha
Where ?????? has at least one non-zero bit, so it's a 32-bit number >=0 and <=0x00FFFFFFFF
Remember we choose an alpha range that makes the high byte always zero

即_mm256_sub_epi32(RGB_inrange, _mm_set1_epi32(1))。我们只需要每个 32 位元素的高字节来获得我们想要的 alpha 值，因为我们使用字节混合将它与源 RGB 值合并。对于 AVX512，这避免了 VPMOVM2D zmm1, k1 指令将比较结果转换回 0/-1 的向量，或者（更昂贵）将每个掩码位与 3 个零交错以将其用于字节混合。

这个sub 而不是cmp 即使对于 AVX2 也有一点优势：sub_epi32 在 Skylake 上的更多端口上运行（p0/p1/p5 与 p0/p1 的 pcmpgt/ pcmpeq）。在所有其他 CPU 上，向量整数加/减在与向量整数比较相同的端口上运行。 (Agner Fog's instruction tables)。

此外，如果您在带有 AVX512 的 CPU 上使用 -march=native 编译 _mm256_cmpeq_epi32()，或者启用 AVX512 然后编译正常的 AVX2 内部函数，一些编译器会愚蠢地使用 AVX512 compare-into-mask，然后展开回向量而不是仅仅使用 VEX 编码的vpcmpeqd。因此，我们使用sub 而不是cmp，即使对于_mm256 内在函数版本也是如此，因为我已经花时间弄清楚它并表明它在为常规AVX2 编译的正常情况下至少同样有效。 （虽然_mm256_setzero_si256() 比set1(1) 便宜；vpxor 可以廉价地将寄存器归零而不是加载常量，但是这种设置发生在循环之外。）

#include <immintrin.h>

#ifdef __AVX2__
// inclusive min and max
__m256i  setAlphaFromRangeCheck_AVX2(__m256i src, __m256i mins, __m256i max_minus_min)
{
    __m256i tmp = _mm256_sub_epi8(src, mins);   // out-of-range wraps to a high signed value

    // (x-min) <= (max-min)  equivalent to:
    // (x-min) - (max-min) saturates to zero
    __m256i RGB_inrange = _mm256_subs_epu8(tmp, max_minus_min);
    // 0x00000000 for in-range pixels, 0x00?????? (some higher value) otherwise

    // this has minor advantages over compare against zero, see full comments on Godbolt    
    __m256i new_alpha = _mm256_sub_epi32(RGB_inrange, _mm256_set1_epi32(1));
    // 0x00000000 - 1  = 0xFFFFFFFF
    // 0x00?????? - 1  = 0x00??????    high byte = new alpha value

    const __m256i RGB_mask = _mm256_set1_epi32(0x00FFFFFF);  // blend mask
    // without AVX512, the only byte-granularity blend is a 2-uop variable-blend with a control register
    // On Ryzen, it's only 1c latency, so probably 1 uop that can only run on one port.  (1c throughput).
    // For 256-bit, that's 2 uops of course.
    __m256i alpha_replaced = _mm256_blendv_epi8(new_alpha, src, RGB_mask);  // RGB from src, 0/FF from new_alpha

    return alpha_replaced;
}
#endif  // __AVX2__

为此函数设置向量参数并使用_mm256_load_si256 / _mm256_store_si256 循环遍历您的数组。 （如果不能保证对齐，则加载u/storeu。）

这个compiles very efficiently (Godbolt Compiler explorer) 带有 gcc、clang 和 MSVC。 （Godbolt 上的 AVX2 版本很好，AVX512 和 SSE 版本仍然是一团糟，还不是所有的技巧都适用于它们。）

;; MSVC's inner loop from a caller that loops over an array with it:
;; see the Godbolt link
$LL4@:
    vmovdqu ymm3, YMMWORD PTR [rdx+rax*4]
    vpsubb   ymm0, ymm3, ymm7
    vpsubusb ymm1, ymm0, ymm6
    vpsubd   ymm2, ymm1, ymm5
    vpblendvb ymm3, ymm2, ymm3, ymm4
    vmovdqu YMMWORD PTR [rcx+rax*4], ymm3
    add      eax, 8
    cmp      eax, r8d
    jb       SHORT $LL4@

因此，MSVC 在内联后设法提升了常量设置。我们从 gcc/clang 中得到了类似的循环。

循环有 4 个向量 ALU 指令，其中一个需要 2 个微指令。总共 5 个向量 ALU 微指令。但是 Haswell/Skylake 上的总融合域微指令 = 9，没有展开，所以幸运的是，它可以每 2.25 个时钟周期以 32 个字节（1 个向量）运行。使用 L1d 或 L2 高速缓存中的热数据可能接近于实际实现这一目标，但 L3 或内存将成为瓶颈。展开后，它可能会限制 L2 缓存带宽。

一个 AVX512 版本（也包含在 Godbolt 链接中），只需要 1 uop 来混合，并且每个周期可以在向量中运行得更快，因此使用 512 字节向量的速度是两倍以上。

Answer 2

这是使此功能与 SSE 指令一起使用的一种可能方法。我使用 SSE 而不是 AVX，因为我想让答案简单。一旦您了解了解决方案的工作原理，使用 AVX 内部函数重写函数应该不是什么大问题。

编辑：请注意，我的方法与PeterCordes 的方法非常相似，但他的代码应该更快，因为他使用的是 AVX。如果你想用 AVX 内部函数重写下面的函数，请将 step 值更改为 8。

void CopyImageBitsWithAlphaRGBA(
  unsigned char *dest,
  const unsigned char *src, int w, int stride, int h,
  unsigned char minred, unsigned char mingre, unsigned char minblu,
  unsigned char maxred, unsigned char maxgre, unsigned char maxblu)
{
  char low = 0x80; // -128
  char high = 0x7f; // 127
  char mnr = *(char*)(&minred) - low;
  char mng = *(char*)(&mingre) - low;
  char mnb = *(char*)(&minblu) - low;
  int32_t lowest = mnr | (mng << 8) | (mnb << 16) | (low << 24);

  char mxr = *(char*)(&maxred) - low;
  char mxg = *(char*)(&maxgre) - low;
  char mxb = *(char*)(&maxblu) - low;
  int32_t highest = mxr | (mxg << 8) | (mxb << 16) | (high << 24);

  // SSE
  int step = 4;
  int sse_width = (w / step)*step;

  for (int y = 0; y < h; ++y)
  {
    for (int x = 0; x < w; x += step)
    {
      if (x == sse_width)
      {
        x = w - step;
      }

      int ptr_offset = y * stride + x;
      const unsigned char* src_ptr = src + ptr_offset;
      unsigned char* dst_ptr = dest + ptr_offset;

      __m128i loaded = _mm_loadu_si128((__m128i*)src_ptr);

      // subtract 128 from every 8-bit int
      __m128i subtracted = _mm_sub_epi8(loaded, _mm_set1_epi8(low));

      // greater than top limit? 
      __m128i masks_hi = _mm_cmpgt_epi8(subtracted, _mm_set1_epi32(highest));

     // lower that bottom limit?
     __m128i masks_lo = _mm_cmplt_epi8(subtracted, _mm_set1_epi32(lowest));

     // perform OR operation on both masks
     __m128i combined = _mm_or_si128(masks_hi, masks_lo);

     // are 32-bit integers equal to zero?
     __m128i eqzer = _mm_cmpeq_epi32(combined, _mm_setzero_si128());

     __m128i shifted = _mm_slli_epi32(eqzer, 24);

    // EDIT: fixed a bug:
     __m128 alpha_unmasked = _mm_and_si128(loaded, _mm_set1_epi32(0x00ffffff));

     __m128i combined = _mm_or_si128(alpha_unmasked, shifted);

     _mm_storeu_si128((__m128i*)dst_ptr, combined);
    }
  }
}

编辑：正如 @PeterCordes 在 cmets 中所述，代码包含一个现已修复的错误。

Answer 3

基于@PeterCordes 解决方案，但将 shift+compare 替换为饱和减法和加法：

// mins_compl shall be [255-minR, 255-minG, 255-minB, 0]
// maxs       shall be [maxR, maxG, maxB, 0]
__m256i  setAlphaFromRangeCheck(__m256i src, __m256i mins_compl, __m256i maxs)
{
    __m256i in_lo = _mm256_adds_epu8(src, mins_compl); // is 255 iff src+mins_coml>=255, i.e. src>=mins
    __m256i in_hi = _mm256_subs_epu8(src, maxs);       // is 0 iff src - maxs <= 0, i.e., src <= maxs

    __m256i inbounds_components = _mm256_andnot_si256(in_hi, in_lo);
    // per-component mask, 0xff, iff (mins<=src && src<=maxs).
    // alpha-channel is always (~src & src) == 0

    // Use a 32-bit element compare to check that all 3 components are in-range
    __m256i RGB_mask = _mm256_set1_epi32(0x00FFFFFF);
    __m256i inbounds = _mm256_cmpeq_epi32(inbounds_components, RGB_mask);

    __m256i new_alpha = _mm256_slli_epi32(inbounds, 24);
    // alternatively _mm256_andnot_si256(RGB_mask, inbounds) ?

    // byte blends (vpblendvb) are at least 2 uops, and Haswell requires port5
    // instead clear alpha and then OR in the new alpha (0 or 0xFF)
    __m256i alphacleared = _mm256_and_si256(src, RGB_mask);   // off the critical path
    __m256i new_alpha_applied = _mm256_or_si256(alphacleared, new_alpha);

    return new_alpha_applied;
}

这节省了vpxor（无需修改src）和一个vpand（alpha 通道自动为0——我猜彼得的解决方案也可以通过相应地选择边界来实现） .

Godbolt-Link，显然，gcc 和 clang 都不认为将RGB_mask 用于这两种用途是值得的......

使用 SSE2 变体进行简单测试：https://wandbox.org/permlink/eVzFHljxfTX5HDcq（您可以玩弄源和边界）