什么系列的内在函数将完成这个 paeth 预测代码？

Question

我有一个在数组上运行的 Paeth 预测函数：

std::array<std::uint8_t,4> birunji::paeth_prediction
    (const std::array<std::uint8_t,4>& a,
     const std::array<std::uint8_t,4>& b,
     const std::array<std::uint8_t,4>& c)
{
    std::array<std::int16_t,4> pa;
    std::array<std::int16_t,4> pb;
    std::array<std::int16_t,4> pc;

    std::array<std::uint8_t,4> results;

    for(std::size_t i = 0; i < 4; ++i)
    {
        pa[i] = b[i] - c[i];
        pb[i] = a[i] - c[i];
        pc[i] = pa[i] + pb[i];

        pa[i] = std::abs(pa[i]);
        pb[i] = std::abs(pb[i]);
        pc[i] = std::abs(pc[i]);

        if(pa[i] <= pb[i] && pa[i] <= pc[i])
            results[i] = a[i];
        else if(pb[i] <= pc[i])
            results[i] = b[i];
        else
            results[i] = c[i];
    }

    return results;
}

我正在尝试手动使用内在函数来向量化代码（出于学习目的）。

__m128i birunji::paeth_prediction(const __m128i& a,
                                  const __m128i& b,
                                  const __m128i& c)
{
    __m128i pa = _mm_sub_epi16(b, c);
    __m128i pb = _mm_sub_epi16(a, c);
    __m128i pc = _mm_add_epi16(pa, pb);
    
    pa = _mm_abs_epi16(pa);
    pb = _mm_abs_epi16(pb);
    pc = _mm_abs_epi16(pc);

    __m128i pa_le_pb = _mm_cmpgt_epi16(pb, pa);
    __m128i pa_le_pc = _mm_cmpgt_epi16(pc, pa);
    __m128i pb_le_pc = _mm_cmpgt_epi16(pc, pb);

    return
    _mm_and_si128(_mm_and_si128(pa_le_pb, pa_le_pc),
                  _mm_and_si128(_mm_and_si128(pb_le_pc,b),a));
}

我遇到的问题是条件语句。我如何成功地将这些向量化？我不确定我上面的尝试是否正确。

Answer 1

_mm_cmpgt_epi16 可用于比较。请注意 _mm_cmpgt_epi16(a, b) = !(a <= b)，但是 _mm_cmpgt_epi16(b, a) != (a <= b)，因为它不是大于或等于比较，而是严格的大于比较。所以面具倒过来了，但这在这种情况下同样有用，不需要显式倒置。

这个函数不应该return一个条件本身，它应该select来自a和b和c根据条件。如果 SSE4.1 可用，_mm_blendv_epi8 可用于实现 selection。例如（未测试）：

__m128i paeth_prediction(__m128i a, __m128i b, __m128i c)
{
    __m128i pa = _mm_sub_epi16(b, c);
    __m128i pb = _mm_sub_epi16(a, c);
    __m128i pc = _mm_add_epi16(pa, pb);
    
    pa = _mm_abs_epi16(pa);
    pb = _mm_abs_epi16(pb);
    pc = _mm_abs_epi16(pc);

    __m128i not_pa_le_pb = _mm_cmpgt_epi16(pa, pb);
    __m128i not_pa_le_pc = _mm_cmpgt_epi16(pa, pc);
    __m128i not_pb_le_pc = _mm_cmpgt_epi16(pb, pc);
    __m128i not_take_a = _mm_or_si128(not_pa_le_pb, not_pa_le_pc);
    __m128i t = _mm_blendv_epi8(b, c, not_pb_le_pc);
    return _mm_blendv_epi8(a, t, not_take_a);
}

最后两行实现的逻辑如下：

if PB is not less-than-or-equal-to PC, take C, otherwise take B.
if PA is not less-than-or-equal-to PB or PA is not less-than-or-equal-to PC, take the result from the previous step, otherwise take A.

如果没有 SSE4.1，混合可以使用 AND/ANDNOT/OR 实现。

我已经更改了函数的签名，因此它按值获取向量，不需要通过 const 引用获取它们（向量很容易复制）并且可以增加间接的开销，尽管这样的开销很可能如果函数最终被编译器内联，则被删除。

作为变体，_mm_min_epi16 可用于实现部分逻辑：

__m128i paeth_prediction(__m128i a, __m128i b, __m128i c)
{
    __m128i pa = _mm_sub_epi16(b, c);
    __m128i pb = _mm_sub_epi16(a, c);
    __m128i pc = _mm_add_epi16(pa, pb);
    
    pa = _mm_abs_epi16(pa);
    pb = _mm_abs_epi16(pb);
    pc = _mm_abs_epi16(pc);

    __m128i not_pb_le_pc = _mm_cmpgt_epi16(pb, pc);
    __m128i take_a = _mm_cmpeq_epi16(pa, _mm_min_epi16(pa, _mm_min_epi16(pb, pc)));
    __m128i t = _mm_blendv_epi8(b, c, not_pb_le_pc);
    return _mm_blendv_epi8(t, a, take_a);
}

因为条件pa <= pb && pa <= pc等价于pa == min(pa, pb, pc).

生成的汇编代码看起来好一点，但我没有以任何方式测试它，包括性能。

Answer 2

您可以通过完全避免任何转换为 int16_t 来简化您的计算。首先，请注意当且仅当 a<=c<=b 或 b<=c<=a 时 pa<=pc 和 pb<=pc 为真。如果 c 小于或等于两者，则返回 max(a,b)；如果 c 大于或等于，则返回 min(a,b)。

所以我们可以首先使用min和max操作“排序”a，b，

A = min(a,b)
B = max(a,b)

剩下三种可能的情况：

A<=B<=c  --> A
c<=A<=B  --> B
A< c< B  --> c

这意味着在 C++ 代码中

std::array<std::uint8_t,4> birunji::paeth_prediction
    (const std::array<std::uint8_t,4>& a,
     const std::array<std::uint8_t,4>& b,
     const std::array<std::uint8_t,4>& c)
{
    std::array<std::uint8_t,4> results;

    for(std::size_t i = 0; i < 4; ++i)
    {
        uint8_t A = std::min(a[i],b[i]);
        uint8_t B = std::max(a[i],b[i]);
        if     (B<=c[i]) results[i] = A;
        else if(c[i]<=A) results[i] = B;
        else             results[i] = c[i];
    }

    return results;
}

不幸的是，没有无符号 SIMD 比较（在 AVX-512 之前），但我们可以使用 (x<=y) == (max(x,y)==y)（或进行饱和减法并与零进行比较）来模拟它。

可能的（未经测试的）SIMD 实现（这也适用于任意多个元素——但您可以只加载最低 32 位中的四个元素并忽略结果的其余部分）：

__m128i paeth_prediction(__m128i a, __m128i b, __m128i c)
{
    __m128i A = _mm_min_epu8(a, b);
    __m128i B = _mm_max_epu8(a, b);

    __m128i A_greater_equal_c = _mm_cmpeq_epi8(_mm_max_epu8(A, c), A);
    __m128i B_less_equal_c    = _mm_cmpeq_epi8(_mm_min_epu8(B, c), B);

    // if you don't have SSE 4.1, this can be done using bitwise and/or operations:
    __m128i t = _mm_blendv_epi8(b, c, A_greater_equal_c);
    return _mm_blendv_epi8(a, t, B_less_equal_c);
}

什么系列的内在函数将完成这个 paeth 预测代码？

What series of intrinsics will complete this paeth prediction code?

c++

sse

intrinsics