ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/clrtab.c
Revision: 2.5
Committed: Tue Feb 23 15:11:59 1993 UTC (31 years, 2 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.4: +9 -3 lines
Log Message:
bug fix -- overflow problem in mktabent()

File Contents

# User Rev Content
1 greg 2.4 /* Copyright (c) 1993 Regents of the University of California */
2 greg 2.1
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 greg 2.2 /* define CLOSEST to get closest colors */
35 greg 2.4 #ifndef CLOSEST
36     #ifdef SPEED
37     #if SPEED > 8
38     #define CLOSEST 1 /* this step takes a little longer */
39     #endif
40     #endif
41     #endif
42 greg 2.1
43    
44     new_histo() /* clear our histogram */
45     {
46     bzero((char *)histo, sizeof(histo));
47     }
48    
49    
50     cnt_pixel(col) /* add pixel to our histogram */
51     register BYTE col[];
52     {
53     cndx(col)++;
54     }
55    
56    
57     cnt_colrs(cs, n) /* add a scanline to our histogram */
58     register COLR *cs;
59     register int n;
60     {
61     while (n-- > 0) {
62     cndx(cs[0])++;
63     cs++;
64     }
65     }
66    
67    
68     new_clrtab(ncolors) /* make new color table using ncolors */
69     int ncolors;
70     {
71     if (ncolors < 1)
72     return(0);
73     if (ncolors > 256)
74     ncolors = 256;
75     /* partition color space */
76     cut(CLRCUBE, 0, ncolors);
77 greg 2.2 #ifdef CLOSEST
78     closest(ncolors); /* ensure colors picked are closest */
79     #endif
80 greg 2.1 /* return new color table size */
81     return(ncolors);
82     }
83    
84    
85     int
86     map_pixel(col) /* get pixel for color */
87     register BYTE col[];
88     {
89     return(cndx(col));
90     }
91    
92    
93     map_colrs(bs, cs, n) /* convert a scanline to color index values */
94     register BYTE *bs;
95     register COLR *cs;
96     register int n;
97     {
98     while (n-- > 0) {
99     *bs++ = cndx(cs[0]);
100     cs++;
101     }
102     }
103    
104    
105     dith_colrs(bs, cs, n) /* convert scanline to dithered index values */
106     register BYTE *bs;
107     register COLR *cs;
108     int n;
109     {
110     static short (*cerr)[3];
111     static int N = 0;
112     int err[3], errp[3];
113     register int x, i;
114    
115     if (n != N) { /* get error propogation array */
116     if (N)
117     cerr = (short (*)[3])realloc((char *)cerr,
118     3*n*sizeof(short));
119     else
120     cerr = (short (*)[3])malloc(3*n*sizeof(short));
121     if (cerr == NULL) {
122     N = 0;
123     map_colrs(bs, cs, n);
124     return;
125     }
126     N = n;
127     bzero((char *)cerr, 3*N*sizeof(short));
128     }
129     err[0] = err[1] = err[2] = 0;
130     for (x = 0; x < n; x++) {
131     for (i = 0; i < 3; i++) { /* dither value */
132     errp[i] = err[i];
133     err[i] += cerr[x][i];
134     #ifdef MAXERR
135     if (err[i] > MAXERR) err[i] = MAXERR;
136     else if (err[i] < -MAXERR) err[i] = -MAXERR;
137     #endif
138     err[i] += cs[x][i];
139     if (err[i] < 0) err[i] = 0;
140     else if (err[i] > 255) err[i] = 255;
141     }
142     bs[x] = cndx(err);
143     for (i = 0; i < 3; i++) { /* propagate error */
144     err[i] -= clrtab[bs[x]][i];
145     err[i] /= 3;
146     cerr[x][i] = err[i] + errp[i];
147     }
148     }
149     }
150    
151    
152     static
153     cut(box, c0, c1) /* partition color space */
154     register int box[3][2];
155     int c0, c1;
156     {
157     register int branch;
158     int kb[3][2];
159    
160     if (c1-c0 <= 1) { /* assign pixel */
161     mktabent(c0, box);
162     return;
163     }
164     /* split box */
165     branch = split(box);
166     bcopy((char *)box, (char *)kb, sizeof(kb));
167     /* do left (lesser) branch */
168     kb[prim(branch)][1] = part(branch);
169     cut(kb, c0, (c0+c1)>>1);
170     /* do right branch */
171     kb[prim(branch)][0] = part(branch);
172     kb[prim(branch)][1] = box[prim(branch)][1];
173     cut(kb, (c0+c1)>>1, c1);
174     }
175    
176    
177     static int
178     split(box) /* find median cut for box */
179     register int box[3][2];
180     {
181     #define c0 r
182     register int r, g, b;
183     int pri;
184     int t[HMAX], med;
185     /* find dominant axis */
186     pri = RED;
187     if (box[GRN][1]-box[GRN][0] > box[pri][1]-box[pri][0])
188     pri = GRN;
189     if (box[BLU][1]-box[BLU][0] > box[pri][1]-box[pri][0])
190     pri = BLU;
191     /* sum histogram over box */
192     med = 0;
193     switch (pri) {
194     case RED:
195     for (r = box[RED][0]; r < box[RED][1]; r++) {
196     t[r] = 0;
197     for (g = box[GRN][0]; g < box[GRN][1]; g++)
198     for (b = box[BLU][0]; b < box[BLU][1]; b++)
199     t[r] += histo[r][g][b];
200     med += t[r];
201     }
202     break;
203     case GRN:
204     for (g = box[GRN][0]; g < box[GRN][1]; g++) {
205     t[g] = 0;
206     for (b = box[BLU][0]; b < box[BLU][1]; b++)
207     for (r = box[RED][0]; r < box[RED][1]; r++)
208     t[g] += histo[r][g][b];
209     med += t[g];
210     }
211     break;
212     case BLU:
213     for (b = box[BLU][0]; b < box[BLU][1]; b++) {
214     t[b] = 0;
215     for (r = box[RED][0]; r < box[RED][1]; r++)
216     for (g = box[GRN][0]; g < box[GRN][1]; g++)
217     t[b] += histo[r][g][b];
218     med += t[b];
219     }
220     break;
221     }
222     if (med < MINSAMP) /* if too sparse, split at midpoint */
223     return(set_branch(pri,(box[pri][0]+box[pri][1])>>1));
224     /* find median position */
225     med >>= 1;
226     for (c0 = box[pri][0]; med > 0; c0++)
227     med -= t[c0];
228     if (c0 > (box[pri][0]+box[pri][1])>>1) /* if past the midpoint */
229     c0--; /* part left of median */
230     return(set_branch(pri,c0));
231     #undef c0
232     }
233    
234    
235     static
236     mktabent(p, box) /* compute average color for box and assign */
237     int p;
238     register int box[3][2];
239     {
240 greg 2.5 unsigned long sum[3];
241     unsigned r, g, n;
242     register unsigned b, c;
243 greg 2.1 /* sum pixels in box */
244     n = 0;
245     sum[RED] = sum[GRN] = sum[BLU] = 0;
246     for (r = box[RED][0]; r < box[RED][1]; r++)
247     for (g = box[GRN][0]; g < box[GRN][1]; g++)
248     for (b = box[BLU][0]; b < box[BLU][1]; b++) {
249     if (c = histo[r][g][b]) {
250     n += c;
251     sum[RED] += (long)c*r;
252     sum[GRN] += (long)c*g;
253     sum[BLU] += (long)c*b;
254     }
255     histo[r][g][b] = p; /* assign pixel */
256     }
257 greg 2.5 if (n >= (1<<23)/HMAX) { /* avoid overflow */
258     sum[RED] /= n;
259     sum[GRN] /= n;
260     sum[BLU] /= n;
261     n = 1;
262     }
263 greg 2.1 if (n) { /* compute average */
264     clrtab[p][RED] = sum[RED]*256/NRED/n;
265     clrtab[p][GRN] = sum[GRN]*256/NGRN/n;
266     clrtab[p][BLU] = sum[BLU]*256/NBLU/n;
267     } else { /* empty box -- use midpoint */
268     clrtab[p][RED] = (box[RED][0]+box[RED][1])*256/NRED/2;
269     clrtab[p][GRN] = (box[GRN][0]+box[GRN][1])*256/NGRN/2;
270     clrtab[p][BLU] = (box[BLU][0]+box[BLU][1])*256/NBLU/2;
271     }
272     }
273 greg 2.2
274    
275     #ifdef CLOSEST
276     #define NBSIZ 32
277     static
278     closest(n) /* make sure we have the closest colors */
279     int n;
280     {
281     BYTE *neigh[256];
282     register int r, g, b;
283     #define i r
284     /* get space for neighbor lists */
285     for (i = 0; i < n; i++) {
286     if ((neigh[i] = (BYTE *)malloc(NBSIZ)) == NULL) {
287     while (i--)
288     free(neigh[i]);
289     return; /* ENOMEM -- abandon effort */
290     }
291     neigh[i][0] = i; /* identity is terminator */
292     }
293     /* make neighbor lists */
294 greg 2.3 for (r = 0; r < NRED; r++)
295     for (g = 0; g < NGRN; g++)
296     for (b = 0; b < NBLU; b++) {
297     if (r < NRED-1 && histo[r][g][b] != histo[r+1][g][b])
298 greg 2.2 addneigh(neigh, histo[r][g][b], histo[r+1][g][b]);
299 greg 2.3 if (g < NGRN-1 && histo[r][g][b] != histo[r][g+1][b])
300 greg 2.2 addneigh(neigh, histo[r][g][b], histo[r][g+1][b]);
301 greg 2.3 if (b < NBLU-1 && histo[r][g][b] != histo[r][g][b+1])
302 greg 2.2 addneigh(neigh, histo[r][g][b], histo[r][g][b+1]);
303     }
304     /* assign closest values */
305     for (r = 0; r < NRED; r++)
306     for (g = 0; g < NGRN; g++)
307     for (b = 0; b < NBLU; b++)
308     setclosest(neigh, r, g, b);
309     /* free neighbor lists */
310     for (i = 0; i < n; i++)
311     free(neigh[i]);
312     #undef i
313     }
314    
315    
316     static
317     addneigh(nl, i, j) /* i and j are neighbors; add them to list */
318     register BYTE *nl[];
319     register int i;
320     int j;
321     {
322     int nc;
323     char *nnl;
324     register int t;
325    
326     for (nc = 0; nc < 2; nc++) { /* do both neighbors */
327     for (t = 0; nl[i][t] != i; t++)
328     if (nl[i][t] == j)
329     break; /* in list already */
330     if (nl[i][t] == i) { /* add to list */
331     nl[i][t++] = j;
332     if (t % NBSIZ == 0) { /* enlarge list */
333     if ((nnl = realloc(nl[i], t+NBSIZ)) == NULL)
334     t--;
335     else
336     nl[i] = (BYTE *)nnl;
337     }
338     nl[i][t] = i; /* terminator */
339     }
340     t = i; i = j; j = t; /* swap and do it again */
341     }
342     }
343    
344    
345     static unsigned
346     dist(col, r, g, b) /* find distance from clrtab entry to r,g,b */
347     register BYTE col[3];
348     int r, g, b;
349     {
350 greg 2.4 register int tmp;
351 greg 2.2 register unsigned sum;
352    
353     tmp = col[RED]*NRED/256 - r;
354     sum = tmp*tmp;
355     tmp = col[GRN]*NGRN/256 - g;
356     sum += tmp*tmp;
357     tmp = col[BLU]*NBLU/256 - b;
358     sum += tmp*tmp;
359     return(sum);
360     }
361    
362    
363     static
364     setclosest(nl, r, g, b) /* find index closest to color and assign */
365     BYTE *nl[];
366     int r, g, b;
367     {
368     int ident;
369     unsigned min;
370     register unsigned d;
371     register BYTE *p;
372     /* get starting value */
373     min = dist(clrtab[ident=histo[r][g][b]], r, g, b);
374     /* find minimum */
375     for (p = nl[ident]; *p != ident; p++)
376     if ((d = dist(clrtab[*p], r, g, b)) < min) {
377     min = d;
378     histo[r][g][b] = *p;
379     }
380     }
381     #endif