home *** CD-ROM | disk | FTP | other *** search
/ High Voltage Shareware / high1.zip / high1 / DIR2 / DVPG30FS.ZIP / JQUANT2.C < prev    next >
C/C++ Source or Header  |  1993-11-28  |  46KB  |  1,213 lines
  1. /*
  2.  * jquant2.c
  3.  *
  4.  * Copyright (C) 1991, 1992, 1993, Thomas G. Lane.
  5.  * This file is part of the Independent JPEG Group's software.
  6.  * For conditions of distribution and use, see the accompanying README file.
  7.  *
  8.  * This file contains 2-pass color quantization (color mapping) routines.
  9.  * These routines are invoked via the methods color_quant_prescan,
  10.  * color_quant_doit, and color_quant_init/term.
  11.  */
  12.  
  13. #include "jinclude.h"
  14.  
  15. #ifdef QUANT_2PASS_SUPPORTED
  16.  
  17. #include "viewdef.h"
  18. extern int more_defaults;
  19.  
  20. /*
  21.  * This module implements the well-known Heckbert paradigm for color
  22.  * quantization.  Most of the ideas used here can be traced back to
  23.  * Heckbert's seminal paper
  24.  *   Heckbert, Paul.  "Color Image Quantization for Frame Buffer Display",
  25.  *   Proc. SIGGRAPH '82, Computer Graphics v.16 #3 (July 1982), pp 297-304.
  26.  *
  27.  * In the first pass over the image, we accumulate a histogram showing the
  28.  * usage count of each possible color.  (To keep the histogram to a reasonable
  29.  * size, we reduce the precision of the input; typical practice is to retain
  30.  * 5 or 6 bits per color, so that 8 or 4 different input values are counted
  31.  * in the same histogram cell.)  Next, the color-selection step begins with a
  32.  * box representing the whole color space, and repeatedly splits the "largest"
  33.  * remaining box until we have as many boxes as desired colors.  Then the mean
  34.  * color in each remaining box becomes one of the possible output colors.
  35.  * The second pass over the image maps each input pixel to the closest output
  36.  * color (optionally after applying a Floyd-Steinberg dithering correction).
  37.  * This mapping is logically trivial, but making it go fast enough requires
  38.  * considerable care.
  39.  *
  40.  * Heckbert-style quantizers vary a good deal in their policies for choosing
  41.  * the "largest" box and deciding where to cut it.  The particular policies
  42.  * used here have proved out well in experimental comparisons, but better ones
  43.  * may yet be found.
  44.  *
  45.  * The most significant difference between this quantizer and others is that
  46.  * this one is intended to operate in YCbCr colorspace, rather than RGB space
  47.  * as is usually done.  Actually we work in scaled YCbCr colorspace, where
  48.  * Y distances are inflated by a factor of 2 relative to Cb or Cr distances.
  49.  * The empirical evidence is that distances in this space correspond to
  50.  * perceptual color differences more closely than do distances in RGB space;
  51.  * and working in this space is inexpensive within a JPEG decompressor, since
  52.  * the input data is already in YCbCr form.  (We could transform to an even
  53.  * more perceptually linear space such as Lab or Luv, but that is very slow
  54.  * and doesn't yield much better results than scaled YCbCr.)
  55.  */
  56.  
  57. #define Y_SCALE 2        /* scale Y distances up by this much */
  58.  
  59. #define MAXNUMCOLORS  (MAXJSAMPLE+1) /* maximum size of colormap */
  60.  
  61.  
  62. /*
  63.  * First we have the histogram data structure and routines for creating it.
  64.  *
  65.  * For work in YCbCr space, it is useful to keep more precision for Y than
  66.  * for Cb or Cr.  We recommend keeping 6 bits for Y and 5 bits each for Cb/Cr.
  67.  * If you have plenty of memory and cycles, 6 bits all around gives marginally
  68.  * better results; if you are short of memory, 5 bits all around will save
  69.  * some space but degrade the results.
  70.  * To maintain a fully accurate histogram, we'd need to allocate a "long"
  71.  * (preferably unsigned long) for each cell.  In practice this is overkill;
  72.  * we can get by with 16 bits per cell.  Few of the cell counts will overflow,
  73.  * and clamping those that do overflow to the maximum value will give close-
  74.  * enough results.  This reduces the recommended histogram size from 256Kb
  75.  * to 128Kb, which is a useful savings on PC-class machines.
  76.  * (In the second pass the histogram space is re-used for pixel mapping data;
  77.  * in that capacity, each cell must be able to store zero to the number of
  78.  * desired colors.  16 bits/cell is plenty for that too.)
  79.  * Since the JPEG code is intended to run in small memory model on 80x86
  80.  * machines, we can't just allocate the histogram in one chunk.  Instead
  81.  * of a true 3-D array, we use a row of pointers to 2-D arrays.  Each
  82.  * pointer corresponds to a Y value (typically 2^6 = 64 pointers) and
  83.  * each 2-D array has 2^5^2 = 1024 or 2^6^2 = 4096 entries.  Note that
  84.  * on 80x86 machines, the pointer row is in near memory but the actual
  85.  * arrays are in far memory (same arrangement as we use for image arrays).
  86.  */
  87.  
  88. #ifndef HIST_Y_BITS        /* so you can override from Makefile */
  89. #define HIST_Y_BITS  6        /* bits of precision in Y histogram */
  90. #endif
  91. #ifndef HIST_C_BITS        /* so you can override from Makefile */
  92. #define HIST_C_BITS  5        /* bits of precision in Cb/Cr histogram */
  93. #endif
  94.  
  95. #define HIST_Y_ELEMS  (1<<HIST_Y_BITS) /* # of elements along histogram axes */
  96. #define HIST_C_ELEMS  (1<<HIST_C_BITS)
  97.  
  98. /* These are the amounts to shift an input value to get a histogram index.
  99.  * For a combination 8/12 bit implementation, would need variables here...
  100.  */
  101.  
  102. #define Y_SHIFT  (BITS_IN_JSAMPLE-HIST_Y_BITS)
  103. #define C_SHIFT  (BITS_IN_JSAMPLE-HIST_C_BITS)
  104.  
  105.  
  106. typedef UINT16 histcell;    /* histogram cell; MUST be an unsigned type */
  107.  
  108. typedef histcell FAR * histptr;    /* for pointers to histogram cells */
  109.  
  110. typedef histcell hist1d[HIST_C_ELEMS]; /* typedefs for the array */
  111. typedef hist1d FAR * hist2d;    /* type for the Y-level pointers */
  112. typedef hist2d * hist3d;    /* type for top-level pointer */
  113.  
  114. static hist3d histogram;    /* pointer to the histogram */
  115.  
  116.  
  117. /*
  118.  * Prescan some rows of pixels.
  119.  * In this module the prescan simply updates the histogram, which has been
  120.  * initialized to zeroes by color_quant_init.
  121.  * Note: workspace is probably not useful for this routine, but it is passed
  122.  * anyway to allow some code sharing within the pipeline controller.
  123.  */
  124.  
  125. METHODDEF void
  126. color_quant_prescan (decompress_info_ptr cinfo, int num_rows,
  127.               JSAMPIMAGE image_data, JSAMPARRAY workspace)
  128. {
  129.   register JSAMPROW ptr0, ptr1, ptr2;
  130.   register histptr histp;
  131.   register int c0, c1, c2;
  132.   int row;
  133.   long col;
  134.   long width = cinfo->image_width;
  135.  
  136.   for (row = 0; row < num_rows; row++) {
  137.      ptr0 = image_data[0][row];
  138.      ptr1 = image_data[1][row];
  139.      ptr2 = image_data[2][row];
  140.      for (col = width; col > 0; col--) {
  141.         /* get pixel value and index into the histogram */
  142.         c0 = GETJSAMPLE(*ptr0++) >> Y_SHIFT;
  143.         c1 = GETJSAMPLE(*ptr1++) >> C_SHIFT;
  144.       c2 = GETJSAMPLE(*ptr2++) >> C_SHIFT;
  145.       histp = & histogram[c0][c1][c2];
  146.       /* increment, check for overflow and undo increment if so. */
  147.       /* We assume unsigned representation here! */
  148.       if (++(*histp) == 0)
  149.     (*histp)--;
  150.     }
  151.   }
  152. }
  153.  
  154.  
  155.  
  156. /*
  157.  * Now we have the really interesting routines: selection of a colormap
  158.  * given the completed histogram.
  159.  * These routines work with a list of "boxes", each representing a rectangular
  160.  * subset of the input color space (to histogram precision).
  161.  */
  162.  
  163. typedef struct {
  164.     /* The bounds of the box (inclusive); expressed as histogram indexes */
  165.     int c0min, c0max;
  166.     int c1min, c1max;
  167.     int c2min, c2max;
  168.     /* The number of nonzero histogram cells within this box */
  169.     long colorcount;
  170.       } box;
  171. typedef box * boxptr;
  172.  
  173. static boxptr boxlist;        /* array with room for desired # of boxes */
  174. static int numboxes;        /* number of boxes currently in boxlist */
  175.  
  176. static JSAMPARRAY my_colormap;    /* the finished colormap (in YCbCr space) */
  177.  
  178.  
  179. LOCAL boxptr
  180. find_biggest_color_pop (void)
  181. /* Find the splittable box with the largest color population */
  182. /* Returns NULL if no splittable boxes remain */
  183. {
  184.   register boxptr boxp;
  185.   register int i;
  186.   register long max = 0;
  187.   boxptr which = NULL;
  188.  
  189.   for (i = 0, boxp = boxlist; i < numboxes; i++, boxp++) {
  190.     if (boxp->colorcount > max) {
  191.       if (boxp->c0max > boxp->c0min || boxp->c1max > boxp->c1min ||
  192.       boxp->c2max > boxp->c2min) {
  193.     which = boxp;
  194.     max = boxp->colorcount;
  195.       }
  196.     }
  197.   }
  198.   return which;
  199. }
  200.  
  201.  
  202. LOCAL boxptr
  203. find_biggest_volume (void)
  204. /* Find the splittable box with the largest (scaled) volume */
  205. /* Returns NULL if no splittable boxes remain */
  206. {
  207.   register boxptr boxp;
  208.   register int i;
  209.   register INT32 max = 0;
  210.   register INT32 norm, c0,c1,c2;
  211.   boxptr which = NULL;
  212.   
  213.   /* We use 2-norm rather than real volume here.
  214.    * Some care is needed since the differences are expressed in
  215.    * histogram-cell units; if HIST_Y_BITS != HIST_C_BITS, we have to
  216.    * adjust the scaling to get the proper scaled-YCbCr-space distance.
  217.    * This code won't work right if HIST_Y_BITS < HIST_C_BITS,
  218.    * but that shouldn't ever be true.
  219.    * Note norm > 0 iff box is splittable, so need not check separately.
  220.    */
  221.   
  222.   for (i = 0, boxp = boxlist; i < numboxes; i++, boxp++) {
  223.     c0 = (boxp->c0max - boxp->c0min) * Y_SCALE;
  224.     c1 = (boxp->c1max - boxp->c1min) << (HIST_Y_BITS-HIST_C_BITS);
  225.     c2 = (boxp->c2max - boxp->c2min) << (HIST_Y_BITS-HIST_C_BITS);
  226.     norm = c0*c0 + c1*c1 + c2*c2;
  227.     if (norm > max) {
  228.       which = boxp;
  229.       max = norm;
  230.     }
  231.   }
  232.   return which;
  233. }
  234.  
  235.  
  236. LOCAL void
  237. update_box (boxptr boxp)
  238. /* Shrink the min/max bounds of a box to enclose only nonzero elements, */
  239. /* and recompute its population */
  240. {
  241.   histptr histp;
  242.   int c0,c1,c2;
  243.   int c0min,c0max,c1min,c1max,c2min,c2max;
  244.   long ccount;
  245.   
  246.   c0min = boxp->c0min;  c0max = boxp->c0max;
  247.   c1min = boxp->c1min;  c1max = boxp->c1max;
  248.   c2min = boxp->c2min;  c2max = boxp->c2max;
  249.   
  250.   if (c0max > c0min)
  251.     for (c0 = c0min; c0 <= c0max; c0++)
  252.       for (c1 = c1min; c1 <= c1max; c1++) {
  253.     histp = & histogram[c0][c1][c2min];
  254.     for (c2 = c2min; c2 <= c2max; c2++)
  255.       if (*histp++ != 0) {
  256.         boxp->c0min = c0min = c0;
  257.         goto have_c0min;
  258.       }
  259.       }
  260.  have_c0min:
  261.   if (c0max > c0min)
  262.     for (c0 = c0max; c0 >= c0min; c0--)
  263.       for (c1 = c1min; c1 <= c1max; c1++) {
  264.     histp = & histogram[c0][c1][c2min];
  265.     for (c2 = c2min; c2 <= c2max; c2++)
  266.       if (*histp++ != 0) {
  267.         boxp->c0max = c0max = c0;
  268.         goto have_c0max;
  269.       }
  270.       }
  271.  have_c0max:
  272.   if (c1max > c1min)
  273.     for (c1 = c1min; c1 <= c1max; c1++)
  274.       for (c0 = c0min; c0 <= c0max; c0++) {
  275.     histp = & histogram[c0][c1][c2min];
  276.     for (c2 = c2min; c2 <= c2max; c2++)
  277.       if (*histp++ != 0) {
  278.         boxp->c1min = c1min = c1;
  279.         goto have_c1min;
  280.       }
  281.       }
  282.  have_c1min:
  283.   if (c1max > c1min)
  284.     for (c1 = c1max; c1 >= c1min; c1--)
  285.       for (c0 = c0min; c0 <= c0max; c0++) {
  286.     histp = & histogram[c0][c1][c2min];
  287.     for (c2 = c2min; c2 <= c2max; c2++)
  288.       if (*histp++ != 0) {
  289.         boxp->c1max = c1max = c1;
  290.         goto have_c1max;
  291.       }
  292.       }
  293.  have_c1max:
  294.   if (c2max > c2min)
  295.     for (c2 = c2min; c2 <= c2max; c2++)
  296.       for (c0 = c0min; c0 <= c0max; c0++) {
  297.     histp = & histogram[c0][c1min][c2];
  298.     for (c1 = c1min; c1 <= c1max; c1++, histp += HIST_C_ELEMS)
  299.       if (*histp != 0) {
  300.         boxp->c2min = c2min = c2;
  301.         goto have_c2min;
  302.       }
  303.       }
  304.  have_c2min:
  305.   if (c2max > c2min)
  306.     for (c2 = c2max; c2 >= c2min; c2--)
  307.       for (c0 = c0min; c0 <= c0max; c0++) {
  308.     histp = & histogram[c0][c1min][c2];
  309.     for (c1 = c1min; c1 <= c1max; c1++, histp += HIST_C_ELEMS)
  310.       if (*histp != 0) {
  311.         boxp->c2max = c2max = c2;
  312.         goto have_c2max;
  313.       }
  314.       }
  315.  have_c2max:
  316.   
  317.   /* Now scan remaining volume of box and compute population */
  318.   ccount = 0;
  319.   for (c0 = c0min; c0 <= c0max; c0++)
  320.      for (c1 = c1min; c1 <= c1max; c1++) {
  321.       histp = & histogram[c0][c1][c2min];
  322.       for (c2 = c2min; c2 <= c2max; c2++, histp++)
  323.     if (*histp != 0) {
  324.       ccount++;
  325.     }
  326.     }
  327.   boxp->colorcount = ccount;
  328. }
  329.  
  330.  
  331. LOCAL void
  332. median_cut (int desired_colors)
  333. /* Repeatedly select and split the largest box until we have enough boxes */
  334. {
  335.   int n,lb;
  336.   int c0,c1,c2,cmax;
  337.   register boxptr b1,b2;
  338.  
  339.   while (numboxes < desired_colors) {
  340.     /* Select box to split */
  341.     /* Current algorithm: by population for first half, then by volume */
  342.     if (numboxes*2 <= desired_colors) {
  343.       b1 = find_biggest_color_pop();
  344.     } else {
  345.       b1 = find_biggest_volume();
  346.     }
  347.     if (b1 == NULL)        /* no splittable boxes left! */
  348.       break;
  349.     b2 = &boxlist[numboxes];    /* where new box will go */
  350.     /* Copy the color bounds to the new box. */
  351.     b2->c0max = b1->c0max; b2->c1max = b1->c1max; b2->c2max = b1->c2max;
  352.     b2->c0min = b1->c0min; b2->c1min = b1->c1min; b2->c2min = b1->c2min;
  353.     /* Choose which axis to split the box on.
  354.      * Current algorithm: longest scaled axis.
  355.      * See notes in find_biggest_volume about scaling...
  356.      */
  357.     c0 = (b1->c0max - b1->c0min) * Y_SCALE;
  358.     c1 = (b1->c1max - b1->c1min) << (HIST_Y_BITS-HIST_C_BITS);
  359.     c2 = (b1->c2max - b1->c2min) << (HIST_Y_BITS-HIST_C_BITS);
  360.     cmax = c0; n = 0;
  361.     if (c1 > cmax) { cmax = c1; n = 1; }
  362.     if (c2 > cmax) { n = 2; }
  363.     /* Choose split point along selected axis, and update box bounds.
  364.       * Current algorithm: split at halfway point.
  365.      * (Since the box has been shrunk to minimum volume,
  366.      * any split will produce two nonempty subboxes.)
  367.      * Note that lb value is max for lower box, so must be < old max.
  368.      */
  369.     switch (n) {
  370.     case 0:
  371.       lb = (b1->c0max + b1->c0min) / 2;
  372.       b1->c0max = lb;
  373.       b2->c0min = lb+1;
  374.       break;
  375.      case 1:
  376.       lb = (b1->c1max + b1->c1min) / 2;
  377.       b1->c1max = lb;
  378.       b2->c1min = lb+1;
  379.       break;
  380.     case 2:
  381.       lb = (b1->c2max + b1->c2min) / 2;
  382.       b1->c2max = lb;
  383.       b2->c2min = lb+1;
  384.       break;
  385.     }
  386.     /* Update stats for boxes */
  387.     update_box(b1);
  388.     update_box(b2);
  389.     numboxes++;
  390.   }
  391. }
  392.  
  393.  
  394. LOCAL void
  395. compute_color (boxptr boxp, int icolor)
  396. /* Compute representative color for a box, put it in my_colormap[icolor] */
  397. {
  398.   /* Current algorithm: mean weighted by pixels (not colors) */
  399.   /* Note it is important to get the rounding correct! */
  400.   histptr histp;
  401.   int c0,c1,c2;
  402.   int c0min,c0max,c1min,c1max,c2min,c2max;
  403.   long count;
  404.   long total = 0;
  405.   long c0total = 0;
  406.   long c1total = 0;
  407.   long c2total = 0;
  408.  
  409.   c0min = boxp->c0min;  c0max = boxp->c0max;
  410.   c1min = boxp->c1min;  c1max = boxp->c1max;
  411.   c2min = boxp->c2min;  c2max = boxp->c2max;
  412.   
  413.   for (c0 = c0min; c0 <= c0max; c0++)
  414.     for (c1 = c1min; c1 <= c1max; c1++) {
  415.       histp = & histogram[c0][c1][c2min];
  416.       for (c2 = c2min; c2 <= c2max; c2++) {
  417.     if ((count = *histp++) != 0) {
  418.       total += count;
  419.       c0total += ((c0 << Y_SHIFT) + ((1<<Y_SHIFT)>>1)) * count;
  420.       c1total += ((c1 << C_SHIFT) + ((1<<C_SHIFT)>>1)) * count;
  421.       c2total += ((c2 << C_SHIFT) + ((1<<C_SHIFT)>>1)) * count;
  422.     }
  423.       }
  424.     }
  425.   
  426.   my_colormap[0][icolor] = (JSAMPLE) ((c0total + (total>>1)) / total);
  427.   my_colormap[1][icolor] = (JSAMPLE) ((c1total + (total>>1)) / total);
  428.   my_colormap[2][icolor] = (JSAMPLE) ((c2total + (total>>1)) / total);
  429. }
  430.  
  431.  
  432. LOCAL void
  433. remap_colormap (decompress_info_ptr cinfo)
  434. /* Remap the internal colormap to the output colorspace */
  435. {
  436.   /* This requires a little trickery since color_convert expects to
  437.    * deal with 3-D arrays (a 2-D sample array for each component).
  438.    * We must promote the colormaps into one-row 3-D arrays.
  439.    */
  440.   short ci;
  441.   JSAMPARRAY input_hack[3];
  442.   JSAMPARRAY output_hack[10];    /* assume no more than 10 output components */
  443.  
  444.   for (ci = 0; ci < 3; ci++)
  445.     input_hack[ci] = &(my_colormap[ci]);
  446.   for (ci = 0; ci < cinfo->color_out_comps; ci++)
  447.     output_hack[ci] = &(cinfo->colormap[ci]);
  448.  
  449.   (*cinfo->methods->color_convert) (cinfo, 1,
  450.                     (long) cinfo->actual_number_of_colors,
  451.                     input_hack, output_hack);
  452. }
  453.  
  454.  
  455. LOCAL void
  456. select_colors (decompress_info_ptr cinfo)
  457. /* Master routine for color selection */
  458. {
  459.   int desired = cinfo->desired_number_of_colors;
  460.   int i;
  461.  
  462.   /* Allocate workspace for box list */
  463.   boxlist = (boxptr) (*cinfo->emethods->alloc_small) (desired * SIZEOF(box));
  464.   /* Initialize one box containing whole space */
  465.   numboxes = 1;
  466.   boxlist[0].c0min = 0;
  467.   boxlist[0].c0max = MAXJSAMPLE >> Y_SHIFT;
  468.   boxlist[0].c1min = 0;
  469.   boxlist[0].c1max = MAXJSAMPLE >> C_SHIFT;
  470.   boxlist[0].c2min = 0;
  471.   boxlist[0].c2max = MAXJSAMPLE >> C_SHIFT;
  472.   /* Shrink it to actually-used volume and set its statistics */
  473.   update_box(& boxlist[0]);
  474.   /* Perform median-cut to produce final box list */
  475.   median_cut(desired);
  476.   /* Compute the representative color for each box, fill my_colormap[] */
  477.   for (i = 0; i < numboxes; i++)
  478.     compute_color(& boxlist[i], i);
  479.   cinfo->actual_number_of_colors = numboxes;
  480.   /* Produce an output colormap in the desired output colorspace */
  481.   remap_colormap(cinfo);
  482.   TRACEMS1(cinfo->emethods, 1, "Selected %d colors for quantization",
  483.        numboxes);
  484.   /* Done with the box list */
  485.   (*cinfo->emethods->free_small) ((void *) boxlist);
  486. }
  487.  
  488.  
  489. /*
  490.  * These routines are concerned with the time-critical task of mapping input
  491.  * colors to the nearest color in the selected colormap.
  492.  *
  493.  * We re-use the histogram space as an "inverse color map", essentially a
  494.  * cache for the results of nearest-color searches.  All colors within a
  495.  * histogram cell will be mapped to the same colormap entry, namely the one
  496.  * closest to the cell's center.  This may not be quite the closest entry to
  497.  * the actual input color, but it's almost as good.  A zero in the cache
  498.  * indicates we haven't found the nearest color for that cell yet; the array
  499.  * is cleared to zeroes before starting the mapping pass.  When we find the
  500.  * nearest color for a cell, its colormap index plus one is recorded in the
  501.  * cache for future use.  The pass2 scanning routines call fill_inverse_cmap
  502.  * when they need to use an unfilled entry in the cache.
  503.  *
  504.  * Our method of efficiently finding nearest colors is based on the "locally
  505.  * sorted search" idea described by Heckbert and on the incremental distance
  506.  * calculation described by Spencer W. Thomas in chapter III.1 of Graphics
  507.  * Gems II (James Arvo, ed.  Academic Press, 1991).  Thomas points out that
  508.  * the distances from a given colormap entry to each cell of the histogram can
  509.  * be computed quickly using an incremental method: the differences between
  510.  * distances to adjacent cells themselves differ by a constant.  This allows a
  511.  * fairly fast implementation of the "brute force" approach of computing the
  512.  * distance from every colormap entry to every histogram cell.  Unfortunately,
  513.  * it needs a work array to hold the best-distance-so-far for each histogram
  514.  * cell (because the inner loop has to be over cells, not colormap entries).
  515.  * The work array elements have to be INT32s, so the work array would need
  516.  * 256Kb at our recommended precision.  This is not feasible in DOS machines.
  517.  * Another disadvantage of the brute force approach is that it computes
  518.  * distances to every cell of the cubical histogram.  When working with YCbCr
  519.  * input, only about a quarter of the cube represents realizable colors, so
  520.  * many of the cells will never be used and filling them is wasted effort.
  521.  *
  522.  * To get around these problems, we apply Thomas' method to compute the
  523.  * nearest colors for only the cells within a small subbox of the histogram.
  524.  * The work array need be only as big as the subbox, so the memory usage
  525.  * problem is solved.  A subbox is processed only when some cell in it is
  526.  * referenced by the pass2 routines, so we will never bother with cells far
  527.  * outside the realizable color volume.  An additional advantage of this
  528.  * approach is that we can apply Heckbert's locality criterion to quickly
  529.  * eliminate colormap entries that are far away from the subbox; typically
  530.  * three-fourths of the colormap entries are rejected by Heckbert's criterion,
  531.  * and we need not compute their distances to individual cells in the subbox.
  532.  * The speed of this approach is heavily influenced by the subbox size: too
  533.  * small means too much overhead, too big loses because Heckbert's criterion
  534.  * can't eliminate as many colormap entries.  Empirically the best subbox
  535.  * size seems to be about 1/512th of the histogram (1/8th in each direction).
  536.  *
  537.  * Thomas' article also describes a refined method which is asymptotically
  538.  * faster than the brute-force method, but it is also far more complex and
  539.  * cannot efficiently be applied to small subboxes.  It is therefore not
  540.  * useful for programs intended to be portable to DOS machines.  On machines
  541.  * with plenty of memory, filling the whole histogram in one shot with Thomas'
  542.  * refined method might be faster than the present code --- but then again,
  543.  * it might not be any faster, and it's certainly more complicated.
  544.  */
  545.  
  546.  
  547. #ifndef BOX_Y_LOG        /* so you can override from Makefile */
  548. #define BOX_Y_LOG  (HIST_Y_BITS-3) /* log2(hist cells in update box, Y axis) */
  549. #endif
  550. #ifndef BOX_C_LOG        /* so you can override from Makefile */
  551. #define BOX_C_LOG  (HIST_C_BITS-3) /* log2(hist cells in update box, C axes) */
  552. #endif
  553.  
  554. #define BOX_Y_ELEMS  (1<<BOX_Y_LOG) /* # of hist cells in update box */
  555. #define BOX_C_ELEMS  (1<<BOX_C_LOG)
  556.  
  557. #define BOX_Y_SHIFT  (Y_SHIFT + BOX_Y_LOG)
  558. #define BOX_C_SHIFT  (C_SHIFT + BOX_C_LOG)
  559.  
  560.  
  561. /*
  562.  * The next three routines implement inverse colormap filling.  They could
  563.  * all be folded into one big routine, but splitting them up this way saves
  564.  * some stack space (the mindist[] and bestdist[] arrays need not coexist)
  565.  * and may allow some compilers to produce better code by registerizing more
  566.  * inner-loop variables.
  567.  */
  568.  
  569. LOCAL int
  570. find_nearby_colors (decompress_info_ptr cinfo, int minc0, int minc1, int minc2,
  571.             JSAMPLE colorlist[])
  572. /* Locate the colormap entries close enough to an update box to be candidates
  573.  * for the nearest entry to some cell(s) in the update box.  The update box
  574.  * is specified by the center coordinates of its first cell.  The number of
  575.  * candidate colormap entries is returned, and their colormap indexes are
  576.  * placed in colorlist[].
  577.  * This routine uses Heckbert's "locally sorted search" criterion to select
  578.  * the colors that need further consideration.
  579.  */
  580. {
  581.   int numcolors = cinfo->actual_number_of_colors;
  582.   int maxc0, maxc1, maxc2;
  583.   int centerc0, centerc1, centerc2;
  584.   int i, x, ncolors;
  585.   INT32 minmaxdist, min_dist, max_dist, tdist;
  586.   INT32 mindist[MAXNUMCOLORS];    /* min distance to colormap entry i */
  587.  
  588.   /* Compute true coordinates of update box's upper corner and center.
  589.    * Actually we compute the coordinates of the center of the upper-corner
  590.    * histogram cell, which are the upper bounds of the volume we care about.
  591.    * Note that since ">>" rounds down, the "center" values may be closer to
  592.    * min than to max; hence comparisons to them must be "<=", not "<".
  593.    */
  594.   maxc0 = minc0 + ((1 << BOX_Y_SHIFT) - (1 << Y_SHIFT));
  595.   centerc0 = (minc0 + maxc0) >> 1;
  596.   maxc1 = minc1 + ((1 << BOX_C_SHIFT) - (1 << C_SHIFT));
  597.   centerc1 = (minc1 + maxc1) >> 1;
  598.   maxc2 = minc2 + ((1 << BOX_C_SHIFT) - (1 << C_SHIFT));
  599.   centerc2 = (minc2 + maxc2) >> 1;
  600.  
  601.   /* For each color in colormap, find:
  602.    *  1. its minimum squared-distance to any point in the update box
  603.    *     (zero if color is within update box);
  604.    *  2. its maximum squared-distance to any point in the update box.
  605.    * Both of these can be found by considering only the corners of the box.
  606.    * We save the minimum distance for each color in mindist[];
  607.    * only the smallest maximum distance is of interest.
  608.    * Note we have to scale Y to get correct distance in scaled space.
  609.    */
  610.   minmaxdist = 0x7FFFFFFFL;
  611.  
  612.   for (i = 0; i < numcolors; i++) {
  613.     /* We compute the squared-c0-distance term, then add in the other two. */
  614.     x = GETJSAMPLE(my_colormap[0][i]);
  615.     if (x < minc0) {
  616.       tdist = (x - minc0) * Y_SCALE;
  617.       min_dist = tdist*tdist;
  618.       tdist = (x - maxc0) * Y_SCALE;
  619.       max_dist = tdist*tdist;
  620.     } else if (x > maxc0) {
  621.       tdist = (x - maxc0) * Y_SCALE;
  622.       min_dist = tdist*tdist;
  623.       tdist = (x - minc0) * Y_SCALE;
  624.       max_dist = tdist*tdist;
  625.     } else {
  626.       /* within cell range so no contribution to min_dist */
  627.       min_dist = 0;
  628.         if (x <= centerc0) {
  629.     tdist = (x - maxc0) * Y_SCALE;
  630.     max_dist = tdist*tdist;
  631.       } else {
  632.     tdist = (x - minc0) * Y_SCALE;
  633.     max_dist = tdist*tdist;
  634.       }
  635.     }
  636.  
  637.     x = GETJSAMPLE(my_colormap[1][i]);
  638.     if (x < minc1) {
  639.         tdist = x - minc1;
  640.       min_dist += tdist*tdist;
  641.       tdist = x - maxc1;
  642.       max_dist += tdist*tdist;
  643.     } else if (x > maxc1) {
  644.       tdist = x - maxc1;
  645.       min_dist += tdist*tdist;
  646.       tdist = x - minc1;
  647.       max_dist += tdist*tdist;
  648.     } else {
  649.       /* within cell range so no contribution to min_dist */
  650.       if (x <= centerc1) {
  651.     tdist = x - maxc1;
  652.     max_dist += tdist*tdist;
  653.       } else {
  654.     tdist = x - minc1;
  655.     max_dist += tdist*tdist;
  656.       }
  657.     }
  658.  
  659.     x = GETJSAMPLE(my_colormap[2][i]);
  660.     if (x < minc2) {
  661.       tdist = x - minc2;
  662.       min_dist += tdist*tdist;
  663.       tdist = x - maxc2;
  664.       max_dist += tdist*tdist;
  665.     } else if (x > maxc2) {
  666.       tdist = x - maxc2;
  667.       min_dist += tdist*tdist;
  668.       tdist = x - minc2;
  669.       max_dist += tdist*tdist;
  670.     } else {
  671.       /* within cell range so no contribution to min_dist */
  672.         if (x <= centerc2) {
  673.     tdist = x - maxc2;
  674.     max_dist += tdist*tdist;
  675.       } else {
  676.     tdist = x - minc2;
  677.     max_dist += tdist*tdist;
  678.       }
  679.     }
  680.  
  681.     mindist[i] = min_dist;    /* save away the results */
  682.     if (max_dist < minmaxdist)
  683.         minmaxdist = max_dist;
  684.   }
  685.  
  686.   /* Now we know that no cell in the update box is more than minmaxdist
  687.    * away from some colormap entry.  Therefore, only colors that are
  688.    * within minmaxdist of some part of the box need be considered.
  689.    */
  690.   ncolors = 0;
  691.   for (i = 0; i < numcolors; i++) {
  692.     if (mindist[i] <= minmaxdist)
  693.       colorlist[ncolors++] = (JSAMPLE) i;
  694.   }
  695.   return ncolors;
  696. }
  697.  
  698.  
  699. LOCAL void
  700. find_best_colors (decompress_info_ptr cinfo, int minc0, int minc1, int minc2,
  701.           int numcolors, JSAMPLE colorlist[], JSAMPLE bestcolor[])
  702. /* Find the closest colormap entry for each cell in the update box,
  703.  * given the list of candidate colors prepared by find_nearby_colors.
  704.  * Return the indexes of the closest entries in the bestcolor[] array.
  705.  * This routine uses Thomas' incremental distance calculation method to
  706.  * find the distance from a colormap entry to successive cells in the box.
  707.  */
  708. {
  709.   int ic0, ic1, ic2;
  710.   int i, icolor;
  711.   register INT32 * bptr;    /* pointer into bestdist[] array */
  712.   JSAMPLE * cptr;        /* pointer into bestcolor[] array */
  713.   INT32 dist0, dist1;        /* initial distance values */
  714.   register INT32 dist2;        /* current distance in inner loop */
  715.   INT32 xx0, xx1;        /* distance increments */
  716.   register INT32 xx2;
  717.   INT32 inc0, inc1, inc2;    /* initial values for increments */
  718.   /* This array holds the distance to the nearest-so-far color for each cell */
  719.   INT32 bestdist[BOX_Y_ELEMS * BOX_C_ELEMS * BOX_C_ELEMS];
  720.  
  721.   /* Initialize best-distance for each cell of the update box */
  722.   bptr = bestdist;
  723.   for (i = BOX_Y_ELEMS*BOX_C_ELEMS*BOX_C_ELEMS-1; i >= 0; i--)
  724.     *bptr++ = 0x7FFFFFFFL;
  725.   
  726.   /* For each color selected by find_nearby_colors,
  727.     * compute its distance to the center of each cell in the box.
  728.    * If that's less than best-so-far, update best distance and color number.
  729.    * Note we have to scale Y to get correct distance in scaled space.
  730.    */
  731.   
  732.   /* Nominal steps between cell centers ("x" in Thomas article) */
  733. #define STEP_Y  ((1 << Y_SHIFT) * Y_SCALE)
  734. #define STEP_C  (1 << C_SHIFT)
  735.   
  736.   for (i = 0; i < numcolors; i++) {
  737.     icolor = GETJSAMPLE(colorlist[i]);
  738.     /* Compute (square of) distance from minc0/c1/c2 to this color */
  739.     inc0 = (minc0 - (int) GETJSAMPLE(my_colormap[0][icolor])) * Y_SCALE;
  740.     dist0 = inc0*inc0;
  741.     inc1 = minc1 - (int) GETJSAMPLE(my_colormap[1][icolor]);
  742.     dist0 += inc1*inc1;
  743.     inc2 = minc2 - (int) GETJSAMPLE(my_colormap[2][icolor]);
  744.     dist0 += inc2*inc2;
  745.     /* Form the initial difference increments */
  746.     inc0 = inc0 * (2 * STEP_Y) + STEP_Y * STEP_Y;
  747.     inc1 = inc1 * (2 * STEP_C) + STEP_C * STEP_C;
  748.     inc2 = inc2 * (2 * STEP_C) + STEP_C * STEP_C;
  749.     /* Now loop over all cells in box, updating distance per Thomas method */
  750.     bptr = bestdist;
  751.     cptr = bestcolor;
  752.     xx0 = inc0;
  753.     for (ic0 = BOX_Y_ELEMS-1; ic0 >= 0; ic0--) {
  754.       dist1 = dist0;
  755.       xx1 = inc1;
  756.       for (ic1 = BOX_C_ELEMS-1; ic1 >= 0; ic1--) {
  757.     dist2 = dist1;
  758.     xx2 = inc2;
  759.     for (ic2 = BOX_C_ELEMS-1; ic2 >= 0; ic2--) {
  760.       if (dist2 < *bptr) {
  761.         *bptr = dist2;
  762.         *cptr = (JSAMPLE) icolor;
  763.       }
  764.       dist2 += xx2;
  765.       xx2 += 2 * STEP_C * STEP_C;
  766.       bptr++;
  767.       cptr++;
  768.     }
  769.     dist1 += xx1;
  770.     xx1 += 2 * STEP_C * STEP_C;
  771.         }
  772.       dist0 += xx0;
  773.       xx0 += 2 * STEP_Y * STEP_Y;
  774.     }
  775.   }
  776. }
  777.  
  778.  
  779. LOCAL void
  780. fill_inverse_cmap (decompress_info_ptr cinfo, int c0, int c1, int c2)
  781. /* Fill the inverse-colormap entries in the update box that contains */
  782. /* histogram cell c0/c1/c2.  (Only that one cell MUST be filled, but */
  783. /* we can fill as many others as we wish.) */
  784. {
  785.   int minc0, minc1, minc2;    /* lower left corner of update box */
  786.   int ic0, ic1, ic2;
  787.   register JSAMPLE * cptr;    /* pointer into bestcolor[] array */
  788.   register histptr cachep;    /* pointer into main cache array */
  789.   /* This array lists the candidate colormap indexes. */
  790.   JSAMPLE colorlist[MAXNUMCOLORS];
  791.   int numcolors;        /* number of candidate colors */
  792.   /* This array holds the actually closest colormap index for each cell. */
  793.   JSAMPLE bestcolor[BOX_Y_ELEMS * BOX_C_ELEMS * BOX_C_ELEMS];
  794.  
  795.   /* Convert cell coordinates to update box ID */
  796.   c0 >>= BOX_Y_LOG;
  797.   c1 >>= BOX_C_LOG;
  798.   c2 >>= BOX_C_LOG;
  799.  
  800.   /* Compute true coordinates of update box's origin corner.
  801.    * Actually we compute the coordinates of the center of the corner
  802.    * histogram cell, which are the lower bounds of the volume we care about.
  803.    */
  804.   minc0 = (c0 << BOX_Y_SHIFT) + ((1 << Y_SHIFT) >> 1);
  805.   minc1 = (c1 << BOX_C_SHIFT) + ((1 << C_SHIFT) >> 1);
  806.   minc2 = (c2 << BOX_C_SHIFT) + ((1 << C_SHIFT) >> 1);
  807.   
  808.   /* Determine which colormap entries are close enough to be candidates
  809.    * for the nearest entry to some cell in the update box.
  810.    */
  811.   numcolors = find_nearby_colors(cinfo, minc0, minc1, minc2, colorlist);
  812.  
  813.   /* Determine the actually nearest colors. */
  814.   find_best_colors(cinfo, minc0, minc1, minc2, numcolors, colorlist,
  815.             bestcolor);
  816.  
  817.   /* Save the best color numbers (plus 1) in the main cache array */
  818.   c0 <<= BOX_Y_LOG;        /* convert ID back to base cell indexes */
  819.   c1 <<= BOX_C_LOG;
  820.   c2 <<= BOX_C_LOG;
  821.   cptr = bestcolor;
  822.   for (ic0 = 0; ic0 < BOX_Y_ELEMS; ic0++) {
  823.     for (ic1 = 0; ic1 < BOX_C_ELEMS; ic1++) {
  824.       cachep = & histogram[c0+ic0][c1+ic1][c2];
  825.       for (ic2 = 0; ic2 < BOX_C_ELEMS; ic2++) {
  826.     *cachep++ = (histcell) (GETJSAMPLE(*cptr++) + 1);
  827.       }
  828.     }
  829.   }
  830. }
  831.  
  832.  
  833. /*
  834.  * These routines perform second-pass scanning of the image: map each pixel to
  835.  * the proper colormap index, and output the indexes to the output file.
  836.  *
  837.  * output_workspace is a one-component array of pixel dimensions at least
  838.  * as large as the input image strip; it can be used to hold the converted
  839.  * pixels' colormap indexes.
  840.  */
  841.  
  842. METHODDEF void
  843. pass2_nodither (decompress_info_ptr cinfo, int num_rows,
  844.         JSAMPIMAGE image_data, JSAMPARRAY output_workspace)
  845. /* This version performs no dithering */
  846. {
  847.   register JSAMPROW ptr0, ptr1, ptr2, outptr;
  848.   register histptr cachep;
  849.   register int c0, c1, c2;
  850.   int row;
  851.   long col;
  852.   long width = cinfo->image_width;
  853.  
  854.   /* Convert data to colormap indexes, which we save in output_workspace */
  855.   for (row = 0; row < num_rows; row++) {
  856.     ptr0 = image_data[0][row];
  857.     ptr1 = image_data[1][row];
  858.     ptr2 = image_data[2][row];
  859.      outptr = output_workspace[row];
  860.     for (col = width; col > 0; col--) {
  861.       /* get pixel value and index into the cache */
  862.       c0 = GETJSAMPLE(*ptr0++) >> Y_SHIFT;
  863.       c1 = GETJSAMPLE(*ptr1++) >> C_SHIFT;
  864.       c2 = GETJSAMPLE(*ptr2++) >> C_SHIFT;
  865.       cachep = & histogram[c0][c1][c2];
  866.       /* If we have not seen this color before, find nearest colormap entry */
  867.       /* and update the cache */
  868.       if (*cachep == 0)
  869.     fill_inverse_cmap(cinfo, c0,c1,c2);
  870.       /* Now emit the colormap index for this cell */
  871.       *outptr++ = (JSAMPLE) (*cachep - 1);
  872.     }
  873.   }
  874.   /* Emit converted rows to the output file */
  875.   (*cinfo->methods->put_pixel_rows) (cinfo, num_rows, &output_workspace);
  876. }
  877.  
  878.  
  879. /* Declarations for Floyd-Steinberg dithering.
  880.  *
  881.  * Errors are accumulated into the array fserrors[], at a resolution of
  882.  * 1/16th of a pixel count.  The error at a given pixel is propagated
  883.  * to its not-yet-processed neighbors using the standard F-S fractions,
  884.  *        ...    (here)    7/16
  885.  *        3/16    5/16    1/16
  886.  * We work left-to-right on even rows, right-to-left on odd rows.
  887.  *
  888.  * We can get away with a single array (holding one row's worth of errors)
  889.  * by using it to store the current row's errors at pixel columns not yet
  890.  * processed, but the next row's errors at columns already processed.  We
  891.  * need only a few extra variables to hold the errors immediately around the
  892.  * current column.  (If we are lucky, those variables are in registers, but
  893.  * even if not, they're probably cheaper to access than array elements are.)
  894.  *
  895.  * The fserrors[] array has (#columns + 2) entries; the extra entry at
  896.  * each end saves us from special-casing the first and last pixels.
  897.  * Each entry is three values long, one value for each color component.
  898.  *
  899.  * Note: on a wide image, we might not have enough room in a PC's near data
  900.  * segment to hold the error array; so it is allocated with alloc_medium.
  901.  */
  902.  
  903. #ifdef EIGHT_BIT_SAMPLES
  904. typedef INT16 FSERROR;        /* 16 bits should be enough */
  905. typedef int LOCFSERROR;        /* use 'int' for calculation temps */
  906. #else
  907. typedef INT32 FSERROR;        /* may need more than 16 bits */
  908. typedef INT32 LOCFSERROR;    /* be sure calculation temps are big enough */
  909. #endif
  910.  
  911. typedef FSERROR FAR *FSERRPTR;    /* pointer to error array (in FAR storage!) */
  912.  
  913. static FSERRPTR fserrors;    /* accumulated errors */
  914. static boolean on_odd_row;    /* flag to remember which row we are on */
  915.  
  916.  
  917. METHODDEF void
  918. pass2_dither (decompress_info_ptr cinfo, int num_rows,
  919.           JSAMPIMAGE image_data, JSAMPARRAY output_workspace)
  920. /* This version performs Floyd-Steinberg dithering */
  921. {
  922.   register LOCFSERROR cur0, cur1, cur2;    /* current error or pixel value */
  923.   LOCFSERROR belowerr0, belowerr1, belowerr2; /* error for pixel below cur */
  924.   LOCFSERROR bpreverr0, bpreverr1, bpreverr2; /* error for below/prev col */
  925.   register FSERRPTR errorptr;    /* => fserrors[] at column before current */
  926.   JSAMPROW ptr0, ptr1, ptr2;    /* => current input pixel */
  927.   JSAMPROW outptr;        /* => current output pixel */
  928.   histptr cachep;
  929.   int dir;            /* +1 or -1 depending on direction */
  930.   int dir3;            /* 3*dir, for advancing errorptr */
  931.   int row;
  932.   long col;
  933.   long width = cinfo->image_width;
  934.   JSAMPLE *range_limit = cinfo->sample_range_limit;
  935.   JSAMPROW colormap0 = my_colormap[0];
  936.   JSAMPROW colormap1 = my_colormap[1];
  937.   JSAMPROW colormap2 = my_colormap[2];
  938.   SHIFT_TEMPS
  939.  
  940.   /* Convert data to colormap indexes, which we save in output_workspace */
  941.   for (row = 0; row < num_rows; row++) {
  942.     ptr0 = image_data[0][row];
  943.     ptr1 = image_data[1][row];
  944.     ptr2 = image_data[2][row];
  945.     outptr = output_workspace[row];
  946.     if (on_odd_row) {
  947.         /* work right to left in this row */
  948.       ptr0 += width - 1;    /* so point to rightmost pixel */
  949.       ptr1 += width - 1;
  950.       ptr2 += width - 1;
  951.       outptr += width - 1;
  952.       dir = -1;
  953.       dir3 = -3;
  954.       errorptr = fserrors + (width+1)*3; /* point to entry after last column */
  955.       on_odd_row = FALSE;    /* flip for next time */
  956.     } else {
  957.       /* work left to right in this row */
  958.       dir = 1;
  959.       dir3 = 3;
  960.       errorptr = fserrors;    /* point to entry before first real column */
  961.       on_odd_row = TRUE;    /* flip for next time */
  962.     }
  963.     /* Preset error values: no error propagated to first pixel from left */
  964.     cur0 = cur1 = cur2 = 0;
  965.     /* and no error propagated to row below yet */
  966.     belowerr0 = belowerr1 = belowerr2 = 0;
  967.      bpreverr0 = bpreverr1 = bpreverr2 = 0;
  968.  
  969.      for (col = width; col > 0; col--) {
  970.         /* curN holds the error propagated from the previous pixel on the
  971.          * current line.  Add the error propagated from the previous line
  972.          * to form the complete error correction term for this pixel, and
  973.          * round the error term (which is expressed * 16) to an integer.
  974.          * RIGHT_SHIFT rounds towards minus infinity, so adding 8 is correct
  975.          * for either sign of the error value.
  976.          * Note: errorptr points to *previous* column's array entry.
  977.          */
  978.         cur0 = RIGHT_SHIFT(cur0 + errorptr[dir3+0] + 8, 4);
  979.         cur1 = RIGHT_SHIFT(cur1 + errorptr[dir3+1] + 8, 4);
  980.         cur2 = RIGHT_SHIFT(cur2 + errorptr[dir3+2] + 8, 4);
  981.  
  982. #define CLAMP 32
  983.          if (cur0 > CLAMP) cur0 = CLAMP;
  984.          else{ if (cur0 < -CLAMP) cur0 = -CLAMP;}
  985.  
  986.          if (cur1 > CLAMP) cur1 = CLAMP;
  987.          else{ if (cur1 < -CLAMP) cur1 = -CLAMP;}
  988.  
  989.          if (cur2 > CLAMP) cur2 = CLAMP;
  990.          else if (cur2 < -CLAMP) cur2 = -CLAMP;
  991.  
  992.         /* Form pixel value + error, and range-limit to 0..MAXJSAMPLE.
  993.          * The maximum error is +- MAXJSAMPLE; this sets the required size
  994.          * of the range_limit array.
  995.          */
  996.  
  997.         cur0 += GETJSAMPLE(*ptr0);
  998.         cur1 += GETJSAMPLE(*ptr1);
  999.       cur2 += GETJSAMPLE(*ptr2);
  1000.       cur0 = GETJSAMPLE(range_limit[cur0]);
  1001.       cur1 = GETJSAMPLE(range_limit[cur1]);
  1002.       cur2 = GETJSAMPLE(range_limit[cur2]);
  1003.       /* Index into the cache with adjusted pixel value */
  1004.       cachep = & histogram[cur0 >> Y_SHIFT][cur1 >> C_SHIFT][cur2 >> C_SHIFT];
  1005.       /* If we have not seen this color before, find nearest colormap */
  1006.       /* entry and update the cache */
  1007.       if (*cachep == 0)
  1008.     fill_inverse_cmap(cinfo, cur0>>Y_SHIFT, cur1>>C_SHIFT, cur2>>C_SHIFT);
  1009.       /* Now emit the colormap index for this cell */
  1010.       { register int pixcode = *cachep - 1;
  1011.     *outptr = (JSAMPLE) pixcode;
  1012.     /* Compute representation error for this pixel */
  1013.     cur0 -= GETJSAMPLE(colormap0[pixcode]);
  1014.     cur1 -= GETJSAMPLE(colormap1[pixcode]);
  1015.     cur2 -= GETJSAMPLE(colormap2[pixcode]);
  1016.       }
  1017.       /* Compute error fractions to be propagated to adjacent pixels.
  1018.        * Add these into the running sums, and simultaneously shift the
  1019.        * next-line error sums left by 1 column.
  1020.        */
  1021.       { register LOCFSERROR bnexterr, delta;
  1022.  
  1023.     bnexterr = cur0;    /* Process component 0 */
  1024.     delta = cur0 * 2;
  1025.     cur0 += delta;        /* form error * 3 */
  1026.     errorptr[0] = (FSERROR) (bpreverr0 + cur0);
  1027.     cur0 += delta;        /* form error * 5 */
  1028.     bpreverr0 = belowerr0 + cur0;
  1029.     belowerr0 = bnexterr;
  1030.     cur0 += delta;        /* form error * 7 */
  1031.     bnexterr = cur1;    /* Process component 1 */
  1032.     delta = cur1 * 2;
  1033.     cur1 += delta;        /* form error * 3 */
  1034.     errorptr[1] = (FSERROR) (bpreverr1 + cur1);
  1035.     cur1 += delta;        /* form error * 5 */
  1036.     bpreverr1 = belowerr1 + cur1;
  1037.     belowerr1 = bnexterr;
  1038.     cur1 += delta;        /* form error * 7 */
  1039.     bnexterr = cur2;    /* Process component 2 */
  1040.     delta = cur2 * 2;
  1041.     cur2 += delta;        /* form error * 3 */
  1042.     errorptr[2] = (FSERROR) (bpreverr2 + cur2);
  1043.     cur2 += delta;        /* form error * 5 */
  1044.     bpreverr2 = belowerr2 + cur2;
  1045.     belowerr2 = bnexterr;
  1046.     cur2 += delta;        /* form error * 7 */
  1047.       }
  1048.       /* At this point curN contains the 7/16 error value to be propagated
  1049.        * to the next pixel on the current line, and all the errors for the
  1050.        * next line have been shifted over.  We are therefore ready to move on.
  1051.        */
  1052.       ptr0 += dir;        /* Advance pixel pointers to next column */
  1053.       ptr1 += dir;
  1054.       ptr2 += dir;
  1055.       outptr += dir;
  1056.       errorptr += dir3;        /* advance errorptr to current column */
  1057.     }
  1058.     /* Post-loop cleanup: we must unload the final error values into the
  1059.      * final fserrors[] entry.  Note we need not unload belowerrN because
  1060.      * it is for the dummy column before or after the actual array.
  1061.      */
  1062.     errorptr[0] = (FSERROR) bpreverr0; /* unload prev errs into array */
  1063.     errorptr[1] = (FSERROR) bpreverr1;
  1064.     errorptr[2] = (FSERROR) bpreverr2;
  1065.   }
  1066.   /* Emit converted rows to the output file */
  1067.   (*cinfo->methods->put_pixel_rows) (cinfo, num_rows, &output_workspace);
  1068. }
  1069.  
  1070.  
  1071. /*
  1072.  * Initialize for two-pass color quantization.
  1073.  */
  1074.  
  1075. METHODDEF void
  1076. color_quant_init (decompress_info_ptr cinfo)
  1077. {
  1078.   int i;
  1079.  
  1080.   /* Lower bound on # of colors ... somewhat arbitrary as long as > 0 */
  1081.   if (cinfo->desired_number_of_colors < 8)
  1082.     ERREXIT(cinfo->emethods, "Cannot request less than 8 quantized colors");
  1083.   /* Make sure colormap indexes can be represented by JSAMPLEs */
  1084.   if (cinfo->desired_number_of_colors > MAXNUMCOLORS)
  1085.     ERREXIT1(cinfo->emethods, "Cannot request more than %d quantized colors",
  1086.          MAXNUMCOLORS);
  1087.  
  1088.   /* Allocate and zero the histogram */
  1089.   histogram = (hist3d) (*cinfo->emethods->alloc_small)
  1090.                 (HIST_Y_ELEMS * SIZEOF(hist2d));
  1091.   for (i = 0; i < HIST_Y_ELEMS; i++) {
  1092.     histogram[i] = (hist2d) (*cinfo->emethods->alloc_medium)
  1093.                 (HIST_C_ELEMS*HIST_C_ELEMS * SIZEOF(histcell));
  1094.     jzero_far((void FAR *) histogram[i],
  1095.           HIST_C_ELEMS*HIST_C_ELEMS * SIZEOF(histcell));
  1096.   }
  1097.  
  1098.   /* Allocate storage for the internal and external colormaps. */
  1099.   /* We do this now since it is FAR storage and may affect the memory */
  1100.   /* manager's space calculations. */
  1101.   my_colormap = (*cinfo->emethods->alloc_small_sarray)
  1102.             ((long) cinfo->desired_number_of_colors,
  1103.              (long) 3);
  1104.   cinfo->colormap = (*cinfo->emethods->alloc_small_sarray)
  1105.             ((long) cinfo->desired_number_of_colors,
  1106.              (long) cinfo->color_out_comps);
  1107.  
  1108.   /* Allocate Floyd-Steinberg workspace if necessary */
  1109.   /* This isn't needed until pass 2, but again it is FAR storage. */
  1110.   if (cinfo->use_dithering) {
  1111.     size_t arraysize = (size_t) ((cinfo->image_width + 2L) *
  1112.                  (3 * SIZEOF(FSERROR)));
  1113.  
  1114.     fserrors = (FSERRPTR) (*cinfo->emethods->alloc_medium) (arraysize);
  1115.     /* Initialize the propagated errors to zero. */
  1116.     jzero_far((void FAR *) fserrors, arraysize);
  1117.     on_odd_row = FALSE;
  1118.   }
  1119.  
  1120.   /* Indicate number of passes needed, excluding the prescan pass. */
  1121.   cinfo->total_passes++;    /* I always use one pass */
  1122. }
  1123.  
  1124.  
  1125. /*
  1126.  * Perform two-pass quantization: rescan the image data and output the
  1127.  * converted data via put_color_map and put_pixel_rows.
  1128.  * The source_method is a routine that can scan the image data; it can
  1129.  * be called as many times as desired.  The processing routine called by
  1130.  * source_method has the same interface as color_quantize does in the
  1131.  * one-pass case, except it must call put_pixel_rows itself.  (This allows
  1132.  * me to use multiple passes in which earlier passes don't output anything.)
  1133.  */
  1134.  
  1135. METHODDEF void
  1136. color_quant_doit (decompress_info_ptr cinfo, quantize_caller_ptr source_method)
  1137. {
  1138.   int i;
  1139.  
  1140.   /* Select the representative colors */
  1141.   select_colors(cinfo);
  1142.   /* Pass the external colormap to the output module. */
  1143.   /* NB: the output module may continue to use the colormap until shutdown. */
  1144.   (*cinfo->methods->put_color_map) (cinfo, cinfo->actual_number_of_colors,
  1145.                     cinfo->colormap);
  1146.   /* Re-zero the histogram so pass 2 can use it as nearest-color cache */
  1147.   for (i = 0; i < HIST_Y_ELEMS; i++) {
  1148.     jzero_far((void FAR *) histogram[i],
  1149.           HIST_C_ELEMS*HIST_C_ELEMS * SIZEOF(histcell));
  1150.   }
  1151.   /* Perform pass 2 */
  1152.   if (cinfo->use_dithering)
  1153.     (*source_method) (cinfo, pass2_dither);
  1154.   else
  1155.     (*source_method) (cinfo, pass2_nodither);
  1156. }
  1157.  
  1158.  
  1159. /*
  1160.  * Finish up at the end of the file.
  1161.  */
  1162.  
  1163. METHODDEF void
  1164. color_quant_term (decompress_info_ptr cinfo)
  1165. {
  1166.   /* no work (we let free_all release the histogram/cache and colormaps) */
  1167.   /* Note that we *mustn't* free the external colormap before free_all, */
  1168.   /* since output module may use it! */
  1169. }
  1170.  
  1171.  
  1172. /*
  1173.  * Map some rows of pixels to the output colormapped representation.
  1174.  * Not used in two-pass case.
  1175.  */
  1176.  
  1177. METHODDEF void
  1178. color_quantize (decompress_info_ptr cinfo, int num_rows,
  1179.         JSAMPIMAGE input_data, JSAMPARRAY output_data)
  1180. {
  1181.   ERREXIT(cinfo->emethods, "Should not get here!");
  1182. }
  1183.  
  1184.  
  1185. /*
  1186.  * The method selection routine for 2-pass color quantization.
  1187.  */
  1188.  
  1189. GLOBAL void
  1190. jsel2quantize (decompress_info_ptr cinfo)
  1191. {
  1192.   if (cinfo->two_pass_quantize) {
  1193.      /* Make sure jdmaster didn't give me a case I can't handle */
  1194.      if (cinfo->num_components != 3 || cinfo->jpeg_color_space != CS_YCbCr)
  1195.         ERREXIT(cinfo->emethods, "2-pass quantization only handles YCbCr input");
  1196.      cinfo->methods->color_quant_init = color_quant_init;
  1197.      cinfo->methods->color_quant_prescan = color_quant_prescan;
  1198.      cinfo->methods->color_quant_doit = color_quant_doit;
  1199.     cinfo->methods->color_quant_term = color_quant_term;
  1200.     cinfo->methods->color_quantize = color_quantize;
  1201.     /* Quantized grayscale output is normally done by jquant1.c (which will do
  1202.      * a much better job of it).  But if the program is configured with only
  1203.      * 2-pass quantization, then I have to do the job.  In this situation,
  1204.      * jseldcolor's clearing of the Cb/Cr component_needed flags is incorrect,
  1205.      * because I will look at those components before conversion.
  1206.      */
  1207.     cinfo->cur_comp_info[1]->component_needed = TRUE;
  1208.      cinfo->cur_comp_info[2]->component_needed = TRUE;
  1209.   }
  1210. }
  1211.  
  1212. #endif /* QUANT_2PASS_SUPPORTED */
  1213.