ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/clrtab.c
Revision: 2.7
Committed: Thu Jun 10 15:45:18 1993 UTC (30 years, 10 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.6: +3 -2 lines
Log Message:
another 16-bit integer fix in mktabent()

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;
242 unsigned long n;
243 register unsigned b, c;
244 /* sum pixels in box */
245 n = 0;
246 sum[RED] = sum[GRN] = sum[BLU] = 0;
247 for (r = box[RED][0]; r < box[RED][1]; r++)
248 for (g = box[GRN][0]; g < box[GRN][1]; g++)
249 for (b = box[BLU][0]; b < box[BLU][1]; b++) {
250 if (c = histo[r][g][b]) {
251 n += c;
252 sum[RED] += (long)c*r;
253 sum[GRN] += (long)c*g;
254 sum[BLU] += (long)c*b;
255 }
256 histo[r][g][b] = p; /* assign pixel */
257 }
258 if (n >= (1L<<23)/HMAX) { /* avoid overflow */
259 sum[RED] /= n;
260 sum[GRN] /= n;
261 sum[BLU] /= n;
262 n = 1;
263 }
264 if (n) { /* compute average */
265 clrtab[p][RED] = sum[RED]*256/NRED/n;
266 clrtab[p][GRN] = sum[GRN]*256/NGRN/n;
267 clrtab[p][BLU] = sum[BLU]*256/NBLU/n;
268 } else { /* empty box -- use midpoint */
269 clrtab[p][RED] = (box[RED][0]+box[RED][1])*256/NRED/2;
270 clrtab[p][GRN] = (box[GRN][0]+box[GRN][1])*256/NGRN/2;
271 clrtab[p][BLU] = (box[BLU][0]+box[BLU][1])*256/NBLU/2;
272 }
273 }
274
275
276 #ifdef CLOSEST
277 #define NBSIZ 32
278 static
279 closest(n) /* make sure we have the closest colors */
280 int n;
281 {
282 BYTE *neigh[256];
283 register int r, g, b;
284 #define i r
285 /* get space for neighbor lists */
286 for (i = 0; i < n; i++) {
287 if ((neigh[i] = (BYTE *)malloc(NBSIZ)) == NULL) {
288 while (i--)
289 free(neigh[i]);
290 return; /* ENOMEM -- abandon effort */
291 }
292 neigh[i][0] = i; /* identity is terminator */
293 }
294 /* make neighbor lists */
295 for (r = 0; r < NRED; r++)
296 for (g = 0; g < NGRN; g++)
297 for (b = 0; b < NBLU; b++) {
298 if (r < NRED-1 && histo[r][g][b] != histo[r+1][g][b])
299 addneigh(neigh, histo[r][g][b], histo[r+1][g][b]);
300 if (g < NGRN-1 && histo[r][g][b] != histo[r][g+1][b])
301 addneigh(neigh, histo[r][g][b], histo[r][g+1][b]);
302 if (b < NBLU-1 && histo[r][g][b] != histo[r][g][b+1])
303 addneigh(neigh, histo[r][g][b], histo[r][g][b+1]);
304 }
305 /* assign closest values */
306 for (r = 0; r < NRED; r++)
307 for (g = 0; g < NGRN; g++)
308 for (b = 0; b < NBLU; b++)
309 setclosest(neigh, r, g, b);
310 /* free neighbor lists */
311 for (i = 0; i < n; i++)
312 free(neigh[i]);
313 #undef i
314 }
315
316
317 static
318 addneigh(nl, i, j) /* i and j are neighbors; add them to list */
319 register BYTE *nl[];
320 register int i;
321 int j;
322 {
323 int nc;
324 char *nnl;
325 register int t;
326
327 for (nc = 0; nc < 2; nc++) { /* do both neighbors */
328 for (t = 0; nl[i][t] != i; t++)
329 if (nl[i][t] == j)
330 break; /* in list already */
331 if (nl[i][t] == i) { /* add to list */
332 nl[i][t++] = j;
333 if (t % NBSIZ == 0) { /* enlarge list */
334 if ((nnl = realloc(nl[i], t+NBSIZ)) == NULL)
335 t--;
336 else
337 nl[i] = (BYTE *)nnl;
338 }
339 nl[i][t] = i; /* terminator */
340 }
341 t = i; i = j; j = t; /* swap and do it again */
342 }
343 }
344
345
346 static unsigned
347 dist(col, r, g, b) /* find distance from clrtab entry to r,g,b */
348 register BYTE col[3];
349 int r, g, b;
350 {
351 register int tmp;
352 register unsigned sum;
353
354 tmp = col[RED]*NRED/256 - r;
355 sum = tmp*tmp;
356 tmp = col[GRN]*NGRN/256 - g;
357 sum += tmp*tmp;
358 tmp = col[BLU]*NBLU/256 - b;
359 sum += tmp*tmp;
360 return(sum);
361 }
362
363
364 static
365 setclosest(nl, r, g, b) /* find index closest to color and assign */
366 BYTE *nl[];
367 int r, g, b;
368 {
369 int ident;
370 unsigned min;
371 register unsigned d;
372 register BYTE *p;
373 /* get starting value */
374 min = dist(clrtab[ident=histo[r][g][b]], r, g, b);
375 /* find minimum */
376 for (p = nl[ident]; *p != ident; p++)
377 if ((d = dist(clrtab[*p], r, g, b)) < min) {
378 min = d;
379 histo[r][g][b] = *p;
380 }
381 }
382 #endif