src/px/clrtab.c

/* Copyright (c) 1992 Regents of the University of California */

#ifndef lint
static char SCCSid[] = "$SunId$ LBL";
#endif

/*
 * Simple median-cut color quantization based on colortab.c
 */

#include "standard.h"

#include "color.h"
                                /* histogram resolution */
#define NRED            36
#define NGRN            48
#define NBLU            24
#define HMAX            NGRN
                                /* minimum box count for adaptive partition */
#define MINSAMP         7
                                /* color partition */
#define set_branch(p,c) ((c)<<2|(p))
#define part(cn)        ((cn)>>2)
#define prim(cn)        ((cn)&3)
                                /* our color table (global) */
BYTE    clrtab[256][3];
                                /* histogram of colors / color assignments */
static unsigned histo[NRED][NGRN][NBLU];
#define cndx(c)         histo[((c)[RED]*NRED)>>8][((c)[GRN]*NGRN)>>8][((c)[BLU]*NBLU)>>8]
                                /* initial color cube boundary */
static int      CLRCUBE[3][2] = {0,NRED,0,NGRN,0,NBLU};
                                /* maximum propagated error during dithering */
#define MAXERR          20


new_histo()             /* clear our histogram */
{
        bzero((char *)histo, sizeof(histo));
}


cnt_pixel(col)          /* add pixel to our histogram */
register BYTE   col[];
{
        cndx(col)++;
}


cnt_colrs(cs, n)                /* add a scanline to our histogram */
register COLR   *cs;
register int    n;
{
        while (n-- > 0) {
                cndx(cs[0])++;
                cs++;
        }
}


new_clrtab(ncolors)             /* make new color table using ncolors */
int     ncolors;
{
        if (ncolors < 1)
                return(0);
        if (ncolors > 256)
                ncolors = 256;
                                /* partition color space */
        cut(CLRCUBE, 0, ncolors);
                                /* return new color table size */
        return(ncolors);
}


int
map_pixel(col)                  /* get pixel for color */
register BYTE   col[];
{
    return(cndx(col));
}


map_colrs(bs, cs, n)            /* convert a scanline to color index values */
register BYTE   *bs;
register COLR   *cs;
register int    n;
{
        while (n-- > 0) {
                *bs++ = cndx(cs[0]);
                cs++;
        }
}


dith_colrs(bs, cs, n)           /* convert scanline to dithered index values */
register BYTE   *bs;
register COLR   *cs;
int     n;
{
        static short    (*cerr)[3];
        static int      N = 0;
        int     err[3], errp[3];
        register int    x, i;

        if (n != N) {           /* get error propogation array */
                if (N)
                        cerr = (short (*)[3])realloc((char *)cerr,
                                        3*n*sizeof(short));
                else
                        cerr = (short (*)[3])malloc(3*n*sizeof(short));
                if (cerr == NULL) {
                        N = 0;
                        map_colrs(bs, cs, n);
                        return;
                }
                N = n;
                bzero((char *)cerr, 3*N*sizeof(short));
        }
        err[0] = err[1] = err[2] = 0;
        for (x = 0; x < n; x++) {
                for (i = 0; i < 3; i++) {       /* dither value */
                        errp[i] = err[i];
                        err[i] += cerr[x][i];
#ifdef MAXERR
                        if (err[i] > MAXERR) err[i] = MAXERR;
                        else if (err[i] < -MAXERR) err[i] = -MAXERR;
#endif
                        err[i] += cs[x][i];
                        if (err[i] < 0) err[i] = 0;
                        else if (err[i] > 255) err[i] = 255;
                }
                bs[x] = cndx(err);
                for (i = 0; i < 3; i++) {       /* propagate error */
                        err[i] -= clrtab[bs[x]][i];
                        err[i] /= 3;
                        cerr[x][i] = err[i] + errp[i];
                }
        }
}


static
cut(box, c0, c1)                        /* partition color space */
register int    box[3][2];
int     c0, c1;
{
        register int    branch;
        int     kb[3][2];
        
        if (c1-c0 <= 1) {               /* assign pixel */
                mktabent(c0, box);
                return;
        }
                                        /* split box */
        branch = split(box);
        bcopy((char *)box, (char *)kb, sizeof(kb));
                                                /* do left (lesser) branch */
        kb[prim(branch)][1] = part(branch);
        cut(kb, c0, (c0+c1)>>1);
                                                /* do right branch */
        kb[prim(branch)][0] = part(branch);
        kb[prim(branch)][1] = box[prim(branch)][1];
        cut(kb, (c0+c1)>>1, c1);
}


static int
split(box)                              /* find median cut for box */
register int    box[3][2];
{
#define c0      r
        register int    r, g, b;
        int     pri;
        int     t[HMAX], med;
                                        /* find dominant axis */
        pri = RED;
        if (box[GRN][1]-box[GRN][0] > box[pri][1]-box[pri][0])
                pri = GRN;
        if (box[BLU][1]-box[BLU][0] > box[pri][1]-box[pri][0])
                pri = BLU;
                                        /* sum histogram over box */
        med = 0;
        switch (pri) {
        case RED:
                for (r = box[RED][0]; r < box[RED][1]; r++) {
                        t[r] = 0;
                        for (g = box[GRN][0]; g < box[GRN][1]; g++)
                                for (b = box[BLU][0]; b < box[BLU][1]; b++)
                                        t[r] += histo[r][g][b];
                        med += t[r];
                }
                break;
        case GRN:
                for (g = box[GRN][0]; g < box[GRN][1]; g++) {
                        t[g] = 0;
                        for (b = box[BLU][0]; b < box[BLU][1]; b++)
                                for (r = box[RED][0]; r < box[RED][1]; r++)
                                        t[g] += histo[r][g][b];
                        med += t[g];
                }
                break;
        case BLU:
                for (b = box[BLU][0]; b < box[BLU][1]; b++) {
                        t[b] = 0;
                        for (r = box[RED][0]; r < box[RED][1]; r++)
                                for (g = box[GRN][0]; g < box[GRN][1]; g++)
                                        t[b] += histo[r][g][b];
                        med += t[b];
                }
                break;
        }
        if (med < MINSAMP)              /* if too sparse, split at midpoint */
                return(set_branch(pri,(box[pri][0]+box[pri][1])>>1));
                                        /* find median position */
        med >>= 1;
        for (c0 = box[pri][0]; med > 0; c0++)
                med -= t[c0];
        if (c0 > (box[pri][0]+box[pri][1])>>1)  /* if past the midpoint */
                c0--;                           /* part left of median */
        return(set_branch(pri,c0));
#undef c0
}


static
mktabent(p, box)        /* compute average color for box and assign */
int     p;
register int    box[3][2];
{
        long    sum[3];
        int     r, g, n;
        register int    b, c;
                                                /* sum pixels in box */
        n = 0;
        sum[RED] = sum[GRN] = sum[BLU] = 0;
        for (r = box[RED][0]; r < box[RED][1]; r++)
            for (g = box[GRN][0]; g < box[GRN][1]; g++)
                for (b = box[BLU][0]; b < box[BLU][1]; b++) {
                    if (c = histo[r][g][b]) {
                        n += c;
                        sum[RED] += (long)c*r;
                        sum[GRN] += (long)c*g;
                        sum[BLU] += (long)c*b;
                    }
                    histo[r][g][b] = p;         /* assign pixel */
                }
        if (n) {                                /* compute average */
                clrtab[p][RED] = sum[RED]*256/NRED/n;
                clrtab[p][GRN] = sum[GRN]*256/NGRN/n;
                clrtab[p][BLU] = sum[BLU]*256/NBLU/n;
        } else {                                /* empty box -- use midpoint */
                clrtab[p][RED] = (box[RED][0]+box[RED][1])*256/NRED/2;
                clrtab[p][GRN] = (box[GRN][0]+box[GRN][1])*256/NGRN/2;
                clrtab[p][BLU] = (box[BLU][0]+box[BLU][1])*256/NBLU/2;
        }
}
Revision:	2.1
Committed:	Mon Oct 12 12:58:39 1992 UTC (33 years, 1 month ago) by greg
Content type:	text/plain
Branch:	MAIN
Log Message:	Initial revision
#	User	Rev	Content
1	greg	2.1	/* Copyright (c) 1992 Regents of the University of California */
2
3			#ifndef lint
4			static char SCCSid[] = "$SunId$ LBL";
5			#endif
6
7			/*
8			* Simple median-cut color quantization based on colortab.c
9			*/
10
11			#include "standard.h"
12
13			#include "color.h"
14			/* histogram resolution */
15			#define NRED 36
16			#define NGRN 48
17			#define NBLU 24
18			#define HMAX NGRN
19			/* minimum box count for adaptive partition */
20			#define MINSAMP 7
21			/* color partition */
22			#define set_branch(p,c) ((c)<<2\|(p))
23			#define part(cn) ((cn)>>2)
24			#define prim(cn) ((cn)&3)
25			/* our color table (global) */
26			BYTE clrtab[256][3];
27			/* histogram of colors / color assignments */
28			static unsigned histo[NRED][NGRN][NBLU];
29			#define cndx(c) histo[((c)[RED]NRED)>>8][((c)[GRN]NGRN)>>8][((c)[BLU]*NBLU)>>8]
30			/* initial color cube boundary */
31			static int CLRCUBE[3][2] = {0,NRED,0,NGRN,0,NBLU};
32			/* maximum propagated error during dithering */
33			#define MAXERR 20
34
35
36			new_histo() /* clear our histogram */
37			{
38			bzero((char *)histo, sizeof(histo));
39			}
40
41
42			cnt_pixel(col) /* add pixel to our histogram */
43			register BYTE col[];
44			{
45			cndx(col)++;
46			}
47
48
49			cnt_colrs(cs, n) /* add a scanline to our histogram */
50			register COLR *cs;
51			register int n;
52			{
53			while (n-- > 0) {
54			cndx(cs[0])++;
55			cs++;
56			}
57			}
58
59
60			new_clrtab(ncolors) /* make new color table using ncolors */
61			int ncolors;
62			{
63			if (ncolors < 1)
64			return(0);
65			if (ncolors > 256)
66			ncolors = 256;
67			/* partition color space */
68			cut(CLRCUBE, 0, ncolors);
69			/* return new color table size */
70			return(ncolors);
71			}
72
73
74			int
75			map_pixel(col) /* get pixel for color */
76			register BYTE col[];
77			{
78			return(cndx(col));
79			}
80
81
82			map_colrs(bs, cs, n) /* convert a scanline to color index values */
83			register BYTE *bs;
84			register COLR *cs;
85			register int n;
86			{
87			while (n-- > 0) {
88			*bs++ = cndx(cs[0]);
89			cs++;
90			}
91			}
92
93
94			dith_colrs(bs, cs, n) /* convert scanline to dithered index values */
95			register BYTE *bs;
96			register COLR *cs;
97			int n;
98			{
99			static short (*cerr)[3];
100			static int N = 0;
101			int err[3], errp[3];
102			register int x, i;
103
104			if (n != N) { /* get error propogation array */
105			if (N)
106			cerr = (short ()[3])realloc((char )cerr,
107			3nsizeof(short));
108			else
109			cerr = (short ()[3])malloc(3n*sizeof(short));
110			if (cerr == NULL) {
111			N = 0;
112			map_colrs(bs, cs, n);
113			return;
114			}
115			N = n;
116			bzero((char )cerr, 3N*sizeof(short));
117			}
118			err[0] = err[1] = err[2] = 0;
119			for (x = 0; x < n; x++) {
120			for (i = 0; i < 3; i++) { /* dither value */
121			errp[i] = err[i];
122			err[i] += cerr[x][i];
123			#ifdef MAXERR
124			if (err[i] > MAXERR) err[i] = MAXERR;
125			else if (err[i] < -MAXERR) err[i] = -MAXERR;
126			#endif
127			err[i] += cs[x][i];
128			if (err[i] < 0) err[i] = 0;
129			else if (err[i] > 255) err[i] = 255;
130			}
131			bs[x] = cndx(err);
132			for (i = 0; i < 3; i++) { /* propagate error */
133			err[i] -= clrtab[bs[x]][i];
134			err[i] /= 3;
135			cerr[x][i] = err[i] + errp[i];
136			}
137			}
138			}
139
140
141			static
142			cut(box, c0, c1) /* partition color space */
143			register int box[3][2];
144			int c0, c1;
145			{
146			register int branch;
147			int kb[3][2];
148
149			if (c1-c0 <= 1) { /* assign pixel */
150			mktabent(c0, box);
151			return;
152			}
153			/* split box */
154			branch = split(box);
155			bcopy((char )box, (char )kb, sizeof(kb));
156			/* do left (lesser) branch */
157			kb[prim(branch)][1] = part(branch);
158			cut(kb, c0, (c0+c1)>>1);
159			/* do right branch */
160			kb[prim(branch)][0] = part(branch);
161			kb[prim(branch)][1] = box[prim(branch)][1];
162			cut(kb, (c0+c1)>>1, c1);
163			}
164
165
166			static int
167			split(box) /* find median cut for box */
168			register int box[3][2];
169			{
170			#define c0 r
171			register int r, g, b;
172			int pri;
173			int t[HMAX], med;
174			/* find dominant axis */
175			pri = RED;
176			if (box[GRN][1]-box[GRN][0] > box[pri][1]-box[pri][0])
177			pri = GRN;
178			if (box[BLU][1]-box[BLU][0] > box[pri][1]-box[pri][0])
179			pri = BLU;
180			/* sum histogram over box */
181			med = 0;
182			switch (pri) {
183			case RED:
184			for (r = box[RED][0]; r < box[RED][1]; r++) {
185			t[r] = 0;
186			for (g = box[GRN][0]; g < box[GRN][1]; g++)
187			for (b = box[BLU][0]; b < box[BLU][1]; b++)
188			t[r] += histo[r][g][b];
189			med += t[r];
190			}
191			break;
192			case GRN:
193			for (g = box[GRN][0]; g < box[GRN][1]; g++) {
194			t[g] = 0;
195			for (b = box[BLU][0]; b < box[BLU][1]; b++)
196			for (r = box[RED][0]; r < box[RED][1]; r++)
197			t[g] += histo[r][g][b];
198			med += t[g];
199			}
200			break;
201			case BLU:
202			for (b = box[BLU][0]; b < box[BLU][1]; b++) {
203			t[b] = 0;
204			for (r = box[RED][0]; r < box[RED][1]; r++)
205			for (g = box[GRN][0]; g < box[GRN][1]; g++)
206			t[b] += histo[r][g][b];
207			med += t[b];
208			}
209			break;
210			}
211			if (med < MINSAMP) /* if too sparse, split at midpoint */
212			return(set_branch(pri,(box[pri][0]+box[pri][1])>>1));
213			/* find median position */
214			med >>= 1;
215			for (c0 = box[pri][0]; med > 0; c0++)
216			med -= t[c0];
217			if (c0 > (box[pri][0]+box[pri][1])>>1) /* if past the midpoint */
218			c0--; /* part left of median */
219			return(set_branch(pri,c0));
220			#undef c0
221			}
222
223
224			static
225			mktabent(p, box) /* compute average color for box and assign */
226			int p;
227			register int box[3][2];
228			{
229			long sum[3];
230			int r, g, n;
231			register int b, c;
232			/* sum pixels in box */
233			n = 0;
234			sum[RED] = sum[GRN] = sum[BLU] = 0;
235			for (r = box[RED][0]; r < box[RED][1]; r++)
236			for (g = box[GRN][0]; g < box[GRN][1]; g++)
237			for (b = box[BLU][0]; b < box[BLU][1]; b++) {
238			if (c = histo[r][g][b]) {
239			n += c;
240			sum[RED] += (long)c*r;
241			sum[GRN] += (long)c*g;
242			sum[BLU] += (long)c*b;
243			}
244			histo[r][g][b] = p; /* assign pixel */
245			}
246			if (n) { /* compute average */
247			clrtab[p][RED] = sum[RED]*256/NRED/n;
248			clrtab[p][GRN] = sum[GRN]*256/NGRN/n;
249			clrtab[p][BLU] = sum[BLU]*256/NBLU/n;
250			} else { /* empty box -- use midpoint */
251			clrtab[p][RED] = (box[RED][0]+box[RED][1])*256/NRED/2;
252			clrtab[p][GRN] = (box[GRN][0]+box[GRN][1])*256/NGRN/2;
253			clrtab[p][BLU] = (box[BLU][0]+box[BLU][1])*256/NBLU/2;
254			}
255			}