[comp.graphics] Variance-based Color Quantization Program

cek@notecnirp.Princeton.EDU (Craig Kolb) (01/13/89)

Due to the number of requests I've received (40+), I've decided to simply
post the code.  Again, "colorvar" currently only runs on the (color) Sun and
Iris 4D workstations.  Porting to different environments should prove to be
relatively easy.

Note that "colorvar" is copyrighted in a GNU/Utah Raster-like fashion.
Feel free to copy and distribute the program as long as the copyright
information remains intact and you don't make any money, directly or
indirectly, in the process.

Several people have asked about the availability of the Utah Raster Toolkit.
I've enclosed "BLURB.UTAH", the "blurb" file from the latest distribution
of the toolkit, which gives such information.

Let me know of any braindamage you may find as well as any useful
modifications/ports you make.  I'd also be interested in hearing how
it compares (perceptually) with other color quantization programs.

Happy quantizing,
Craig Kolb
cek@princeton.edu
craig@lom1.math.yale.edu

#! /bin/sh
# This is a shell archive, meaning:
# 1. Remove everything above the #! /bin/sh line.
# 2. Save the resulting text in a file.
# 3. Execute the file with /bin/sh (not csh) to create:
#	BLURB.UTAH
#	Copyright
#	Makefile
#	README
#	args.c
#	boxes.c
#	color.c
#	colorvar.1
#	display.c
#	main.c
#	setup.c
#	var.h
# This archive created: Thu Jan 12 10:44:44 1989
export PATH; PATH=/bin:/usr/bin:$PATH
echo shar: "extracting 'BLURB.UTAH'" '(6402 characters)'
if test -f 'BLURB.UTAH'
then
	echo shar: "will not over-write existing file 'BLURB.UTAH'"
else
sed 's/^	X//' << \SHAR_EOF > 'BLURB.UTAH'
	X
	X		       THE UTAH RASTER TOOLKIT
	X
	XThe Utah Raster toolkit is a collection of programs and C routines for
	Xdealing with raster images commonly encountered in computer graphics.  It
	Xprovides the following major functions:
	X
	X	* A device and system independent image format for storing images
	X	  and information about them.  Called the RLE format, it uses
	X	  run length encoding to reduce storage space for most images.
	X
	X	* A library of C routines for reading, writing and manipulating
	X	  images stored in the RLE format.
	X
	X	* A collections of programs for manipulating and displaying RLE
	X	  images.
	X
	X
	XThe Format:
	X
	X  The device independent RLE file has two parts, the header, which stores
	X  information about the image (size, position, channel information,
	X  color maps, comments, etc), and the actual image data in a run length
	X  encoded format.  The RLE format often requires about a third of the
	X  available space for most "image synthesis" style images.  If the image
	X  does not compress well, the RLE format stores it as straight pixel data
	X  with little extra overhead.  The format has been developed over the past
	X  five years at Utah.
	X
	XThe Library:
	X
	X  C routines are provided for setting up and reading the image header,
	X  and for reading and writing the image a scanline at a time.  Images can
	X  be read or written using two different methods.  Using the "row" method,
	X  the library performs the RLE encoding and decoding.  With the "raw" method,
	X  scanlines are constructed directly with RLE opcodes.  Additional routines
	X  are available for generating dither matrices (e.g., for display programs
	X  running on devices with less than 24 bits of color).
	X
	XThe Tools:
	X  applymap   - Apply color map values to pixel values.
	X  avg4 	     - Downfilter an image by 1/4, generating a matte channel if one
	X	       didn't previously exist
	X  comp 	     - Digital image compositor.  Provides the operations over, atop,
	X	       in, out, xor, plus, minus and diff on two images.
	X  crop 	     - Crop an image.
	X  dvi2rle    - Convert TeX output into anti-aliased images.
	X  fant 	     - Rotate and/or scale in image by an arbitrary (float) value.
	X  mcut       - Quantize an image from 24 to eight bits using the median cut
	X	       algorithm.
	X  mergechan  - Merge several channels from different files into a single
	X	       RLE file.
	X  pyrmask    - Blend images using Gaussian pyrimids.
	X  repos      - Change the position in the RLE header.
	X  rleClock   - Generate an image of a clock.
	X  rleaddcom  - Add a comment to the RLE file's header.
	X  rlebg      - Generate a solid or variable background.
	X  rlebox     - Find the actual non-background area of an image.
	X  rleflip    - Rotate an image by 90/180 degree increments.
	X  rlehdr     - Dump the contents of the RLE header in human readable form.
	X  rlehisto   - Generate the histogram of an RLE file.
	X  rleldmap   - Load a color map into an RLE file from a variety of sources.
	X  rlemandl   - Generate Mandlebrot sets as RLE files.
	X  rlenoise   - Adds random noise to an image.	
	X  rlepatch   - Overlay several smaller images over a larger one.
	X  rlescale   - Generates gray scale and color scale RLE files.
	X  rlesetbg   - Set the background color stored in the RLE header.
	X  rlesplit   - Split a file containing several images into several files.
	X  rleswap    - Swap, copy or delete channels in an RLE file.
	X  rletops    - Convert an RLE image to PostScript (graylevel).
	X  rlezoom    - Enlarge an image with pixel replication.
	X  smush      - Perform a simple Gaussian filter on an image.
	X  to8 	     - Convert a 24 bit RGB image to an eight bit dithered one.
	X  tobw 	     - Convert 24 bits to 8 bits black and white.
	X  unexp      - Convert an "exponential" image to a displayable one.
	X  unslice    - Quickly assemble an image from several horizontal strips
	X
	X  Format conversion programs are provided for:
	X   	- Simple pixel streams (color & B&W)
	X	- Targa image format
	X	- Cubicomp image format
	X	- PostScript
	X	- MacPaint
	X	- Sun rasterfiles
	X	- Wastatch paint systems
	X
	X  Display programs are provided for:
	X  getap	     - Apollo workstations
	X  getbob     - HP Series 300 ("bobocat") running Windows 9000
	X  getcx3d    - Chromatics CX1500 display
	X  getfb      - BRL "libfb" displays
	X  getgmr     - Grinnell GMR-27 (remember those?)
	X  getX	     - Workstations running the X window system
	X  getX11     - Workstations running X11
	X  getOrion   - Orion displays
	X  getren     - HP 98721 "Rennasance" display
	X  getsun     - Suns running Suntools
	X  getmac     - Macintosh.
	X  getmex     - Iris running Mex
	X  getqcr     - Photograph images with the Matrix QCR-Z camera.
	X  getiris    - Iris in raw 24 bit mode.
	X             - [Note display programs for a particular device are
	X		simple to add]
	X
	X  All the tools are designed to pipe together, so they can be used as 
	X  filters on images much like the standard Unix tools filter text.
	X
	XPlus:
	X
	X  The raster toolkit also includes Unix man pages for the library and
	X  commands, some sample images, and additional documentation.
	X
	XSystem Requirements:
	X
	X  We have successfully ported the Raster Toolkit to a number of Unix
	X  systems, including 4.2/4.3bsd (Vax, Sun, etc), Apollo Domain/IX, HP
	X  Series 3000, SGI Iris, Gould UTX.  Display programs are included for
	X  several devices.  Creating display programs for additional devices is
	X  a straightforward task.
	X
	XDistribution:
	X
	X  For ARPAnet sites, the toolkit may be obtained via anonymous FTP to the
	X  site cs.utah.edu, in the file pub/toolkit-2.0.tar (or, if you cannot FTP
	X  that large a file at once, in pub/toolkit-2.0.tar.1, pub/toolkit-2.0.tar.2
	X  and pub/toolkit-2.0.tar.3).  Sites not on the ARPAnet can obtain the Raster
	X  Toolkit on a 9-track, 1600 bpi tar format tape by sending check or
	X  money order for $200.00, payable to the Department of Computer Science,
	X  to:
	X
	X	Attn: Utah Raster Toolkit, Loretta Looser
	X	Department of Computer Science
	X	University of Utah
	X 	Salt Lake City, UT, 84112
	X
	X  Courtesy Mike Muuss at BRL, the Raster Toolkit is also included as
	X  contributed software in the BRL-CAD distribution.
	X
	X  [Note: because of the size of the distribution, we can not distribute
	X  it via mail or UUCP]
	X
	X  Although the Raster Toolkit software is copyrighted, it may be freely 
	X  re-distributed on a "GNU-like" basis.
	X
	XFor further technical information on the Raster Toolkit, send mail
	Xto:
	X	toolkit-request@cs.utah.edu		(ARPA)
	X	{ihnp4,decvax}!utah-cs!toolkit-request	(UUCP)
	X
	X
SHAR_EOF
if test 6402 -ne "`wc -c < 'BLURB.UTAH'`"
then
	echo shar: "error transmitting 'BLURB.UTAH'" '(should have been 6402 characters)'
fi
fi
echo shar: "extracting 'Copyright'" '(216 characters)'
if test -f 'Copyright'
then
	echo shar: "will not over-write existing file 'Copyright'"
else
sed 's/^	X//' << \SHAR_EOF > 'Copyright'
	XCopyright (C) 1988 Craig Kolb
	XThis program may be freely copied, modified, and distributed provided that
	Xit is given free-of-charge for strictly non-commercial purposes and that
	Xthis copyright notice remains intact.
SHAR_EOF
if test 216 -ne "`wc -c < 'Copyright'`"
then
	echo shar: "error transmitting 'Copyright'" '(should have been 216 characters)'
fi
fi
echo shar: "extracting 'Makefile'" '(773 characters)'
if test -f 'Makefile'
then
	echo shar: "will not over-write existing file 'Makefile'"
else
sed 's/^	X//' << \SHAR_EOF > 'Makefile'
	X#
	X# Makefile for variance-based color quantization program
	X#
	X# Craig Kolb 12/88
	X#
	X# Location of Utah Raster include directory and library.
	XRLEINC = /usr/u/utah/include
	XRLELIB = /usr/u/utah/lib/librle.a
	X#
	X# Libraries.
	X# For Iris:		LIBS = $(RLELIB) -lm -lgl -lmalloc
	X# For Sun:		LIBS = $(RLELIB) -lm -lsuntool -lsunwindow -lpixrect
	X#
	XLIBS = $(RLELIB) -lm -lgl -lmalloc
	X#
	X# C-compiler flags.
	X# For Iris 4D:		CFLAGS = -I$(RLEINC) -DIRIS4D
	X# For Iris 4D/GT:	CFLAGS = -I$(RLEINC) -DIRIS4DGT
	X# For Sun (no fpa):	CFLAGS = -I$(RLEINC) -DSUN -f68881
	X#
	XCFLAGS = -O -I$(RLEINC) -DIRIS4D
	X#
	X# End of configuration section.
	X#
	X
	XOBJS = main.o color.o setup.o boxes.o args.o display.o
	X
	Xcolorvar:	$(OBJS)
	X	cc $(CFLAGS) -o colorvar $(OBJS) $(LIBS)
	X
	X$(OBJS): var.h
	X
	Xclean:
	X	/bin/rm -f *.o core
SHAR_EOF
if test 773 -ne "`wc -c < 'Makefile'`"
then
	echo shar: "error transmitting 'Makefile'" '(should have been 773 characters)'
fi
fi
echo shar: "extracting 'README'" '(2026 characters)'
if test -f 'README'
then
	echo shar: "will not over-write existing file 'README'"
else
sed 's/^	X//' << \SHAR_EOF > 'README'
	X"colorvar" performs variance-based color quantization on a given Utah Raster
	Xformat image and displays the results on one one of several types
	Xof framebuffers. Support in this version is provided for Sun and Iris
	Xworkstations.  See the manual page for command-line syntax and the like.
	X
	XThe program is based upon an algorithm presented in "Variance-based Color
	XImage Quantization For Frame Buffer Display" by Wan, Prusinkiewicz and Wong,
	Xwhich is as yet unpublished.  However, earlier results were published in ACM
	XTransactions on Mathematical Software, Vol. 14, #2 (June, 1988) pp. 153-162.
	X
	XTo compile, edit the Makefile to match the needs of your hardware and to
	Xset directory names.
	X
	XNote that this program is, alas, not public-domain.  Although I have no
	Xdoubt that any person with half a brain and copy of the appropriate paper
	Xcould whip together something similar in a day or so, I have no desire to
	Xsee somebody else make money (laugh) off of this.
	X
	XBe warned that dithering has not been debugged, or even thought about at any
	Xgreat length.  Also, the code has not been speed-optimized.  In addition, it
	Xonly works on 24-bit images.  We never deal with anything but 24-bit images,
	Xso I've never worried about it.
	X
	XPorting Notes:
	XThis program has not been compiled on a strictly System V-based machine.
	XCalls to "bzero" will have to be changed to the appropriate thing.
	X
	XWriting versions of display_image() for different framebuffers should be
	Xrelatively easy.  See the version for the IRIS in display.c for a
	Xnext-to-generic example of how to do this.
	X
	XThe prime target for optimization is the getneighbors()/makenearest()
	Xcombination.   The scheme used here for closest-neighboring mapping is
	Xfar from optimal.
	X
	XPlease let me know if (when) you find bugs/braindamage and if you make
	Xany useful modifications.  I'd also appreciate feedback on the dithering
	Xcode -- I haven't had the time to do this properly.
	X
	XCraig Kolb  1/9/89
	XYale University Math Dept.
	Xcek@princeton.edu or
	Xcraig@lom1.math.yale.edu
	X(203) 432-7053
SHAR_EOF
if test 2026 -ne "`wc -c < 'README'`"
then
	echo shar: "error transmitting 'README'" '(should have been 2026 characters)'
fi
fi
echo shar: "extracting 'args.c'" '(1974 characters)'
if test -f 'args.c'
then
	echo shar: "will not over-write existing file 'args.c'"
else
sed 's/^	X//' << \SHAR_EOF > 'args.c'
	X/*
	X * args.c
	X *
	X * Craig Kolb 10/88
	X *
	X * Copyright (C) 1988 Craig Kolb
	X * This program may be freely copied, modified, and distributed provided that
	X * it is given free-of-charge for strictly non-commercial purposes and that
	X * this copyright notice remains intact.
	X */
	X#include <stdio.h>
	X#include "var.h"
	X
	Xchar *Progname;		/* argv[0] */
	X
	X/*
	X * Parse command-line arguments.
	X */
	Xparse_args(argc, argv)
	Xint argc;
	Xchar **argv;
	X{
	X	extern FILE *fimage;
	X	double atof();
	X
	X	OutColors = DEFAULT_OUTCOLORS;
	X	disp_gamma = DEFAULT_GAMMA;
	X	Verbose = FALSE;
	X
	X	Progname = argv[0];
	X
	X	argv++; argc--;
	X
	X	while(argc) {
	X		if (*argv[0] != '-')
	X			break;
	X		switch(argv[0][1]) {
	X#ifdef SUN
	X			case 'd':
	X				Dither = TRUE;
	X				FastColor = TRUE;
	X				break;
	X#endif
	X			case 'c':
	X				OutColors = atoi(argv[1]);
	X#ifdef SUN
	X				if (OutColors > 256) {
	X					fprintf(stderr,"A maximum of 256 output colors may be used.\n");
	X					exit(1);
	X				}
	X#endif
	X				if (OutColors < 1) {
	X					fprintf(stderr,"%d colors?\n", OutColors);
	X					exit(1);
	X				}
	X				argv++; argc--;
	X				break;
	X			/*
	X			 * Toggle FastColor.
	X			 */
	X			case 'f':
	X				FastColor = FastColor ? FALSE : TRUE;
	X				break;
	X			case 'g':
	X				disp_gamma = atof(argv[1]);
	X				argv++; argc--;
	X				break;
	X			case 'v':
	X				Verbose = TRUE;
	X				break;
	X			default:
	X				usage();
	X				exit(1);
	X				break;
	X		}
	X		argv++; argc--;
	X	}
	X
	X	if (argc == 1) {
	X		filename = argv[0];
	X		fimage = fopen(filename, "r");
	X		if (fimage == (FILE *)NULL) {
	X			fprintf(stderr,"Cannot open %s.\n",argv[0]);
	X			exit(0);
	X		}
	X	} else
	X		fimage = stdin;
	X
	X	disp_gamma = 1. / disp_gamma;
	X}
	X
	Xusage()
	X{
	X#ifdef SUN
	X	fprintf(stderr,"usage: %s [-c colors] [-d] [-f] [-g gamma] [-v] [filename]\n",
	X			Progname);
	X#else
	X	fprintf(stderr,"usage: %s [-c colors] [-f] [-g gamme] [-v] [filename]\n", Progname);
	X#endif
	X	fprintf(stderr,"If no gamma is specified a value of %lf is assumed.\n",
	X			DEFAULT_GAMMA);
	X	fprintf(stderr,"By default, the image is quantized to %d colors.\n",
	X			DEFAULT_OUTCOLORS);
	X}
	X
SHAR_EOF
if test 1974 -ne "`wc -c < 'args.c'`"
then
	echo shar: "error transmitting 'args.c'" '(should have been 1974 characters)'
fi
fi
echo shar: "extracting 'boxes.c'" '(6646 characters)'
if test -f 'boxes.c'
then
	echo shar: "will not over-write existing file 'boxes.c'"
else
sed 's/^	X//' << \SHAR_EOF > 'boxes.c'
	X/*
	X * boxes.c
	X *
	X * Routines to perform box-cutting.
	X *
	X * Craig Kolb 10/88
	X *
	X * Copyright (C) 1988 Craig Kolb
	X * This program may be freely copied, modified, and distributed provided that
	X * it is given free-of-charge for strictly non-commercial purposes and that
	X * this copyright notice remains intact.
	X */
	X#include <math.h>
	X#include "var.h"
	X
	X/*
	X * Initialize first color box and perform box-cutting.
	X */
	Xcompute_boxes()
	X{
	X	int curbox;
	X
	X	/*
	X	 * The initial box encompasses the entire color cube.
	X	 * Note that the span of a color box is [low, high).
	X	 */
	X	Boxes[0].low[RED] = Boxes[0].low[GREEN] = Boxes[0].low[BLUE] = 0;
	X	Boxes[0].high[RED] = Boxes[0].high[GREEN] = Boxes[0].high[BLUE] = 33;
	X	Boxes[0].weight = NPixels;
	X
	X	/*
	X	 * Gather statistics on the initial box.
	X	 */
	X	get_stats(&Boxes[0]);
	X
	X	/*
	X	 * Cut away...
	X	 */
	X	for (curbox = 0; curbox < OutColors -1; curbox++) {
	X		/*
	X		 * Cut the box with the largest variance,
	X		 * creating two new boxes, one of which
	X		 * replaces the box which was cut while the other
	X		 * is placed at the end of the list of color boxes.
	X		 */
	X		if (makecut(&Boxes[GreatestVariance(curbox)],
	X			&Boxes[curbox +1]) == FALSE) {
	X				/*
	X				 * Couldn't cut box with greatest variance,
	X				 * thus we cannot make any more cuts.
	X				 */
	X				OutColors = 1 + curbox;
	X				return;
	X			}
	X	}
	X
	X}
	X
	X/*
	X * Find the box with the greatest variance.  Restrict search
	X * to those boxes numbered 0 - n.
	X */
	XGreatestVariance(n)
	Xint n;
	X{
	X	double max;
	X	int whichbox, i;
	X
	X	max = -1;
	X	for (i = 0; i <= n; i++) {
	X		if (Boxes[i].weightedvar > max) {
	X			max = Boxes[i].weightedvar;
	X			whichbox = i;
	X		}
	X	}
	X	return whichbox;
	X}
	X
	X/*
	X * Compute the projected variances of a given box.
	X */
	Xget_stats(box)
	XBox *box;
	X{
	X	register int b, g, r;
	X	int i, color, myfreq;
	X	double var[3];
	X
	X	if(box->weight == 0) {
	X		box->weightedvar = 0;
	X		return;
	X	}
	X	/*
	X	 * Zero out projected frequency arrays.
	X	 */
	X	bzero(box->freq[0], 33*sizeof(int));
	X	bzero(box->freq[1], 33*sizeof(int));
	X	bzero(box->freq[2], 33*sizeof(int));
	X
	X	/*
	X	 * Compute projected frequencies.
	X	 */
	X	for (r = box->low[RED] ; r < box->high[RED] ; r++) {
	X		for (g = box->low[GREEN] ; g < box->high[GREEN] ; g++) {
	X			for (b = box->low[BLUE] ; b < box->high[BLUE] ; b++) {
	X				if ((myfreq = frequency[r][g][b]) == 0)
	X					continue;
	X				box->freq[RED][r] += myfreq;
	X				box->freq[GREEN][g] += myfreq;
	X				box->freq[BLUE][b] += myfreq;
	X			}
	X		}
	X	}
	X	/*
	X	 * Compute weighted variances.
	X	 */
	X	for (color = 0; color < 3; color++) {
	X		box->mean[color] = var[color] = 0.;
	X		for (i = box->low[color] ; i < box->high[color]; i++) {
	X			box->mean[color] += (double)(i * box->freq[color][i]);
	X			var[color] += (double)(i*i * box->freq[color][i]);
	X		}
	X	}
	X	for (color = 0; color < 3; color++) {
	X		box->mean[color] /= box->weight;
	X		var[color] -= box->mean[color] * box->mean[color]
	X						* box->weight;
	X	}
	X
	X	box->weightedvar = (var[RED]+var[GREEN]+var[BLUE]) / NPixels;
	X}
	X
	X/*
	X * Compute optimal cut point for a given box.
	X */
	Xmakecut(box, newbox)
	XBox *box, *newbox;
	X{
	X	int i;
	X	double totalvar[3];
	X	Box newboxes[3][2];
	X
	X	if (box->weightedvar == 0. || box->weight == 0)
	X		return FALSE;	/* Can't cut this box. */
	X	/*
	X	 * For each axis, find the sum of the weighted variances of the
	X	 * two boxes which would result by cutting along that axis. The
	X	 * actual cut is made along the axis which minimizes this sum.
	X	 */
	X
	X	/*
	X 	 * Perform cut along each axis and sum weighted variances.
	X	 */
	X	for (i = 0; i < 3; i++) {
	X		find_cutpoint(box, i, &newboxes[i][0], &newboxes[i][1]);
	X		totalvar[i] = newboxes[i][0].weightedvar +
	X				newboxes[i][1].weightedvar;
	X	}
	X
	X	/*
	X	 * Find which cut minimizes the sum of the weighted variances.
	X	 * The current box is replaced by one of the two new
	X	 * boxes, while "newbox" is set equal to the other.
	X	 */
	X	if (totalvar[RED] <= totalvar[GREEN] &&
	X	    totalvar[RED] <= totalvar[BLUE]) {
	X		*box = newboxes[RED][0];
	X		*newbox = newboxes[RED][1];
	X	}
	X	else if (totalvar[GREEN] <= totalvar[RED] &&
	X		 totalvar[GREEN] <= totalvar[BLUE]) {
	X		*box = newboxes[GREEN][0];
	X		*newbox = newboxes[GREEN][1];
	X	} else {
	X		*box = newboxes[BLUE][0];
	X		*newbox = newboxes[BLUE][1];
	X	}
	X
	X	return TRUE;	/* Cut successfully. */
	X}
	X
	X/*
	X * Find "optimal" cutpoint along an axis of the given box.
	X * Create and calculate variances of the two new boxes which
	X * result from this cut.  "color" indicates the axis along which to cut.
	X */
	Xfind_cutpoint(box, color, newbox1, newbox2)
	XBox *box, *newbox1, *newbox2;
	Xint color;
	X{
	X	double u, v, max;
	X	int i, maxindex, minindex, optweight, curweight, cutpoint;
	X
	X	if (box->low[color] + 1 == box->high[color]) {
	X		/*
	X		 * Axis has zero length, so we can't cut it.
	X		 */
	X		newbox1->weightedvar = newbox2->weightedvar = HUGE;
	X		return;
	X	}
	X	/*
	X	 * minindex and maxindex represent the range in which we
	X	 * must search to find the optimal cutpoint.
	X	 */
	X	minindex = (int)((box->low[color] + box->mean[color]) / 2.);
	X	maxindex = (int)((box->mean[color] + box->high[color]) / 2.);
	X
	X	/*
	X	 * If the full weight of the box is contained in
	X	 * box->freq[color][minindex], then we know we should cut
	X	 * at minindex.  Thus, we initialize the cutpoint to be there.
	X	 */
	X	cutpoint = minindex;
	X	optweight = box->weight;
	X
	X	curweight = 0.;
	X	/*
	X	 * Sum that (projected) weight of the box that occurs at
	X	 * points less than minindex.
	X	 */
	X	for (i = box->low[color] ; i < minindex ; i++)
	X		curweight += box->freq[color][i];
	X
	X	/*
	X	 * The mess below is (16) from Wan, Prusinkiewicz and Wong.
	X	 * curweight is the weight of the first interval
	X	 * (box->weight - curweight) is the weight of the second
	X	 * u is the mean of the projected distribution of the
	X	 * 	first interval
	X	 */
	X	u = 0.;
	X	max = -1;
	X	for (i = minindex; i <= maxindex ; i++) {
	X		curweight += box->freq[color][i];
	X		if (curweight == box->weight)
	X			break;
	X		u += (double)i * (double)box->freq[color][i] /
	X					(double)box->weight;
	X		v = ((double)curweight / (double)(box->weight-curweight)) *
	X				(box->mean[color]-u)*(box->mean[color]-u);
	X		/*
	X		 * We cut at the point where v is maximized.
	X		 */
	X		if (v > max) {
	X			max = v;
	X			cutpoint = i;
	X			optweight = curweight;
	X		}
	X	}
	X	/*
	X	 * The cutpoint should represent the place where
	X	 * the second box begins, so we must increment the cutpoint.  
	X	 */
	X	cutpoint++;
	X	/*
	X	 * The two new boxes are nearly identical to the old...
	X	 */
	X	*newbox1 = *newbox2 = *box;
	X	/*
	X	 * newbox1 is the box which falls "to the left" of
	X	 * the cutpoint (the lower interval).
	X	 */
	X	newbox1->weight = optweight;
	X	newbox2->weight -= optweight;
	X	newbox1->high[color] = cutpoint;
	X	newbox2->low[color] = cutpoint;
	X	/*
	X	 * Compute variances of the new boxes.
	X	 */
	X	get_stats(newbox1);
	X	get_stats(newbox2);
	X}
SHAR_EOF
if test 6646 -ne "`wc -c < 'boxes.c'`"
then
	echo shar: "error transmitting 'boxes.c'" '(should have been 6646 characters)'
fi
fi
echo shar: "extracting 'color.c'" '(6760 characters)'
if test -f 'color.c'
then
	echo shar: "will not over-write existing file 'color.c'"
else
sed 's/^	X//' << \SHAR_EOF > 'color.c'
	X/*
	X * color.c
	X *
	X * Routines to compute colormap and NearestColor arrays given an
	X * array of color boxes.
	X *
	X * Craig Kolb 10/88
	X *
	X * Copyright (C) 1988 Craig Kolb
	X * This program may be freely copied, modified, and distributed provided that
	X * it is given free-of-charge for strictly non-commercial purposes and that
	X * this copyright notice remains intact.
	X */
	X#include <math.h>
	X#include <stdio.h>
	X#include "var.h"
	X
	X/*
	X * Form colormap and NearestColor arrays.
	X */
	Xfind_colors()
	X{
	X	register int i;
	X	int num, *neighbors;
	X	Box *box;
	X
	X	neighbors = (int *)malloc(OutColors * sizeof(int));
	X
	X	/*
	X	 * Form color map.  Centroid of ith box is the ith entry.
	X	 */
	X	for (i = 0, box = Boxes; i < OutColors; i++, box++) {
	X		colormap[i][RED] = box->mean[RED] + 0.5;
	X		colormap[i][GREEN] = box->mean[GREEN] + 0.5;
	X		colormap[i][BLUE] = box->mean[BLUE] + 0.5;
	X	}
	X	/*
	X	 * Form map of representative (nearest) colors.
	X	 */
	X	for (i = 0; i < OutColors; i++) {
	X		if (FastColor)
	X			makefastnearest(i);
	X		else {
	X			/*
	X			 * Create list of candidate neighbors and
	X			 * find closest representative for each
	X			 * color in the box.
	X			 */
	X			num = getneighbors(i, neighbors);
	X			makenearest(i, num, neighbors);
	X		}
	X	}
	X	free((char *)neighbors);
	X}
	X
	X/*
	X * In order to minimize our search for 'best representative', we form the
	X * 'neighbors' array.  This array contains the number of the boxes whose
	X * centroids *might* be used as a representative for some color in the
	X * current box.  We need only consider those boxes whose centroids are closer
	X * to one or more of the current box's corners than is the centroid of the
	X * current box. 'Closeness' is measured by Euclidean distance.
	X */
	Xgetneighbors(box, neighbors)
	Xint box;
	Xint *neighbors;
	X{
	X	register int i, j;
	X	int dist, rdist, gdist, bdist, maxdist;
	X	Box *bp;
	X	double realdist, LowR, LowG, LowB, HighR, HighG, HighB;
	X
	X	bp = &Boxes[box];
	X
	X	/*
	X	 * Find the maximum distance between the centroid (representative
	X	 * color) of the box and any point in the box.  Such a point
	X	 * must be one of the corners, so we simply check those.
	X	 */
	X	rdist = bp->low[RED] - colormap[box][RED];
	X	gdist = bp->low[GREEN] - colormap[box][GREEN];
	X	bdist = bp->low[BLUE] - colormap[box][BLUE];
	X	maxdist = rdist*rdist + gdist+gdist + bdist*bdist;
	X
	X	rdist = bp->high[RED] - colormap[box][RED];
	X	gdist = bp->low[GREEN] - colormap[box][GREEN];
	X	bdist = bp->low[BLUE] - colormap[box][BLUE];
	X	dist = rdist*rdist + gdist+gdist + bdist*bdist;
	X	if (dist > maxdist)
	X		maxdist = dist;
	X
	X	rdist = bp->low[RED] - colormap[box][RED];
	X	gdist = bp->high[GREEN] - colormap[box][GREEN];
	X	bdist = bp->low[BLUE] - colormap[box][BLUE];
	X	dist = rdist*rdist + gdist+gdist + bdist*bdist;
	X	if (dist > maxdist)
	X		maxdist = dist;
	X
	X	rdist = bp->low[RED] - colormap[box][RED];
	X	gdist = bp->low[GREEN] - colormap[box][GREEN];
	X	bdist = bp->high[BLUE] - colormap[box][BLUE];
	X	dist = rdist*rdist + gdist+gdist + bdist*bdist;
	X	if (dist > maxdist)
	X		maxdist = dist;
	X
	X	rdist = bp->high[RED] - colormap[box][RED];
	X	gdist = bp->high[GREEN] - colormap[box][GREEN];
	X	bdist = bp->low[BLUE] - colormap[box][BLUE];
	X	dist = rdist*rdist + gdist+gdist + bdist*bdist;
	X	if (dist > maxdist)
	X		maxdist = dist;
	X
	X	rdist = bp->high[RED] - colormap[box][RED];
	X	gdist = bp->low[GREEN] - colormap[box][GREEN];
	X	bdist = bp->high[BLUE] - colormap[box][BLUE];
	X	dist = rdist*rdist + gdist+gdist + bdist*bdist;
	X	if (dist > maxdist)
	X		maxdist = dist;
	X
	X	rdist = bp->low[RED] - colormap[box][RED];
	X	gdist = bp->high[GREEN] - colormap[box][GREEN];
	X	bdist = bp->high[BLUE] - colormap[box][BLUE];
	X	dist = rdist*rdist + gdist+gdist + bdist*bdist;
	X	if (dist > maxdist)
	X		maxdist = dist;
	X
	X	rdist = bp->high[RED] - colormap[box][RED];
	X	gdist = bp->high[GREEN] - colormap[box][GREEN];
	X	bdist = bp->high[BLUE] - colormap[box][BLUE];
	X	dist = rdist*rdist + gdist+gdist + bdist*bdist;
	X	if (dist > maxdist)
	X		maxdist = dist;
	X
	X	realdist = sqrt((double)maxdist);
	X
	X	/*
	X	 * Loop over all colors in the colormap, the ith entry of which
	X	 * corresponds to the ith box.
	X	 *
	X	 * If the centroid of a box is as close to any corner of the
	X	 * current box as is the centroid of the current box, add that
	X	 * box to the list of "neighbors" of the current box.
	X	 */
	X	HighR = (double)bp->high[RED] + realdist;
	X	HighG = (double)bp->high[GREEN] + realdist;
	X	HighB = (double)bp->high[BLUE] + realdist;
	X	LowR = (double)bp->low[RED] - realdist;
	X	LowG = (double)bp->low[GREEN] - realdist;
	X	LowB = (double)bp->low[BLUE] - realdist;
	X	for (i = j = 0, bp = Boxes; i < OutColors; i++, bp++) {
	X		if (LowR <= bp->mean[RED] && HighR >= bp->mean[RED] &&
	X		    LowG <= bp->mean[GREEN] && HighG >= bp->mean[GREEN] &&
	X		    LowB <= bp->mean[BLUE] && HighB >= bp->mean[BLUE])
	X			neighbors[j++] = i;
	X	}
	X
	X	return j;	/* Return the number of neighbors found. */
	X
	X}
	X
	X/*
	X * Assign representative colors to every pixel in a given box through
	X * the construction of the NearestColor array.  For each color in the
	X * given box, we look at the list of neighbors passed to find the
	X * one whose centroid is closest to the current color.
	X */
	Xmakenearest(boxnum, nneighbors, neighbors)
	Xint boxnum;
	Xint nneighbors, *neighbors;
	X{
	X	register int n, b, g, r;
	X	int rdist, gdist, bdist, dist, mindist, which;
	X	Box *box;
	X
	X	box = &Boxes[boxnum];
	X
	X	for (r = box->low[RED]; r < box->high[RED]; r++) {
	X		for (g = box->low[GREEN]; g < box->high[GREEN]; g++) {
	X			for (b = box->low[BLUE]; b < box->high[BLUE]; b++) {
	X#ifdef SUN
	X				/*
	X				 * If we're dithering, we have to have
	X				 * a representative for every possible color.
	X				 */
	X				if (!Dither && frequency[r][g][b] == 0)
	X					continue;
	X#else
	X				if (frequency[r][g][b] == 0)
	X					continue;
	X#endif
	X				mindist = HUGE;
	X				/*
	X				 * Find the colormap entry which is
	X				 * closest to the current color.
	X				 */
	X				for (n = 0; n < nneighbors; n++) {
	X					rdist = r - (int)colormap[neighbors[n]][RED];
	X					gdist = g - (int)colormap[neighbors[n]][GREEN];
	X					bdist = b - (int)colormap[neighbors[n]][BLUE];
	X					dist = rdist*rdist + gdist*gdist + bdist*bdist;
	X					if (dist < mindist) {
	X						mindist = dist;
	X						which = neighbors[n];
	X					}
	X				}
	X				/*
	X				 * The colormap entry closest to this
	X				 * color is used as a representative.
	X				 */
	X				NearestColor[r][g][b] = which;
	X			}
	X		}
	X	}
	X}
	X
	X/*
	X * Assign representative colors to every pixel in a given box through
	X * the construction of the NearestColor array.  Here, we simply assign the
	X * centroid of a box to be the representative for every point in the box.
	X */
	Xmakefastnearest(boxnum)
	Xint boxnum;
	X{
	X	register int b, g, r;
	X	Box *box;
	X
	X	box = &Boxes[boxnum];
	X
	X	for (r = box->low[RED]; r < box->high[RED]; r++) {
	X		for (g = box->low[GREEN]; g < box->high[GREEN]; g++) {
	X			for (b = box->low[BLUE]; b < box->high[BLUE]; b++)
	X				NearestColor[r][g][b] = boxnum;
	X		}
	X	}
	X}
SHAR_EOF
if test 6760 -ne "`wc -c < 'color.c'`"
then
	echo shar: "error transmitting 'color.c'" '(should have been 6760 characters)'
fi
fi
echo shar: "extracting 'colorvar.1'" '(2668 characters)'
if test -f 'colorvar.1'
then
	echo shar: "will not over-write existing file 'colorvar.1'"
else
sed 's/^	X//' << \SHAR_EOF > 'colorvar.1'
	X.\" Manual page for colorvar, -man format.
	X.\" Initial version December 22, 1988 Craig Kolb
	X.\"
	X.TH COLORVAR 1G "December 22, 1988"
	X.UC 4
	X.SH NAME
	Xcolorvar \- a variance-based color quantization program
	X.SH
	X.B colorvar
	X[
	X.B options
	X] [
	X.B filename
	X]
	X.SH DESCRIPTION
	X.I Colorvar
	Xreads a Utah Raster format image file and displays a quantized color
	Ximage on one of several types of framebuffers.  If no filename is
	Xgiven the image is read from the standard input.
	X.I Colorvar
	Xquantizes colors in the original 24-bit image by constructing a map
	Xof representative colors
	Xand then mapping each pixel in the input image
	Xto an entry in the colormap.  The program uses a variance-based quantization
	Xalgorithm to construct the map of representative colors.
	X.SH OPTIONS
	X.TP
	X\fB\-c\fI output_colors\fR
	XSpecifies the number of colors in the final image.  The default is 256.
	X.TP
	X.B \-d
	X(Sun workstation only.)  Perform dithering on the final image.  This option
	Ximplies
	X.I \-f.
	X.TP
	X.B \-f
	XToggles fast full-color to colormap entry mapping.  Rather than finding the
	Xclosest representative for a given color, this option forces the use of
	Xthe centroid of the color box in which a color is contained
	Xas the representative.  The resulting error is usually tolerable.
	XNote that
	X.I "\-d \-f"
	Xwill cause both dithering and the more precise method of
	Xrepresentative-finding to be used.
	X.TP
	X\fB\-g\fI gamma\fR
	XSpecifies the gamma value of the display monitor.  The default is 1.6 for
	Xthe Sun workstation and 1.0 for the Iris workstation.
	X.TP
	X.B \-v
	XThis "verbose" flag causes the total number of distinct 15-bit colors to be
	Xprinted to the standard error before constructing the colormap.  (\fIColorvar\fR
	X"prequantizes" 24-bit images to 15 bits in order to speed up the algorithm and
	Xreduce storage space.)
	X.SH CAVEATS
	X.PP
	XThe dithering option does not work properly.  Or at least not to my
	Xsatisfaction.
	X.PP
	X.I Colorvar
	Xcurrently only works with 24-bit (8 each of red, green and blue) images.
	X.PP
	XA maximum of 256 colormap entries can be used on the Sun workstation.
	X.PP
	XThere is currently no option to change where in the Iris' lookup table
	Xthe colormap is written.  Running more than one
	X.I colorvar
	Xat a time on the Iris will almost surely not work properly, as each
	Xinvocation will clobber the colormap of the last invocation.
	X.SH REFERENCES
	X.PP
	XWan, Wong and Prusinkiewicz, \fIAn Algorithm for Multidimensional Data
	XClustering\fR.  ACM Transactions on Mathematical Software, Vol. 14,
	XNo. 2, June 1988. pp. 153-162
	X.PP
	XWan, Prusinkiewicz and Wong, \fIVariance-based Color Image Quantization For
	XFrame Buffer Display\fR
	X.SH AUTHOR
	XCraig Kolb (cek@princeton.edu)
	X.SH "SEE ALSO"
	Xrle(5)
SHAR_EOF
if test 2668 -ne "`wc -c < 'colorvar.1'`"
then
	echo shar: "error transmitting 'colorvar.1'" '(should have been 2668 characters)'
fi
fi
echo shar: "extracting 'display.c'" '(7187 characters)'
if test -f 'display.c'
then
	echo shar: "will not over-write existing file 'display.c'"
else
sed 's/^	X//' << \SHAR_EOF > 'display.c'
	X/*
	X * display.c
	X *
	X * Routines to display quantized image on Sun and Iris workstations.
	X *
	X * Craig Kolb 11/88
	X * Copyright (C) 1988 Craig Kolb
	X * This program may be freely copied, modified, and distributed provided that
	X * it is given free-of-charge for strictly non-commercial purposes and that
	X * this copyright notice remains intact.
	X */
	X#include <stdio.h>
	X#include <math.h>
	X#include "var.h"
	X#include "svfb_global.h"
	X
	Xextern FILE *fimage;
	X
	Xunsigned char Gamma();
	X
	X#ifdef SUN
	X
	X#include <suntool/sunview.h>
	X#include <suntool/canvas.h>
	X
	X#define COLORNAME	"var_colormap"	/* Name of new colormap. */
	X#define FACT		(255./32.)
	X
	X/*
	X * Display image on a color Sun workstation.
	X */
	Xdisplay_image()
	X{
	X	unsigned char *outscan, redmap[256], greenmap[256], bluemap[256];
	X	float scaledred[256], scaledgreen[256], scaledblue[256];
	X	short *curr, *curg, *curb, *nextr, *nextg, *nextb;
	X	int x, y, r, g, b, i;
	X	float RErr, GErr, BErr;
	X	char *cptmp, *rindex(), colorname[256];
	X	Pixwin *pw;
	X	Canvas canvas;
	X	Frame frame;
	X	struct pixrect *pm;
	X
	X	/*
	X	 * Find basename of image file, use that & resolution to title window.
	X	 */
	X	if ((cptmp = rindex(filename, '/')) != (char *)NULL)
	X		filename = cptmp + 1;
	X
	X	sprintf(filename, "%s (%d x %d) %d colors", filename, Xres, Yres, OutColors);
	X
	X	/*
	X	 * Create a new frame and canvas of the proper size.
	X	 */
	X	frame = window_create(NULL, FRAME, 
	X				FRAME_SHOW_LABEL, TRUE,
	X				FRAME_LABEL, filename,
	X				FRAME_NO_CONFIRM, TRUE, 0);
	X	canvas = window_create(frame, CANVAS,
	X				WIN_HEIGHT, Yres,
	X				WIN_WIDTH, Xres, 0);
	X	window_fit(frame);
	X
	X	/*
	X	 * Perform gamma correction on the color map.
	X	 * In addition, precompute 'scaled' color values which
	X	 * are used often when dithering.
	X	 */
	X	for (x = 0; x < OutColors; x++) {
	X		redmap[x] = Gamma((int)colormap[x][0]);
	X		greenmap[x] = Gamma((int)colormap[x][1]);
	X		bluemap[x] = Gamma((int)colormap[x][2]);
	X		if (Dither) {
	X			scaledred[x] = FACT * colormap[x][0];
	X			scaledgreen[x] = FACT * colormap[x][1];
	X			scaledblue[x] = FACT * colormap[x][2];
	X		}
	X	}
	X	/*
	X	 * Zero out end of color map.
	X	 */
	X	for (x = OutColors ; x < 256; x++)
	X		redmap[x] = greenmap[x] = bluemap[x] = 0;
	X
	X	/*
	X	 * Give canvas pixwin new colormap and colormap name.
	X	 */
	X	pw = canvas_pixwin(canvas);
	X	sprintf(colorname, "%s%d",COLORNAME, getpid());
	X	pw_setcmsname(pw, colorname);
	X	pw_putcolormap(pw, 0, 256, redmap, greenmap, bluemap);
	X	outscan = (unsigned char *)malloc(Xres * sizeof(unsigned char));
	X	/*
	X	 * Replace each RGB value in the image with its representative
	X	 * color.
	X	 * The 'cur' arrays contain the red, green, and blue values of
	X	 * the current scanline of the original image.
	X	 * The 'next' arrays contain the R, G and B values of the next
	X	 * scanline of the original image.
	X	 */
	X	nextr = Image[0][0];
	X	nextg = Image[0][1];
	X	nextb = Image[0][2];
	X	for (y = 0; y < Yres -1; y++) {
	X		curr = nextr;
	X		curg = nextg;
	X		curb = nextb;
	X		nextr = Image[y+1][RED];
	X		nextg = Image[y+1][GREEN];
	X		nextb = Image[y+1][BLUE];
	X		for (x = 0; x < Xres; x++) {
	X			if (!Dither) {
	X				/*
	X				 * Image has already been quantized to
	X				 * [0,32].
	X				 */
	X				outscan[x] = NearestColor[curr[x]][curg[x]][curb[x]];
	X				continue;
	X			}
	X			/*
	X			 * Quantize to [0,32].
	X			 */
	X			r = Quantize(curr[x]);
	X			g = Quantize(curg[x]);
	X			b = Quantize(curb[x]); 
	X			if (r > 32) r = 32;
	X			else if (r < 0) r = 0;
	X			if (g > 32) g = 32;
	X			else if (g < 0) g = 0;
	X			if (b > 32) b = 32;
	X			else if (b < 0) b = 0;
	X
	X			/*
	X			 * Find representative.
	X			 */
	X			outscan[x] = NearestColor[r][g][b];
	X
	X			i = (int)outscan[x];
	X
	X			/*
	X			 * Compute red, green and blue error.
	X			 */
	X			RErr = (float)curr[x] - scaledred[i];
	X			GErr = (float)curg[x] - scaledgreen[i];
	X			BErr = (float)curb[x] - scaledblue[i];
	X
	X			/*
	X			 * Add 3/8 error upwards and to the right,
	X			 * 1/4 error diagonally upwards.
	X			 */
	X			if (RErr) {
	X				curr[x+1] += (short)(.375 * RErr);
	X				nextr[x] += (short)(.375 * RErr);
	X				nextr[x+1] += (short)(.25 * RErr);
	X			}
	X			if (GErr) {
	X				curg[x+1] += (short)(.375 * GErr);
	X				nextg[x] += (short)(.375 * GErr);
	X				nextg[x+1] += (short)(.25 * GErr);
	X			}
	X			if (BErr) {
	X				curb[x+1] += (short)(.375 * BErr);
	X				nextb[x] += (short)(.375 * BErr);
	X				nextb[x+1] += (short)(.25 * BErr);
	X			}
	X		}
	X		/*
	X		 * Write scanline to the window.
	X		 */
	X		pm = mem_point(Xres, 1, 8, outscan);
	X		pw_rop(pw, 0, Yres-y-1, Xres, 1, PIX_SRC, pm, 0, 0);
	X		free((char *)pm);
	X		free((char *)curr);
	X		free((char *)curg);
	X		free((char *)curb);
	X	}
	X	/*
	X	 * Output final scanline.
	X	 */
	X	for (x = 0; x < Xres; x++) {
	X		if (!Dither)
	X			outscan[x] = NearestColor[nextr[x]][nextg[x]][nextb[x]];
	X		else {
	X			r = Quantize(nextr[x]);
	X			g = Quantize(nextg[x]);
	X			b = Quantize(nextb[x]);
	X			if (r > 32) r = 32;
	X			else if (r < 0) r = 0;
	X			if (g > 32) g = 32;
	X			else if (g < 0) g = 0;
	X			if (b > 32) b = 32;
	X			else if (b < 0) b = 0;
	X			outscan[x] = NearestColor[r][g][b];
	X		}
	X	}
	X	pm = mem_point(Xres, 1, 8, outscan);
	X	pw_rop(pw, 0, 0, Xres, 1, PIX_SRC, pm, 0, 0);
	X	free((char *)pm);
	X	free((char *)nextr);
	X	free((char *)nextg);
	X	free((char *)nextb);
	X	window_main_loop(frame);
	X}
	X
	X#endif
	X
	X#ifdef IRIS4D || ifdef IRIS4DGT
	X
	X#include <device.h>
	X
	X/*
	X * MAP_OFFSET is the location in device colormap to begin writing
	X * image colormap.
	X * Note that this means that you can only have one copy of this
	X * program running at a time and have correct results, as each
	X * time the program is run the colormap is overwritten.
	X */
	X#define MAP_OFFSET	256
	X
	Xunsigned short *scanline;			/* Output scanline */
	X
	X/*
	X * Display quantized image on an IRIS4D
	X */
	Xdisplay_image()
	X{
	X	int x;
	X	char *cptmp, *rindex();
	X	short val;
	X
	X	/*
	X	 * Find basename of image file.
	X	 */
	X	if ((cptmp = rindex(filename, '/')) != (char *)NULL)
	X		filename = cptmp + 1;
	X
	X	prefsize(Xres, Yres);
	X	winopen(filename);
	X	sprintf(filename, "%s (%d x %d) %d colors", filename, Xres, Yres, OutColors);
	X	wintitle(filename);
	X	qdevice(REDRAW);
	X	qdevice(WINQUIT);
	X	unqdevice(INPUTCHANGE);
	X	qenter(REDRAW);
	X	color(0);
	X	clear();
	X	/*
	X	 * Load gamma-corrected colormap.
	X	 */
	X	for (x = 0; x < OutColors; x++)
	X		mapcolor(x+MAP_OFFSET, 	(short)Gamma((int)colormap[x][RED]), 
	X				(short)Gamma((int)colormap[x][GREEN]),
	X				(short)Gamma((int)colormap[x][BLUE]));
	X
	X	scanline = (unsigned short *)malloc(Xres * sizeof(unsigned short));
	X	while (1) {
	X		switch (qread(&val)) {
	X			case REDRAW:
	X				reshapeviewport();
	X				loadimage();
	X				break;
	X			case WINQUIT:
	X				exit(0);
	X				break;
	X		}
	X	}
	X}
	X
	X/*
	X * Draw image on screen.
	X */
	Xloadimage()
	X{
	X	int x, y, r, g, b;
	X
	X	for (y = 0; y < Yres ; y++) {
	X		for (x = 0; x < Xres; x ++) {
	X			r = (int)Image[y][RED][x];
	X			g = (int)Image[y][GREEN][x];
	X			b = (int)Image[y][BLUE][x];
	X			/*
	X			 * Use current color's representative.
	X			 */
	X			scanline[x] = (unsigned short)NearestColor[r][g][b] + MAP_OFFSET;
	X		}
	X#ifdef IRIS4DGT
	X		rectwrite(0, y, Xres, y, scanline);
	X#else
	X		cmov2i(0, y);
	X		writepixels(Xres, scanline);
	X#endif
	X	}
	X}
	X#endif
	X
	X/*
	X * Perform gamma correction on values ranging from 0 to 32, returning
	X * values between 0 and 255.
	X */
	Xunsigned char
	XGamma(val)
	Xint val;
	X{
	X	val = (int)(0.5 + 255. * pow((double)val/ 32., disp_gamma));
	X
	X	if (val < 0)
	X		val = 0;
	X	else if (val > 255)
	X		val = 255;
	X	return (unsigned char)val;
	X}
	X
SHAR_EOF
if test 7187 -ne "`wc -c < 'display.c'`"
then
	echo shar: "error transmitting 'display.c'" '(should have been 7187 characters)'
fi
fi
echo shar: "extracting 'main.c'" '(1525 characters)'
if test -f 'main.c'
then
	echo shar: "will not over-write existing file 'main.c'"
else
sed 's/^	X//' << \SHAR_EOF > 'main.c'
	X/*
	X * main.c
	X *
	X * Craig Kolb 10/88
	X *
	X * Copyright (C) 1988 Craig Kolb
	X * This program may be freely copied, modified, and distributed provided that
	X * it is given free-of-charge for strictly non-commercial purposes and that
	X * this copyright notice remains intact.
	X */
	X#include <stdio.h>
	X#include <math.h>
	X#include "svfb_global.h"
	X#include "var.h"
	X
	X/*
	X * Since the sun has a maximum colormap size of 256, we can
	X * stores all our indicies into it in unsigned chars.
	X */
	X#ifdef SUN
	Xunsigned char	NearestColor[33][33][33];
	X#else
	Xshort		NearestColor[33][33][33];
	X#endif
	X
	XBox		*Boxes;				/* Array of color boxes. */
	Xshort		***Image;			/* Original image */
	Xunsigned char	**colormap;			/* Representative colors */
	Xint		frequency[33][33][33],		/* Frequency of each color. */
	X		OutColors,			/* # of entries in colormap */
	X		NPixels,			/* # of pixels in image */
	X		Xres, Yres,			/* Resolution of image */
	X		FastColor,			/* Simple rep. finding. */
	X		Verbose;			/* Print # of orig. colors */
	Xchar		*filename = DEFAULT_FILENAME;	/* Name of input image. */
	XFILE		*fimage;			/* Input image file */
	Xdouble		disp_gamma;			/* Gamma of framebuffer. */
	X#ifdef SUN
	Xint		Dither;				/* Dithering flag */
	X#endif
	X
	Xmain(argc, argv)
	Xint argc;
	Xchar **argv;
	X{	
	X	parse_args(argc, argv);
	X
	X	/*
	X	 * Read image & gather statistics.
	X	 */
	X	setup();
	X
	X	/*
	X 	 * Compute color boxes.
	X	 */
	X	compute_boxes();
	X
	X	/*
	X 	 * Compute color map and NearestColor array.
	X	 */
	X	find_colors();
	X	/*
	X	 * Display the image on the framebuffer.
	X	 */
	X	display_image();
	X
	X	exit(0);
	X}
SHAR_EOF
if test 1525 -ne "`wc -c < 'main.c'`"
then
	echo shar: "error transmitting 'main.c'" '(should have been 1525 characters)'
fi
fi
echo shar: "extracting 'setup.c'" '(3368 characters)'
if test -f 'setup.c'
then
	echo shar: "will not over-write existing file 'setup.c'"
else
sed 's/^	X//' << \SHAR_EOF > 'setup.c'
	X/*
	X * setup.c
	X *
	X * Read Utah Raster image & gather initial statistics.
	X *
	X * Craig Kolb 10/88
	X *
	X * Copyright (C) 1988 Craig Kolb
	X * This program may be freely copied, modified, and distributed provided that
	X * it is given free-of-charge for strictly non-commercial purposes and that
	X * this copyright notice remains intact.
	X */
	X#include <stdio.h>
	X#include "svfb_global.h"
	X#include "var.h"
	X
	X/*
	X * Read image and gather statistics.
	X */
	Xsetup()
	X{
	X	register int i, j, k;
	X	int ColorsUsed, red, green, blue, y;
	X	unsigned char *inline[3];
	X	extern FILE *fimage;
	X	extern char *Progname;
	X
	X	/*
	X 	 * Allocate space for boxes.
	X	 */
	X	Boxes = (Box *)malloc(sizeof(Box) * OutColors);
	X
	X	/*
	X	 * Allocate space for colormap.
	X	 */
	X	colormap = (unsigned char **)malloc(OutColors * sizeof(unsigned char *));
	X	for (i = 0; i < OutColors; i++)
	X		colormap[i] = (unsigned char *)malloc(3*sizeof(unsigned char));
	X	/*
	X	 * Read Utah Raster header, calculate image resolution.
	X	 */
	X	sv_globals.svfb_fd = fimage;
	X
	X	if (rle_get_setup(&sv_globals) < 0) {
	X		fprintf(stderr,"%s: Error reading RLE header.\n", Progname);
	X		exit(1);
	X	}
	X
	X	SV_CLR_BIT(sv_globals, SV_ALPHA);
	X	Xres = sv_globals.sv_xmax - sv_globals.sv_xmin + 1;
	X	Yres = sv_globals.sv_ymax - sv_globals.sv_ymin + 1;
	X	sv_globals.sv_xmax -= sv_globals.sv_xmin;
	X	sv_globals.sv_xmin = 0;
	X
	X	/*
	X	 * Allocate input scan lines.
	X	 */
	X	inline[0] = (unsigned char *)malloc(Xres * sizeof(unsigned char));
	X	inline[1] = (unsigned char *)malloc(Xres * sizeof(unsigned char));
	X	inline[2] = (unsigned char *)malloc(Xres * sizeof(unsigned char));
	X	Image = (short ***)malloc(Yres * sizeof(short **));
	X	/*
	X	 * Read image one line at a time.
	X	 */
	X	for (y = 0; y < Yres; y++) {
	X		Image[y] = (short **)malloc(3 * sizeof(short *));
	X		Image[y][RED] = (short *)malloc(Xres * sizeof(short));
	X		Image[y][GREEN] = (short *)malloc(Xres * sizeof(short));
	X		Image[y][BLUE] = (short *)malloc(Xres * sizeof(short));
	X		rle_getrow(&sv_globals, inline);
	X		/*
	X		 * Quantize each of RGB from 0-255 to 0-32.  Increment
	X		 * the frequency count for this RGB triplet.
	X		 */
	X		for (i = 0; i < Xres; i++) {
	X			/*
	X			 * Divide by 8, rouding up, to reduce the
	X			 * size of the space we must search for
	X			 * cut-points.
	X			 * Rounding up ensures that we have a full range
	X			 * of color intensties in the final image and
	X			 * necessitates arrays of size 33 (values from
	X			 * 0 to 32 inclusive).
	X			 */
	X			red = Quantize(inline[RED][i]);
	X			green = Quantize(inline[GREEN][i]);
	X			blue =  Quantize(inline[BLUE][i]);
	X
	X			frequency[red][green][blue]++;
	X#ifdef SUN
	X			if (Dither) {
	X				/*
	X				 * Store original pixel values.
	X				 */
	X				Image[y][RED][i] = (short)inline[RED][i];
	X				Image[y][GREEN][i] = (short)inline[GREEN][i];
	X				Image[y][BLUE][i] = (short)inline[BLUE][i];
	X			} else
	X#endif
	X			{
	X				/*
	X				 * Store prequantized values.
	X				 */
	X				Image[y][RED][i] = (short)red;
	X				Image[y][GREEN][i] = (short)green;
	X				Image[y][BLUE][i] = (short)blue;
	X			}
	X		}
	X	}
	X
	X	fclose(fimage);
	X	free(inline[0]);
	X	free(inline[1]);
	X	free(inline[2]);
	X	NPixels = Xres * Yres;
	X
	X	if (Verbose) {
	X		/*
	X		 * Report the number of colors before cutting.
	X		 */
	X		ColorsUsed = 0;
	X		for (i = 0;i < 33; i++) {
	X			for (j = 0; j < 33; j++)
	X				for (k = 0; k < 33; k++)
	X					if (frequency[i][j][k] != 0)
	X						ColorsUsed++;
	X		}
	X		fprintf(stderr,"Before cutting, there are %d distinct colors.\n",ColorsUsed);
	X	}
	X}
SHAR_EOF
if test 3368 -ne "`wc -c < 'setup.c'`"
then
	echo shar: "error transmitting 'setup.c'" '(should have been 3368 characters)'
fi
fi
echo shar: "extracting 'var.h'" '(1336 characters)'
if test -f 'var.h'
then
	echo shar: "will not over-write existing file 'var.h'"
else
sed 's/^	X//' << \SHAR_EOF > 'var.h'
	X/*
	X * var.h
	X *
	X * Header file for variance-based color quantization program.
	X *
	X * Craig Kolb 10/88
	X *
	X * Copyright (C) 1988 Craig Kolb
	X * This program may be freely copied, modified, and distributed provided that
	X * it is given free-of-charge for strictly non-commercial purposes and that
	X * this copyright notice remains intact.
	X */
	X#define DEFAULT_OUTCOLORS	256		/* Default output map size */
	X#define DEFAULT_FILENAME	"stdin"		/* stdin filename */
	X
	X#define TRUE			(0==0)
	X#define FALSE			(1==0)
	X
	X#ifdef	SUN
	X#	define DEFAULT_GAMMA	1.6
	X#else
	X#	define DEFAULT_GAMMA	1.0
	X#endif	SUN
	X
	X#define RED	0
	X#define GREEN	1
	X#define BLUE	2
	X
	X/*
	X * Quantize(i) return (int)(i / 8 + 0.5) bu uses integer arithmetic.
	X */
	X#define Quantize(i)		((i >> 3) + ((i & 4) >> 2))
	X
	Xtypedef struct {
	X	double	weightedvar;		/* Weighted variance of box. */
	X	int weight;			/* Weight of box. */
	X	double	mean[3];		/* Centroid of box. */
	X	int	freq[3][33];		/* Projected frequencies. */
	X	int	low[3], high[3];	/* Box bounds. */
	X} Box;
	X
	X#ifdef SUN
	Xextern unsigned char	NearestColor[33][33][33];
	X#else
	Xextern short		NearestColor[33][33][33];
	X#endif
	X
	Xextern Box		*Boxes;
	Xextern unsigned char	**colormap;
	Xextern int		frequency[33][33][33], OutColors, NPixels,
	X			Xres, Yres, Dither, FastColor, Verbose;
	Xextern double		disp_gamma; 
	Xextern char		*filename;
	Xextern short		***Image;
SHAR_EOF
if test 1336 -ne "`wc -c < 'var.h'`"
then
	echo shar: "error transmitting 'var.h'" '(should have been 1336 characters)'
fi
fi
exit 0
#	End of shell archive