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