#! /usr/bin/gawk -f
# Last edited on 2007-01-10 03:39:29 by stolfi

BEGIN {
  abort = -1;
  usage = ( \
    "complement-hue \\\n" \
    "  [ -v debug={BOOL} ] \\\n" \
    "  [ -v showlum={BOOL} ] \\\n" \
    "  < INFILE.txt \\\n" \
    "  > OUTFILE.txt" \
  ); 

  # Reads a list of colors, at most one per line. Outputs the 
  # colors with same luminosity Y, complementary hue, and
  # same relative saturation, also one per line.
  #
  # Each color is given as a triple "{R} {G} {B}" followed 
  # by an optional denominator "/ {D}".
  #
  # The program ignores '#'-comments, blank lines, commas, and 
  # parentheses of types "[]{}()<>".
  
  if (debug == "") { debug = 0; }
  if (showlum == "") { showlum = 0; }
}

(abort >= 0) { exit abort; }

/^[ ]*([\#]|$)/ { next; }

// {
  # Remove '#'-comments, if any:
  gsub(/[\#].*$/, "", $0);
  
  # Remove delimiters:
  gsub(/[][(){}<>,\011\015]/, " ", $0);
  
  # Parse RGB components:
  if ((NF != 3) && (NF != 5)) { data_error(("bad color = \"" $0 "\"\n")); }

  R = check_value($1);
  G = check_value($2); 
  B = check_value($3);
  
  # Parse and apply denominator {D}, if any, and save its text form {txtD}: 
  if ((NF >= 5) && ($4 == "/"))
    { txtD = $5;
      D = check_value($5);
      R /= D; G /= D; B /= D;
    }
  else
    { txtD = ""; D = 1; }
  
  # Clip to unit cube:
  R = clip_value(R);
  G = clip_value(G);
  B = clip_value(B);
  if (debug) { debug_color("input color", R, G, B); }
  
  # Compute luminosity {Y} (European TV standard formula):
  Y = lum(R,G,B);
  if (debug) { debug_color("equiv gray", Y, Y, Y); }
  if (debug) { debug_scalar("lum of equiv gray", lum(Y,Y,Y)); }
  
  if ((Y < 0.0001) || (Y > 0.9999))
    { # Not worth computing the exact answer: 
      nR = Y; nG = Y; nB = Y; 
    }
  else
    { # Compute the chroma vector:
      cR = R - Y;
      cG = G - Y;
      cB = B - Y;
      if (debug) { debug_color("chroma vector", cR, cG, cB); }
      if (debug) { debug_scalar("lum of chroma vec", lum(cR,cG,cB)); }
      if (debug) { debug_color("supplementary color", Y-cR, Y-cG, Y-cB); }
      if (debug) { debug_scalar("lum of supp color", lum(Y-cR,Y-cG,Y-cB)); }

      # Find the relative saturation {S} of {cR,cR,cG}:
      S = sat(R,G,B)
      if (debug) { debug_scalar("rel saturation (S)", "%7.5f", S); }
      
      assert((S >= 0) && (S <= 1)); 
      
      # Find the relative saturation T of {-cR,-cG,-cB}:
      T = sat(Y-cR,Y-cG,Y-cB);
      if (debug) { debug_scalar("cmp saturation (T)", "%7.5f", T); }
      
      assert(T >= 0);
      
      # Compute the scale factor {F} for {-cR,-cG,-cB}:
      if (T == 0) { assert(S < 0.000001); F = 0; } else { F = S/T; }
      
      # Apply the reversed chroma, scaled, to {Y,Y,Y}:
      nR = Y - F*cR;
      nG = Y - F*cG;
      nB = Y - F*cB;
      if (debug) { debug_color("result", nR, nG, nB); }
      if (debug) { debug_scalar("lum of result", lum(nR,nG,nB)); }
      if (debug) { debug_color("chroma of result", Y-nR, Y-nG, Y-nB); }
    }      

  # Just in case, clip to unit cube:
  nR = clip_value(nR);
  nG = clip_value(nG);
  nB = clip_value(nB);

  # Output result:
  if ((txtD == "") || (D <= 1.0))
    { printf "%5.3f %5.3f %5.3f", nR, nG, nB; }
  else if (D >= 1000)
    { printf "%5d %5d %5d / %s", round(nR*D), round(nG*D), round(nB*D), txtD; }
  else if (D >= 100)
    { printf "%5.1f %5.1f %5.1f / %s", nR*D, nG*D, nB*D, txtD; }
  else if (D >= 10)
    { printf "%5.2f %5.2f %5.2f / %s", nR*D, nG*D, nB*D, txtD; }
  else
    { printf "%5.3f %5.3f %5.3f / %s", nR*D, nG*D, nB*D, txtD; }
    
  # Output luminosity if requested:
  if (showlum) { printf " # Y = %6.4f", Y; }
  
  printf "\n";
}
  
function round(x  )
{
  # Rounds {x} to integer:
  return int(x + (x < 0 ? -0.5 : +0.5));
}
  
function lum(R,G,B)
{
  # Luminosity of {R,G,B} (European TV formula):
  return 0.298911*R + 0.586611*G + 0.114478*B;
}

function sat(R,G,B,  Y,S,sR,sG,sB)
{
  Y = lum(R,G,B);
  S = 0;
  sR = sat1(R-Y,Y); if (sR > S) { S = sR; }
  sG = sat1(G-Y,Y); if (sG > S) { S = sG; }
  sB = sat1(B-Y,Y); if (sB > S) { S = sB; }
  return S;
}

function sat1(C,Y  )
{
  # Saturation of {C} relative to the range {[-Y _ 1-Y]}
  return (C < 0 ? -C/Y : C/(1-Y));
}

function check_value(x,  v)
{
  # Checks format of {x}.
  if (x !~ /^[ ]*[-+]?[0-9]*([.][0-9]*|[0-9])[ ]*$/)
    { data_error(("bad number = \"" x "\"\n")); }
  v = x + 0.000;
  return v;
}

function clip_value(V,  eps)
{
  # Returns V clipped to range [0_1]; warns if outside.
  eps = 1.0e-5; # Fudge factor to allow for roundoff errors.
  if ((V < -eps) || (V > 1+eps))
    { data_warning(("value " V " not in unit cube")); }
  if (V < 0) { V = 0; } else if (V > 1) { V = 1; }
  return V;
}

function assert(cond)
{
  if (! cond) { prog_error("bug"); }
}

function data_warning(msg)
{
  printf "%s:%s: warning - %s\n", FILENAME, FNR, msg > "/dev/stderr";
}

function data_error(msg)
{
  printf "%s:%s: ** %s\n", FILENAME, FNR, msg > "/dev/stderr";
  abort = 1;
  exit abort;
}

function prog_error(msg)
{
  printf "** %s\n", msg > "/dev/stderr";
  abort = 1;
  exit abort;
}

function debug_color(lab,R,G,B)
{
  printf "%-20s = ( %+6.3f %+6.3f %+6.3f )\n", lab, R, G, B > "/dev/stderr";
}

function debug_scalar(lab,fmt,V)
{
  printf "%-20s = ", lab > "/dev/stderr";
  printf fmt, V > "/dev/stderr";
  printf "\n" > "/dev/stderr";
}