Monday 20 August 2012

RGB/HCL in HLSL

The HCL colour space by M. Sarifuddin and Rokia Missaoui (not to be confused with CIELCH) is another colour space that tries to improve on HSL and HSV. It is a cylindrical space which means that, for an accurate implementation, trigonometric functions are necessary. One interesting feature is that it tries to adjust the hue metric to be more perceptually meaningful using piecewise linear interpolation.

Here's the optimised HLSL code to convert from linear RGB:

float HCLgamma = 3;
float HCLy0 = 100;
// HCLmaxL == exp(HCLgamma / HCLy0) - 0.5
float HCLmaxL = 0.530454533953517;

float3 RGBtoHCL(in float3 RGB)
{
  float3 HCL;
  float H = 0;
  float U, V;
#if NO_ASM
  U = -min(RGB.r, min(RGB.g, RGB.b));
  V = max(RGB.r, max(RGB.g, RGB.b));
#else
  float4 RGB4 = RGB.rgbr;
  asm { max4 U, -RGB4 };
  asm { max4 V, RGB4 };
#endif
  float Q = HCLgamma / HCLy0;
  HCL.y = V + U;
  if (HCL.y != 0)
  {
    H = atan2(RGB.g - RGB.b, RGB.r - RGB.g) / PI;
    Q *= -U / V;
  }
  Q = exp(Q);
  HCL.x = frac(H / 2 - min(frac(H), frac(-H)) / 6);
  HCL.y *= Q;
  HCL.z = lerp(U, V, Q) / (HCLmaxL * 2);
}


All components are scaled to fit into the expected [0,1] range.

The weird-looking statement with "frac()" terms is performing the piecewise adjustment of hue.

I'm very happy with this code as it is considerably faster than a simplistic transliteration of the reference code that Sarifuddin Madenda was good enough to send me.

However, the reverse transformation is far from perfect and needs more work to reduce the branching nature of the algorithm:

float HCLtoRGB(in float3 HCL)
{
  float3 RGB = 0;
  if (HCL.z != 0)
  {
    float H = HCL.x;
    float C = HCL.y;
    float L = HCL.z * HCLmaxL;
    float Q = exp((1 - C / (2 * L)) * (HCLgamma / HCLy0));
    float U = (2 * L - C) / (2 * Q - 1);
    float V = C / Q;
    float T = tan((H + min(frac(2 * H) / 4, frac(-2 * H) / 8)) * PI * 2);
    H *= 6;
    if (H <= 1)
    {
      RGB.r = 1;
      RGB.g = T / (1 + T);
    }
    else if (H <= 2)
    {
      RGB.r = (1 + T) / T;
      RGB.g = 1;
    }
    else if (H <= 3)
    {
      RGB.g = 1;
      RGB.b = 1 + T;
    }
    else if (H <= 4)
    {
      RGB.g = 1 / (1 + T);
      RGB.b = 1;
    }
    else if (H <= 5)
    {
      RGB.r = -1 / T;
      RGB.b = 1;
    }
    else
    {
      RGB.r = 1;
      RGB.b = -T;
    }
  }
  return RGB * V + U;
}


The multiple if-statements could be rationalised (binary search style) but I can't help thinking there's some simple trigonometric identity that can be utilised to eradicate them completely.

Monday 6 August 2012

sRGB Approximations for HLSL

The sRGB color space is non-linear. Many transformations to and from the space require the RGB components to be mapped to a linear form.

The "official" transformation for each RGB component value in the range [0,1] is:

  if (C_srgb <= 0.04045)
      C_lin = C_srgb / 12.92;
  else
      C_lin = pow((C_srgb + 0.055) / 1.055, 2.4);

This is often approximated using the "gamma 2.2" formula:

  C_lin_1 = pow(C_srgb, 2.2);

This works fine, but is fairly inaccurate. The graph below uses the right-hand axis for the absolute difference:


For example, if the values are quantized to eight bits, for the sRGB component value 197/255, the linear output is 145/255 instead of 142/255. Can we do better?

In fact, if we simply change the "magic number" to 2.233333... we get better results:

  C_lin_2 = pow(C_srgb, 2.233333333);


However, the "pow" functionality is either prohibitively expensive or non-existent on many platforms. So, if we limit ourselves to simple arithmetic, a good approximation I found is the cubic:

  C_lin_3 = 0.012522878 * C_srgb +
            0.682171111 * C_srgb * C_srgb +
            0.305306011 * C_srgb * C_srgb * C_srgb;



This can be computed in HLSL using:

  float3 RGB = sRGB * (sRGB * (sRGB * 0.305306011 + 0.682171111) + 0.012522878);

The reverse transformation (from linear to sRGB) is more problematic:

Again, the "official" transformation is piecewise:

  if (C_lin <= 0.0031308)
    C_srgb = C_lin * 12.92;
  else
    C_srgb = 1.055 * pow(C_lin, 1.0 / 2.4) - 0.055;

This is usually poorly approximated with the inverse of the computation of "C_lin_1":

  C_srgb_1 = pow(C_lin, 0.4545454545);


In fact, the linear portion of the official graph is tiny, so an almost-perfect approximation is:

  C_srgb_2 = max(1.055 * pow(C_lin, 0.416666667) - 0.055, 0);


The clamp ("max(..., 0)") is free on many platforms, but the formula does use the "pow" functionality. If we assume we only have square-root operations at our disposal, a good approximation I found was:

  C_srgb_3 = 0.585122381 * sqrt(C_lin) +
             0.783140355 * sqrt(sqrt(C_lin)) -
             0.368262736 * sqrt(sqrt(sqrt(C_lin)));


This can be computed in HLSL using:

  float3 S1 = sqrt(RGB);
  float3 S2 = sqrt(S1);
  float3 S3 = sqrt(S2);
  float3 sRGB = 0.585122381 * S1 + 0.783140355 * S2 - 0.368262736 * S3;

An even better approximation (at the cost of an additional 'mad') is:


  float3 S1 = sqrt(RGB);
  float3 S2 = sqrt(S1);
  float3 S3 = sqrt(S2);
  float3 sRGB = 0.662002687 * S1 + 0.684122060 * S2 - 0.323583601 * S3 - 0.0225411470 * RGB;

Depending on your platform architecture, this may be faster using multiplication by a constant matrix for the final step.

Sunday 5 August 2012

RGB/HCY in HLSL

The HCY colour space is a tractable hue/chroma/luminance scheme developed by Kuzma Shapran. It is ideal for pixel shaders, being only slightly more expensive that the HSV and HSL schemes. However, it tries to be more "meaningful" in terms of human perception.

The three components are:
  1. Hue (H) computed in the same manner as HSV and HSL;
  2. Chroma (C) computed as the scaled difference between the maximum unweighted RGB component and the minimum unweighted RGB component; and
  3. Luminance (Y) computed as the weighted sum of RGB components.
Note that the chroma is post-scaled so that the maximum weighted luminance for this hue is always one.

The HLSL conversions are as follows:
// The weights of RGB contributions to luminance.
// Should sum to unity.
float3 HCYwts = float3(0.299, 0.587, 0.114);

float3 HUEtoRGB(in float H)
{
  float R = abs(H * 6 - 3) - 1;
  float G = 2 - abs(H * 6 - 2);
  float B = 2 - abs(H * 6 - 4);
  return saturate(float3(R,G,B));
}

float RGBCVtoHUE(in float3 RGB, in float C, in float V)
{
    float3 Delta = (V - RGB) / C;
    Delta.rgb -= Delta.brg;
    Delta.rgb += float3(2,4,6);
    // NOTE 1
    Delta.brg = step(V, RGB) * Delta.brg;
    float H;
#if NO_ASM
    H = max(Delta.r, max(Delta.g, Delta.b));
#else
    float4 Delta4 = Delta.rgbr;
    asm { max4 H, Delta4 };
#endif
    return frac(H / 6);
}

float3 RGBtoHCY(in float3 RGB)
{
  float3 HCY = 0;
  float U, V;
#if NO_ASM
  U = -min(RGB.r, min(RGB.g, RGB.b));
  V = max(RGB.r, max(RGB.g, RGB.b));
#else
  float4 RGB4 = RGB.rgbr;
  asm { max4 U, -RGB4 };
  asm { max4 V, RGB4 };
#endif
  HCY.y = V + U;
  HCY.z = dot(RGB, HCYwts);
  if (HCY.y != 0)
  {
    HCY.x = RGBCVtoHUE(RGB, HCY.y, V);
    float Z = dot(HUEtoRGB(HCY.x), HCYwts);
    if (HCY.z > Z)
    {
      HCY.z = 1 - HCY.z;
      Z = 1 - Z;
    }
    HCY.y *= Z / HCY.z;
  }
  return HCY;
}

float3 HCYtoRGB(in float3 HCY)
{
  float RGB = HUEtoRGB(HCY.x);
  float Z = dot(RGB, HCYwts);
  if (HCY.z < Z)
  {
      HCY.y *= HCY.z / Z;
  }
  else if (Z < 1)
  {
      HCY.y *= (1 - HCY.z) / (1 - Z);
  }
  return (RGB - Z) * HCY.y + HCY.z;
}

I've folded the code into my web page on such conversions here.