makesRGB.c 11 KB

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