Home | History | Annotate | Download | only in tools 
      1 #!/bin/sh
      2 #
      3 # intgamma.sh
      4 #
      5 # Last changed in libpng 1.6.0 [February 14, 2013]
      6 #
      7 # COPYRIGHT: Written by John Cunningham Bowler, 2013.
      8 # To the extent possible under law, the author has waived all copyright and
      9 # related or neighboring rights to this work.  This work is published from:
     10 # United States.
     11 #
     12 # Shell script to generate png.c 8-bit and 16-bit log tables (see the code in
     13 # png.c for details).
     14 #
     15 # This script uses the "bc" arbitrary precision calculator to calculate 32-bit
     16 # fixed point values of logarithms appropriate to finding the log of an 8-bit
     17 # (0..255) value and a similar table for the exponent calculation.
     18 #
     19 # "bc" must be on the path when the script is executed, and the math library
     20 # (-lm) must be available
     21 #
     22 # function to print out a list of numbers as integers; the function truncates
     23 # the integers which must be one-per-line
     24 function print(){
     25    awk 'BEGIN{
     26       str = ""
     27    }
     28    {
     29       sub("\\.[0-9]*$", "")
     30       if ($0 == "")
     31          $0 = "0"
     32 
     33       if (str == "")
     34          t = "   " $0 "U"
     35       else
     36          t = str ", " $0 "U"
     37 
     38       if (length(t) >= 80) {
     39          print str ","
     40          str = "   " $0 "U"
     41       } else
     42          str = t
     43    }
     44    END{
     45       print str
     46    }'
     47 }
     48 #
     49 # The logarithm table.
     50 cat <<END
     51 /* 8-bit log table: png_8bit_l2[128]
     52  * This is a table of -log(value/255)/log(2) for 'value' in the range 128 to
     53  * 255, so it's the base 2 logarithm of a normalized 8-bit floating point
     54  * mantissa.  The numbers are 32-bit fractions.
     55  */
     56 static const png_uint_32
     57 png_8bit_l2[128] =
     58 {
     59 END
     60 #
     61 bc -lqws <<END | print
     62 f=65536*65536/l(2)
     63 for (i=128;i<256;++i) { .5 - l(i/255)*f; }
     64 END
     65 echo '};'
     66 echo
     67 #
     68 # The exponent table.
     69 cat <<END
     70 /* The 'exp()' case must invert the above, taking a 20-bit fixed point
     71  * logarithmic value and returning a 16 or 8-bit number as appropriate.  In
     72  * each case only the low 16 bits are relevant - the fraction - since the
     73  * integer bits (the top 4) simply determine a shift.
     74  *
     75  * The worst case is the 16-bit distinction between 65535 and 65534; this
     76  * requires perhaps spurious accuracy in the decoding of the logarithm to
     77  * distinguish log2(65535/65534.5) - 10^-5 or 17 bits.  There is little chance
     78  * of getting this accuracy in practice.
     79  *
     80  * To deal with this the following exp() function works out the exponent of the
     81  * frational part of the logarithm by using an accurate 32-bit value from the
     82  * top four fractional bits then multiplying in the remaining bits.
     83  */
     84 static const png_uint_32
     85 png_32bit_exp[16] =
     86 {
     87 END
     88 #
     89 bc -lqws <<END | print
     90 f=l(2)/16
     91 for (i=0;i<16;++i) {
     92    x = .5 + e(-i*f)*2^32;
     93    if (x >= 2^32) x = 2^32-1;
     94    x;
     95 }
     96 END
     97 echo '};'
     98 echo
     99 #
    100 # And the table of adjustment values.
    101 cat <<END
    102 /* Adjustment table; provided to explain the numbers in the code below. */
    103 #if 0
    104 END
    105 bc -lqws <<END | awk '{ printf "%5d %s\n", 12-NR, $0 }'
    106 for (i=11;i>=0;--i){
    107    (1 - e(-(2^i)/65536*l(2))) * 2^(32-i)
    108 }
    109 END
    110 echo '#endif'
    111