Clamp with desaturation (partial code)

clamp_with_desaturation.hh

// RGB to YUV conversion matrix
#define LUMA_COEFFICIENTS      0.29900,  0.58700,  0.11400
#define CHROMA1_COEFFICIENTS  -0.14713, -0.28886,  0.43600
#define CHROMA2_COEFFICIENTS   0.61500, -0.51499, -0.10001
// YUV to RGB conversion matrix
#define R_YUV_COEFFICIENTS     1,  0.00000,  1.13983
#define G_YUV_COEFFICIENTS     1, -0.39465, -0.58060
#define B_YUV_COEFFICIENTS     1,  2.03211,  0.00000

template<typename VecType>
unsigned ClampWithDesaturationAndDithering(VecType rgb, unsigned x,unsigned y)
{
    constexpr float bayer4x4_f[16] =
    {
         0/16.f, 8/16.f, 1/16.f, 9/16.f,
        12/16.f, 4/16.f,13/16.f, 5/16.f,
         3/16.f,11/16.f, 2/16.f,10/16.f,
        15/16.f, 7/16.f,14/16.f, 6/16.f
    };

    float l = rgb.DotProduct(VecType( LUMA_COEFFICIENTS, nullptr ));
    //float l = rgb.HorizontalSum() / 3.f; // faster but less accurate
    if(l >= 255.f)
        return 0xFFFFFF;
    else if(__builtin_expect(l <= 0.f, false))
        return 0x000000;
    else
    {
        float s = 1.f;
        for(unsigned n=0; n<3; ++n)
            if(rgb[n] > 255.f)
                s = std::min(s, (l-255.f) / (l-rgb[n]));
            else if(__builtin_expect(rgb[n] < 0.f, false))
                s = std::min(s, l         / (l-rgb[n]));
        if(s != 1.f)
            rgb = (rgb - l) * s + l;
    }
    auto r = rgb + bayer4x4_f[(y%4)*4+(x%4)];
    return (unsigned(r[0]) << 16) + (unsigned(r[1]) << 8) + (unsigned(r[2]));
}