home *** CD-ROM | disk | FTP | other *** search
/ Practical Algorithms for Image Analysis / Practical Algorithms for Image Analysis.iso / TARFILE.GZ / tarfile / libtiff / tools / tiffmedian.c < prev    next >
Encoding:
C/C++ Source or Header  |  1999-09-11  |  22.6 KB  |  888 lines
  1. /* $Header: /usr/people/sam/tiff/tools/RCS/tiffmedian.c,v 1.18 1995/10/10 00:35:22 sam Rel $ */
  2.  
  3. /*
  4.  * Apply median cut on an image.
  5.  *
  6.  * tiffmedian [-c n] [-f] input output
  7.  *     -C n        - set colortable size.  Default is 256.
  8.  *     -f        - use Floyd-Steinberg dithering.
  9.  *     -c lzw        - compress output with LZW
  10.  *     -c none        - use no compression on output
  11.  *     -c packbits    - use packbits compression on output
  12.  *     -r n        - create output with n rows/strip of data
  13.  * (by default the compression scheme and rows/strip are taken
  14.  *  from the input file)
  15.  *
  16.  * Notes:
  17.  *
  18.  * [1] Floyd-Steinberg dither:
  19.  *  I should point out that the actual fractions we used were, assuming
  20.  *  you are at X, moving left to right:
  21.  *
  22.  *            X     7/16
  23.  *         3/16   5/16  1/16    
  24.  *
  25.  *  Note that the error goes to four neighbors, not three.  I think this
  26.  *  will probably do better (at least for black and white) than the
  27.  *  3/8-3/8-1/4 distribution, at the cost of greater processing.  I have
  28.  *  seen the 3/8-3/8-1/4 distribution described as "our" algorithm before,
  29.  *  but I have no idea who the credit really belongs to.
  30.  
  31.  *  Also, I should add that if you do zig-zag scanning (see my immediately
  32.  *  previous message), it is sufficient (but not quite as good) to send
  33.  *  half the error one pixel ahead (e.g. to the right on lines you scan
  34.  *  left to right), and half one pixel straight down.  Again, this is for
  35.  *  black and white;  I've not tried it with color.
  36.  *  -- 
  37.  *                        Lou Steinberg
  38.  *
  39.  * [2] Color Image Quantization for Frame Buffer Display, Paul Heckbert,
  40.  *    Siggraph '82 proceedings, pp. 297-307
  41.  */
  42.  
  43. #include <stdio.h>
  44. #include <stdlib.h>
  45. #include <string.h>
  46.  
  47. #include "tiffio.h"
  48.  
  49. #define    MAX_CMAP_SIZE    256
  50.  
  51. #define    streq(a,b)    (strcmp(a,b) == 0)
  52. #define    strneq(a,b,n)    (strncmp(a,b,n) == 0)
  53.  
  54. #define    COLOR_DEPTH    8
  55. #define    MAX_COLOR    256
  56.  
  57. #define    B_DEPTH        5        /* # bits/pixel to use */
  58. #define    B_LEN        (1L<<B_DEPTH)
  59.  
  60. #define    C_DEPTH        2
  61. #define    C_LEN        (1L<<C_DEPTH)    /* # cells/color to use */
  62.  
  63. #define    COLOR_SHIFT    (COLOR_DEPTH-B_DEPTH)
  64.  
  65. typedef    struct colorbox {
  66.     struct    colorbox *next, *prev;
  67.     int    rmin, rmax;
  68.     int    gmin, gmax;
  69.     int    bmin, bmax;
  70.     int    total;
  71. } Colorbox;
  72.  
  73. typedef struct {
  74.     int    num_ents;
  75.     int    entries[MAX_CMAP_SIZE][2];
  76. } C_cell;
  77.  
  78. uint16    rm[MAX_CMAP_SIZE], gm[MAX_CMAP_SIZE], bm[MAX_CMAP_SIZE];
  79. int    bytes_per_pixel;
  80. int    num_colors;
  81. int    histogram[B_LEN][B_LEN][B_LEN];
  82. Colorbox *freeboxes;
  83. Colorbox *usedboxes;
  84. C_cell    **ColorCells;
  85. TIFF    *in, *out;
  86. uint32    rowsperstrip = (uint32) -1;
  87. uint16    compression = (uint16) -1;
  88. uint16    bitspersample = 1;
  89. uint16    samplesperpixel;
  90. uint32    imagewidth;
  91. uint32    imagelength;
  92. uint16    predictor = 0;
  93.  
  94. static    void get_histogram(TIFF*, Colorbox*);
  95. static    void splitbox(Colorbox*);
  96. static    void shrinkbox(Colorbox*);
  97. static    void map_colortable(void);
  98. static    void quant(TIFF*, TIFF*);
  99. static    void quant_fsdither(TIFF*, TIFF*);
  100. static    Colorbox* largest_box(void);
  101.  
  102. static    void usage(void);
  103. static    int processCompressOptions(char*);
  104.  
  105. #define    CopyField(tag, v) \
  106.     if (TIFFGetField(in, tag, &v)) TIFFSetField(out, tag, v)
  107.  
  108. int
  109. main(int argc, char* argv[])
  110. {
  111.     int i, dither = 0;
  112.     uint16 shortv, config, photometric;
  113.     Colorbox *box_list, *ptr;
  114.     float floatv;
  115.     uint32 longv;
  116.     int c;
  117.     extern int optind;
  118.     extern char* optarg;
  119.  
  120.     num_colors = MAX_CMAP_SIZE;
  121.     while ((c = getopt(argc, argv, "c:C:r:f")) != -1)
  122.         switch (c) {
  123.         case 'c':        /* compression scheme */
  124.             if (!processCompressOptions(optarg))
  125.                 usage();
  126.             break;
  127.         case 'C':        /* set colormap size */
  128.             num_colors = atoi(optarg);
  129.             if (num_colors > MAX_CMAP_SIZE) {
  130.                 fprintf(stderr,
  131.                    "-c: colormap too big, max %d\n",
  132.                    MAX_CMAP_SIZE);
  133.                 usage();
  134.             }
  135.             break;
  136.         case 'f':        /* dither */
  137.             dither = 1;
  138.             break;
  139.         case 'r':        /* rows/strip */
  140.             rowsperstrip = atoi(optarg);
  141.             break;
  142.         case '?':
  143.             usage();
  144.             /*NOTREACHED*/
  145.         }
  146.     if (argc - optind != 2)
  147.         usage();
  148.     in = TIFFOpen(argv[optind], "r");
  149.     if (in == NULL)
  150.         return (-1);
  151.     TIFFGetField(in, TIFFTAG_IMAGEWIDTH, &imagewidth);
  152.     TIFFGetField(in, TIFFTAG_IMAGELENGTH, &imagelength);
  153.     TIFFGetField(in, TIFFTAG_BITSPERSAMPLE, &bitspersample);
  154.     TIFFGetField(in, TIFFTAG_SAMPLESPERPIXEL, &samplesperpixel);
  155.     if (bitspersample != 8 && bitspersample != 16) {
  156.         fprintf(stderr, "%s: Image must have at least 8-bits/sample\n",
  157.             argv[optind]);
  158.         return (-3);
  159.     }
  160.     if (!TIFFGetField(in, TIFFTAG_PHOTOMETRIC, &photometric) ||
  161.         photometric != PHOTOMETRIC_RGB || samplesperpixel < 3) {
  162.         fprintf(stderr, "%s: Image must have RGB data\n", argv[optind]);
  163.         return (-4);
  164.     }
  165.     TIFFGetField(in, TIFFTAG_PLANARCONFIG, &config);
  166.     if (config != PLANARCONFIG_CONTIG) {
  167.         fprintf(stderr, "%s: Can only handle contiguous data packing\n",
  168.             argv[optind]);
  169.         return (-5);
  170.     }
  171.  
  172.     /*
  173.      * STEP 1:  create empty boxes
  174.      */
  175.     usedboxes = NULL;
  176.     box_list = freeboxes = (Colorbox *)_TIFFmalloc(num_colors*sizeof (Colorbox));
  177.     freeboxes[0].next = &freeboxes[1];
  178.     freeboxes[0].prev = NULL;
  179.     for (i = 1; i < num_colors-1; ++i) {
  180.         freeboxes[i].next = &freeboxes[i+1];
  181.         freeboxes[i].prev = &freeboxes[i-1];
  182.     }
  183.     freeboxes[num_colors-1].next = NULL;
  184.     freeboxes[num_colors-1].prev = &freeboxes[num_colors-2];
  185.  
  186.     /*
  187.      * STEP 2: get histogram, initialize first box
  188.      */
  189.     ptr = freeboxes;
  190.     freeboxes = ptr->next;
  191.     if (freeboxes)
  192.         freeboxes->prev = NULL;
  193.     ptr->next = usedboxes;
  194.     usedboxes = ptr;
  195.     if (ptr->next)
  196.         ptr->next->prev = ptr;
  197.     get_histogram(in, ptr);
  198.  
  199.     /*
  200.      * STEP 3: continually subdivide boxes until no more free
  201.      * boxes remain or until all colors assigned.
  202.      */
  203.     while (freeboxes != NULL) {
  204.         ptr = largest_box();
  205.         if (ptr != NULL)
  206.             splitbox(ptr);
  207.         else
  208.             freeboxes = NULL;
  209.     }
  210.  
  211.     /*
  212.      * STEP 4: assign colors to all boxes
  213.      */
  214.     for (i = 0, ptr = usedboxes; ptr != NULL; ++i, ptr = ptr->next) {
  215.         rm[i] = ((ptr->rmin + ptr->rmax) << COLOR_SHIFT) / 2;
  216.         gm[i] = ((ptr->gmin + ptr->gmax) << COLOR_SHIFT) / 2;
  217.         bm[i] = ((ptr->bmin + ptr->bmax) << COLOR_SHIFT) / 2;
  218.     }
  219.  
  220.     /* We're done with the boxes now */
  221.     _TIFFfree(box_list);
  222.     freeboxes = usedboxes = NULL;
  223.  
  224.     /*
  225.      * STEP 5: scan histogram and map all values to closest color
  226.      */
  227.     /* 5a: create cell list as described in Heckbert[2] */
  228.     ColorCells = (C_cell **)_TIFFmalloc(C_LEN*C_LEN*C_LEN*sizeof (C_cell*));
  229.     _TIFFmemset(ColorCells, 0, C_LEN*C_LEN*C_LEN*sizeof (C_cell*));
  230.     /* 5b: create mapping from truncated pixel space to color
  231.        table entries */
  232.     map_colortable();
  233.  
  234.     /*
  235.      * STEP 6: scan image, match input values to table entries
  236.      */
  237.     out = TIFFOpen(argv[optind+1], "w");
  238.     if (out == NULL)
  239.         return (-2);
  240.  
  241.     CopyField(TIFFTAG_SUBFILETYPE, longv);
  242.     CopyField(TIFFTAG_IMAGEWIDTH, longv);
  243.     TIFFSetField(out, TIFFTAG_BITSPERSAMPLE, (short)COLOR_DEPTH);
  244.     if (compression != (uint16)-1) {
  245.         TIFFSetField(out, TIFFTAG_COMPRESSION, compression);
  246.         switch (compression) {
  247.         case COMPRESSION_LZW:
  248.         case COMPRESSION_DEFLATE:
  249.             if (predictor != 0)
  250.                 TIFFSetField(out, TIFFTAG_PREDICTOR, predictor);
  251.             break;
  252.         }
  253.     } else
  254.         CopyField(TIFFTAG_COMPRESSION, compression);
  255.     TIFFSetField(out, TIFFTAG_PHOTOMETRIC, (short)PHOTOMETRIC_PALETTE);
  256.     CopyField(TIFFTAG_ORIENTATION, shortv);
  257.     TIFFSetField(out, TIFFTAG_SAMPLESPERPIXEL, (short)1);
  258.     CopyField(TIFFTAG_PLANARCONFIG, shortv);
  259.     TIFFSetField(out, TIFFTAG_ROWSPERSTRIP,
  260.         TIFFDefaultStripSize(out, rowsperstrip));
  261.     CopyField(TIFFTAG_MINSAMPLEVALUE, shortv);
  262.     CopyField(TIFFTAG_MAXSAMPLEVALUE, shortv);
  263.     CopyField(TIFFTAG_RESOLUTIONUNIT, shortv);
  264.     CopyField(TIFFTAG_XRESOLUTION, floatv);
  265.     CopyField(TIFFTAG_YRESOLUTION, floatv);
  266.     CopyField(TIFFTAG_XPOSITION, floatv);
  267.     CopyField(TIFFTAG_YPOSITION, floatv);
  268.  
  269.     if (dither)
  270.         quant_fsdither(in, out);
  271.     else
  272.         quant(in, out);
  273.     /*
  274.      * Scale colormap to TIFF-required 16-bit values.
  275.      */
  276. #define    SCALE(x)    (((x)*((1L<<16)-1))/255)
  277.     for (i = 0; i < MAX_CMAP_SIZE; ++i) {
  278.         rm[i] = SCALE(rm[i]);
  279.         gm[i] = SCALE(gm[i]);
  280.         bm[i] = SCALE(bm[i]);
  281.     }
  282.     TIFFSetField(out, TIFFTAG_COLORMAP, rm, gm, bm);
  283.     (void) TIFFClose(out);
  284.     return (0);
  285. }
  286.  
  287. static int
  288. processCompressOptions(char* opt)
  289. {
  290.     if (streq(opt, "none"))
  291.         compression = COMPRESSION_NONE;
  292.     else if (streq(opt, "packbits"))
  293.         compression = COMPRESSION_PACKBITS;
  294.     else if (strneq(opt, "lzw", 3)) {
  295.         char* cp = strchr(opt, ':');
  296.         if (cp)
  297.             predictor = atoi(cp+1);
  298.         compression = COMPRESSION_LZW;
  299.     } else if (strneq(opt, "zip", 3)) {
  300.         char* cp = strchr(opt, ':');
  301.         if (cp)
  302.             predictor = atoi(cp+1);
  303.         compression = COMPRESSION_DEFLATE;
  304.     } else
  305.         return (0);
  306.     return (1);
  307. }
  308.  
  309. char* stuff[] = {
  310. "usage: tiffmedian [options] input.tif output.tif",
  311. "where options are:",
  312. " -r #        make each strip have no more than # rows",
  313. " -C #        create a colormap with # entries",
  314. " -f        use Floyd-Steinberg dithering",
  315. " -c lzw[:opts]    compress output with Lempel-Ziv & Welch encoding",
  316. " -c zip[:opts]    compress output with deflate encoding",
  317. " -c packbits    compress output with packbits encoding",
  318. " -c none    use no compression algorithm on output",
  319. "",
  320. "LZW and deflate options:",
  321. " #        set predictor value",
  322. "For example, -c lzw:2 to get LZW-encoded data with horizontal differencing",
  323. NULL
  324. };
  325.  
  326. static void
  327. usage(void)
  328. {
  329.     char buf[BUFSIZ];
  330.     int i;
  331.  
  332.     setbuf(stderr, buf);
  333.     for (i = 0; stuff[i] != NULL; i++)
  334.         fprintf(stderr, "%s\n", stuff[i]);
  335.     exit(-1);
  336. }
  337.  
  338. static void
  339. get_histogram(TIFF* in, Colorbox* box)
  340. {
  341.     register unsigned char *inptr;
  342.     register int red, green, blue;
  343.     register uint32 j, i;
  344.     unsigned char *inputline;
  345.  
  346.     inputline = (unsigned char *)_TIFFmalloc(TIFFScanlineSize(in));
  347.     if (inputline == NULL) {
  348.         fprintf(stderr, "No space for scanline buffer\n");
  349.         exit(-1);
  350.     }
  351.     box->rmin = box->gmin = box->bmin = 999;
  352.     box->rmax = box->gmax = box->bmax = -1;
  353.     box->total = imagewidth * imagelength;
  354.  
  355.     { register int *ptr = &histogram[0][0][0];
  356.       for (i = B_LEN*B_LEN*B_LEN; i-- > 0;)
  357.         *ptr++ = 0;
  358.     }
  359.     for (i = 0; i < imagelength; i++) {
  360.         if (TIFFReadScanline(in, inputline, i, 0) <= 0)
  361.             break;
  362.         inptr = inputline;
  363.         for (j = imagewidth; j-- > 0;) {
  364.             red = *inptr++ >> COLOR_SHIFT;
  365.             green = *inptr++ >> COLOR_SHIFT;
  366.             blue = *inptr++ >> COLOR_SHIFT;
  367.             if (red < box->rmin)
  368.                 box->rmin = red;
  369.                 if (red > box->rmax)
  370.                 box->rmax = red;
  371.                 if (green < box->gmin)
  372.                 box->gmin = green;
  373.                 if (green > box->gmax)
  374.                 box->gmax = green;
  375.                 if (blue < box->bmin)
  376.                 box->bmin = blue;
  377.                 if (blue > box->bmax)
  378.                 box->bmax = blue;
  379.                 histogram[red][green][blue]++;
  380.         }
  381.     }
  382.     _TIFFfree(inputline);
  383. }
  384.  
  385. static Colorbox *
  386. largest_box(void)
  387. {
  388.     register Colorbox *p, *b;
  389.     register int size;
  390.  
  391.     b = NULL;
  392.     size = -1;
  393.     for (p = usedboxes; p != NULL; p = p->next)
  394.         if ((p->rmax > p->rmin || p->gmax > p->gmin ||
  395.             p->bmax > p->bmin) &&  p->total > size)
  396.                 size = (b = p)->total;
  397.     return (b);
  398. }
  399.  
  400. static void
  401. splitbox(Colorbox* ptr)
  402. {
  403.     int        hist2[B_LEN];
  404.     int        first, last;
  405.     register Colorbox    *new;
  406.     register int    *iptr, *histp;
  407.     register int    i, j;
  408.     register int    ir,ig,ib;
  409.     register int sum, sum1, sum2;
  410.     enum { RED, GREEN, BLUE } axis;
  411.  
  412.     /*
  413.      * See which axis is the largest, do a histogram along that
  414.      * axis.  Split at median point.  Contract both new boxes to
  415.      * fit points and return
  416.      */
  417.     i = ptr->rmax - ptr->rmin;
  418.     if (i >= ptr->gmax - ptr->gmin  && i >= ptr->bmax - ptr->bmin)
  419.         axis = RED;
  420.     else if (ptr->gmax - ptr->gmin >= ptr->bmax - ptr->bmin)
  421.         axis = GREEN;
  422.     else
  423.         axis = BLUE;
  424.     /* get histogram along longest axis */
  425.     switch (axis) {
  426.     case RED:
  427.         histp = &hist2[ptr->rmin];
  428.             for (ir = ptr->rmin; ir <= ptr->rmax; ++ir) {
  429.             *histp = 0;
  430.             for (ig = ptr->gmin; ig <= ptr->gmax; ++ig) {
  431.                 iptr = &histogram[ir][ig][ptr->bmin];
  432.                 for (ib = ptr->bmin; ib <= ptr->bmax; ++ib)
  433.                     *histp += *iptr++;
  434.             }
  435.             histp++;
  436.             }
  437.             first = ptr->rmin;
  438.         last = ptr->rmax;
  439.             break;
  440.     case GREEN:
  441.             histp = &hist2[ptr->gmin];
  442.             for (ig = ptr->gmin; ig <= ptr->gmax; ++ig) {
  443.             *histp = 0;
  444.             for (ir = ptr->rmin; ir <= ptr->rmax; ++ir) {
  445.                 iptr = &histogram[ir][ig][ptr->bmin];
  446.                 for (ib = ptr->bmin; ib <= ptr->bmax; ++ib)
  447.                     *histp += *iptr++;
  448.             }
  449.             histp++;
  450.             }
  451.             first = ptr->gmin;
  452.         last = ptr->gmax;
  453.             break;
  454.     case BLUE:
  455.             histp = &hist2[ptr->bmin];
  456.             for (ib = ptr->bmin; ib <= ptr->bmax; ++ib) {
  457.             *histp = 0;
  458.             for (ir = ptr->rmin; ir <= ptr->rmax; ++ir) {
  459.                 iptr = &histogram[ir][ptr->gmin][ib];
  460.                 for (ig = ptr->gmin; ig <= ptr->gmax; ++ig) {
  461.                     *histp += *iptr;
  462.                     iptr += B_LEN;
  463.                 }
  464.             }
  465.             histp++;
  466.             }
  467.             first = ptr->bmin;
  468.         last = ptr->bmax;
  469.             break;
  470.     }
  471.     /* find median point */
  472.     sum2 = ptr->total / 2;
  473.     histp = &hist2[first];
  474.     sum = 0;
  475.     for (i = first; i <= last && (sum += *histp++) < sum2; ++i)
  476.         ;
  477.     if (i == first)
  478.         i++;
  479.  
  480.     /* Create new box, re-allocate points */
  481.     new = freeboxes;
  482.     freeboxes = new->next;
  483.     if (freeboxes)
  484.         freeboxes->prev = NULL;
  485.     if (usedboxes)
  486.         usedboxes->prev = new;
  487.     new->next = usedboxes;
  488.     usedboxes = new;
  489.  
  490.     histp = &hist2[first];
  491.     for (sum1 = 0, j = first; j < i; j++)
  492.         sum1 += *histp++;
  493.     for (sum2 = 0, j = i; j <= last; j++)
  494.         sum2 += *histp++;
  495.     new->total = sum1;
  496.     ptr->total = sum2;
  497.  
  498.     new->rmin = ptr->rmin;
  499.     new->rmax = ptr->rmax;
  500.     new->gmin = ptr->gmin;
  501.     new->gmax = ptr->gmax;
  502.     new->bmin = ptr->bmin;
  503.     new->bmax = ptr->bmax;
  504.     switch (axis) {
  505.     case RED:
  506.         new->rmax = i-1;
  507.             ptr->rmin = i;
  508.             break;
  509.     case GREEN:
  510.             new->gmax = i-1;
  511.             ptr->gmin = i;
  512.             break;
  513.     case BLUE:
  514.             new->bmax = i-1;
  515.             ptr->bmin = i;
  516.             break;
  517.     }
  518.     shrinkbox(new);
  519.     shrinkbox(ptr);
  520. }
  521.  
  522. static void
  523. shrinkbox(Colorbox* box)
  524. {
  525.     register int *histp, ir, ig, ib;
  526.  
  527.     if (box->rmax > box->rmin) {
  528.         for (ir = box->rmin; ir <= box->rmax; ++ir)
  529.             for (ig = box->gmin; ig <= box->gmax; ++ig) {
  530.                 histp = &histogram[ir][ig][box->bmin];
  531.                     for (ib = box->bmin; ib <= box->bmax; ++ib)
  532.                     if (*histp++ != 0) {
  533.                         box->rmin = ir;
  534.                         goto have_rmin;
  535.                     }
  536.             }
  537.     have_rmin:
  538.         if (box->rmax > box->rmin)
  539.             for (ir = box->rmax; ir >= box->rmin; --ir)
  540.                 for (ig = box->gmin; ig <= box->gmax; ++ig) {
  541.                     histp = &histogram[ir][ig][box->bmin];
  542.                     ib = box->bmin;
  543.                     for (; ib <= box->bmax; ++ib)
  544.                         if (*histp++ != 0) {
  545.                             box->rmax = ir;
  546.                             goto have_rmax;
  547.                         }
  548.                     }
  549.     }
  550. have_rmax:
  551.     if (box->gmax > box->gmin) {
  552.         for (ig = box->gmin; ig <= box->gmax; ++ig)
  553.             for (ir = box->rmin; ir <= box->rmax; ++ir) {
  554.                 histp = &histogram[ir][ig][box->bmin];
  555.                     for (ib = box->bmin; ib <= box->bmax; ++ib)
  556.                 if (*histp++ != 0) {
  557.                     box->gmin = ig;
  558.                     goto have_gmin;
  559.                 }
  560.             }
  561.     have_gmin:
  562.         if (box->gmax > box->gmin)
  563.             for (ig = box->gmax; ig >= box->gmin; --ig)
  564.                 for (ir = box->rmin; ir <= box->rmax; ++ir) {
  565.                     histp = &histogram[ir][ig][box->bmin];
  566.                     ib = box->bmin;
  567.                     for (; ib <= box->bmax; ++ib)
  568.                         if (*histp++ != 0) {
  569.                             box->gmax = ig;
  570.                             goto have_gmax;
  571.                         }
  572.                     }
  573.     }
  574. have_gmax:
  575.     if (box->bmax > box->bmin) {
  576.         for (ib = box->bmin; ib <= box->bmax; ++ib)
  577.             for (ir = box->rmin; ir <= box->rmax; ++ir) {
  578.                 histp = &histogram[ir][box->gmin][ib];
  579.                     for (ig = box->gmin; ig <= box->gmax; ++ig) {
  580.                     if (*histp != 0) {
  581.                         box->bmin = ib;
  582.                         goto have_bmin;
  583.                     }
  584.                     histp += B_LEN;
  585.                     }
  586.                 }
  587.     have_bmin:
  588.         if (box->bmax > box->bmin)
  589.             for (ib = box->bmax; ib >= box->bmin; --ib)
  590.                 for (ir = box->rmin; ir <= box->rmax; ++ir) {
  591.                     histp = &histogram[ir][box->gmin][ib];
  592.                     ig = box->gmin;
  593.                     for (; ig <= box->gmax; ++ig) {
  594.                         if (*histp != 0) {
  595.                             box->bmax = ib;
  596.                             goto have_bmax;
  597.                         }
  598.                         histp += B_LEN;
  599.                     }
  600.                     }
  601.     }
  602. have_bmax:
  603.     ;
  604. }
  605.  
  606. static C_cell *
  607. create_colorcell(int red, int green, int blue)
  608. {
  609.     register int ir, ig, ib, i;
  610.     register C_cell *ptr;
  611.     int mindist, next_n;
  612.     register int tmp, dist, n;
  613.  
  614.     ir = red >> (COLOR_DEPTH-C_DEPTH);
  615.     ig = green >> (COLOR_DEPTH-C_DEPTH);
  616.     ib = blue >> (COLOR_DEPTH-C_DEPTH);
  617.     ptr = (C_cell *)_TIFFmalloc(sizeof (C_cell));
  618.     *(ColorCells + ir*C_LEN*C_LEN + ig*C_LEN + ib) = ptr;
  619.     ptr->num_ents = 0;
  620.  
  621.     /*
  622.      * Step 1: find all colors inside this cell, while we're at
  623.      *       it, find distance of centermost point to furthest corner
  624.      */
  625.     mindist = 99999999;
  626.     for (i = 0; i < num_colors; ++i) {
  627.         if (rm[i]>>(COLOR_DEPTH-C_DEPTH) != ir  ||
  628.             gm[i]>>(COLOR_DEPTH-C_DEPTH) != ig  ||
  629.             bm[i]>>(COLOR_DEPTH-C_DEPTH) != ib)
  630.             continue;
  631.         ptr->entries[ptr->num_ents][0] = i;
  632.         ptr->entries[ptr->num_ents][1] = 0;
  633.         ++ptr->num_ents;
  634.             tmp = rm[i] - red;
  635.             if (tmp < (MAX_COLOR/C_LEN/2))
  636.             tmp = MAX_COLOR/C_LEN-1 - tmp;
  637.             dist = tmp*tmp;
  638.             tmp = gm[i] - green;
  639.             if (tmp < (MAX_COLOR/C_LEN/2))
  640.             tmp = MAX_COLOR/C_LEN-1 - tmp;
  641.             dist += tmp*tmp;
  642.             tmp = bm[i] - blue;
  643.             if (tmp < (MAX_COLOR/C_LEN/2))
  644.             tmp = MAX_COLOR/C_LEN-1 - tmp;
  645.             dist += tmp*tmp;
  646.             if (dist < mindist)
  647.             mindist = dist;
  648.     }
  649.  
  650.     /*
  651.      * Step 3: find all points within that distance to cell.
  652.      */
  653.     for (i = 0; i < num_colors; ++i) {
  654.         if (rm[i] >> (COLOR_DEPTH-C_DEPTH) == ir  &&
  655.             gm[i] >> (COLOR_DEPTH-C_DEPTH) == ig  &&
  656.             bm[i] >> (COLOR_DEPTH-C_DEPTH) == ib)
  657.             continue;
  658.         dist = 0;
  659.             if ((tmp = red - rm[i]) > 0 ||
  660.             (tmp = rm[i] - (red + MAX_COLOR/C_LEN-1)) > 0 )
  661.             dist += tmp*tmp;
  662.             if ((tmp = green - gm[i]) > 0 ||
  663.             (tmp = gm[i] - (green + MAX_COLOR/C_LEN-1)) > 0 )
  664.             dist += tmp*tmp;
  665.             if ((tmp = blue - bm[i]) > 0 ||
  666.             (tmp = bm[i] - (blue + MAX_COLOR/C_LEN-1)) > 0 )
  667.             dist += tmp*tmp;
  668.             if (dist < mindist) {
  669.             ptr->entries[ptr->num_ents][0] = i;
  670.             ptr->entries[ptr->num_ents][1] = dist;
  671.             ++ptr->num_ents;
  672.             }
  673.     }
  674.  
  675.     /*
  676.      * Sort color cells by distance, use cheap exchange sort
  677.      */
  678.     for (n = ptr->num_ents - 1; n > 0; n = next_n) {
  679.         next_n = 0;
  680.         for (i = 0; i < n; ++i)
  681.             if (ptr->entries[i][1] > ptr->entries[i+1][1]) {
  682.                 tmp = ptr->entries[i][0];
  683.                 ptr->entries[i][0] = ptr->entries[i+1][0];
  684.                 ptr->entries[i+1][0] = tmp;
  685.                 tmp = ptr->entries[i][1];
  686.                 ptr->entries[i][1] = ptr->entries[i+1][1];
  687.                 ptr->entries[i+1][1] = tmp;
  688.                 next_n = i;
  689.                 }
  690.     }
  691.     return (ptr);
  692. }
  693.  
  694. static void
  695. map_colortable(void)
  696. {
  697.     register int *histp = &histogram[0][0][0];
  698.     register C_cell *cell;
  699.     register int j, tmp, d2, dist;
  700.     int ir, ig, ib, i;
  701.  
  702.     for (ir = 0; ir < B_LEN; ++ir)
  703.         for (ig = 0; ig < B_LEN; ++ig)
  704.             for (ib = 0; ib < B_LEN; ++ib, histp++) {
  705.                 if (*histp == 0) {
  706.                     *histp = -1;
  707.                     continue;
  708.                 }
  709.                 cell = *(ColorCells +
  710.                     (((ir>>(B_DEPTH-C_DEPTH)) << C_DEPTH*2) +
  711.                     ((ig>>(B_DEPTH-C_DEPTH)) << C_DEPTH) +
  712.                     (ib>>(B_DEPTH-C_DEPTH))));
  713.                 if (cell == NULL )
  714.                     cell = create_colorcell(
  715.                         ir << COLOR_SHIFT,
  716.                         ig << COLOR_SHIFT,
  717.                         ib << COLOR_SHIFT);
  718.                 dist = 9999999;
  719.                 for (i = 0; i < cell->num_ents &&
  720.                     dist > cell->entries[i][1]; ++i) {
  721.                     j = cell->entries[i][0];
  722.                     d2 = rm[j] - (ir << COLOR_SHIFT);
  723.                     d2 *= d2;
  724.                     tmp = gm[j] - (ig << COLOR_SHIFT);
  725.                     d2 += tmp*tmp;
  726.                     tmp = bm[j] - (ib << COLOR_SHIFT);
  727.                     d2 += tmp*tmp;
  728.                     if (d2 < dist) {
  729.                         dist = d2;
  730.                         *histp = j;
  731.                     }
  732.                 }
  733.             }
  734. }
  735.  
  736. /*
  737.  * straight quantization.  Each pixel is mapped to the colors
  738.  * closest to it.  Color values are rounded to the nearest color
  739.  * table entry.
  740.  */
  741. static void
  742. quant(TIFF* in, TIFF* out)
  743. {
  744.     unsigned char    *outline, *inputline;
  745.     register unsigned char    *outptr, *inptr;
  746.     register uint32 i, j;
  747.     register int red, green, blue;
  748.  
  749.     inputline = (unsigned char *)_TIFFmalloc(TIFFScanlineSize(in));
  750.     outline = (unsigned char *)_TIFFmalloc(imagewidth);
  751.     for (i = 0; i < imagelength; i++) {
  752.         if (TIFFReadScanline(in, inputline, i, 0) <= 0)
  753.             break;
  754.         inptr = inputline;
  755.         outptr = outline;
  756.         for (j = 0; j < imagewidth; j++) {
  757.             red = *inptr++ >> COLOR_SHIFT;
  758.             green = *inptr++ >> COLOR_SHIFT;
  759.             blue = *inptr++ >> COLOR_SHIFT;
  760.             *outptr++ = histogram[red][green][blue];
  761.         }
  762.         if (TIFFWriteScanline(out, outline, i, 0) < 0)
  763.             break;
  764.     }
  765.     _TIFFfree(inputline);
  766.     _TIFFfree(outline);
  767. }
  768.  
  769. #define    SWAP(type,a,b)    { type p; p = a; a = b; b = p; }
  770.  
  771. #define    GetInputLine(tif, row, bad)                \
  772.     if (TIFFReadScanline(tif, inputline, row, 0) <= 0)    \
  773.         bad;                        \
  774.     inptr = inputline;                    \
  775.     nextptr = nextline;                    \
  776.     for (j = 0; j < imagewidth; ++j) {            \
  777.         *nextptr++ = *inptr++;                \
  778.         *nextptr++ = *inptr++;                \
  779.         *nextptr++ = *inptr++;                \
  780.     }
  781. #define    GetComponent(raw, cshift, c)                \
  782.     cshift = raw;                        \
  783.     if (cshift < 0)                        \
  784.         cshift = 0;                    \
  785.     else if (cshift >= MAX_COLOR)                \
  786.         cshift = MAX_COLOR-1;                \
  787.     c = cshift;                        \
  788.     cshift >>= COLOR_SHIFT;
  789.  
  790. static void
  791. quant_fsdither(TIFF* in, TIFF* out)
  792. {
  793.     unsigned char *outline, *inputline, *inptr;
  794.     short *thisline, *nextline;
  795.     register unsigned char    *outptr;
  796.     register short *thisptr, *nextptr;
  797.     register uint32 i, j;
  798.     uint32 imax, jmax;
  799.     int lastline, lastpixel;
  800.  
  801.     imax = imagelength - 1;
  802.     jmax = imagewidth - 1;
  803.     inputline = (unsigned char *)_TIFFmalloc(TIFFScanlineSize(in));
  804.     thisline = (short *)_TIFFmalloc(imagewidth * 3 * sizeof (short));
  805.     nextline = (short *)_TIFFmalloc(imagewidth * 3 * sizeof (short));
  806.     outline = (unsigned char *) _TIFFmalloc(TIFFScanlineSize(out));
  807.  
  808.     GetInputLine(in, 0, goto bad);        /* get first line */
  809.     for (i = 1; i < imagelength; ++i) {
  810.         SWAP(short *, thisline, nextline);
  811.         lastline = (i == imax);
  812.         GetInputLine(in, i, break);
  813.         thisptr = thisline;
  814.         nextptr = nextline;
  815.         outptr = outline;
  816.         for (j = 0; j < imagewidth; ++j) {
  817.             int red, green, blue;
  818.             register int oval, r2, g2, b2;
  819.  
  820.             lastpixel = (j == jmax);
  821.             GetComponent(*thisptr++, r2, red);
  822.             GetComponent(*thisptr++, g2, green);
  823.             GetComponent(*thisptr++, b2, blue);
  824.             oval = histogram[r2][g2][b2];
  825.             if (oval == -1) {
  826.                 int ci;
  827.                 register int cj, tmp, d2, dist;
  828.                 register C_cell    *cell;
  829.  
  830.                 cell = *(ColorCells +
  831.                     (((r2>>(B_DEPTH-C_DEPTH)) << C_DEPTH*2) +
  832.                     ((g2>>(B_DEPTH-C_DEPTH)) << C_DEPTH ) +
  833.                     (b2>>(B_DEPTH-C_DEPTH))));
  834.                 if (cell == NULL)
  835.                     cell = create_colorcell(red,
  836.                         green, blue);
  837.                 dist = 9999999;
  838.                 for (ci = 0; ci < cell->num_ents && dist > cell->entries[ci][1]; ++ci) {
  839.                     cj = cell->entries[ci][0];
  840.                     d2 = (rm[cj] >> COLOR_SHIFT) - r2;
  841.                     d2 *= d2;
  842.                     tmp = (gm[cj] >> COLOR_SHIFT) - g2;
  843.                     d2 += tmp*tmp;
  844.                     tmp = (bm[cj] >> COLOR_SHIFT) - b2;
  845.                     d2 += tmp*tmp;
  846.                     if (d2 < dist) {
  847.                         dist = d2;
  848.                         oval = cj;
  849.                     }
  850.                 }
  851.                 histogram[r2][g2][b2] = oval;
  852.             }
  853.             *outptr++ = oval;
  854.             red -= rm[oval];
  855.             green -= gm[oval];
  856.             blue -= bm[oval];
  857.             if (!lastpixel) {
  858.                 thisptr[0] += blue * 7 / 16;
  859.                 thisptr[1] += green * 7 / 16;
  860.                 thisptr[2] += red * 7 / 16;
  861.             }
  862.             if (!lastline) {
  863.                 if (j != 0) {
  864.                     nextptr[-3] += blue * 3 / 16;
  865.                     nextptr[-2] += green * 3 / 16;
  866.                     nextptr[-1] += red * 3 / 16;
  867.                 }
  868.                 nextptr[0] += blue * 5 / 16;
  869.                 nextptr[1] += green * 5 / 16;
  870.                 nextptr[2] += red * 5 / 16;
  871.                 if (!lastpixel) {
  872.                     nextptr[3] += blue / 16;
  873.                         nextptr[4] += green / 16;
  874.                         nextptr[5] += red / 16;
  875.                 }
  876.                 nextptr += 3;
  877.             }
  878.         }
  879.         if (TIFFWriteScanline(out, outline, i-1, 0) < 0)
  880.             break;
  881.     }
  882. bad:
  883.     _TIFFfree(inputline);
  884.     _TIFFfree(thisline);
  885.     _TIFFfree(nextline);
  886.     _TIFFfree(outline);
  887. }
  888.