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