ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/clrtab.c
Revision: 2.19
Committed: Fri May 20 02:06:39 2011 UTC (12 years, 11 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R4, rad5R2, rad4R2P2, rad5R0, rad5R1, rad4R2, rad4R1, rad4R2P1, rad5R3, HEAD
Changes since 2.18: +17 -17 lines
Log Message:
Changed every instance of BYTE to uby8 to avoid conflicts

File Contents

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