home *** CD-ROM | disk | FTP | other *** search
/ Amiga ISO Collection / AmigaUtilCD1.iso / GFX / viewer / flick-1.4.lha / src / median.c < prev    next >
Encoding:
C/C++ Source or Header  |  1995-01-03  |  7.3 KB  |  234 lines
  1. /*****************************************************************************
  2.  
  3.         Flick FLI-format Animation Viewer v1.2          19 Feb 1994
  4.         --------------------------------------
  5.  
  6.  
  7. This program plays FLI/FLC-format bitmapped animation files on any ECS
  8. or AGA Amiga running OS2.04 or higher.  FLI/FLC-format files are
  9. produced by Autodesk Animator and Autodesk 3D Studio on a PC, as well
  10. as by other programs.
  11.  
  12. The files in this archive may be distributed anywhere provided they are
  13. unmodified and are not sold for profit.
  14.  
  15. Ownership and copyright of all files remains with the author:
  16.  
  17.     Peter McGavin, 86 Totara Crescent, Lower Hutt, New Zealand.
  18.     e-mail: peterm@maths.grace.cri.nz
  19.  
  20. *****************************************************************************/
  21.  
  22. #include "includes.h"
  23.  
  24. #define NBITS 4            /* for each r, g, b */
  25. #define NSIDE (1<<NBITS)
  26. #define NPALETTE 256      /* # of palette entries */
  27.  
  28. typedef enum rgb_type {R, G, B};
  29.  
  30. struct box_type {                  /* Defines a parallelopiped */
  31.   UBYTE start[3];
  32.   UBYTE end[3];
  33. };
  34.  
  35. struct node_type {
  36.   struct box_type box;        /* corners of the rectangular box */
  37.   UWORD n;            /* total number of pixels in box */
  38.   enum rgb_type long_primary;    /* longest direction (red, green or blue) */
  39.   ULONG badness;        /* spread of pixels within box */
  40.   UWORD primaryhist[3][NSIDE];    /* 3 histograms (one for each of R,G,B) */
  41.   UBYTE palette[3];        /* centre of gravity */
  42. };
  43.  
  44.  
  45. static UWORD hist[NSIDE][NSIDE][NSIDE];    /* 3D histogram (cube) */
  46.  
  47. static void update_entry (struct node_type *node)
  48. /* given a rectangular box of given dimensions (in node), calculate and
  49.    fill in all the other parts of the node. */
  50. {
  51.   UWORD r, g, b, i, mean;
  52.   ULONG c, wsum, bad;
  53.   LONG signed_diff;
  54.   enum rgb_type primary;
  55.  
  56.   /* build histogram for each primary */
  57.   memset (node->primaryhist, 0, 3 * NSIDE * sizeof(UWORD));
  58.   node->n = 0;
  59.   for (r = node->box.start[R]; r <= node->box.end[R]; r++)
  60.     for (g = node->box.start[G]; g <= node->box.end[G]; g++)
  61.       for (b = node->box.start[B]; b <= node->box.end[B]; b++) {
  62.         c = hist[r][g][b];
  63.         node->primaryhist[R][r] += c;
  64.         node->primaryhist[G][g] += c;
  65.         node->primaryhist[B][b] += c;
  66.         node->n += c;
  67.       }
  68.  
  69.   if (node->n == 0)
  70.     die ("Empty box!  This should never happen.");
  71.  
  72.   /* shrink the box */
  73.   for (primary = R; primary <= B; primary++) {
  74.     i = node->box.start[primary];
  75.     while (node->primaryhist[primary][i] == 0)
  76.       i++;
  77.     node->box.start[primary] = i;
  78.     i = node->box.end[primary];
  79.     while (node->primaryhist[primary][i] == 0)
  80.       i--;
  81.     node->box.end[primary] = i;
  82.   }
  83.  
  84.   /* compute the badness & choose the primary with the greatest badness */
  85.   node->badness = 0;
  86.   for (primary = R; primary <= B; primary++) {
  87.     wsum = 0;
  88.     for (i = node->box.start[primary]; i <= node->box.end[primary]; i++)
  89.       wsum += (node->primaryhist[primary][i] * (ULONG)i);
  90.     mean = ((wsum << 1) + node->n) / (node->n << 1); /* round(wsum/node->n) */
  91.     node->palette[primary] = (UBYTE)(mean << 4);
  92.     bad = 0;
  93.     for (i = node->box.start[primary]; i <= node->box.end[primary]; i++) {
  94.       signed_diff = ((WORD)mean) - ((WORD)(i));
  95.       bad += node->primaryhist[primary][i] * (ULONG)(signed_diff * signed_diff);
  96.     }
  97.     if (bad >= node->badness) {
  98.       node->badness = bad;
  99.       node->long_primary = primary;
  100.     }
  101.   }
  102. }
  103.  
  104.  
  105. static void cut (struct node_type *node0, struct node_type *node1)
  106. /* cut node0 into 2 pieces in the node0->long_primary direction and store
  107.    the 2 pieces back into node0 and node1 */
  108. {
  109.   ULONG sum, goal;
  110.   UWORD split;
  111.  
  112.   if (node0->badness == 0)
  113.     die ("Badness == 0!  This should never happen");
  114.  
  115.   /* find the median */
  116.   goal = node0->n >> 1;    /* how many we want on each side of the cut */
  117.   sum = 0;
  118.   split = node0->box.start[node0->long_primary];
  119.   while (sum <= goal)
  120.     sum += node0->primaryhist[node0->long_primary][split++];
  121.  
  122.   /* go back a bit if necessary to get a better balance */
  123.   if ((sum - goal) >
  124.       (goal - (sum - node0->primaryhist[node0->long_primary][split - 1])))
  125.     sum -= node0->primaryhist[node0->long_primary][--split];
  126.  
  127.   /* copy box */
  128.   node1->box = node0->box;
  129.  
  130.   /* set new dimensions of boxes */
  131.   node0->box.end[node0->long_primary] = split - 1;
  132.   node1->box.start[node0->long_primary] = split;
  133.  
  134.   /* update the nodes */
  135.   update_entry (node0);
  136.   update_entry (node1);
  137. }
  138.  
  139.  
  140. static struct node_type node[32];          /* node list */
  141.  
  142. void median_cut (UBYTE colourtable[NPALETTE][3],
  143.                  UWORD *viewcolourtable,
  144.                  UBYTE *xlate,
  145.                  enum mode_type mode)
  146. {
  147.   UWORD idx, this_idx, worst_idx, c, r, g, b, best_idx, max_nodes;
  148.   WORD dr, dg, db;
  149.   ULONG max_badness, max_n, dist, best_dist;
  150.   enum rgb_type primary;
  151.  
  152.   switch (mode) {
  153.     case MODE_COLOUR4:
  154.       max_nodes = 16;
  155.       break;
  156.     case MODE_EHB:
  157.       max_nodes = 32;
  158.       break;
  159.     default:
  160.       die ("Wrong mode in median_cut!");
  161.   }
  162.  
  163.   /* build histogram from colourtable (unlike conventional median split
  164.      algorithm which uses pixels from image, but that would be too slow) */
  165.   memset (hist, 0, NSIDE * NSIDE * NSIDE * sizeof(UWORD));
  166.   for (c = 0; c < NPALETTE; c++) {
  167.     r = colourtable[c][R] >> 4;
  168.     g = colourtable[c][G] >> 4;
  169.     b = colourtable[c][B] >> 4;
  170.     if (mode == MODE_EHB && r < 8 && g < 8 && b < 8)    /* EHB if possible */
  171.       hist[r << 1][g << 1][b << 1]++;
  172.     else
  173.       hist[r][g][b]++;
  174.   }
  175.  
  176.   /* set up an initial box containing the entire colour cube */
  177.   for (primary = R; primary <= B; primary++) {
  178.     node[0].box.start[primary] = 0;
  179.     node[0].box.end[primary] = NSIDE - 1;
  180.   }
  181.   update_entry (&node[0]);
  182.  
  183.   /* locate and subdivide the worst box until there are not enough
  184.      palette entries for more boxes or all the boxes are too small to
  185.      subdivide further */
  186.   for (this_idx = 1; this_idx < max_nodes; this_idx++) {
  187.     max_badness = 0;
  188.     max_n = 0;
  189.     for (idx = 0; idx < this_idx; idx++)
  190.       if (node[idx].badness > max_badness ||
  191.           (node[idx].badness == max_badness && node[idx].n > max_n)) {
  192.         max_badness = node[idx].badness;
  193.         max_n = node[idx].n;
  194.         worst_idx = idx;
  195.       }
  196.     if (max_badness == 0)
  197.       break;
  198.     cut (&node[worst_idx], &node[this_idx]);
  199.   }
  200.  
  201.   /* fill the viewcolourtable */
  202.   for (idx = 0; idx < this_idx; idx++)
  203.     viewcolourtable[idx] = ((node[idx].palette[R] >> 4) << 8) |
  204.                            ((node[idx].palette[G] >> 4) << 4) |
  205.                             (node[idx].palette[B] >> 4);
  206.  
  207.   /* calculate pixel translation table */
  208.   for (c = 0; c < NPALETTE; c++) {
  209.     r = colourtable[c][R];
  210.     g = colourtable[c][G];
  211.     b = colourtable[c][B];
  212.     best_dist = 32760;
  213.     for (idx = 0; idx < this_idx; idx++) {
  214.       dr = ((WORD)r) - (WORD)node[idx].palette[R];
  215.       dg = ((WORD)g) - (WORD)node[idx].palette[G];
  216.       db = ((WORD)b) - (WORD)node[idx].palette[B];
  217.       if ((dist = dr * dr + dg * dg + db * db) < best_dist) {
  218.         best_dist = dist;
  219.         best_idx = idx;
  220.       }
  221.       if (mode == MODE_EHB) {
  222.         dr = ((WORD)r) - (WORD)(node[idx].palette[R] >> 1);
  223.         dg = ((WORD)g) - (WORD)(node[idx].palette[G] >> 1);
  224.         db = ((WORD)b) - (WORD)(node[idx].palette[B] >> 1);
  225.         if ((dist = dr * dr + dg * dg + db * db) < best_dist) {
  226.           best_dist = dist;
  227.           best_idx = idx + 32;    /* EHB colour entry */
  228.         }
  229.       }
  230.     }
  231.     xlate[c] = best_idx;
  232.   }
  233. }
  234.