contrib/tools/makesRGB.c

   1 /* makesRGB.c -- build sRGB-to-linear and linear-to-sRGB conversion tables
   2  *
   3  * COPYRIGHT: Written by John Cunningham Bowler, 2013.
   4  * To the extent possible under law, the author has waived all copyright and
   5  * related or neighboring rights to this work.  This work is published from:
   6  * United States.
   7  *
   8  * Make a table to convert 8-bit sRGB encoding values into the closest 16-bit
   9  * linear value.
  10  *
  11  * Make two tables to take a linear value scaled to 255*65535 and return an
  12  * approximation to the 8-bit sRGB encoded value.  Calculate the error in these
  13  * tables and display it.
  14  */
  15
  16 #define _C99_SOURCE 1
  17 #include <stdio.h>
  18 #include <math.h>
  19 #include <stdlib.h>
  20
  21 /* pngpriv.h includes the definition of 'PNG_sRGB_FROM_LINEAR' which is required
  22  * to verify the actual code.
  23  */
  24 #include "../../pngpriv.h"
  25
  26 #include "sRGB.h"
  27
  28 /* The tables are declared 'const' in pngpriv.h, so this redefines the tables to
  29  * be used.
  30  */
  31 #define png_sRGB_table sRGB_table
  32 #define png_sRGB_base sRGB_base
  33 #define png_sRGB_delta sRGB_delta
  34
  35 static png_uint_16 png_sRGB_table[256];
  36 static png_uint_16 png_sRGB_base[512];
  37 static png_byte png_sRGB_delta[512];
  38
  39 static const unsigned int max_input = 255*65535;
  40
  41 double
  42 fsRGB(double l)
  43 {
  44    return sRGB_from_linear(l/max_input);
  45 }
  46
  47 double
  48 sRGB(unsigned int i)
  49 {
  50    return fsRGB(i);
  51 }
  52
  53 double
  54 finvsRGB(unsigned int i)
  55 {
  56    return 65535 * linear_from_sRGB(i/255.);
  57 }
  58
  59 png_uint_16
  60 invsRGB(unsigned int i)
  61 {
  62    unsigned int x = nearbyint(finvsRGB(i));
  63
  64    if (x > 65535)
  65    {
  66       fprintf(stderr, "invsRGB(%u) overflows to %u\n", i, x);
  67       exit(1);
  68    }
  69
  70    return (png_uint_16)x;
  71 }
  72
  73 int
  74 main(int argc, char **argv)
  75 {
  76    unsigned int i, i16, ibase;
  77    double min_error = 0;
  78    double max_error = 0;
  79    double min_error16 = 0;
  80    double max_error16 = 0;
  81    double adjust;
  82    double adjust_lo = 0.4, adjust_hi = 0.6, adjust_mid = 0.5;
  83    unsigned int ec_lo = 0, ec_hi = 0, ec_mid = 0;
  84    unsigned int error_count = 0;
  85    unsigned int error_count16 = 0;
  86    int test_only = 0;
  87
  88    if (argc > 1)
  89       test_only = strcmp("--test", argv[1]) == 0;
  90
  91    /* Initialize the encoding table first. */
  92    for (i=0; i<256; ++i)
  93    {
  94       png_sRGB_table[i] = invsRGB(i);
  95    }
  96
  97    /* Now work out the decoding tables (this is where the error comes in because
  98     * there are 512 set points and 512 straight lines between them.)
  99     */
 100    for (;;)
 101    {
 102       if (ec_lo == 0)
 103          adjust = adjust_lo;
 104
 105       else if (ec_hi == 0)
 106          adjust = adjust_hi;
 107
 108       else if (ec_mid == 0)
 109          adjust = adjust_mid;
 110
 111       else if (ec_mid < ec_hi)
 112          adjust = (adjust_mid + adjust_hi)/2;
 113
 114       else if (ec_mid < ec_lo)
 115          adjust = (adjust_mid + adjust_lo)/2;
 116
 117       else
 118       {
 119          fprintf(stderr, "not reached: %u .. %u .. %u\n", ec_lo, ec_mid, ec_hi);
 120          exit(1);
 121       }
 122
 123       /* Calculate the table using the current 'adjust' */
 124       for (i=0; i<=511; ++i)
 125       {
 126          double lo = 255 * sRGB(i << 15);
 127          double hi = 255 * sRGB((i+1) << 15);
 128          unsigned int calc;
 129
 130          calc = nearbyint((lo+adjust) * 256);
 131          if (calc > 65535)
 132          {
 133             fprintf(stderr, "table[%d][0]: overflow %08x (%d)\n", i, calc,
 134                calc);
 135             exit(1);
 136          }
 137          png_sRGB_base[i] = calc;
 138
 139          calc = nearbyint((hi-lo) * 32);
 140          if (calc > 255)
 141          {
 142             fprintf(stderr, "table[%d][1]: overflow %08x (%d)\n", i, calc,
 143                calc);
 144             exit(1);
 145          }
 146          png_sRGB_delta[i] = calc;
 147       }
 148
 149       /* Check the 16-bit linear values alone: */
 150       error_count16 = 0;
 151       for (i16=0; i16 <= 65535; ++i16)
 152       {
 153          unsigned int i = 255*i16;
 154          unsigned int iexact = nearbyint(255*sRGB(i));
 155          unsigned int icalc = PNG_sRGB_FROM_LINEAR(i);
 156
 157          if (icalc != iexact)
 158             ++error_count16;
 159       }
 160
 161       /* Now try changing the adjustment. */
 162       if (ec_lo == 0)
 163          ec_lo = error_count16;
 164
 165       else if (ec_hi == 0)
 166          ec_hi = error_count16;
 167
 168       else if (ec_mid == 0)
 169       {
 170          ec_mid = error_count16;
 171          printf("/* initial error counts: %u .. %u .. %u */\n", ec_lo, ec_mid,
 172             ec_hi);
 173       }
 174
 175       else if (error_count16 < ec_mid)
 176       {
 177          printf("/* adjust (mid ): %f: %u -> %u */\n", adjust, ec_mid,
 178             error_count16);
 179          ec_mid = error_count16;
 180          adjust_mid = adjust;
 181       }
 182
 183       else if (adjust < adjust_mid && error_count16 < ec_lo)
 184       {
 185          printf("/* adjust (low ): %f: %u -> %u */\n", adjust, ec_lo,
 186             error_count16);
 187          ec_lo = error_count16;
 188          adjust_lo = adjust;
 189       }
 190
 191       else if (adjust > adjust_mid && error_count16 < ec_hi)
 192       {
 193          printf("/* adjust (high): %f: %u -> %u */\n", adjust, ec_hi,
 194             error_count16);
 195          ec_hi = error_count16;
 196          adjust_hi = adjust;
 197       }
 198
 199       else
 200       {
 201          adjust = adjust_mid;
 202          printf("/* adjust: %f: %u */\n", adjust, ec_mid);
 203          break;
 204       }
 205    }
 206
 207    /* For each entry in the table try to adjust it to minimize the error count
 208     * in that entry.  Each entry corresponds to 128 input values.
 209     */
 210    for (ibase=0; ibase<65536; ibase+=128)
 211    {
 212       png_uint_16 base = png_sRGB_base[ibase >> 7], trybase = base, ob=base;
 213       png_byte delta = png_sRGB_delta[ibase >> 7], trydelta = delta, od=delta;
 214       unsigned int ecbase = 0, eco;
 215
 216       for (;;)
 217       {
 218          png_sRGB_base[ibase >> 7] = trybase;
 219          png_sRGB_delta[ibase >> 7] = trydelta;
 220
 221          /* Check the 16-bit linear values alone: */
 222          error_count16 = 0;
 223          for (i16=ibase; i16 < ibase+128; ++i16)
 224          {
 225             unsigned int i = 255*i16;
 226             unsigned int iexact = nearbyint(255*sRGB(i));
 227             unsigned int icalc = PNG_sRGB_FROM_LINEAR(i);
 228
 229             if (icalc != iexact)
 230                ++error_count16;
 231          }
 232
 233          if (error_count16 == 0)
 234             break;
 235
 236          if (ecbase == 0)
 237          {
 238             eco = ecbase = error_count16;
 239             ++trybase; /* First test */
 240          }
 241
 242          else if (error_count16 < ecbase)
 243          {
 244             if (trybase > base)
 245             {
 246                base = trybase;
 247                ++trybase;
 248             }
 249             else if (trybase < base)
 250             {
 251                base = trybase;
 252                --trybase;
 253             }
 254             else if (trydelta > delta)
 255             {
 256                delta = trydelta;
 257                ++trydelta;
 258             }
 259             else if (trydelta < delta)
 260             {
 261                delta = trydelta;
 262                --trydelta;
 263             }
 264             else
 265             {
 266                fprintf(stderr, "makesRGB: impossible\n");
 267                exit(1);
 268             }
 269             ecbase = error_count16;
 270          }
 271
 272          else
 273          {
 274             if (trybase > base)
 275                trybase = base-1;
 276             else if (trybase < base)
 277             {
 278                trybase = base;
 279                ++trydelta;
 280             }
 281             else if (trydelta > delta)
 282                trydelta = delta-1;
 283             else if (trydelta < delta)
 284                break; /* end of tests */
 285          }
 286       }
 287
 288       png_sRGB_base[ibase >> 7] = base;
 289       png_sRGB_delta[ibase >> 7] = delta;
 290       if (base != ob || delta != od)
 291       {
 292          printf("/* table[%u]={%u,%u} -> {%u,%u} %u -> %u errors */\n",
 293             ibase>>7, ob, od, base, delta, eco, ecbase);
 294       }
 295       else if (0)
 296          printf("/* table[%u]={%u,%u} %u errors */\n", ibase>>7, ob, od,
 297             ecbase);
 298    }
 299
 300    /* Only do the full (slow) test at the end: */
 301    min_error = -.4999;
 302    max_error = .4999;
 303    error_count = 0;
 304
 305    for (i=0; i <= max_input; ++i)
 306    {
 307       unsigned int iexact = nearbyint(255*sRGB(i));
 308       unsigned int icalc = PNG_sRGB_FROM_LINEAR(i);
 309
 310       if (icalc != iexact)
 311       {
 312          double err = 255*sRGB(i) - icalc;
 313
 314          if (err > (max_error+.001) || err < (min_error-.001))
 315          {
 316             printf(
 317                "/* 0x%08x: exact: %3d, got: %3d [tables: %08x, %08x] (%f) */\n",
 318                i, iexact, icalc, png_sRGB_base[i>>15],
 319                png_sRGB_delta[i>>15], err);
 320          }
 321
 322          ++error_count;
 323          if (err > max_error)
 324             max_error = err;
 325          else if (err < min_error)
 326             min_error = err;
 327       }
 328    }
 329
 330    /* Re-check the 16-bit cases too, including the warning if there is an error
 331     * bigger than 1.
 332     */
 333    error_count16 = 0;
 334    max_error16 = 0;
 335    min_error16 = 0;
 336    for (i16=0; i16 <= 65535; ++i16)
 337    {
 338       unsigned int i = 255*i16;
 339       unsigned int iexact = nearbyint(255*sRGB(i));
 340       unsigned int icalc = PNG_sRGB_FROM_LINEAR(i);
 341
 342       if (icalc != iexact)
 343       {
 344          double err = 255*sRGB(i) - icalc;
 345
 346          ++error_count16;
 347          if (err > max_error16)
 348             max_error16 = err;
 349          else if (err < min_error16)
 350             min_error16 = err;
 351
 352          if (abs(icalc - iexact) > 1)
 353             printf(
 354                "/* 0x%04x: exact: %3d, got: %3d [tables: %08x, %08x] (%f) */\n",
 355                i16, iexact, icalc, png_sRGB_base[i>>15],
 356                png_sRGB_delta[i>>15], err);
 357       }
 358    }
 359
 360    /* Check the round trip for each 8-bit sRGB value. */
 361    for (i16=0; i16 <= 255; ++i16)
 362    {
 363       unsigned int i = 255 * png_sRGB_table[i16];
 364       unsigned int iexact = nearbyint(255*sRGB(i));
 365       unsigned int icalc = PNG_sRGB_FROM_LINEAR(i);
 366
 367       if (i16 != iexact)
 368       {
 369          fprintf(stderr, "8-bit rounding error: %d -> %d\n", i16, iexact);
 370          exit(1);
 371       }
 372
 373       if (icalc != i16)
 374       {
 375          double finv = finvsRGB(i16);
 376
 377          printf("/* 8-bit roundtrip error: %d -> %f -> %d(%f) */\n",
 378             i16, finv, icalc, fsRGB(255*finv));
 379       }
 380    }
 381
 382
 383    printf("/* error: %g - %g, %u (%g%%) of readings inexact */\n",
 384       min_error, max_error, error_count, (100.*error_count)/max_input);
 385    printf("/* 16-bit error: %g - %g, %u (%g%%) of readings inexact */\n",
 386       min_error16, max_error16, error_count16, (100.*error_count16)/65535);
 387
 388    if (!test_only)
 389    {
 390       printf("const png_uint_16 png_sRGB_table[256] =\n{\n   ");
 391       for (i=0; i<255; )
 392       {
 393          do
 394          {
 395             printf("%d,", png_sRGB_table[i++]);
 396          }
 397          while ((i & 0x7) != 0 && i<255);
 398          if (i<255) printf("\n   ");
 399       }
 400       printf("%d\n};\n\n", png_sRGB_table[i]);
 401
 402
 403       printf("const png_uint_16 png_sRGB_base[512] =\n{\n   ");
 404       for (i=0; i<511; )
 405       {
 406          do
 407          {
 408             printf("%d,", png_sRGB_base[i++]);
 409          }
 410          while ((i & 0x7) != 0 && i<511);
 411          if (i<511) printf("\n   ");
 412       }
 413       printf("%d\n};\n\n", png_sRGB_base[i]);
 414
 415       printf("const png_byte png_sRGB_delta[512] =\n{\n   ");
 416       for (i=0; i<511; )
 417       {
 418          do
 419          {
 420             printf("%d,", png_sRGB_delta[i++]);
 421          }
 422          while ((i & 0xf) != 0 && i<511);
 423          if (i<511) printf("\n   ");
 424       }
 425       printf("%d\n};\n\n", png_sRGB_delta[i]);
 426    }
 427
 428    return 0;
 429 }