tools/tiffmedian.c

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