[libpng17] Moved READ_GAMMA to pngrtran.c. This makes everything in pngrtrans.c

depend on READ_TRANSFORMS and moves all the transform gamma (READ_GAMMA)
code from elsewhere to png_tran.c.  There are no code changes.
There is one remaining use of the gamma (16-bit) code in the simplified
API in pngread.c but that is because of a long-standing libpng bug,
namely that the gamma corrected palette is no produced by
png_read_update_info (John Bowler).

1 files changed, 3 insertions, 910 deletions

diff --git a/png.c b/png.c
index f1462dd1a..8642abd76 100644
--- a/png.c
+++ b/png.c
@@ -689,13 +689,13 @@ png_get_copyright(png_const_structrp png_ptr)
 #else
 #  ifdef __STDC__
    return PNG_STRING_NEWLINE \
-     "libpng version 1.7.0beta63 - June 4, 2015" PNG_STRING_NEWLINE \
+     "libpng version 1.7.0beta63 - June 5, 2015" PNG_STRING_NEWLINE \
      "Copyright (c) 1998-2015 Glenn Randers-Pehrson" PNG_STRING_NEWLINE \
      "Copyright (c) 1996-1997 Andreas Dilger" PNG_STRING_NEWLINE \
      "Copyright (c) 1995-1996 Guy Eric Schalnat, Group 42, Inc." \
      PNG_STRING_NEWLINE;
 #  else
-      return "libpng version 1.7.0beta63 - June 4, 2015\
+      return "libpng version 1.7.0beta63 - June 5, 2015\
       Copyright (c) 1998-2015 Glenn Randers-Pehrson\
       Copyright (c) 1996-1997 Andreas Dilger\
       Copyright (c) 1995-1996 Guy Eric Schalnat, Group 42, Inc.";
@@ -3258,7 +3258,7 @@ png_muldiv(png_fixed_point_p res, png_fixed_point a, png_int_32 times,
 
    return 0;
 }
-#endif /* READ_GAMMA || INCH_CONVERSIONS */
+#endif /* GAMMA || INCH_CONVERSIONS */
 
 #ifdef PNG_GAMMA_SUPPORTED /* more fixed point functions for gamma */
 /* Calculate a reciprocal, return 0 on div-by-zero or overflow. */
@@ -3291,913 +3291,6 @@ png_gamma_significant(png_fixed_point gamma_val)
 }
 #endif
 
-#ifdef PNG_READ_GAMMA_SUPPORTED
-#ifndef PNG_FLOATING_ARITHMETIC_SUPPORTED
-/* A local convenience routine. */
-static png_fixed_point
-png_product2(png_fixed_point a, png_fixed_point b)
-{
-   /* The required result is 1/a * 1/b; the following preserves accuracy. */
-   png_fixed_point res;
-
-   if (png_muldiv(&res, a, b, 100000) != 0)
-      return res;
-
-   return 0; /* overflow */
-}
-#endif /* FLOATING_ARITHMETIC */
-
-/* The inverse of the above. */
-png_fixed_point
-png_reciprocal2(png_fixed_point a, png_fixed_point b)
-{
-   /* The required result is 1/a * 1/b; the following preserves accuracy. */
-#ifdef PNG_FLOATING_ARITHMETIC_SUPPORTED
-   if (a != 0 && b != 0)
-   {
-      double r = 1E15/a;
-      r /= b;
-      r = floor(r+.5);
-
-      if (r <= 2147483647. && r >= -2147483648.)
-         return (png_fixed_point)r;
-   }
-#else
-   /* This may overflow because the range of png_fixed_point isn't
-    * symmetric, but this API is only used for the product of file and
-    * screen gamma so it doesn't matter that the smallest number it can
-    * produce is 1/21474, not 1/100000
-    */
-   png_fixed_point res = png_product2(a, b);
-
-   if (res != 0)
-      return png_reciprocal(res);
-#endif
-
-   return 0; /* overflow */
-}
-#endif /* READ_GAMMA */
-
-#ifdef PNG_READ_GAMMA_SUPPORTED /* gamma table code */
-#ifndef PNG_FLOATING_ARITHMETIC_SUPPORTED
-/* Fixed point gamma.
- *
- * The code to calculate the tables used below can be found in the shell script
- * contrib/tools/intgamma.sh
- *
- * To calculate gamma this code implements fast log() and exp() calls using only
- * fixed point arithmetic.  This code has sufficient precision for either 8-bit
- * or 16-bit sample values.
- *
- * The tables used here were calculated using simple 'bc' programs, but C double
- * precision floating point arithmetic would work fine.
- *
- * 8-bit log table
- *   This is a table of -log(value/255)/log(2) for 'value' in the range 128 to
- *   255, so it's the base 2 logarithm of a normalized 8-bit floating point
- *   mantissa.  The numbers are 32-bit fractions.
- */
-static const png_uint_32
-png_8bit_l2[128] =
-{
-   4270715492U, 4222494797U, 4174646467U, 4127164793U, 4080044201U, 4033279239U,
-   3986864580U, 3940795015U, 3895065449U, 3849670902U, 3804606499U, 3759867474U,
-   3715449162U, 3671346997U, 3627556511U, 3584073329U, 3540893168U, 3498011834U,
-   3455425220U, 3413129301U, 3371120137U, 3329393864U, 3287946700U, 3246774933U,
-   3205874930U, 3165243125U, 3124876025U, 3084770202U, 3044922296U, 3005329011U,
-   2965987113U, 2926893432U, 2888044853U, 2849438323U, 2811070844U, 2772939474U,
-   2735041326U, 2697373562U, 2659933400U, 2622718104U, 2585724991U, 2548951424U,
-   2512394810U, 2476052606U, 2439922311U, 2404001468U, 2368287663U, 2332778523U,
-   2297471715U, 2262364947U, 2227455964U, 2192742551U, 2158222529U, 2123893754U,
-   2089754119U, 2055801552U, 2022034013U, 1988449497U, 1955046031U, 1921821672U,
-   1888774511U, 1855902668U, 1823204291U, 1790677560U, 1758320682U, 1726131893U,
-   1694109454U, 1662251657U, 1630556815U, 1599023271U, 1567649391U, 1536433567U,
-   1505374214U, 1474469770U, 1443718700U, 1413119487U, 1382670639U, 1352370686U,
-   1322218179U, 1292211689U, 1262349810U, 1232631153U, 1203054352U, 1173618059U,
-   1144320946U, 1115161701U, 1086139034U, 1057251672U, 1028498358U, 999877854U,
-   971388940U, 943030410U, 914801076U, 886699767U, 858725327U, 830876614U,
-   803152505U, 775551890U, 748073672U, 720716771U, 693480120U, 666362667U,
-   639363374U, 612481215U, 585715177U, 559064263U, 532527486U, 506103872U,
-   479792461U, 453592303U, 427502463U, 401522014U, 375650043U, 349885648U,
-   324227938U, 298676034U, 273229066U, 247886176U, 222646516U, 197509248U,
-   172473545U, 147538590U, 122703574U, 97967701U, 73330182U, 48790236U,
-   24347096U, 0U
-
-#if 0
-   /* The following are the values for 16-bit tables - these work fine for the
-    * 8-bit conversions but produce very slightly larger errors in the 16-bit
-    * log (about 1.2 as opposed to 0.7 absolute error in the final value).  To
-    * use these all the shifts below must be adjusted appropriately.
-    */
-   65166, 64430, 63700, 62976, 62257, 61543, 60835, 60132, 59434, 58741, 58054,
-   57371, 56693, 56020, 55352, 54689, 54030, 53375, 52726, 52080, 51439, 50803,
-   50170, 49542, 48918, 48298, 47682, 47070, 46462, 45858, 45257, 44661, 44068,
-   43479, 42894, 42312, 41733, 41159, 40587, 40020, 39455, 38894, 38336, 37782,
-   37230, 36682, 36137, 35595, 35057, 34521, 33988, 33459, 32932, 32408, 31887,
-   31369, 30854, 30341, 29832, 29325, 28820, 28319, 27820, 27324, 26830, 26339,
-   25850, 25364, 24880, 24399, 23920, 23444, 22970, 22499, 22029, 21562, 21098,
-   20636, 20175, 19718, 19262, 18808, 18357, 17908, 17461, 17016, 16573, 16132,
-   15694, 15257, 14822, 14390, 13959, 13530, 13103, 12678, 12255, 11834, 11415,
-   10997, 10582, 10168, 9756, 9346, 8937, 8531, 8126, 7723, 7321, 6921, 6523,
-   6127, 5732, 5339, 4947, 4557, 4169, 3782, 3397, 3014, 2632, 2251, 1872, 1495,
-   1119, 744, 372
-#endif
-};
-
-static png_int_32
-png_log8bit(unsigned int x)
-{
-   unsigned int lg2 = 0;
-   /* Each time 'x' is multiplied by 2, 1 must be subtracted off the final log,
-    * because the log is actually negate that means adding 1.  The final
-    * returned value thus has the range 0 (for 255 input) to 7.994 (for 1
-    * input), return -1 for the overflow (log 0) case, - so the result is
-    * always at most 19 bits.
-    */
-   if ((x &= 0xff) == 0)
-      return -1;
-
-   if ((x & 0xf0) == 0)
-      lg2  = 4, x <<= 4;
-
-   if ((x & 0xc0) == 0)
-      lg2 += 2, x <<= 2;
-
-   if ((x & 0x80) == 0)
-      lg2 += 1, x <<= 1;
-
-   /* result is at most 19 bits, so this cast is safe: */
-   return (png_int_32)((lg2 << 16) + ((png_8bit_l2[x-128]+32768)>>16));
-}
-
-/* The above gives exact (to 16 binary places) log2 values for 8-bit images,
- * for 16-bit images we use the most significant 8 bits of the 16-bit value to
- * get an approximation then multiply the approximation by a correction factor
- * determined by the remaining up to 8 bits.  This requires an additional step
- * in the 16-bit case.
- *
- * We want log2(value/65535), we have log2(v'/255), where:
- *
- *    value = v' * 256 + v''
- *          = v' * f
- *
- * So f is value/v', which is equal to (256+v''/v') since v' is in the range 128
- * to 255 and v'' is in the range 0 to 255 f will be in the range 256 to less
- * than 258.  The final factor also needs to correct for the fact that our 8-bit
- * value is scaled by 255, whereas the 16-bit values must be scaled by 65535.
- *
- * This gives a final formula using a calculated value 'x' which is value/v' and
- * scaling by 65536 to match the above table:
- *
- *   log2(x/257) * 65536
- *
- * Since these numbers are so close to '1' we can use simple linear
- * interpolation between the two end values 256/257 (result -368.61) and 258/257
- * (result 367.179).  The values used below are scaled by a further 64 to give
- * 16-bit precision in the interpolation:
- *
- * Start (256): -23591
- * Zero  (257):      0
- * End   (258):  23499
- */
-#ifdef PNG_16BIT_SUPPORTED
-static png_int_32
-png_log16bit(png_uint_32 x)
-{
-   unsigned int lg2 = 0;
-
-   /* As above, but now the input has 16 bits. */
-   if ((x &= 0xffff) == 0)
-      return -1;
-
-   if ((x & 0xff00) == 0)
-      lg2  = 8, x <<= 8;
-
-   if ((x & 0xf000) == 0)
-      lg2 += 4, x <<= 4;
-
-   if ((x & 0xc000) == 0)
-      lg2 += 2, x <<= 2;
-
-   if ((x & 0x8000) == 0)
-      lg2 += 1, x <<= 1;
-
-   /* Calculate the base logarithm from the top 8 bits as a 28-bit fractional
-    * value.
-    */
-   lg2 <<= 28;
-   lg2 += (png_8bit_l2[(x>>8)-128]+8) >> 4;
-
-   /* Now we need to interpolate the factor, this requires a division by the top
-    * 8 bits.  Do this with maximum precision.
-    */
-   x = ((x << 16) + (x >> 9)) / (x >> 8);
-
-   /* Since we divided by the top 8 bits of 'x' there will be a '1' at 1<<24,
-    * the value at 1<<16 (ignoring this) will be 0 or 1; this gives us exactly
-    * 16 bits to interpolate to get the low bits of the result.  Round the
-    * answer.  Note that the end point values are scaled by 64 to retain overall
-    * precision and that 'lg2' is current scaled by an extra 12 bits, so adjust
-    * the overall scaling by 6-12.  Round at every step.
-    */
-   x -= 1U << 24;
-
-   if (x <= 65536U) /* <= '257' */
-      lg2 += ((23591U * (65536U-x)) + (1U << (16+6-12-1))) >> (16+6-12);
-
-   else
-      lg2 -= ((23499U * (x-65536U)) + (1U << (16+6-12-1))) >> (16+6-12);
-
-   /* Safe, because the result can't have more than 20 bits: */
-   return (png_int_32)((lg2 + 2048) >> 12);
-}
-#endif /* 16BIT */
-
-/* The 'exp()' case must invert the above, taking a 20-bit fixed point
- * logarithmic value and returning a 16 or 8-bit number as appropriate.  In
- * each case only the low 16 bits are relevant - the fraction - since the
- * integer bits (the top 4) simply determine a shift.
- *
- * The worst case is the 16-bit distinction between 65535 and 65534. This
- * requires perhaps spurious accuracy in the decoding of the logarithm to
- * distinguish log2(65535/65534.5) - 10^-5 or 17 bits.  There is little chance
- * of getting this accuracy in practice.
- *
- * To deal with this the following exp() function works out the exponent of the
- * frational part of the logarithm by using an accurate 32-bit value from the
- * top four fractional bits then multiplying in the remaining bits.
- */
-static const png_uint_32
-png_32bit_exp[16] =
-{
-   /* NOTE: the first entry is deliberately set to the maximum 32-bit value. */
-   4294967295U, 4112874773U, 3938502376U, 3771522796U, 3611622603U, 3458501653U,
-   3311872529U, 3171459999U, 3037000500U, 2908241642U, 2784941738U, 2666869345U,
-   2553802834U, 2445529972U, 2341847524U, 2242560872U
-};
-
-/* Adjustment table; provided to explain the numbers in the code below. */
-#if 0
-for (i=11;i>=0;--i){ print i, " ", (1 - e(-(2^i)/65536*l(2))) * 2^(32-i), "\n"}
-   11 44937.64284865548751208448
-   10 45180.98734845585101160448
-    9 45303.31936980687359311872
-    8 45364.65110595323018870784
-    7 45395.35850361789624614912
-    6 45410.72259715102037508096
-    5 45418.40724413220722311168
-    4 45422.25021786898173001728
-    3 45424.17186732298419044352
-    2 45425.13273269940811464704
-    1 45425.61317555035558641664
-    0 45425.85339951654943850496
-#endif
-
-static png_uint_32
-png_exp(png_fixed_point x)
-{
-   if (x > 0 && x <= 0xfffff) /* Else overflow or zero (underflow) */
-   {
-      /* Obtain a 4-bit approximation */
-      png_uint_32 e = png_32bit_exp[(x >> 12) & 0xf];
-
-      /* Incorporate the low 12 bits - these decrease the returned value by
-       * multiplying by a number less than 1 if the bit is set.  The multiplier
-       * is determined by the above table and the shift. Notice that the values
-       * converge on 45426 and this is used to allow linear interpolation of the
-       * low bits.
-       */
-      if (x & 0x800)
-         e -= (((e >> 16) * 44938U) +  16U) >> 5;
-
-      if (x & 0x400)
-         e -= (((e >> 16) * 45181U) +  32U) >> 6;
-
-      if (x & 0x200)
-         e -= (((e >> 16) * 45303U) +  64U) >> 7;
-
-      if (x & 0x100)
-         e -= (((e >> 16) * 45365U) + 128U) >> 8;
-
-      if (x & 0x080)
-         e -= (((e >> 16) * 45395U) + 256U) >> 9;
-
-      if (x & 0x040)
-         e -= (((e >> 16) * 45410U) + 512U) >> 10;
-
-      /* And handle the low 6 bits in a single block. */
-      e -= (((e >> 16) * 355U * (x & 0x3fU)) + 256U) >> 9;
-
-      /* Handle the upper bits of x. */
-      e >>= x >> 16;
-      return e;
-   }
-
-   /* Check for overflow */
-   if (x <= 0)
-      return png_32bit_exp[0];
-
-   /* Else underflow */
-   return 0;
-}
-
-static png_byte
-png_exp8bit(png_const_structrp png_ptr, png_fixed_point lg2)
-{
-   /* Get a 32-bit value: */
-   png_uint_32 x = png_exp(lg2);
-
-   /* Convert the 32-bit value to 0..255 by multiplying by 256-1. Note that the
-    * second, rounding, step can't overflow because of the first, subtraction,
-    * step.
-    */
-   x -= x >> 8;
-   return png_check_byte(png_ptr, (x + 0x7fffffU) >> 24);
-   PNG_UNUSEDRC(png_ptr)
-}
-
-static png_uint_16
-png_exp16bit(png_const_structrp png_ptr, png_fixed_point lg2)
-{
-   /* Get a 32-bit value: */
-   png_uint_32 x = png_exp(lg2);
-
-   /* Convert the 32-bit value to 0..65535 by multiplying by 65536-1: */
-   x -= x >> 16;
-   return png_check_u16(png_ptr, (x + 32767U) >> 16);
-   PNG_UNUSEDRC(png_ptr)
-}
-#endif /* FLOATING_ARITHMETIC */
-
-png_byte
-png_gamma_8bit_correct(png_const_structrp png_ptr, png_uint_32 value,
-   png_fixed_point gamma_val)
-{
-   if (value > 0 && value < 255)
-   {
-#     ifdef PNG_FLOATING_ARITHMETIC_SUPPORTED
-         /* 'value' is unsigned, ANSI-C90 requires the compiler to correctly
-          * convert this to a floating point value.  This includes values that
-          * would overflow if 'value' were to be converted to 'int'.
-          *
-          * Apparently GCC, however, does an intermediate conversion to (int)
-          * on some (ARM) but not all (x86) platforms, possibly because of
-          * hardware FP limitations.  (E.g. if the hardware conversion always
-          * assumes the integer register contains a signed value.)  This results
-          * in ANSI-C undefined behavior for large values.
-          *
-          * Other implementations on the same machine might actually be ANSI-C90
-          * conformant and therefore compile spurious extra code for the large
-          * values.
-          *
-          * We can be reasonably sure that an unsigned to float conversion
-          * won't be faster than an int to float one.  Therefore this code
-          * assumes responsibility for the undefined behavior, which it knows
-          * can't happen because of the check above.
-          *
-          * Note the argument to this routine is an (unsigned int) because, on
-          * 16-bit platforms, it is assigned a value which might be out of
-          * range for an (int); that would result in undefined behavior in the
-          * caller if the *argument* ('value') were to be declared (int).
-          */
-         double r = floor(255*pow((int)/*SAFE*/value/255.,gamma_val*.00001)+.5);
-         if (r >= 0 && r <= 255)
-             return (png_byte)/*SAFE*/r;
-#     else
-         png_int_32 lg2 = png_log8bit(value);
-         png_fixed_point res;
-
-         if (png_muldiv(&res, gamma_val, lg2, PNG_FP_1) != 0)
-            return png_exp8bit(png_ptr, res);
-#     endif
-
-      /* Overflow. */
-      handled("8-bit gamma overflow");
-      return 0;
-   }
-
-   return png_check_byte(png_ptr, value);
-   PNG_UNUSED(png_ptr) /* Only used in non-release builds */
-}
-
-png_uint_16
-png_gamma_16bit_correct(png_const_structrp png_ptr, png_uint_32 value,
-   png_fixed_point gamma_val)
-{
-   if (value > 0 && value < 65535)
-   {
-#     ifdef PNG_FLOATING_ARITHMETIC_SUPPORTED
-         /* The same (unsigned int)->(double) constraints apply here as above,
-          * however in this case the (unsigned int) to (int) conversion can
-          * overflow on an ANSI-C90 compliant system so the cast needs to ensure
-          * that this is not possible.
-          */
-         double r = floor(65535.*pow((png_int_32)value/65535.,
-                     gamma_val*.00001)+.5);
-         if (r >= 0 && r <= 65535)
-             return (png_uint_16)/*SAFE*/r;
-#     else
-         png_int_32 lg2 = png_log16bit(value);
-         png_fixed_point res;
-
-         if (png_muldiv(&res, gamma_val, lg2, PNG_FP_1) != 0)
-            return png_exp16bit(png_ptr, res);
-#     endif
-
-      /* Overflow. */
-      handled("16-bit gamma overflow");
-      return 0;
-   }
-
-   return png_check_u16(png_ptr, value);
-   PNG_UNUSED(png_ptr) /* Only used in non-release builds */
-}
-
-#define PNG_GAMMA_TABLE_8       0 /* 8-bit entries in png_byte */
-#define PNG_GAMMA_TABLE_8_IN_16 1 /* 8-bit entries * 257 in png_uint_16 */
-#define PNG_GAMMA_TABLE_16      2 /* 16-bit entries in png_uint_16 */
-
-typedef struct
-{
-   png_fixed_point gamma;
-   png_uint_32     mult;
-   unsigned int    add;
-   unsigned int    shift;  /* input value is (i * mult + add) >> shift */
-   int             output; /* One of the above values */
-   int             adjust; /* Divide or multiple output by 257 */
-   png_voidp       table;  /* Lookup table */
-}  gamma_table_data;
-
-static unsigned int
-write_gamma_table_entry(png_const_structrp png_ptr,
-   const gamma_table_data *data, png_uint_32 i)
-   /* Calculate and write a single entry into table[i], the value of the entry
-    * written is returned.
-    */
-{
-   png_uint_32 in = (i * data->mult + data->add) >> data->shift;
-   unsigned int out;
-   
-   /* If the output is TABLE_8 with no adjust, or the output is not with an
-    * adjust, use 8-bit correction.
-    */
-   if ((data->output == PNG_GAMMA_TABLE_8) != (data->adjust != 0))
-   {
-      out = png_gamma_8bit_correct(png_ptr, in, data->gamma);
-
-      if (data->adjust != 0)
-         out *= 257U;
-   }
-
-   else /* 16-bit correction */
-   {
-      out = png_gamma_16bit_correct(png_ptr, in, data->gamma);
-
-      if (data->adjust != 0)
-         out = PNG_DIV257(out);
-   }
-
-   PNG_UNUSEDRC(png_ptr)
-   if (data->output == PNG_GAMMA_TABLE_8)
-      png_upcast(png_bytep, data->table)[i] = png_check_byte(png_ptr, out);
-
-   else
-      png_upcast(png_uint_16p, data->table)[i] = png_check_u16(png_ptr, out);
-
-   return out;
-}
-
-static void
-write_gamma_table(png_const_structrp png_ptr, const gamma_table_data *data,
-   png_uint_32 lo, unsigned int loval, png_uint_32 hi, unsigned int hival)
-   /* Fill in gamma table entries between lo and hi, exclusive.  The entries at
-    * table[lo] and table[hi] have already been written, the intervening entries
-    * are written.
-    */
-{
-   if (hi > lo+1) /* Else nothing to fill in */
-   {
-      if (hival == loval)
-      {
-         /* All intervening entries must be the same. */
-         if (data->output == PNG_GAMMA_TABLE_8)
-         {
-            png_bytep table8 = png_voidcast(png_bytep, data->table);
-
-            while (++lo < hi)
-               table8[lo] = png_check_byte(png_ptr, loval);
-         }
-
-         else
-         {
-            png_uint_16p table16 = png_voidcast(png_uint_16p, data->table);
-
-            while (++lo < hi)
-               table16[lo] = png_check_u16(png_ptr, loval);
-         }
-      }
-
-      else
-      {
-         png_uint_32 mid = (lo+hi) >> 1;
-         unsigned int midval = write_gamma_table_entry(png_ptr, data, mid);
-
-         /* The algorithm used is to divide the entries to be written in half
-          * and fill in the middle.  For all practical tables with significant
-          * gamma this will result in a performance gain because the expensive
-          * gamma correction arithmetic is avoided for some entries.
-          */
-         write_gamma_table(png_ptr, data, lo, loval, mid, midval);
-         write_gamma_table(png_ptr, data, mid, midval, hi, hival);
-      }
-   }
-}
-
-static void *
-png_build_gamma_table(png_structrp png_ptr, png_fixed_point gamma_val,
-   int output/*as above*/, int input_depth, int use_shift)
-   /* Build a gamma lookup table to encode input_depth bit input values.
-    * The table will have 2^input_depth entries plus an extra one if use_shift
-    * is specified.  With shift the table is accessed:
-    *
-    *    table[(original-value + rounding) >> shift]
-    *
-    * And an extra entry exists to accomodate overflow of original-value on
-    * rounding.  If use_shift is not specified the table is accessed with an
-    * input_depth bit value and the original values must have been correctly
-    * scaled to this range (not using a shift!)
-    *
-    * Each table entry contains input-value^gamma_val rounded to the output
-    * precision.  This is 8-bit precision unless output is specified as
-    * PNG_GAMMA_TABLE_16, in which case it is 16-bit precision.  For
-    * PNG_GAMMA_TABLE_8_IN_16 the 8-bit value is scaled to 16-bits by
-    * multiplying by 257.
-    */
-{
-   png_uint_32 size;
-   unsigned int hival;
-   gamma_table_data data;
-
-   /* If use_shift is true or if the input or output is not 8-bit the gamma
-    * correction will use the 16-bit correction code.  This requires a value in
-    * the range 0..65535.  For use_shift the value is simply:
-    *
-    *    input << shift
-    *
-    * For the scaling case the value is:
-    *
-    *    round(input * 65535 / ((1<<input_depth)-1)
-    *
-    * Both these expressions can be rewritten as:
-    *
-    *    (input * mult + add) >> shift;
-    *
-    * With 'mult' and 'add' chosen to minimize the error for all input values
-    * in the range 0..((1<<input_depth)-1).  The following table does this for
-    * the scaling case.  In fact all the answers are except except for the
-    * 13-bit case, where the maximum error (from the exact value) is 0.500183.
-    *
-    * This table can be produced using the code in contrib/tools/scale.c
-    */
-   static const struct
-   {
-      png_uint_32 mult;
-      png_uint_16 add;
-      png_byte    shift;
-   }  multadd65535[16] =
-   {
-      {      65535,     0,  0 }, /* 65535/1 */
-      {      21845,     0,  0 }, /* 65535/3 */
-      {      37449,     0,  2 }, /* 65535/7 */
-      {       4369,     0,  0 }, /* 65535/15 */
-      {      33825,     0,  4 }, /* 65535/31 */
-      {     266301,   121,  8 }, /* 65535/63 */
-      {    1056817,   970, 11 }, /* 65535/127 */
-      {        257,     0,  0 }, /* 65535/255 */
-      {     262653,  1020, 11 }, /* 65535/511 */
-      {    1049585,  8165, 14 }, /* 65535/1023 */
-      {    2098145, 31774, 16 }, /* 65535/2047 */
-      {      65551,  2055, 12 }, /* 65535/4095 */
-      {      65543,  4100, 13 }, /* 65535/8191 ERROR: .5+0.000183128*/
-      {      65539,  8193, 14 }, /* 65535/16383 */
-      {      32769,     0, 14 }, /* 65535/32767 */
-      {          1,     0,  0 }  /* 65535/65535 */
-#  if 0 /* inverse */
-      {      1,         0, 15 }, /* 1/65535 */
-      {      3,     32769, 16 }, /* 3/65535 */
-      {  28673, 134188470, 28 }, /* 7/65535 */
-      {     15,     32775, 16 }, /* 15/65535 */
-      {  31745,  33522654, 26 }, /* 31/65535 */
-      {  64513,  33552693, 26 }, /* 63/65535 */
-      {  65025,  16776620, 25 }, /* 127/65535 */
-      {    255,     32895, 16 }, /* 255/65535 */
-      {  65409,   4194134, 23 }, /* 511/65535 */
-      {  65473,   2097037, 22 }, /* 1023/65535 */
-      {  65505,   1048544, 21 }, /* 2047/65535 */
-      {  65521,    524167, 20 }, /* 4095/65535 */
-      {  65529,    262136, 19 }, /* 8191/65535 */
-      {  65533,    131065, 18 }, /* 16383/65535 */
-      {      1,         0,  1 }, /* 32767/65535 */
-      {      1,         0,  0 }  /* 65535/65535 */
-#  endif
-   };
-
-   /* When both the input and output are 8-bit (i.e. the output is not
-    * PNG_GAMMA_TABLE_16 and the input_depth is <9) the 8-bit gamma correction
-    * code can be used, it is slightly faster.  This requires values scaled to
-    * 255, not 65535:
-    */
-   static const struct
-   {
-      png_uint_16 mult;
-      png_byte    add;
-      png_byte    shift;
-   }  multadd255[8] =
-   {
-      { 255,    0,  0 }, /* 255/1 */
-      {  85,    0,  0 }, /* 255/3 */
-      {  73,    0,  1 }, /* 255/7 */
-      {  17,    0,  0 }, /* 255/15 */
-      { 527,   23,  6 }, /* 255/31 */
-      { 259,   33,  6 }, /* 255/63 */
-      { 129,    0,  6 }, /* 255/127 */
-      {   1,    0,  0 }  /* 255/255 */
-#  if 0 /* inverse */
-      {   1,    0,  7 }, /* 1/255 */
-      {   3,  129,  8 }, /* 3/255 */
-      { 225, 4060, 13 }, /* 7/255 */
-      {  15,  135,  8 }, /* 15/255 */
-      { 249, 1014, 11 }, /* 31/255 */
-      { 253,  505, 10 }, /* 63/255 */
-      {   1,    0,  1 }, /* 127/255 */
-      {   1,    0,  0 }  /* 255/255 */
-#  endif
-   };
-
-   /* Basic table size, increased by one below in the use_shift case where the
-    * input is rounded.
-    */
-   size = 1U << input_depth;
-   data.gamma = gamma_val;
-   data.output = output;
-
-   if (output < PNG_GAMMA_TABLE_16 && input_depth <= 8)
-   {
-      /* The 8-bit correction can only be used if both input and output have no
-       * more than 8 bits of precision.
-       */
-      data.adjust = output > PNG_GAMMA_TABLE_8;
-
-      if (use_shift != 0)
-      {
-         /* The multiplier does the shift: */
-         data.mult = 1U << (8-input_depth);
-         data.add = 0;
-         data.shift = 0;
-         if (input_depth < 8) ++size;
-      }
-
-      else
-      {
-         data.mult = multadd255[input_depth-1].mult;
-         data.add = multadd255[input_depth-1].add;
-         data.shift = multadd255[input_depth-1].shift;
-      }
-   }
-
-   else
-   {
-      /* 16-bit correction is used for cases where input or output require more
-       * than 8 bits.
-       */
-      data.adjust = output == PNG_GAMMA_TABLE_8;
-
-      if (use_shift != 0)
-      {
-         data.mult = 1U << (16-input_depth);
-         data.add = 0;
-         data.shift = 0;
-         if (input_depth < 16) ++size;
-      }
-
-      else
-      {
-         data.mult = multadd65535[input_depth-1].mult;
-         data.add = multadd65535[input_depth-1].add;
-         data.shift = multadd65535[input_depth-1].shift;
-      }
-   }
-
-   if (output == PNG_GAMMA_TABLE_8)
-   {
-      data.table = png_malloc(png_ptr, size * sizeof (png_byte));
-      ((png_bytep)data.table)[0] = 0;
-      hival = ((png_bytep)data.table)[size-1] = 255;
-   }
-
-   else
-   {
-      /* Output is 16 bits, although it may only have 8 bits of precision */
-      data.table = png_malloc(png_ptr, size * sizeof (png_uint_16));
-      ((png_uint_16p)data.table)[0] = 0;
-      hival = ((png_uint_16p)data.table)[size-1] = 65535;
-   }
-
-   if (png_gamma_significant(gamma_val) != 0)
-      write_gamma_table(png_ptr, &data, 0, 0, size-1, hival);
-
-   else /* gamma_val not significant */
-   {
-      if (output == PNG_GAMMA_TABLE_8)
-      {
-         png_uint_32 i;
-         png_bytep table8 = ((png_bytep)data.table);
-
-         if (data.adjust)
-            for (i=1; i<size-1; ++i)
-               table8[i] = png_check_byte(png_ptr,
-                  PNG_DIV257((i * data.mult + data.add) >> data.shift));
-
-         else
-            for (i=1; i<size-1; ++i)
-               table8[i] = png_check_byte(png_ptr,
-                  (i * data.mult + data.add) >> data.shift);
-      }
-
-      else
-      {
-         png_uint_32 i;
-         png_uint_16p table16 = (png_uint_16p)data.table;
-
-         if (data.adjust)
-            for (i=1; i<size-1; ++i)
-               table16[i] = png_check_u16(png_ptr,
-                  ((i * data.mult + data.add) >> data.shift) * 257U);
-
-         else
-            for (i=1; i<size-1; ++i)
-               table16[i] = png_check_u16(png_ptr,
-                  (i * data.mult + data.add) >> data.shift);
-      }
-   }
-
-   return data.table;
-}
-
-/* Used from png_read_destroy and below to release the memory used by the gamma
- * tables.
- */
-void /* PRIVATE */
-png_destroy_gamma_table(png_structrp png_ptr)
-{
-   png_free(png_ptr, png_ptr->gamma_table);
-   png_ptr->gamma_table = NULL;
-
-   png_free(png_ptr, png_ptr->gamma_16_table);
-   png_ptr->gamma_16_table = NULL;
-
-#if defined(PNG_READ_BACKGROUND_SUPPORTED) || \
-   defined(PNG_READ_ALPHA_MODE_SUPPORTED) || \
-   defined(PNG_READ_RGB_TO_GRAY_SUPPORTED)
-   png_free(png_ptr, png_ptr->gamma_from_1);
-   png_ptr->gamma_from_1 = NULL;
-   png_free(png_ptr, png_ptr->gamma_to_1);
-   png_ptr->gamma_to_1 = NULL;
-
-   png_free(png_ptr, png_ptr->gamma_16_from_1);
-   png_ptr->gamma_16_from_1 = NULL;
-   png_free(png_ptr, png_ptr->gamma_16_to_1);
-   png_ptr->gamma_16_to_1 = NULL;
-#endif /* READ_BACKGROUND || READ_ALPHA_MODE || RGB_TO_GRAY */
-}
-
-/* We build the 8- or 16-bit gamma tables here.  Note that for 16-bit
- * tables, we don't make a full table if we are reducing to 8-bit in
- * the future.  Note also how the gamma_16 tables are segmented so that
- * we don't need to allocate > 64K chunks for a full 16-bit table.
- */
-void /* PRIVATE */
-png_build_gamma_tables(png_structrp png_ptr, int bit_depth)
-{
-  png_debug(1, "in png_build_gamma_table");
-
-  /* Remove any existing table; this copes with multiple calls to
-   * png_read_update_info.  The warning is because building the gamma tables
-   * multiple times is a performance hit - it's harmless but the ability to call
-   * png_read_update_info() multiple times is new in 1.5.6 so it seems sensible
-   * to warn if the app introduces such a hit.
-   */
-  if (png_ptr->gamma_table != NULL || png_ptr->gamma_16_table != NULL)
-  {
-    png_warning(png_ptr, "gamma table being rebuilt");
-    png_destroy_gamma_table(png_ptr);
-  }
-
-  if (bit_depth <= 8)
-  {
-     png_ptr->gamma_table = png_voidcast(png_bytep, png_build_gamma_table(
-         png_ptr, png_ptr->screen_gamma > 0 ?
-         png_reciprocal2(png_ptr->colorspace.gamma, png_ptr->screen_gamma) :
-         PNG_FP_1, PNG_GAMMA_TABLE_8, 8/*input depth*/, 0/*scale*/));
-
-#if defined(PNG_READ_BACKGROUND_SUPPORTED) || \
-   defined(PNG_READ_ALPHA_MODE_SUPPORTED) || \
-   defined(PNG_READ_RGB_TO_GRAY_SUPPORTED)
-     if ((png_ptr->transformations & (PNG_COMPOSE | PNG_RGB_TO_GRAY)) != 0)
-     {
-        /* This sets the accuracy of 8-bit composition and the 8-bit RGB to gray
-         * conversion - PNG_MAX_GAMMA_8 (the number of bits in the sixteen bit
-         * value that are considered significant.)
-         */
-        png_ptr->gamma_to_1 = png_voidcast(png_uint_16p, png_build_gamma_table(
-            png_ptr, png_reciprocal(png_ptr->colorspace.gamma),
-            PNG_GAMMA_TABLE_16, 8/*input depth*/, 0/*scale*/));
-
-        png_ptr->gamma_from_1 = png_voidcast(png_bytep, png_build_gamma_table(
-            png_ptr, png_ptr->screen_gamma > 0 ?
-            png_reciprocal(png_ptr->screen_gamma) :
-            png_ptr->colorspace.gamma/* Probably doing rgb_to_gray */,
-            PNG_GAMMA_TABLE_8, PNG_MAX_GAMMA_8/*input depth*/, 1/*shift*/));
-
-        png_ptr->gamma_shift = 16-PNG_MAX_GAMMA_8;
-     }
-#endif /* READ_BACKGROUND || READ_ALPHA_MODE || RGB_TO_GRAY */
-  }
-  else
-  {
-     png_byte shift, sig_bit;
-     int table_type;
-
-#    ifdef PNG_16BIT_SUPPORTED
-         table_type = PNG_GAMMA_TABLE_16;
-#    else
-         table_type = PNG_GAMMA_TABLE_8_IN_16;
-#    endif
-
-     if ((png_ptr->color_type & PNG_COLOR_MASK_COLOR) != 0)
-     {
-        sig_bit = png_ptr->sig_bit.red;
-
-        if (png_ptr->sig_bit.green > sig_bit)
-           sig_bit = png_ptr->sig_bit.green;
-
-        if (png_ptr->sig_bit.blue > sig_bit)
-           sig_bit = png_ptr->sig_bit.blue;
-     }
-     else
-        sig_bit = png_ptr->sig_bit.gray;
-
-     /* shift == insignificant bits */
-     if (sig_bit > 0 && sig_bit < 16U)
-        shift = png_check_byte(png_ptr, 16U - sig_bit);
-
-     else
-        shift = 0; /* keep all 16 bits */
-
-     if ((png_ptr->transformations & (PNG_16_TO_8 | PNG_SCALE_16_TO_8)) != 0)
-     {
-        /* PNG_MAX_GAMMA_8 is the number of bits to keep - effectively
-         * the significant bits in the *input* when the output will
-         * eventually be 8 bits.
-         */
-        if (shift < (16U - PNG_MAX_GAMMA_8))
-           shift = (16U - PNG_MAX_GAMMA_8);
-
-        table_type = PNG_GAMMA_TABLE_8_IN_16;
-     }
-
-     png_ptr->gamma_shift = shift;
-
-     png_ptr->gamma_16_table = png_voidcast(png_uint_16p, png_build_gamma_table(
-         png_ptr, png_ptr->screen_gamma > 0 ? png_reciprocal2(
-         png_ptr->colorspace.gamma, png_ptr->screen_gamma) : PNG_FP_1,
-         table_type, (16-shift)/*input depth*/, 1/*shift*/));
-
-#if defined(PNG_READ_BACKGROUND_SUPPORTED) || \
-   defined(PNG_READ_ALPHA_MODE_SUPPORTED) || \
-   defined(PNG_READ_RGB_TO_GRAY_SUPPORTED)
-     if ((png_ptr->transformations & (PNG_COMPOSE | PNG_RGB_TO_GRAY)) != 0)
-     {
-        png_ptr->gamma_16_to_1 = png_voidcast(png_uint_16p,
-            png_build_gamma_table(png_ptr,
-            png_reciprocal(png_ptr->colorspace.gamma), PNG_GAMMA_TABLE_16,
-            (16-shift)/*input depth*/, 1/*shift*/));
-
-        /* Notice that the '16 from 1' table should be full precision, however
-         * the lookup on this table still uses gamma_shift, so it can't be.
-         * TODO: fix this.
-         */
-        png_ptr->gamma_16_from_1 = png_voidcast(png_uint_16p,
-            png_build_gamma_table(png_ptr, png_ptr->screen_gamma > 0 ?
-            png_reciprocal(png_ptr->screen_gamma) :
-            png_ptr->colorspace.gamma/* Probably doing rgb_to_gray */,
-            PNG_GAMMA_TABLE_16, (16-shift)/*input depth*/, 1/*shift*/));
-     }
-#endif /* READ_BACKGROUND || READ_ALPHA_MODE || RGB_TO_GRAY */
-  }
-}
-#endif /* READ_GAMMA */
-
 /* HARDWARE OPTION SUPPORT */
 #ifdef PNG_SET_OPTION_SUPPORTED
 int PNGAPI