ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/clrtab.c
Revision: 2.6
Committed: Wed May 5 10:18:05 1993 UTC (31 years ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.5: +1 -1 lines
Log Message:
changed int to long in split() (16-bit bug fix)

File Contents

# Content
1 /* Copyright (c) 1993 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 /* define CLOSEST to get closest colors */
35 #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
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 #ifdef CLOSEST
78 closest(ncolors); /* ensure colors picked are closest */
79 #endif
80 /* 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 long 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 unsigned long sum[3];
241 unsigned r, g, n;
242 register unsigned b, c;
243 /* 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 if (n >= (1<<23)/HMAX) { /* avoid overflow */
258 sum[RED] /= n;
259 sum[GRN] /= n;
260 sum[BLU] /= n;
261 n = 1;
262 }
263 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
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 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 addneigh(neigh, histo[r][g][b], histo[r+1][g][b]);
299 if (g < NGRN-1 && histo[r][g][b] != histo[r][g+1][b])
300 addneigh(neigh, histo[r][g][b], histo[r][g+1][b]);
301 if (b < NBLU-1 && histo[r][g][b] != histo[r][g][b+1])
302 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 register int tmp;
351 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